Study of high-temperature uniaxial low cyclic fatigue and creep–fatigue behavior of P92 steel using a unified viscoplastic constitutive model

Abstract

In this study, within the framework of the unified viscoplastic cyclic constitutive theory, a cyclic constitutive model is developed based on the modified Chaboche nonlinear kinematic hardening rule. The proposed model improves the prediction of stress relaxation and relaxation stress decay behavior of P92 steel with cycling through the flow rule of hyperbolic sine and static recovery term coupled with accumulative plastic strain and introduces the plastic strain range memory effect to better describe the cyclic softening characteristics of P92 steel. The proposed model is written into the user subroutine of finite element software ABAQUS. Strain-controlled low-cycle fatigue (LCF) and creep–fatigue interaction behavior simulations are performed for the ferritic–martensitic stainless steel P92 at different temperatures and different mechanical strain ranges, respectively. Moreover, a preliminary prediction of the LCF behavior of P92 steel under non-isothermal conditions is provided. By comparing with the experimental results, it can be found that the proposed modified model can provide positive simulation results.

Graphical abstract

Appendices

Appendix 1

Appendix 1.1: Numerical discretization

Since the experiments were performed under isothermal conditions, the temperature rate and temperature strain rate were zero. First, the backward Euler method was used to discretize the constitutive equations, and the following expressions were obtained for the time interval from the $n{\text{th}}$ to the $n + 1{\text{th}}$ step ($\Delta t_{n + 1} = t_{n + 1} - t_{n}$):

$${{\varvec{\upvarepsilon}}}_{{{\varvec{n}} + 1}} = {{\varvec{\upvarepsilon}}}_{{{\varvec{n}} + 1}}^{e} + {{\varvec{\upvarepsilon}}}_{{{\varvec{n}} + 1}}^{p}$$

(30)

$${{\varvec{\upvarepsilon}}}_{{{\text{n}} + 1}}^{p} = {{\varvec{\upvarepsilon}}}_{{\varvec{n}}}^{p} + \Delta {{\varvec{\upvarepsilon}}}_{{{\varvec{n}} + 1}}^{p}$$

(31)

$${{\varvec{\upsigma}}}_{n + 1} = {\mathbf{D}}^{e} :\left( {{{\varvec{\upvarepsilon}}}_{{{\varvec{n}} + 1}} - {{\varvec{\upvarepsilon}}}_{{{\varvec{n}} + 1}}^{p} } \right)$$

(32)

$$\Delta {{\varvec{\upvarepsilon}}}_{n + 1}^{p} = \sqrt{\frac{3}{2}} \Delta p_{n + 1} {\mathbf{n}}_{n + 1}$$

(33)

$$\Delta p_{n + 1} = \phi \left( {\Delta p,{{\varvec{\upalpha}}},R} \right)\Delta t = \alpha {\text{sinh}}\beta \left( {F_{{y\left( {n + 1} \right)}} } \right)\Delta t$$

(34)

$${\mathbf{n}}_{n + 1} = \sqrt{\frac{3}{2}} \frac{{{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} }}{{J\left( {{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} } \right)}} = \frac{{{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} }}{{\left\| {{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} } \right\|}}$$

(35)

The yield surface is discretized as Eq. (36)

$$F_{{y\left( {n + 1} \right)}} = \sqrt {1.5\left( {{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} } \right):\left( {{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} } \right)} - R_{n + 1} - \sigma_{0}$$

(36)

Next, the isotropic hardening is discretized as follows:

$$dR_{n + 1} = k\left[ { - D\left( {Q_{sa} + R_{n} } \right)d\Delta p_{n + 1} } \right]$$

(37)

where $k = \frac{1}{{1 + D\Delta p_{n + 1} }}$.

The mean stress evolution tensor $Y_{n + 1}^{(i)}$ is discretized to obtain:

$$\Delta {\mathbf{Y}}_{n + 1}^{(i)} = - \alpha_{b}^{(i)} \left( {\frac{3}{2}Y_{st}^{(i)} \frac{{{{\varvec{\upalpha}}}_{n + 1}^{(i)} }}{{J\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)}} + {\dot{\mathbf{Y}}}_{n + 1}^{(i)} } \right)J\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)^{{m^{(i)} }} \Delta t$$

(38)

Simplifying Eq. (38), the following is given:

$${\mathbf{Y}}_{n + 1}^{(i)} = \varphi^{(i)} \left( {{\mathbf{Y}}_{n}^{(i)} - \rho^{(i)} {{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)$$

(39)

where $\varphi^{(i)} = \frac{1}{{1 + \alpha_{b}^{(i)} J\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)^{{m^{(i)} }} \Delta t}}$, $\rho^{(i)} = \alpha_{b}^{(i)} Y_{st}^{(i)} J\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)^{{m^{(i)} - 1}} \Delta t$.

The back stress evolution equation is expressed as Eq. (40):

$$\Delta {{\varvec{\upalpha}}}_{n + 1}^{(i)} = w^{(i)} \left\{ \begin{gathered} \frac{2}{3}C^{(i)} \Delta {{\varvec{\upvarepsilon}}}^{p} - \gamma^{(i)} \delta^{\prime } \Delta p_{n + 1} {{\varvec{\upalpha}}}_{n}^{(i)} - \gamma^{(i)} \left( {1 - \delta^{\prime } } \right)\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} :{\mathbf{n}}_{{{\varvec{n}} + 1}} } \right){\mathbf{n}}_{{{\varvec{n}} + 1}} \Delta p_{n + 1} + \hfill \\ \gamma^{(i)} \varphi^{(i)} \Delta p_{n + 1} \left( {{\mathbf{Y}}_{n}^{(i)} - \rho_{k} {{\varvec{\upalpha}}}_{n}^{(i)} } \right) - b^{(i)} J\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)^{{m^{(i)} - 1}} {{\varvec{\upalpha}}}_{n}^{(i)} \Delta t \hfill \\ \end{gathered} \right\}$$

(40)

where $w^{(i)} = \frac{1}{{1 + \gamma^{(i)} \delta^{\prime}\Delta p_{n + 1} + \gamma^{(i)} \varphi^{(i)} \rho^{(i)} \Delta p_{n + 1} + b^{(i)} J\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)^{{m^{(i)} - 1}} \Delta t}}$.

Appendix 1.2: Implicit stress integration method

The implicit integration method, also known as the radial regression method, is used in this paper. In this method, the stress return mapping uses an elastic prediction-inelastic correction process. Thus, the stress tensor can be written as Eq. (41):

$${{\varvec{\upsigma}}}_{n + 1} = {{\varvec{\upsigma}}}_{n + 1}^{tr} - 2G\Delta p_{n + 1} {\mathbf{n}}_{n + 1}$$

(41)

where G is the shear. The test stress of von Mises is indicated as follows:

$${{\varvec{\upsigma}}}_{e}^{tr} = \left[ {\frac{3}{2}({\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} ):({\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} )} \right]^{\frac{1}{2}}$$

(42)

According to the von Mises test stress, the effect force such as von Mises can be obtained:

$$\sigma_{e} = J({{\varvec{\upsigma}}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} ) = \sigma_{e}^{tr} - 3G\Delta p_{n + 1}$$

(43)

Then, the stress normal can be rewritten as a function of the trial quantities::

$${\mathbf{n}}_{n + 1} = \sqrt{\frac{3}{2}} \frac{{{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} }}{{J\left( {{{\varvec{\upsigma}}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} } \right)}} = \sqrt{\frac{3}{2}} \frac{{{\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} }}{{J({{\varvec{\upsigma}}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} )}}$$

(44)

The viscoplastic flow rule is in hyperbolic sinusoidal form, and then the expression for $\dot{p}_{n + 1}$ is as follows:

$$\dot{p}_{n + 1} = \frac{{\Delta p_{n + 1} }}{\Delta t} = \phi (\Delta p_{n + 1} ,{{\varvec{\upalpha}}}_{n + 1} ,R_{n + 1} ) = \alpha_{0} \sinh \beta (\sigma_{e}^{tr} - 3G\Delta p_{n + 1} - R_{n + 1} - \sigma_{0} )$$

(45)

Newton–Raphson iteration method is used to solve for $\Delta p$:

$$\varphi = \Delta p_{n + 1} - \phi (\Delta p_{n + 1} ,{{\varvec{\upalpha}}}_{n + 1} ,R_{n + 1} ) = 0$$

(46)

$$\varphi + (1 - \frac{\partial \phi }{{\partial \Delta p_{n + 1} }}\Delta t)\Delta p_{n + 1} - \frac{\partial \phi }{{\partial {{\varvec{\upalpha}}}_{n + 1} }}\Delta t:d{{\varvec{\upalpha}}}_{n + 1} - \frac{\partial \phi }{{\partial R_{n + 1} }}\Delta tdR = 0$$

(47)

where

$$\frac{\partial \phi }{{\partial \Delta p_{n + 1} }} = - 3G\alpha_{0} \beta \cosh \beta (\sigma_{e}^{tr} - 3G\Delta p_{n + 1} - R_{n + 1} - \sigma_{0} ) = - 3GZ$$

(48)

$$\frac{\partial \phi }{{\partial R_{n + 1} }} = - \alpha_{0} \beta \cosh \beta (\sigma_{e}^{tr} - 3G\Delta p_{n + 1} - R_{n + 1} - \sigma_{0} ) = - Z$$

(49)

$$\frac{\partial \phi }{{\partial {{\varvec{\upalpha}}}_{n + 1} }} = - \frac{\partial \phi }{{\partial \sigma_{e} }}\frac{{\partial \sigma_{e} }}{{\partial {{\varvec{\upalpha}}}_{n + 1} }} = - \alpha_{0} \beta \cosh \beta (\sigma_{e}^{tr} - 3G\Delta p_{n + 1} - R_{n + 1} - \sigma_{0} )\sqrt{\frac{3}{2}} {\mathbf{n}}_{n + 1} = - \sqrt{\frac{3}{2}} {\mathbf{n}}_{n + 1} Z$$

(50)

And we substitute Eqs. (48)–(50) into Eq. (47) to simplify the formula about $\Delta p$:

$$\Delta p_{n + 1} = \Delta p_{n} + d\Delta p_{n + 1}$$

(51)

$$d\Delta p_{n + 1} = \frac{{\phi \left( {\Delta p,{{\varvec{\upalpha}}}_{n + 1} ,R} \right) - \frac{\Delta p}{{\Delta t}} + \sqrt{\frac{3}{2}} Z{\mathbf{n}}_{n + 1} :\mathop \sum \limits_{1}^{i = 4} \left\{ {w^{(i)} b^{(i)} J\left( {{{\varvec{\upalpha}}}_{n + 1}^{\left( i \right)} } \right)^{{m^{(i)} - 1}} {{\varvec{\upalpha}}}_{n + 1}^{\left( i \right)} \Delta t} \right\}}}{{\frac{1}{\Delta t} + 3GZ + \sum\limits_{1}^{i = 4} {w^{(i)} ZC^{(i)} } - \sqrt{\frac{3}{2}} Z{\mathbf{n}}_{n + 1} :\sum\limits_{1}^{i = 4} {w^{(i)} {{\varvec{\upalpha}}}_{n}^{\left( i \right)} \gamma^{(i)} \delta^{\prime}} - \sqrt{\frac{3}{2}} Z\sum\limits_{1}^{i = 4} {w^{(i)} \gamma^{(i)} \left( {1 - \delta^{\prime}} \right)\left( {{{\varvec{\upalpha}}}_{n + 1}^{\left( i \right)} :{\mathbf{n}}_{n + 1} } \right)} + \sqrt{\frac{3}{2}} Z{\mathbf{n}}_{n + 1} :\sum\limits_{1}^{i = 4} {w^{(i)} \gamma^{(i)} \varphi^{(i)} \left( {{\mathbf{Y}}_{n}^{\left( i \right)} - \rho^{(i)} {{\varvec{\upalpha}}}_{n}^{\left( i \right)} } \right)} - ZkD\left( {Q_{sa} + R_{n} } \right)}}$$

(52)

Appendix 1.3: Consistent material tangent operator

During the finite element implementation, the user needs to provide subroutines containing the constitutive equations and the consistent tangent stiffness matrix to maintain the quadratic convergence of the global Newtonian method. The formulation of the consistent tangent stiffness matrix depends on the chosen integration format and the constitutive equation.

Differentiating the discrete constitutive equations:

$$\delta {{\varvec{\upsigma}}}_{n + 1} = \delta \sigma^{tr} - 2G\delta \Delta p_{n + 1} {\mathbf{n}}_{n + 1}$$

(53)

$$\delta \sigma_{e}^{tr} = \delta \sigma_{e} + 3G\delta \Delta p_{n + 1} = {\mathbf{n}}_{n + 1} :\delta ({\mathbf{s}}_{{^{n + 1} }}^{tr} - {{\varvec{\upalpha}}}_{n + 1} )$$

(54)

The above formula is expressed by the deviation stress tensor:

$$\delta {\mathbf{s}}_{n + 1} = \delta {\mathbf{s}}^{tr} - 2G\delta \Delta p_{n + 1} {\mathbf{n}}_{n + 1}$$

(55)

The partial stress tensor can be written as follows using the stress tensor and the strain tensor:

$$\delta {\mathbf{s}}_{n + 1} = \delta {{\varvec{\upsigma}}}_{n + 1} - K{\mathbf{II}}:\delta {{\varvec{\upvarepsilon}}}_{n + 1}$$

(56)

The test deviation stress tensor can be expressed by the strain tensor:

$$\delta {\mathbf{s}}_{n + 1}^{tr} = 2G\delta {{\varvec{\upvarepsilon}}}_{n + 1} - \frac{2}{3}G(\delta {{\varvec{\upvarepsilon}}}_{n + 1} :{\mathbf{I}}){\mathbf{I}}$$

(57)

where K indicates the bulk modulus and I shows the identity matrix.

The stress normals are expressed as:

$$\sqrt{\frac{3}{2}} \frac{{{\mathbf{s}}_{n + 1} - {{\varvec{\upalpha}}}_{n + 1} }}{{\sigma_{e} }} = \sqrt{\frac{3}{2}} \frac{{{\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} }}{{\sigma_{e}^{tr} }}$$

(58)

Rearranging the above equation:

$${\mathbf{s}}_{n + 1} = {{\varvec{\upalpha}}}_{n + 1} + \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }}({\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} )$$

(59)

$$\delta {\mathbf{s}}_{n + 1} = (1 - \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }})\delta {{\varvec{\upalpha}}}_{n + 1} + \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }}\delta {\mathbf{s}}_{n + 1}^{tr} + \frac{{{\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} }}{{\sigma_{e}^{tr} }}(\delta \sigma_{e} - \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }}\delta \sigma_{e}^{tr} )$$

(60)

Substituting Eqs. (57) and (61) into Eq. (56), we obtain:

$$\delta {{\varvec{\upsigma}}}_{n + 1} = K{\mathbf{II}}:\delta {{\varvec{\upvarepsilon}}}_{n + 1} + \left( {1 - \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }}} \right)\sum\limits_{i = 1}^{4} {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } + \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }}\left( {2G\delta {{\varvec{\upvarepsilon}}}_{n + 1} - \frac{2}{3}G(\delta {{\varvec{\upvarepsilon}}}_{n + 1} :{\mathbf{I}}){\mathbf{I}}} \right) + \frac{{{\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} }}{{\sigma_{e}^{tr} }}\left( {\delta \sigma_{e} - \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }}\delta \sigma_{e}^{tr} } \right)$$

(61)

Differentiating the effective viscoplastic strain increment, we obtain:

$$\delta \Delta p_{n + 1} = (\frac{\partial \phi }{{\partial {{\varvec{\upsigma}}}_{n + 1} }}:\delta {{\varvec{\upsigma}}}_{n + 1} + \frac{\partial \phi }{{\partial {{\varvec{\upalpha}}}_{n + 1} }}:\delta {{\varvec{\upalpha}}}_{n + 1} + \frac{\partial \phi }{{\partial R_{n + 1} }}\delta R_{n + 1} )\Delta t$$

(62)

where

$$\frac{\partial \phi }{{\partial {{\varvec{\upsigma}}}_{n + 1} }} = \frac{\partial \phi }{{\partial \sigma_{e} }}\frac{{\partial \sigma_{e} }}{{\partial {{\varvec{\upsigma}}}_{n + 1} }} = \sqrt{\frac{3}{2}} Z{\mathbf{n}}_{n + 1}$$

(63)

$$\frac{\partial \phi }{{\partial R_{n + 1} }} = - \alpha_{0} \beta \cosh \beta (\sigma_{e}^{{{\text{tr}}}} - 3G\Delta p - R_{n + 1} - \sigma_{0} ) = - Z$$

(64)

$$\frac{\partial \phi }{{\partial {{\varvec{\upalpha}}}_{n + 1} }} = - \frac{\partial \phi }{{\partial \sigma_{e} }}\frac{{\partial \sigma_{e} }}{{\partial {{\varvec{\upalpha}}}_{n + 1} }} = - \alpha_{0} \beta \cosh \beta (\sigma_{e}^{{{\text{tr}}}} - 3G\Delta p - R_{n + 1} - \sigma_{0} )\sqrt{\frac{3}{2}} {\mathbf{n}}_{n + 1} = - \sqrt{\frac{3}{2}} {\mathbf{n}}_{n + 1} Z$$

(65)

Substituting Eqs. (63), (64) and (35) into Eq. (62) and simplifying gives:

$$\delta \Delta p_{n + 1} = \sqrt{\frac{3}{2}} \frac{Z}{{Z_{1} }}{\mathbf{n}}_{n + 1} :\delta {{\varvec{\upsigma}}}_{n + 1}^{tr} + \frac{{Z_{2} }}{{Z_{1} }}$$

(66)

where

$$\begin{gathered} Z_{1} = \frac{1}{\Delta t} + \sqrt 6 GZ + Z\sum\limits_{i = 1}^{4} {C^{(i)} } - \sqrt{\frac{3}{2}} Z{\mathbf{n}}_{n + 1} :\sum\limits_{i = 1}^{4} {\gamma^{(i)} {{\varvec{\upalpha}}}_{n + 1}^{(i)} } \delta{\prime} - \sqrt{\frac{3}{2}} Z\sum\limits_{i = 1}^{4} {\gamma^{(i)} } (1 - \delta{\prime} )({{\varvec{\upalpha}}}_{n + 1}^{(i)} :{\mathbf{n}}_{n + 1} ) \hfill \\ \, + \sqrt{\frac{3}{2}} Z{\mathbf{n}}_{n + 1} :\sum\limits_{i = 1}^{4} {\gamma^{(i)} {\mathbf{Y}}_{n + 1}^{(i)} } - ZD(Q_{sa} + R_{n + 1} ) \hfill \\ \end{gathered}$$

(67)

$$\begin{gathered} Z_{2} = \sqrt{\frac{3}{2}} Z{\mathbf{n}}_{n + 1} :\sum\limits_{i = 1}^{4} {b^{(i)} J\left( {{{\varvec{\upalpha}}}_{{{\text{n}} + 1}}^{(i)} } \right)^{{m^{(i)} - 1}} {{\varvec{\upalpha}}}_{n + 1}^{(i)} \Delta t} \\ \\ \end{gathered}$$

(68)

Substituting $\delta \Delta p$,$\delta \sigma_{e}$ and $\delta \sigma_{e}^{{{\text{tr}}}}$ into Eq. (61) to obtain consistent tangent stiffness (CTS), the definition of the terms in the CTS is described in Eq. (69):

$$\begin{aligned} \delta {{\varvec{\upsigma}}}_{n + 1} = & Y_{1} \delta {{\varvec{\upvarepsilon}}}_{n + 1} + Y_{2} {\mathbf{II}}:\delta {{\varvec{\upvarepsilon}}}_{n + 1} + Y_{3} {\mathbf{n}}_{n + 1} :\delta {{\varvec{\upvarepsilon}}}_{n + 1} - ({\mathbf{n}}_{n + 1} :{\mathbf{I}})Y_{4} {\mathbf{I}}:\delta {{\varvec{\upvarepsilon}}}_{n + 1} - ({\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} )Y_{5} {\mathbf{n}}_{n + 1} :\delta {{\varvec{\upvarepsilon}}}_{n + 1} \\ & \quad + ({\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} )Y_{6} ({\mathbf{n}}_{n + 1} :{\mathbf{I}}){\mathbf{I}}:\delta {{\varvec{\upvarepsilon}}}_{n + 1} + (1 - \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }})\sum\limits_{i = 1}^{4} {\left\{ {\frac{2}{3}C^{(i)} \sqrt{\frac{3}{2}} {\mathbf{n}}_{n + 1} - \gamma^{(i)} \left[ {{{\varvec{\upalpha}}}_{n + 1}^{(i)} \delta{\prime} + (1 - \delta{\prime} )({{\varvec{\upalpha}}}_{n + 1}^{(i)} :{\mathbf{n}}_{n + 1} ){\mathbf{n}}_{n + 1} } \right]} \right.} \\ & \quad \left. { + \gamma^{(i)} {\mathbf{Y}}_{n + 1}^{(i)} } \right\}\frac{{Z_{2} }}{{Z_{1} }}\left( {1 + \frac{{{\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} }}{{\sigma_{e}^{tr} }}} \right) - (1 - \frac{{\sigma_{e} }}{{\sigma_{e}^{tr} }})\sum\limits_{i = 1}^{4} {\left\{ {b^{(i)} J\left( {{{\varvec{\upalpha}}}_{n + 1}^{(i)} } \right)^{{m^{(i)} - 1}} {{\varvec{\upalpha}}}_{n + 1}^{(i)} \Delta t} \right\}} \left( {1 + \frac{{{\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} }}{{\sigma_{e}^{tr} }}} \right) \\ & \quad - \frac{{{\mathbf{s}}_{n + 1}^{tr} - {{\varvec{\upalpha}}}_{n + 1} }}{{\sigma_{e}^{tr} }}3G\frac{{Z_{2} }}{{Z_{1} }} \\ \end{aligned}$$

(69)

From the above equation, the Jacobian matrix is finally obtained as (Table 

Table 4 The definition of the terms in the CTS

4):

$$J = \frac{{\partial \delta {{\varvec{\upsigma}}}_{n + 1} }}{{\partial \delta {{\varvec{\upvarepsilon}}}_{n + 1} }}$$

(70)