Due to the limited testing methods and high costs, researches on the evolution law of residual stress in the post-weld heat treatment (PWHT) process of pressure equipment are mainly by numerical simulation, in which the selection of a reasonable and reliable creep constitutive model is crucial to ensure the accuracy of numerical predictions. It is worth noting that the temperature changes and residual stress variation during the heat treatment process result in the following creep characteristics: the creep occurs under different temperatures during the heat treatment processes of heating, holding, and cooling; the creep occurs under different residual stress level where the residual stress varies during the heat treatment stages. Thus, to accurately simulate the stress release behavior during the PWHT process, the proper creep constitutive model must involve the influences of the temperature and the stress state. Meanwhile, current heat treatment simulation usually adopts simplified qualitative approaches, and the creep constitutive models often ignore the physical mechanisms and characteristics of actual creep processes, such as the Norton, the Norton-Bailey, and the Omega models. As a result, the accuracy of heat treatment simulation is limited. In this paper, a new creep constitutive relation under varying temperature and stress is proposed. The hyperbolic-sinusoidal function is used to represent the complete three stages of creep, which compensate the inaccuracy prediction at the first and the third creep stages. Moreover, the Arrhenius function is introduced to characterize the creep response under varying temperature conditions, which enables the prediction of residual stress evolution throughout the entire heating-holding-cooling process of heat treatment. Then, a typical pressure vessel cylindrical circumferential weld is simulated. The stress release during the heat treatment is calculated, and the efficacy of the proposed model is validated through surface residual stress testing. It is found that the creep behavior is a crucial factor in numerical simulation of heat treatment, and the creep constitutive relation can significantly affect the accuracy of heat treatment simulation. The simulation utilizing the hyperbolic-sinusoidal temperature-dependent creep model is much more accurate compared to the classical simulation. The simulation results are much more consistent to the experimental results. Specifically, the traditional model often underestimates the stress state, the residual stress drops abruptly during the heating stage. The residual stress evolution mechanism is simulated using the new model: The majority of residual stress are released during the heating process which primarily relies on the creep strain transformation (contributes up to 85%). The stress state does not change during the holding stage and slight higher at the cooling stage.