Smart Materials and Structures

PAPER • OPEN ACCESS You may also like

- Macroscopic stress–strain curve, local
A polynomial constitutive model of shape memory strain band behavior and the texture of
NiTi thin sheets
alloys based on kinematic hardening Go Murasawa, Kazuhiro Kitamura, Satoru
Yoneyama et al.

- Mechanical behaviour and functional

To cite this article: Cheng Yang et al 2023 Smart Mater. Struct. 32 075016 properties of porous Ti-45 at. % Ni alloy
produced by self-propagating high-
temperature synthesis
N Resnina, S Belyaev, A Voronkov et al.

- A computational study of shape memory

View the article online for updates and enhancements. effect and pseudoelasticity of NiTi alloy
under uniaxial tension during complete
and partial phase transformation
Rajat Rastogi and S J Pawar

This content was downloaded from IP address 183.83.234.227 on 16/08/2023 at 19:02

 Smart Materials and Structures

Smart Mater. Struct. 32 (2023) 075016 (19pp) https://doi.org/10.1088/1361-665X/acd91a

A polynomial constitutive model of

shape memory alloys based on
kinematic hardening
Cheng Yang1,2, Dayi Zhang1,2, Fabrizio Scarpa3, Qicheng Zhang4,∗, He Sun1,2
and Xiaoyong Zhang1,2
1
School of Energy and Power Engineering, Beihang University, No.37 Xueyuan RD. Haidian Dist.,
Beijing 100191, People’s Republic of China
2
Beijing Key Laboratory of Aero-Engine Structure and Strength, Beijing 100191, People’s Republic of
China
3
Dynamics and Control Research Group, University of Bristol, BS8 1TR Bristol, United Kingdom
4
Faculty of Science and Engineering, Swansea University, Swansea SA1 8EN, United Kingdom

E-mail: [email protected]

Received 16 January 2023, revised 10 May 2023

Accepted for publication 25 May 2023
Published 16 June 2023

Abstract
This paper describes the derivation of a phenomenological model for shape memory alloys
under the framework of classical plasticity theory. The proposed model combines the Souza
constitutive approach with kinematic hardening; the model requires solving only one nonlinear
equation rather than several nonlinear ones, therefore increasing the computational efficiency
and convergence. Moreover, the original Souza model is improved by adding an odd polynomial
function to describe the phase transformation of the shape memory alloys, making it possible to
use a lower number of parameters for the inverse identification of the constitutive properties of
SMAs from simple tensile tests. A tangent stiffness formulation is also derived to simulate the
variation of the elastic modulus during the phase transformation. The tangent stiffness
formulation proposed here extends the one used in classical plasticity and improves the
convergence of the proposed model. The reliability and fidelity of the model described in this
work are benchmarked against experimental data and other models. The numerical results show
that the proposed phenomenological approach can describe well the pseudoelasticity and shape
memory effect of shape memory alloys. The formulation described in this paper can be readily
generalized to finite strains and other formulations based on existing formulations related to
classical plasticity theory.
Keywords: shape memory alloy, pseudoelasticity, constitutive model, plasticity,
kinematic hardening, polynomial phase transformation

(Some figures may appear in colour only in the online journal)

∗
Author to whom any correspondence should be addressed.

Original content from this work may be used under the terms
of the Creative Commons Attribution 4.0 licence. Any fur-
ther distribution of this work must maintain attribution to the author(s) and the
title of the work, journal citation and DOI.

1361-665X/23/075016+19$33.00 Printed in the UK 1 © 2023 The Author(s). Published by IOP Publishing Ltd
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

1. Introduction model developed by Lagoudas and Boyd et al [26–28] is

based on the classical Tanaka model [29]. In this model,
Shape memory alloys (SMAs) have been widely used in the Gibbs free energy function of the SMA is derived with
aerospace [1, 2], instrumentation [3, 4] and biomedical the martensite volume fraction and the phase transformation
applications [5, 6] because of their pseudoelasticity and shape strain as the internal variables. The governing equations of
memory effect. In recent decades, scholar and engineers have the control variables and the thermodynamic driving force
extensively developed constitutive models and prototypes of of martensite are derived using the Clausius-Duhem inequal-
devices involving the mechanical behavior of SMAs [7, 8]. ity. The phase transformation yield condition used in the
Studies have also shown that the pseudoelasticity and shape Lagoudas model states that when the thermodynamic driving
memory effect of SMAs are caused by the phase changes of force reaches the critical value, the SMA will undergo phase
twinned and detwinned martensite and austenite under differ- transformation. The thermodynamic driving force is obtained
ent temperatures and stresses [7, 9–11]. by taking the partial derivative of the Gibbs free energy with
SMA constitutive models can be approximately divided respect to martensite volume fraction, which makes the phys-
into two categories: meso-crystal ones and macroscopic ther- ical meaning of thermodynamic driving force and critical driv-
modynamic phenomenological approaches [12]. We introduce ing force not very intuitive and easy to understand. Lagoudas
first the theoretical fundamentals of the meso-crystal con- and Brinson et al [30–32] also proposed several phase change
stitutive model and its advantages and disadvantages. The hardening functions in exponential and trigonometric forms
meso-crystal constitutive model [13–15] is based on the lat- to better characterize the mechanical behaviors of the SMA
tice deformation of martensite variants. The single-crystal phase transformation. Based on the temperature-stress phase
constitutive relationship is derived by applying the WLR diagram, Popov and Lagoudas [33] introduced the volume
crystallography theory, which is used to predict the crystallo- fraction of twinned martensite and detwinned and austenite
graphic characteristic parameters, with the multi-crystal con- to describe the transformation between different metal phases
stitutive relationship being then derived by combining a self- and complex mechanical behaviors in SMAs under various
consistent method. Meso-crystal models aim to describe the thermal and static loading paths. Zaki and Moumni [34, 35]
microstructural behaviors of SMA grains (i.e. nucleation and developed the Lagoudas model by introducing the martensite
interface movements), which can help to understand the mech- orientation strain, which enables to describe the reorientation
anisms underpinning the mechanics of SMAs. Anand, Wang, and plastic deformation [36] of microstructures in SMAs. On
and Xiao et al [16–19] considered the volume fraction of a the basis of Lagoudas model, Chemisky, Lagoudas, Hartl and
group of martensite variants as the internal variables and intro- Zhang et al have performed several analyses on the tension-
duced an interaction energy matrix to describe the interaction compression asymmetry [37, 38], finite strain [39, 40], large
between martensite variants. The Gibbs free energy function deformation [41], non-proportional loading [39, 42], load rate
can then be obtained by using the meso-mechanical theory, influence [43–45], thermomechanical coupling [46–49], cyc-
resulting in the derivation of an SMA constitutive model. This lic loading [32, 50, 51] and plastic fatigue [52–57] of SMAs,
model can describe well the reorientation [20, 21] of martens- which have extended the assessment of the feasibility of this
ite, detwinning-induced plasticity, and cyclic degeneration model.
[22–24]. However, the material parameters required for the The second class of SMA phenomenological models is
meso-crystal approach need to be extracted from specialistic the Souza theoretical approach. The basic theory underpin-
experiment tests. Moreover, the meso-crystal model involves ning this model and its developments are described as fol-
the use of a significant number of internal variables, leading lows. Based on the use of temperature-stress phase diagram,
to very time-consuming calculations that are not suitable for Souza et al [58] constructed Helmholtz free energy using the
implementations within finite elements frameworks. transformation strain as the internal variable and derived a 3D
The phenomenological constitutive model is more suitable macroscopic constitutive model, in which the phase transform-
to be incorporated with finite element models because of its ation driving force is obtained by taking the partial derivat-
simplicity in predicting the thermomechanical response of ive of the Helmholtz free energy with respect to the trans-
SMA actuators [25]. The phenomenological model only uses formation strain. The phase transformation driving force has
macroscopic state variables to describe the overall mechan- a similar meaning to stress, which makes its physical signi-
ical behavior of the SMA materials, with no need of using ficance intuitive and easy to understand. On this basis, Souza
microscale parameters. This approach is therefore well suited introduced a transformation yield function, which was further
to practical engineering simulations. The material paramet- improved in Souza’s follow-up work [59]. This yield function
ers associated to phenomenological models can are usually describes the occurrence of the phase transformation within
determined using classical experimental measurements, and SMAs. This phase transformation yield function also indic-
the constitutive equations can be also embedded within finite ates that SMAs will undergo phase transformation when the
elements codes. driving force of the phase transformation goes beyond the
Phenomenological models are broadly categorized in two elastic region [60]. Souza et al [59] also derived an evolu-
classes, the Lagoudas and the Souza models. These mod- tion criteria based on the Drucker-Prager inequality in terms
els involve the derivation of the free energy function by of the deviatoric transformation stress and the mean trans-
selecting appropriate internal variables [26]. The fundament- formation stress. The model can describe the shape memory
als about the Lagoudas model are described here first. The effect and pseudoelasticity of SMAs. Arghavani et al [61]

2
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

expressed the transformation strain as the product of a scalar instead. In the fifth section, numerical solutions are provided
and an orientation tensor, which also enabled the Souza model to verify the reliability of the proposed model to accurately
to describe the reorientation phenomenon. Auricchio and Gu describe the complex mechanical behaviors of SMAs, such as
et al [62, 63] have also optimized the numerical integration pseudoelasticity and shape memory effect.
algorithm adopted in the original Souza model to provide a
more stable convergence. Altas et al [64] introduced expo-
nential functions to describe the phase transformation harden- 2. Transform stress–strain curve into kinematic
ing characteristics of SMAs during the forward and reverse hardening curve
phase transformations. Based on the Souza model, Chen et al
[65] considered the plastic deformation after the completion The model described in this paper is established based on
of phase transition and the coupling between the plastic strain the temperature–stress (σ − T) phase diagram, and the mater-
and the transformation strain. Kang and Kan et al [66–68] ial parameters can be obtained by standard uniaxial tensile
first observed ratchetting deformations of super-elastic NiTi and compression tests. Therefore, the parameters of the phe-
alloy under uniaxial stress-controlled cyclic tests and then nomenological model can be measured by simple mechanical
evaluated the dependence of the ratchetting versus the applied tests. The (σ − T) phase diagram contains the key information
stress, stress rates, temperatures and types of loading. The at the beginning and the completion of the SMA phase trans-
Authors indicated that the phenomenon of the transforma- formation. This key information has a strong inherent mech-
tion ratchetting is caused by transformation-induced plasticity anical relation with the stress–strain curve, representing the
(TRIP) and plastic yield. They therefore introduced energy constitutive properties of the SMA. Therefore, the inherent
dissipation into the equations describing the damage evolu- relation between the σ − T phase diagram and the stress–strain
tion and established a predictive multi-axial fatigue life model (σ − ε) curve will be discussed in this Section. How to trans-
for NiTi-SMAs [69]. Kang and Petrini et al [70, 71] intro- form the stress–strain curve of the SMA material into a quasi-
duced the cumulated induced-martensite volume fraction vari- kinematic hardening curve is also presented.
able into the phenomenological model to describe the cyclic The significance of the typical stress-temperature phase
ratchetting behavior. Their model can comprehensively con- diagram of a SMA (figure 1) is described briefly. The metal
sider the evolution of TRIP, tension-compression asymmetry, phases of SMAs can be divided into twinned martensite (M t ),
the dependence of transformation ratchetting on the stress and detwinned martensite (Md ) and austenite (A), which are related
phase transformation hardening. Similarly, Auricchio, Zaki, to temperature and stress. The σ − T phase diagram mainly
Moumni, and Petrini et al have proposed several refinements to includes three types of parameters: the start and the finish of
the Souza model in terms of martensite reorientation [72–75], the complete critical transformation temperature for martens-
tension-compression asymmetry [76], cyclic behavior [77], ite and austenite (T M M A A
s , Tf , T s , T f ), the forward and the reverse
finite strain [72, 78–80], and plastic coupling [81]. critical transformation stresses for martensite and austenite
The Souza model has a clear physical implication, and its (σsM , σfM , σsA , σfA ), and the slope of the critical phase trans-
constitutive equations are similar to those used in classical formation stresses (CM and CA ). R is the radius of the elastic
plastic mechanics theory. The numerical solution of the model region in the σ − T phase diagram. The term τm is the so-called
is however difficult, because a set of nonlinear equations needs Maxwell stress.
to be solved to obtain the transformation strain. To deal with SMAs display different mechanical behaviors at low and
this problem, this paper proposes a macroscopic phenomen- high temperatures. SMAs behaves in a pseudoelastic man-
ological constitutive model for SMAs within the framework ner (PE) at temperatures higher than T Af (see figure 2). When
of classical plasticity theory. This approach is based on the the temperature is lower than TM f , the SMA behaves under
use of the Souza model and on adopting kinematic harden- shape memory effect (SME). There is a strong internal rela-
ing principles. The second section of this paper describes the tion between the σ − T phase diagram and the σ − ε curve. The
internal relation between temperature-stress phase diagrams σ − ε curves of SMAs at different temperatures can be determ-
and stress–strain curves and derives the transformation rela- ined based on the σ − T phase diagram, the elastic modulus of
tion related to the different internal variables. In the third the material at martensite and austenite states, and the max-
section, a 1D SMA macroscopic constitutive mode like a kin- imum phase transformation strain. The internal mechanical
ematic hardening one is firstly derived by using Helmholtz free relations between the phase diagram and stress–strain curves
energy function. Then, an odd polynomial function is intro- are discussed in appendix.
duced to describe the mechanical response during the phase The SMA phenomenological model described in this paper
transformation and an analytical expression related to a con- is derived by referring to the classical plastic mechanics the-
sistent tangent stiffness is given. In the fourth section of this ory. Therefore, the stress–strain curve of the SMA material is
paper, the midpoint integration algorithm is used to estab- transformed into a quasi-kinematic hardening curve. Taking
lish the solution for the numerical integral of the polynomial the pseudoelastic behavior of a shape memory alloy as an
macroscopic constitutive model. When phase transformation example (figure 3), the coordinate point (ετelm , τm ) of the ori-
occurs, the proposed algorithm requires to solve one nonlin- ginal stress–strain curve is shifted to the origin to obtain a
ear equation only to complete the stress correction; this can stress–strain curve which is similar to the one of a kinematic
help to overcome the intrinsic difficulty of the Souza model, hardening model. ετelm is the elastic strain corresponding to the
which requires the solution of a set of nonlinear equations Maxwell stress τm under elastic loading The term εLtr represents

3
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

Figure 1. Typical stress–temperature phase diagram of SMA.

Figure 2. Internal relation between the σ–T phase diagram and the σ–ε curve of SMA ((a) The σ–ε curve of the shape memory effect.
(b) The stress-temperature phase diagram. (c) The σ–ε curve of the pseudoelasticity).

Figure 3. Illustration of the transformations related to the kinematic hardening model. (a) The radius of the elastic
( region
) and the Maxwell
τm
stress are shown in the pseudoelasticity curve of SMA. (b) The pseudoelasticity curve of SMA is shifted from εel , τm to the origin point
to obtain the σ–ε curve of the kinematic hardening model.).

the maximum phase transformation strain. We are calling this equation (1). The back stress amplitudes of the forward and the
model as the quasi-kinematic hardening one. The red dashed reverse transformations are denoted as αamp . To ensure that the
line in the figure is the kinematic hardening track (typical back forward and the reverse back stress are consistent, we assume
stress). The amplitudes of the back stress in the forward and that:
the reverse transformation process can be expressed as:
( )
{ σfM = σsM + CA TAf − T As . (2)
αforward = σfM − σsM
. (1)
αreverse = σsA − σfA
Hence,
( ) maximum back stress amplitude is αmax =
the
During the shape memory state, the stress τm is equal to CA TAf − T As . Based on the above transformation process fol-
zero. In this state, the amplitude of the back stress in the lowing a kinematic hardening model, the function between the
forward transformation process can also be expressed as in stress τm and the temperature is

4
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

This use of the method first requires the value of the specific
free energy (Helmholtz free energy or Gibbs free energy) at
any point of the state of the material, which can be represen-
ted by the state variable. Here, the phenomenological model
of an SMA is derived in terms of Helmholtz free energy. The
selected state variables are the total strain ε, the test temperat-
ure T, the phase transformation strain εtr and the phase trans-
formation hardening parameter Eh . Among them, ε and T are
observable variables that can be directly measured from mech-
anical tests. The parameters εtr and Eh are internal variables.
The term εel is the elastic strain. The constitutive model in this
paper does not consider the effect of the plastic strain during
Figure 4. Radius of the elastic region.
phase transformation. In addition, the thermal strain caused by
the change of temperature is relatively small compared with
τm (T) = β(T − T ∗ ) the phase transformation strain. For the sake of simplicity,
+
(3)
the influence of the thermal strain is ignored in this model.
where the constant β is defined as: β = (CM + CA )/ 2. (•)
+ Therefore, the total strain of a SMA can be decomposed into
represents the positive function, which is defined as (x) =
+ an elastic and a phase change part:
(x + |x|)/ 2. The Maxwell stress is 0 when T ∗ = T. At that
ε = εel + εtr . (5)
particular instant of time, the stress–strain curve of SMA is
the one shown in figure 2(a). The phase transformation of the From a physical point of view, the Helmholtz free energy
SMA is related to the historical loading path. Moreover, the of a representative unit is composed of the elastic potential
proposed model can record the volume fractions of twinned energy Ψ el , the thermal energy Ψ T and the phase transforma-
martensite, detwinned martensite and austenite during the tion energy Ψ tr . Hence, the Helmholtz free energy of a SMA
loading history. Hence, only twinned martensite and austen- can be written as
ite phases undergo forward phase transformation. No reverse [ ( )]
phase transformation occurs at T ∗ = T when the load is applied 1 2 T
Ψ (ε, εtr , T) = E(ε − εtr ) + cv (T − T0 ) − T ln
to SMA. The relationship between the radius of the elastic 2ρ T0
region R and temperature of the quasi-kinematic hardening ˆ
1 1 εtr
model (figure 4) is + τm εtr + f (εtr ) dεtr (6)
ρ ρ 0

 R0 ( ) TM< Ts
M
where E is the elastic modulus of the NiTi SMA. ρ is the dens-
R= R0 + β T − Ts M
Ts ⩽ T < T ∗ (4) ity of NiTi SMA. The temperature T0 is the equilibrium tem-

R1 T∗ < T perature at which the Helmholtz free energy of the austenite
[ ( )]/ and martensite phases are the same. The term cv is the specific
where R0 = σs0 M M
, R1 = σs0 + β TfA − TM
s 2, T ∗ = TM
s + heat at constant volume. The other terms of the Helmholtz free
(R1 − R0 )/ β. energy are:
 2
 Ψ el (ε, εtr ) =[ 2ρ E(ε − εtr ) ( )]
1

3. One-dimensional odd polynomial
Ψ T (T) = cv (T − T0 ) − T ln TT0 . (7)
phenomenological model 
 ´
1 1 εtr
Ψ tr (εtr , T) = ρ τm (T) εtr + ρ 0 f (εtr ) dεtr
In this Section, the driving force of the SMA phase transform-
ation, composed of the stress and the back stress, is derived The equations include the terms E (εtr ) = EA + εεtrL ∆E and
according to the Helmholtz free energy. An odd polynomial ∆E = EM − EA . The term f (εtr ) is the phase transformation
function is constructed to describe the phase transformation stress during the phase transform. The complete differential
process and obtain its tangential modulus. Based on the tan- equation of the Helmholtz free energy function is
gential modulus and the phase transformation yield function, it ∂Ψ ∂Ψ ∂Ψ
can be further derived the evolution law of the transformation Ψ= dε + dεtr + dT
∂ε ∂εtr ∂T
strain and back stress. Finally, the tangent and thermal stiff-
1 1
nesses are deduced by the consistency condition of the yield = E (ε − εtr ) dε + [τm + f (εtr )
function and help to improving the convergence of the numer- ρ ρ
ical solution. ∂Ψ
−E (ε − εtr )] dεtr + dT. (8)
∂T
3.1. Constitutive relations According to the Clausius–Duhem inequality( expressed
)
The local state method is used here to derive the constitutive in terms of Helmholtz free energy (σ : dε − ρ Ψ̇ + sdT −
q·∇T
equation of SMAs based on the Clausius-Duhem inequality. T ⩾ 0), we obtain:

5
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

( ) ( )
∂Ψ ∂Ψ ∂Ψ q · ∇T of the midpoint of the stress–strain curve at the phase trans-
σ−ρ dε − ρ s + dT − ρ dεtr − ⩾0
∂ε ∂T ∂εtr T formation stage. The coefficient C can be set as half of the
(9) maximum back stress, C = αamp /2. In addition, the(coeffi- )
cient A can be obtained combining the starting point ε0el , 0 ,
where q is the heat flux and ∇T is the temperature gradient. L α
the midpoint ( ε2 + ε0el , amp L 0
2 ) and end point (ε + εel , αamp )
The above inequality holds for any values of ε, εtr and T. αamp εL L
Furthermore, it is assumed that the SMA material is strongly of the phase transformation stage, A = ( 2 − B 2 )/( ε2 )m .
dissipative, so the stress and entropy are respectively: Therefore, the tangential modulus of the phase transformation
stage can be obtained by taking the derivative of equation (16)
∂Ψ versus the strain ε:
σ=ρ = E (ε − εtr ) (10)
∂ε ( )m−1
( ) dσ εL
∂Ψ T 1 Et = = mA ε − εel −
0
+ B. (18)
s=− = cv ln − βεtr . (11) dε 2
∂T T0 ρ
Let us take the shape memory effect SME as an example.
Therefore, the phase transformation and the heat dissipa-
The elastic modulus of the SMA remains unchanged during
tion potential energy are respectively:
the phase transformation stage. The stress–strain curves and its
∂Ψ tangent modulus change characteristics during the phase trans-
φ tr = −ρ dεtr ⩾ 0 (12) formation (figure 5). When the exponential of an odd poly-
∂εtr
nomial function is large, the middle segment of the curve is
q · ∇T flatter. Moreover, the minimum tangential modulus of the odd
φT = − ⩾ 0. (13)
T polynomial function is equal to the coefficient B. The layout
of the stress–strain curve described by equation (16) is intuit-
Also, according to equation (12), the driving force of the
ive and simple. At first, only the start and the finish points of
SMA phase transformation is:
the phase transformation stage and the slope of the interme-
∂Ψ diate point need to be determined. Then different indices are
χ = −ρ = E (ε − εtr ) − [τm + f (εtr )] = σ − α. (14) selected according to the flatness of the stress–strain curve of
∂εtr
the SMA phase transformation. All these parameters can be
The right end of equation (14) is involves the use of the obtained by standard uniaxial tension and compression tests.
stress minus the term related to the temperature and phase The resulting odd polynomial function provides a good fit to
transformation strain. By comparing the kinematic hardening the stress–strain curve during transformation.
model with figure 3, α = τm + f (εtr ) can be regarded as the
back stress of the SMAs. Thus, the constitutive relations of
3.3. Evolution of the transformation strain and back stress
the SMA material are established as:
 The evolution law of the SMA constitutive model is derived by
 ε = εel + εtr
referring to the kinematic hardening model of plastic mechan-
σ = E (ε − εtr ) . (15)
 ics theory. The back stress and plastic strain need to be identi-
χ = E (ε − εtr ) − [τm + f (εtr )]
fied when the kinematic hardening model is applied. To obtain
the evolution law of the phase transformation strain, the phase
3.2. Odd polynomial function during phase transformation transformation yield function of the SMA should be estab-
lished first. The SMA model used in this paper is related to
The stress–strain curves during phase transformation are like
one-dimensional structures and the model only involves uni-
odd polynomial functions (for example the segment BC of the
axial tension and compression, which can be represented by
curve in figure 10). Therefore, we introduce an odd polynomial
positive and negative signs. Therefore, the yield function of
function to describe the phase transformation process, and its
the phase transformation can be assumed to follow the form:
form is as follows:
( )m ( ) 1 1
εL εL 2
F (σ, εtr , T) = (σ − τm − α) − R2 .
σ = A ε − ε0el − + B ε − ε0el − +C (16) 2 2
(19)
2 2

where m = 3, 5, 7, · · · . The strain εL is the maximum strain The physical meaning of equation (19) is that when the
generated during the phase transformation stage, which is stress exceeds the interval with the radius R of the center point
equal to the maximum phase transformation strain plus the as τm + α, the SMA will undergo phase transformation to pull
elastic strain during the phase transformation stage: the stress back to the yield surface. For the case of associ-
( ) ated plasticity, the evolution law of the transformation strain
σ0 + αamp σ0 is expressed as:
εL = εLtr + − . (17)
Em Ea

The coefficients A, B and C in equation (16) need to be ∂F

dεtr = dλ = dλ (σ − τm − α) . (20)
determined. The coefficient B can be identified from the slope ∂σ

6
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

Figure 5. Odd polynomial stress–strain function and tangent modulus of the transformation stage. ((a) The curve shape of the odd
polynomial function change with the different exponential m, m = 3,5,7. (b) The tangent modulus of the odd polynomial function change
with the different exponential m, m = 3,5,7).

The evolution law of the back stress requires the harden- to zero when the SMA is converted from twinned martensite
ing modulus of the phase transformation process. However, to detwinned martensite. In this case, the second term at the
the explicit expression of the hardening modulus Eh cannot right end is zero, and the above equation degenerates into the
be obtained directly by using the stress–strain relation of the classical expression of the hardening modulus. Thus, when the
phase transformation mentioned above. Hence, the function phase transformation occurs, the increment of back stress is:
describing the hardening modulus needs to be deduced first
according to the mechanics principles. The total differential dα (εtr ) = Eh dεtr . (24)
form of equation (10) is:
3.4. Consistent tangent stiffness
dσ = E (dε − dεtr ) + dE (ε − εtr ) = E (dε − dεtr ) + dEεel .
(21) The tangent and thermal stiffnesses of a SMA can be deduced
by the total differential form of equation (10) and the con-
( )
where dE = dεtr /εLtr ∆E. In addition, the stress increment sistency condition of the yield function (19). The tangential
can also be obtained by the tangential modulus and the strain stiffness here is the consistent tangent stiffness matrix adop-
increment during the transformation stage: ted in standard implicit commercial finite element codes like
ANSYS or ABAQUS. In the process of phase transformation,
dσ = Et dε. (22) the state variables always remain on the yield surface, and the
yield function is always equal to zero. In other words, the value
By combining equations (21) and (22) and the chain of the yield function does not change with the state variables,
rule, the following relationship of the hardening modulus is i.e. the yield function of phase transformation meets the con-
obtained sistency condition:
dσ dσ dε EEt εel ∆EEt ∂F ∂F ∂F
Eh = = = − . (23) dσ + dεtr + dT = 0 (25)
dεtr dε dεtr E − Et εLtr E − Et ∂σ ∂εtr ∂T

The first term on the right end of equation (23) is the where
classical expression of the hardening modulus used in plastic 


∂F
= − ∂τ
∂F
= − ∂α
∂F
= s − τm − α
mechanics theory when the elastic modulus is constant. The  ∂σ m

second term on the right is the correction caused by the dif- ∂F

= ∂ F ∂α
= −Eh (s − τm − α) . (26)


∂εtr ∂α ∂εtr
ference of the elastic modulus between martensite and aus- 
tenite during phase transformation. Obviously, ∆E is equal
∂F
∂T = ∂ F ∂τm
∂τm ∂ T = −β (s − τm − α)

7
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

 n+1
Equation (21) is arranged and expressed as: 
 σ = σ trail,n+1
 αn+1 = αn
. (34)
dσ = Edε − Qdεtr (27)  εnel+1 = εnel + ∆ε

 n+1
εtr = εntr
where Q = E − ∆εELεel . Substituting equation (27) into
tr ( )
equation (25), the expression of transformation strain incre- 2. F σ trail,n+1 , εntr , T n > 0, the SMA undergoes phase
ment is: transformation. It is therefore necessary to solve the
∂F ∂F
transformation strain increment ∆εtr and then carry out
∂σ Edε ∂ T dT the stress correction.
dεtr = + . (28)
∂σ Q − ∂εtr ∂σ Q − ∂εtr
∂F ∂F ∂F ∂F
(c) If the value of the yield function is positive, the stress cor-
rection will be carried out. It is not necessary to consider
By combining equations (25) and (28), the tangent and the directions of the stress or of the transformation strain
thermal stiffnesses of a SMA can be obtained as: in the 1D model. Hence, the actual stress is:
∂F E ∂F
dσ Q ∂σ E Q ∂εtr EEh σ n+1 = σ n + dσ
DK = = E− = = . [ ]
∂σ Q − ∂εtr Q ∂εtr − ∂σ
dε ∂F ∂F 1 ∂F ∂F Eh + Q
= σ n + (1 − θ) En + θEn+1 (∆ε − ∆εtr )
(29) ∆εtr n+1
+ L ∆Eεel . (35)
−Q ∂∂FT ∂F εtr
dσ ∂T βQ
DT = = = = . (30)
∂σ Q − ∂εtr Q ∂εtr − ∂σ
dT ∂F ∂F 1 ∂F ∂F Eh + Q Set

If EA = EM , Q = E = EM . Equations (29) and (30) can be E = (1 − θ) En + θEn+1 . (36)

simplified as:
{ The increment of back stress can be expressed as:
DK = EEE
h +E
h

. (31) ∆α = Enh+1 ∆εtr (37)

DT = Eβh +EE
where:
4. Numerical integration algorithm framework
En+1 Ent +1 εnel+1 ∆EEnt +1
Enh+1 = − (38)
In this paper, the semi-implicit midpoint integration algorithm En+1 − Etn+1 εLtr En+1 − Ent +1
is used to solve the proposed SMA phenomenological model.
εntr + ∆εtr
The constitutive model here is based on the principles of En+1 = EA + ∆E (39)
the kinematic hardening model, so its numerical integra- εLtr
tion algorithm takes inspiration from the schemes adopted in ( )m−1
plastic kinematic hardening. At the beginning of the current εL
Ent +1 = mA εn + ∆ε − + B. (40)
increment step, the known state variables are σ n , αn , εnel , εntr , 2
∆ε, ∆T. Firstly, it is necessary to solve the plastic multiplier
∆λ and then solve the increments of each variable ∆εtr , ∆α, Therefore, the back stress can be obtained as:
∆σ and obtain each state variable under the current increment
step. The numerical integration framework of the 1D SMA αn+1 = αn + ∆α. (41)
phenomenological model (figure 6) is performed as follows:
In addition
(a) Calculate the trail stress: τmn+1 = β(T n + ∆T − T ∗ ) .
+
(42)
σ trail,n+1 = σ n + En ∆ε. (32) Substituting equations (35)–(42) into equation (19), the
nonlinear equation about ∆εtr can be obtained as follows:
(b) Substitute equation (32) into equation (25) to obtain the
value of the yield function:
( ) 1( )2 1
( ) 1( )2 1 F σ n+1 , εntr+1 , Tn+1 = σ n+1 − τmn+1 − αn+1 − R2
F σ trail,n+1 , εntr , T n
= σ trail,n+1 − τmn − αn − R2 . 2[ 2
2 2 1 n ∆εtr n+1
(33) = σ + E (∆ε − ∆εtr) + L ∆Eεel
2 εtr
( )]2
( ) −τmn+1 − αn + Ehn+1 ∆εtr
1. If F σ trail,n+1 , εntr , T n ⩽ 0, the SMA is in an elastic
state. Therefore, the trail stress is the actual stress, and 1
− R2 = 0. (43)
each state variable is 2

8
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

Figure 6. Algorithm framework of the 1D SMA phenomenological model.

The Newton–Rapson method is used to solve the nonlin-

ear equation (43) iteratively to obtain the transformation ∆εtr
Xn+1 = σ n + E (∆ε − ∆εtr ) + L ∆Eεnel+1
strain increment ∆εtr . The key calculation formulas of the εtr
( n )
iterative process are described as: − τm − α + Eh ∆εtr
n+1 n+ 1
(45)

[
dF θ εnel+1 ( )2
=X n+1
∆E (∆ε − ∆ε tr ) − E + ∆E dEnh+1 ∆E Ent +1 1 ∆EEnt +1
d∆εtr εLtr εLtr =− ( ) + L
( )] d∆εtr εL En+1 − En+1
tr
2
t
εtr En+1 − Ent +1
∆εtr dEnh+1
− L ∆E − ∆εtr + Eh n+ 1 2
εnel+1 (∆E) Ent +1
εtr d∆εtr + [ L( )]2 . (46)
(44) εtr En+1 − Ent +1

9
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

(d) Substitute ∆εtr back to equations (34), (35), (37), and (41) which simplifies the actual situation to some extent. The σ − T
to obtain εntr+1 , σ n+1 , ∆α, αn+1 . phase diagram assumes a linear relationship between the crit-
(e) Finally, obtain the tangent stiffness and thermal stiffness ical stress and temperature. In addition, the critical stress of a
of SMA as: twinned martensite transforming into detwinned martensite is
 n+1 n+1 constant, but the critical stress and temperature are not actually
 Dn+1 = En+1 Eh linear.
K Eh +Qn+1
(47)
 Dn+1 = nβ+Q1 n+1 The numerical prediction of the critical transformation
T E h Qn+1 stress is lower than the experimental data at 30 ◦ C. The reason
( ) is that the transformation forms R phase in SMA wire used in
where Qn+1 = En+1 − εnel+1 /εLtr ∆E. Zaki and Moumni’s test [33], a phenomenon affect the trans-
formation but is not accounted in the proposed model. When
5. Numerical results and discussion the SMA wire is stretched at this temperature, the R phase
of the SMA is generated first, and then it becomes twinned
In this section, the accuracy of the polynomial phenomenolo- martensite. The critical stress of the SMA wire is increased due
gical model is verified for the pseudoelasticity and the shape to this physical process. Therefore, the experimental data are
memory effect in shape memory alloys. The research object larger than the predictions provided by the numerical model.
here is the SMA wire. The one end of the SMA wire is fixed When the temperature is higher than 30.0 ◦ C, the numerical
and the tension is applied to the other end of it. The SMA wire calculation of the critical stress has a small offset around the
is heated or cooled through a temperature control box to con- experimental results. Because the proposed model assumes
trol its temperature. Firstly, the stress–strain curves of the pro- that the critical stress is proportional to the temperature, i.e. the
posed model at different temperatures are simulated and com- stress influence coefficients CM and CA of the martensite and
pared with the predictions provided by the Zaki and Moumni austenite transformations are constant. However, these coeffi-
model and the experimental data of the Zaki and Moumni’s test cients will have uncertain and small changes with the change
[33], in order to verify the pseudoelastic manner of the poly- of temperature, which will lead to the discrepancy of the
nomial phenomenological model. The martensite volume frac- numerical prediction around the experimental results. On the
tion and transformation strain of the SMA are also analyzed whole, the stress–strain curves obtained by the polynomial
during loading and unloading process at different temperat- phenomenological model are however in good agreement with
ures. Secondly, the SMA is heated up and then cooled down to the experimental data obtained by Zaki and Moumni above
simulate the shape memory effect of the proposed model. The 40 ◦ C; this is a confirmation of the overall validity of the
mechanical characteristics of the SMA under complex loading pseudoelastic modeling approach for the SMA.
paths are predicted by the proposed model at last. The stress–strain curves shown in figure 8 describe the
behavior of the model with the temperature. When the ambi-
ent temperature is lower than T As , only the forward phase
5.1. Pseudoelasticity (PE)
transformation occurs during the loading and unloading pro-
The stress–strain curves obtained by Zaki and Moumni’s cess, but no inverse phase transformation happens. Twinned
tensile tests [33] (loading and unloading) of an SMA wire at martensite or the mixture of twinned martensite and austen-
different temperatures are used to verify the validity of the pro- ite is transformed into detwinned martensite. The transform-
posed model. The mechanical material parameters are shown ation strain generated by the loading cannot be recovered
in table 1. during unloading. When the temperature is higher than T Af ,
The integral solver of the polynomial phenomenological the forward phase transformation occurs during the loading
model has been compiled and implemented on a Matlab plat- process, and the austenite changes to detwinned martensite.
form. The program is only solved for a single integral point. During unloading, the inverse phase transformation occurs.
The mechanical material parameters in table 1 are inputed The detwinned martensite changes back to austenite and
into the integral solver. The loading conditions are as fol- the transformation strain disappears completely. However,
lows: at different temperatures, the integral point of the SMA when the ambient temperature is between T As and T Af , such
is first loaded until the forward phase transformation is com- as at 37.0 ◦ C, a forward phase transformation occurs dur-
pleted, and then the force load is unloaded to zero. The res- ing the loading process and only part of the transformation
ulting stress–strain curves are shown in figure 7. Note that the strain undergoes an inverse phase transformation during the
exponent of the polynomial function of the phenomenological unloading process. Therefore, residual transformation strain
model is set to 7. is present after unloading. Figure 7 also shows the critical
As shown in figure 7, the phenomenological model pro- stresses at different temperatures during the forward (inverse)
posed in this paper can reflect more accurately the actual phase transformation. Each critical stress is linearly dependent
mechanical properties of the SMA, compared with the pre- versus the temperature.
dictions provided by the Zaki and Moumni model. However, Figure 8 shows that the proposed model allows to pre-
there are some differences between the calculated results of dict the mechanical properties of SMA actuators within a
the proposed constitutive model and the experimental results. broad range of temperatures. In addition, due to the close rela-
The reason behind these discrepancies is that the polynomial tionship existing between the proposed model and the phase
phenomenological model is based on the σ − T phase diagram, diagram, one can develop from this simulation approach a

10
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

Table 1. Material parameters of SMA wire tensile test in Popov and Lagoudas [33].

Constitutive parameter EA (Mpa) EM (Mpa) CA (Mpa K K−1 ) CM (Mpa K−1 ) T As (K) T Af (K) TM
s (K)

Value 40 624.0 21 364.0 11.5 7.0 311.0 315.0 295.0

−3
Constitutive parameter TM
f (K) σ s (Mpa) h (Mpa) ρ (kg m ) εmax µ m

Value 291.0 190.0 0.0 6500.0 0.04 0.33 7

Figure 7. The stress–strain curves simulated by the proposed model are compared with the predictions provided by the Zaki and Moumni
model and the experimental data of the Zaki and Moumni’s test at different temperatures [33]. ((a) The temperature is 30 ◦ C. (b) The
temperature is 40 ◦ C. (c) The temperature is 50 ◦ C. (d) The temperature is 60 ◦ C).

Figure 8. Stress–strain curves of the proposed phenomenological SMA model at different temperatures.

11
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

Figure 9. Characteristics of the martensite and transformation strains of the SMA at different temperatures. ((a) Different martensite
volume fraction change during loading and unloading process at 10 ◦ C. (b) Different martensite volume fraction change during loading and
unloading process at 20 ◦ C. (c) Different martensite volume fraction change during loading and unloading process at 60 ◦ C. (d) The
transformation strain and the total strain change at different temperatures).

higher fidelity version for practical engineering applications is completed and at this point, ξ t = 0.0, ξ dt = 1.0, while the
by simply determining the key parameters of the phase dia- phase change strain reaches 4%. During the phase transforma-
gram. Hence, the proposed model could be more convenient tion, the increment of the volume fraction of the detwinned
to use for the general design and evaluation of SMA actuators. martensite is equal to the volume fraction decrease of the
The proposed phenomenological model can describe the twinned martensite, dξ dt = −dξ t , i.e. the dissociated twinned
physical characteristics of a SMA in different phase trans- martensite is transformed into twinned martensite. After that
formations. Figure 9 shows the variation of volume fraction of point, the load continues to increase, the volume fraction of
twinned and detwinned martensite during loading and unload- twinned and detwinned martensite do not change and the phase
ing at three typical temperatures (10 ◦ C: twinned martens- change strain remains at 0.04. So does the subsequent unload-
ite, 20 ◦ C: twinned martensite and austenite mixture, and ing process, which is only related to the elastic unloading pro-
60 ◦ C: austenite). Figure 8 also shows the change of the trans- cess of the twinned martensite.
formation strain and of the total strain at the corresponding It must be noticed that the change of volume fraction
temperature. of twinned and detwinned martensite is not linear with the
As shown in figure 9(a), when the ambient temperature change of load during the phase transformation. At the start
is 10 ◦ C (which is lower than TM f −291.0 K), only twinned and finish of the phase transformation, the volume fraction of
martensitic phase exists in the SMA. The twinned martens- twinned and detwinned martensite changes relatively slowly.
ite phase changes to detwinned martensite during loading. The volume fraction changes rapidly in the middle of the phase
The unloading process is elastic unloading of the twinned transformation. The trend of the variation of the martensite
martensite without any phase transformation. At the begin- volume fraction shows the characteristics of slowly increas-
ning of the phase transformation, the volume fraction of ing at first and then rapidly increasing, followed by a soft
the twinned martensite is ξ t = 1.0, and the volume frac- increase. This is because the stresses and strains show an odd
tion of the detwinned martensite is ξ dt = 0. The twinned polynomial function relationship in the phase transformation.
martensite begins to transform into detwinned martensite with At the beginning of phase transformation, the stress increases
increasing load. The volume fraction of the twinned martens- quickly and the total strain increases slowly, making the elastic
ite decreases, while the volume fraction of the detwinned strain the predominant component of the strain. However,
martensite increases. Correspondingly, the transformation the increase of the transformation strain and the martensite
strain shown in figure 9(d) begins to increase. Before the volume fraction is relatively small. The same is true when the
loading step reaches 0.5 (i.e. the stress corresponding to the phase transformation is finished. In the middle of the phase
load almost reaches 275.0 MPa.), the phase transformation transformation, the stress increases very little while the total

12
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

strain increases rapidly. The result is a small proportion of

elastic strain and a large proportion of transformation strain.
Therefore, the volume fraction of the twinned and detwinned
martensite increases rapidly. The change characteristics of the
twinned and detwinned martensite volume fractions reason-
ably reflect the physical characteristics during the phase trans-
formation process.
When the load is zero and the ambient temperature is 20 ◦ C
(i.e. between T M M
s and Tf ), the SMA wires are made of a mix-
ture between twinned martensite and austenite. Figure 9(b)
shows that the volume fraction of the twinned martensite is
about 0.43 without loading; the volume fraction of the twinned
martensite is 0, the austenite volume fraction is about 0.57
and the martensite volume fraction is close to 0.43. The mix-
ture begins to transform into detwinned martensite when the
load exceeds a critical stress. Before the load step reaches
0.5, the phase transformation is completed. At this time,
ξ t = 0, ξ dt = 1.0, and the twinned martensite and austenite
completely become detwinned martensite. The phenomeno-
logical model proposed in this work assumes that twinned
martensite and austenite transform into detwinned martensite Figure 10. Stress–strain curve of a SMA wire under shape memory
in equal proportions during the phase transformation. When effect. The SMA is loaded to 450.0 MPa at 283.15 K and then
the phase transformation is finished, the transformation strain heated up the SMA to 370.0 K and finally cooled down to 283.15 K.
is equal to 4.0%. Thereafter, the transformation strain remains
unchanged. Figure 9(d) also shows that the critical strain at
the phase transformation start at 10 ◦ C is larger than the ana- between austenite and twinned martensite occur in the SMA
logous one at 20 ◦ C. The reason is that the SMA wire is pure wires.
twinned martensite phase at 10 ◦ C, but the mixture between As shown in figure 9, the proposed model can track the
twinned martensite and austenite phase is present at 20 ◦ C. change of the volume fraction of each material phase (twinned
The elastic modulus of the former is less than the one of the and detwinned martensite, and austenite) during the loading
latter. Therefore, the critical strain of the former is larger than history. In other words, the proposed model can represent the
the SMA with twinned martensite and austenite when the crit- dependence of the phase states of the SMA during its loading
ical stresses are equal. history. This can help understand the internal changes of the
When the ambient temperature is higher than 60 ◦ C (i.e. lar- shape memory alloy and provide guidance for the design and
ger than T Af ), only austenite exists within the SMA wire optimization of actuators made with this smart material.
(figure 9(c)). At the initial loading time, ξ t = ξ dt = 0 and
the austenite volume fraction is 1.0. When the load exceeds
5.2. Shape memory effect (SME)
a certain value, the austenite starts to transform into twinned
martensite. Before the load step reaches 1.0, the phase trans- To assess how the proposed model represents the shape
formation is finished and at this point, ξ dt = 1.0, ξ t = 0, the memory effect of the constitutive model, a load correspond-
austenite volume fraction is 0 and the transformation strain is ing to 450.0 MPa was applied at 283.15 K, which is a lower
0.04. The transformation process is only a forward transform- temperature than TM f . Then, keeping the load constant, the
ation from austenite to twinned martensite. From the comple- temperature was increased to 370.0 K and finally lowered to
tion of the phase transformation to the peak loading, it is the 283.15 K.
elastic loading process of the twinned martensite. Load step As shown in the stress–strain curves of figure 10, during
2 represents the unloading process. Elastic unloading of the the initial loading stage, the A-B segment is the elastic load-
twinned martensite occurs first. Then, the twinned martensite ing of the twinned martensite. When the stress exceeds point
transforms into austenite when the stress is lower than the B, the SMA begins to undergo forward phase transformation
initial critical one. This process is a reverse phase transform- until the finish of the phase transformation (point C). The C-
ation. Before the loading step reaches 1.75, the reverse phase D segment is the elastic loading of the twinned martensite.
transformation of the SMA wire is completed. The austenite The slopes of the A-B and C-D points are equal because the
volume fraction is 1.0, ξ dt = 0, ξ t = 0 and the phase trans- elastic modulus of the twinned and detwinned martensite are
formation strain is 0. The following process is represented by equal. After point D, by keeping the load unchanged and the
the elastic unloading during austenite. When the stress drops temperature rising to 370.0 K, the SMA undergoes a reverse
to 0, the total strain of the SMA wire is 0; during loading and phase transformation starting from point E and the total strain
unloading, only the forward and reverse phase transformations decreases accordingly—i.e. the point E moves towards point F.

13
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

Figure 11. Shape memory effect of the SMA wire. ((a) The transformation strain and the total strain change with temperature at 450.0 MPa.
(b): The relationship between the transformation strain and the total strain when the temperature changes.).

Figure 12. The transformation strain and the total strain variation characteristics of the SMA under loading process. The blue solid line and
the green dash line show the change of the stress and temperature at different load steps.

This is because the twinned martensite is transformed into aus- strain is 0. Then the temperature is increased to 370.0 K and
tenite during the reverse phase transformation, and the trans- the strain remains unchanged. During cooling, the strain does
formation strain decreases to 0. Then, the forward phase trans- not change until the temperature drops to 332.1 K. After the
formation occurs when the SMA wire is cooled to 283.15 K. temperature is lower than 332.1 K, the forward phase trans-
The austenite is transformed into twinned martensite again, formation starts and the total and phase transformation strains
and the point F goes to point H through point G. During this increase again until 325.6 K, when the forward phase trans-
process, transformation strain and total strain increase. formation is completed.
Figures 11 and 12 clearly show the shape memory effect Now we analyze the mechanical characteristics of the
of the SMA wire when the temperature is modified. The total SMA under complex loading paths as predicted by the
strain corresponding to point D is 0.061, the transformation model developed in this work. The loading paths 1 and 2 in
strain is 0.04, and the elastic strain is 0.021. The temperature figure 13(a) show that the SMA is loaded to 450.0 MPa at
rises from 283.15 K. The strain remains constant, and no phase the initial temperature of 283.15 K (twinned martensite) and
transformation occurs until the temperature rises to 347.1 K. 323.15 K (austenite), respectively; the load is then kept con-
When the ambient temperature is higher than 347.1 K, the stant and the temperature is increased up to 370.0 K; the SMA
reverse phase transformation starts. The total strain and trans- is finally unloaded to 0 MPa stress at the same temperature.
formation strain decrease slowly at first, then rapidly, and The mechanical behaviors of SMAs under loading path 1 and 2
finally again slowly. When the temperature reaches 351.1 K, are A-B-C-D-E-F-G-H and A-B-C-D-E-F-G-H, respectively.
the total and elastic strains are 1.1%, and the transformation The phase transformation process of loading paths 1 and 2 is

14
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

Figure 13. Mechanical behavior of the SMA under different loading paths. ((a) The SMA is loaded to 450.0 MPa at 283.15 K and
323.15 K, respectively. Then the SMA is heated up to 370.0 K at 450.0 MPa and finally unloaded to 0 MPa stress at 370.0 K. (b) The SMA
is loaded to 450.0 MPa at 283.15 K and 323.15 K, respectively; then the SMA is heated up to 370.0 K and cooled to the initial temperature
(283.15 K and 323.15 K) at 450.0 MPa. Finally, the SMA is unloaded to 0 MPa stress.)

forward, and the critical transformation stress of the twinned and shape memory effect of SMAs themselves. This model is
martensite is lower than that of the austenite. Both loading therefore expected to identify the restoring force and study the
paths lead to a complete transformation into twinned martens- mechanical properties of SMA actuators under complex load-
ite after loading to 450.0 MPa. The reverse phase transforma- ing conditions.
tions of loading paths 1 and 2 in the heating process are also
coincident, i.e. they reproduce a shape memory effect. In this
process, the twinned martensite is reversed into austenite. The 6. Conclusion
slope of the G(g)-H(h) segment is equal to that of the a-b seg-
ment (austenite elastic modulus) and larger than that of the A- This paper has first introduced an internal transforma-
B segment (martensite elastic modulus). After the SMA wire tion relationship between temperature-stress phase diagram
subjected to loading paths 1 and 2 is completely unloaded, the and stress–strain curve of shape memory alloy materials.
total strain also becomes 0. Combined with temperature-stress phase diagram, the stress–
The loading paths 3 and 4 in figure 13(b) are slightly dif- strain curve of a SMA is described by using odd polynomial
ferent from the loading paths 1 and 2. After the temperat- functions, which make the proposed model more intuitive and
ure rises to 370.0 K, the load is not unloaded directly, but simpler to mimic the mechanical behavior of SMA. In this
the SMA is cooled to the initial temperature (283.15 K or way, the mechanical parameters of the proposed model can
323.15 K), and then the load is reduced to 0. During the cool- be measured directly by simple uniaxial tensile tests. A 1D
ing process, loading paths 3 and 4 also experience a forward thermodynamic phenomenological model of the SMA is then
phase transformation process. The forward phase transform- derived and established by referring to the kinematic harden-
ation process of the two loading paths coincides, which is ing model from classical plasticity theory. Due to the differ-
H(h)-I(i) segment. In the loading path 3, the temperature is ence between elastic moduli at martensite and austenite, the
reduced to 283.15 K, and then SMA is unloaded. The unload- elastic modulus of a SMA will change during the phase trans-
ing section is J-K, which is the elastic unloading of twinned formation. Hence, the SMA constitutive models established in
martensite. After complete unloading, the residual transform- the past are approximated by the martensite volume fraction to
ation strain of SMA is 0.04. In the loading path 4, the temper- represent the consistent tangent stiffness matrix or the Jacobi
ature is lowered to 323.15 K and then SMA is unloaded. The matrix of a set of nonlinear equations. In this paper, the expli-
unloading stage is j-k-l-m, which includes elastic unloading cit expression of the equivalent hardening modulus has been
of twinned martensite, reverse phase transformation and aus- derived, and its expression consists of the classical hardening
tenite elastic unloading. After complete unloading, the SMA modulus term plus a correction term due to the change of the
returns to its original length and the total strain is 0. Note that elastic modulus. Furthermore, the exact analytical formula of
loading paths 3 and 4 exhibit two typical mechanical behaviors the consistent tangent stiffness matrix of the SMA phase trans-
of SMA: pseudoelasticity and shape memory effect. formation has been derived in this paper and helps to improv-
According to figure 13 and the analysis above, the proposed ing the convergence of the numerical solution. This makes the
model appears suitable for describing a variety of complex tangent stiffness of classical plasticity theory suitable for more
thermal-mechanical loading conditions. The proposed model general cases.
can effectively predict the mechanical properties of SMAs dur- A semi-implicit midpoint integration algorithm is used
ing their loading history and characterize the pseudoelasticity to establish the numerical integration framework of the

15
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

proposed model. By comparing the numerical results with the Suppose that the maximum phase transformation strain at
experimental data and the model of Zaki and Moumni, it is a specific temperature is εLtr . We denote with σsM and σsA the
proved that the proposed polynomial phenomenological model start transformation stresses σs . The stresses σfM and σfA are
can accurately describe the pseudoelasticity of shape memory the finish transformation stress σf . The total strains and elastic
alloy wires. The numerical results also prove that the proposed strain corresponding to σs and σf are respectively:
model can reliably describe shape memory effects in SMAs. 
Moreover, the phenomenological model proposed in this work ε = εsel = σs / Es
 s
has a strong scalability. For example, if the plastic strain ε fel = σf / Ef . (50)
is considered, the phenomenological model can be directly 
 L
extended to the mixed hardening model of classical plastic εf = εtr + σf / Ef
theory.
There are still currently some limitations in the proposed The stresses σsM , σfM , σsA , σfA , and the corresponding strains
phenomenological model. Further work needs to be done to and elastic strain can be formulated as:
make the model reproduce the 3D pseudo-elasticity and shape  A ( )
 σs = CA ( T − T s )
A

memory effect. An improved model should be also able to 

σ A = CA T − TA
mimic the cyclic degradations present in pseudo-elasticity and f f
( ) (51)
shape memory effects under the cyclic thermo-mechanical 
 σs = σs0 + CM T − TM
M M

 (
s
)
loading conditions, such as the so-called transformation ratch-  M
σf = σfM0 + CM T − TM f
etting and re-orientation ratchetting identified by Kang et al
[66–69].  A /

 εs = εL + σsA EM

 εA = σ A / E A

f f
Data availability statement / (52)
εs = σsM EA

M

 /
 M
All data that support the findings of this study are included εf = εLtr + σfM EM
within the article (and any supplementary files).  As /
 εel = σsA EM


Acknowledgments εAf = σ A / EA

el f
/ . (53)

 el
ε Ms
= σ M
EA

 s
/
This work was financially supported by the National Natural  Mf
Science Foundation of China (Grant Nos. 52175071 and εel = σfM EM
52105083), the Major Projects of Aero-engines and Gas tur-
bines (J2019-I-0008-0008) and the Innovation Centre for ORCID iDs
Advanced Aviation Power (HKCX2020-02-016).
Fabrizio Scarpa  https://orcid.org/0000-0002-5470-4834
Qicheng Zhang  https://orcid.org/0000-0001-7696-6221
Appendix. Internal mechanical relations between
the phase diagram and stress–strain curve
References
From figure 2, it is possible to observe that twinned martens-
[1] Jani J M, Leary M, Subic A and Gibson M A 2014 A review of
ite can be transformed into austenite and vice versa within the shape memory alloy research, applications and
temperature range of T M s −Tf . That is, the SMA has an inter-
M
opportunities Mater. Des. 56 1078–113
val intermediate state between twinned martensite and austen- [2] Leal P B C and Savi M A 2018 Shape memory alloy-based
ite when the stress is lower than σsM and the temperature is mechanism for aeronautical application: theory,
between T M M optimization and experiment Aerosp. Sci. Technol.
s and Tf . Assuming that the forward and reverse
76 155–63
transformations between twinned martensite and austenite are [3] Choi E, Mohammadzadeh B and Kim H S 2019 SMA bending
linearly dependent versus the temperature, the elastic moduli bars as self-centering and damping devices Smart Mater.
corresponding to the beginning critical transformation stresses Struct. 28 025029
(σsM and σsA ) are expressed as [52]: [4] Hu S, Qiu C and Zhu S 2022 Machine learning-driven
performance-based seismic design of hybrid self-centering
 braced frames with SMA braces and viscous dampers Smart

 EM T < TMf
T−TM
Mater. Struct. 31 105024
Es = EM + TM −TM (EA − EM ) Tf ⩽ T < TM
f M (48) [5] Petrini L and Migliavacca F 2011 Biomedical applications of

 s f
s
shape memory alloys J. Metall. 2011 1–15
EA Ts ⩽ T
M
[6] De Bock S, Iannaccone F, De Santis G, De Beule M,
Mortier P, Verhegghe B and Segers P 2012 Our capricious
and vessels: the influence of stent design and vessel geometry
on the mechanics of intracranial aneurysm stent deployment
Ef = EM . (49) J. Biomech. 45 1353–9

16
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

[7] Cisse C, Zaki W and Zineb T B 2016 A review of constitutive [26] Bo Z and Lagoudas D C 1999 Thermomechanical modeling of
models and modeling techniques for shape memory alloys polycrystalline SMAs under cyclic loading, Part I:
Int. J. Plast. 76 244–84 theoretical derivations Int. J. Eng. Sci. 37 1089–140
[8] Jaber M B, Mehrez S and Ghazouani O 2014 A 1D constitutive [27] Leclercq S and Lexcellent C 1996 A general macroscopic
model for shape memory alloy using strain and temperature description of the thermomechanical behavior of shape
as control variables and including martensite reorientation memory alloys J. Mech. Phys. Solids 44 953–80
and asymmetric behaviors Smart Mater. Struct. 23 095026 [28] Lagoudas D C and Bo Z 1999 Thermomechanical modeling of
[9] Zhang X, Wang S, Yan X, Yue D, Sun R and Zhou X 2016 polycrystalline SMAs under cyclic loading, Part II: material
Probabilistic analysis for the functional and structural characterization and experimental results for a stable
fatigue of NiTi wires Mater. Des. 102 213–24 transformation cycle Int. J. Eng. Sci. 37 1141–73
[10] Sahli M L and Necib B 2014 Characterisation and modelling [29] Adeodato A, Vignoli L L, Paiva A, Monteiro L L,
of behaviour of a shape memory alloys J. Adv. Manuf. Pacheco P M, Savi M A 2022 A shape memory alloy
Technol. 70 1847–57 constitutive model with polynomial phase transformation
[11] Wang W, Fang C and Liu J 2016 Large size superelastic SMA kinetics Shape Mem. Superelasticity 8 1–18
bars: heat treatment strategy, mechanical property and [30] Lagoudas D C, Bo Z and Qidwai M A 1996 A unified
seismic application Smart Mater. Struct. 25 075001 thermodynamic constitutive model for SMA and finite
[12] Cisse C, Zaki W and Zineb T B 2016 A review of modeling element analysis of active metal matrix composites Mech.
techniques for advanced effects in shape memory alloy Compos. Mater. Struct. 3 153–79
behavior Smart Mater. Struct. 25 103001 [31] Lagoudas D, Hartl D, Chemisky Y, Machado L and Popov P
[13] Salem A A, Kalidindi S R and Semiatin S L 2005 Strain 2012 Constitutive model for the numerical analysis of phase
hardening due to deformation twinning in α-titanium: transformation in polycrystalline shape memory alloys Int.
constitutive relations and crystal-plasticity modeling Acta J. Plast. 32 155–83
Mater. 53 3495–502 [32] Mehrabi R, Shirani M, Kadkhodaei M and Elahinia M 2015
[14] Salem A A, Kalidindi S R, Doherty R D and Semiatin S L Constitutive modeling of cyclic behavior in shape memory
2006 Strain hardening due to deformation twinning in alloys Int. J. Mech. Sci. 103 181–8
α-titanium: mechanisms Metall. Mater. Trans. A 37 259–68 [33] Popov P and Lagoudas D C 2007 A 3D constitutive model for
[15] Yu C, Kang G, Song D and Kan Q 2015 Effect of martensite shape memory alloys incorporating pseudoelasticity and
reorientation and reorientation-induced plasticity on detwinning of self-accommodated martensite Int. J. Plast.
multiaxial transformation ratchetting of super-elastic NiTi 23 1679–720
shape memory alloy: new consideration in constitutive [34] Zaki W and Moumni Z 2007 A three-dimensional model of
model Int. J. Plast. 67 69–101 the thermomechanical behavior of shape memory alloys J.
[16] Anand L and Gurtin M E 2003 Thermal effects in the Mech. Phys. Solids 55 2455–90
superelasticity of crystalline shape-memory materials J. [35] Zaki W, Zamfir S and Moumni Z 2010 An extension of the
Mech. Phys. Solids 51 1015–58 ZM model for shape memory alloys accounting for plastic
[17] Wang X M, Xu B X and Yue Z F 2008 Micromechanical deformation Mech. Mater. 42 266–74
modelling of the effect of plastic deformation on the [36] Hartl D J, Chatzigeorgiou G and Lagoudas D C 2010
mechanical behaviour in pseudoelastic shape memory Three-dimensional modeling and numerical analysis of
alloys Int. J. Plast. 24 1307–32 rate-dependent irrecoverable deformation in shape memory
[18] Xiao Y, Zeng P and Lei L 2018 Micromechanical modeling on alloys Int. J. Plast. 26 1485–507
thermomechanical coupling of cyclically deformed [37] Paiva A, Savi M A, Braga A M B and Pacheco P M C L 2005
superelastic NiTi shape memory alloy Int. J. Plast. A constitutive model for shape memory alloys considering
107 164–88 tensile–compressive asymmetry and plasticity Int. J. Solids
[19] Xiao Y, Zeng P and Lei L 2019 Micromechanical modelling Struct. 42 3439–57
on thermomechanical coupling of superelastic NiTi alloy [38] Viet N V, Zaki W and Moumni Z 2019 A model for shape
Int. J. Mech. Sci. 153 36–47 memory alloy beams accounting for tensile compressive
[20] Thamburaja P 2005 Constitutive equations for martensitic asymmetry J. Intell. Mater. Syst. Struct. 30 2697–715
reorientation and detwinning in shape-memory alloys J. [39] Xu L, Solomou A, Baxevanis T and Lagoudas D 2021 Finite
Mech. Phys. Solids 53 825–56 strain constitutive modeling for shape memory alloys
[21] Luo J, Xu K, Li C, Li M and Lin Y 2021 The evolution of considering transformation-induced plasticity and two-way
dynamic recrystallization and recrystallization texture shape memory effect Int. J. Solids Struct. 221 42–59
during isothermal compression of NiTi shape memory alloy [40] Zaki W 2012 An efficient implementation for a model of
Mater. Sci. Eng. A 820 141424 martensite reorientation in martensitic shape memory alloys
[22] Ebrahimi P, Arghavani J, Naghdabadi R and McGarry J P 2020 under multiaxial nonproportional loading Int. J. Plast.
On the effect of detwinning-induced plasticity in 37 72–94
compressive cyclic loading of NiTi shape memory alloys [41] Xu L, Baxevanis T and Lagoudas D C 2019 A
Mech. Mater. 148 103451 three-dimensional constitutive model for the martensitic
[23] Gloanec A-L, Bilotta G and Gerland M 2013 Deformation transformation in polycrystalline shape memory alloys
mechanisms in a TiNi shape memory alloy during cyclic under large deformation Smart Mater. Struct. 28 074004
loading Mater. Sci. Eng. A 564 351–8 [42] Chatziathanasiou D, Chemisky Y, Chatzigeorgiou G and
[24] Wang B, Kang G, Yu C, Gu B and Yuan W 2021 Molecular Meraghni F 2016 Modeling of coupled phase
dynamics simulations on one-way shape memory effect of transformation and reorientation in shape memory alloys
nanocrystalline NiTi shape memory alloy and its cyclic under non-proportional thermomechanical loading Int. J.
degeneration Int. J. Mech. Sci. 211 106777 Plast. 82 192–224
[25] Rao Z, Leng J, Yan Z, Tan L and Yan X 2023 A [43] Shakiba S, Yousefi-Koma A and Ayati M 2020 Development of
three-dimensional constitutive model for shape memory a frequency-dependent constitutive model for hysteresis of
alloy considering transformation-induced plasticity, shape memory alloys Proc. Inst. Mech. Eng. L 234 1535–49
two-way shape memory effect, plastic yield and [44] Lexcellent C and Rejzner J 2000 Modeling of the strain rate
tension-compression asymmetry Eur. J. Mech. A 99 104945 effect, creep and relaxation of a Ni-Ti shape memory alloy

17
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

under tension (compression)-torsional proportional [62] Auricchio F and Petrini L 2004 A three-dimensional model
loading in the pseudoelastic range Smart Mater. Struct. describing stress-temperature induced solid phase
9 613 transformations: solution algorithm and boundary value
[45] Roh J-H 2014 Thermomechanical Modeling of shape memory problems Int. J. Numer. Methods Eng. 61 807–36
alloys with rate dependency on the pseudoelastic behavior [63] Gu X, Zaki W, Morin C, Moumni Z and Zhang W 2015 Time
Math. Probl. Eng. 2014 1–9 integration and assessment of a model for shape memory
[46] Zhang X, Yan X, Xie H and Sun R 2014 Modeling evolutions alloys considering multiaxial nonproportional loading cases
of plastic strain, maximum transformation strain and Int. J. Solids Struct. 54 82–99
transformation temperatures in SMA under superelastic [64] Altas E, Khosravi F, Gokkaya H, Arab Maleki V, Akınay Y,
cycling Comput. Mater. Sci. 81 113–22 Ozdemir O, Bayraktar O and Kandas H 2022 Finite
[47] Morin C, Moumni Z and Zaki W 2011 A constitutive model element simulation and experimental investigation on the
for shape memory alloys accounting for thermomechanical effect of temperature on pseudoelastic behavior of
coupling Int. J. Plast. 27 748–67 perforated Ni–Ti shape memory alloy strips Smart Mater.
[48] Morin C, Moumni Z and Zaki W 2011 Thermomechanical Struct. 31 025031
coupling in shape memory alloys under cyclic loadings: [65] Chen B, Peng X, Chen X, Wang J, Wang H and Hu N 2012 A
experimental analysis and constitutive modeling Int. J. three-dimensional model of shape memory alloys under
Plast. 27 1959–80 coupled transformation and plastic deformation Comput.
[49] Zhang Y, Yu C, Zhu Y, Kan Q and Kang G 2022 Mater. Contin. 30 145–76
Thermo-mechanically coupled deformation of pseudoelastic [66] Kang G, Kan Q, Qian L and Liu Y 2009 Ratchetting
NiTi SMA helical spring Int. J. Mech. Sci. 236 107767 deformation of super-elastic and shape-memory NiTi alloys
[50] Shastry V V, Singh G and Ramamurty U 2021 On the fatigue Mech. Mater. 41 139–53
life enhancement due to periodic healing of a NiTi shape [67] Song D, Kang G, Kan Q, Yu C and Zhang C 2014
memory alloy Mater. Sci. Eng. A 815 141272 Non-proportional multiaxial transformation ratchetting of
[51] Yu C, Kang G and Kan Q 2014 A physical mechanism based super-elastic NiTi shape memory alloy: experimental
constitutive model for temperature-dependent observations Mech. Mater. 70 94–105
transformation ratchetting of NiTi shape memory alloy: [68] Zhao T, Kang G, Yu C and Kan Q 2020 Experimental study on
one-dimensional model Mech. Mater. 78 1–10 whole-life one-way shape memory cyclic degradation and
[52] You Y, Zhang Y, Moumni Z, Anlas G and Zhang W 2017 fatigue failure of NiTi shape memory alloy Mater. Today
Effect of the thermomechanical coupling on fatigue crack Commun. 25 101621
propagation in NiTi shape memory alloys Mater. Sci. Eng. [69] Zhao T and Kang G 2022 Experimental study and life
A 685 50–56 prediction on fatigue failure of NiTi shape memory alloy
[53] Chemisky Y, Hartl D J and Meraghni F 2018 under multi-axial one-way shape memory cyclic loadings
Three-dimensional constitutive model for structural and Int. J. Fatigue 155 106609
functional fatigue of shape memory alloy actuators Int. J. [70] Kan Q and Kang G 2010 Constitutive model for uniaxial
Fatigue 112 263–78 transformation ratchetting of super-elastic NiTi shape
[54] Qidwai M A and Lagoudas D C 2000 Numerical memory alloy at room temperature Int. J. Plast.
implementation of a shape memory alloy thermomechanical 26 441–65
constitutive model using return mapping algorithms Int. J. [71] Petrini L and Bertini A 2020 A three-dimensional
Numer. Methods Eng. 47 1123–68 phenomenological model describing cyclic behavior of
[55] Scalet G, Niccoli F, Garion C, Chiggiato P, Maletta C and shape memory alloys Int. J. Plast. 125 348–73
Auricchio F 2019 A three-dimensional phenomenological [72] Auricchio F, Coda A, Reali A and Urbano M 2009 SMA
model for shape memory alloys including two-way shape numerical modeling versus experimental results: parameter
memory effect and plasticity Mech. Mater. 136 103085 identification and model prediction capabilities J. Mater.
[56] Auricchio F, Constantinescu A, Menna C and Scalet G 2016 A Eng. Perform. 18 649–54
shakedown analysis of high cycle fatigue of shape memory [73] Wang J, Moumni Z, Zhang W, Xu Y and Zaki W 2017 A 3D
alloys Int. J. Fatigue 87 112–23 finite-strain-based constitutive model for shape memory
[57] Zhang X, Huang D, Yan X and Zhou X 2016 Modeling alloys accounting for thermomechanical coupling and
functional fatigue of SMA using a more accurate martensite reorientation Smart Mater. Struct. 26 065006
subdivision of martensite volume fractions Mech. Mater. [74] Auricchio F, Bonetti E, Scalet G and Ubertini F 2014
96 12–29 Theoretical and numerical modeling of shape memory
[58] Souza A C, Mamiya E N and Zouain N 1998 alloys accounting for multiple phase transformations and
Three-dimensional model for solids undergoing martensite reorientation Int. J. Plast. 59 30–54
stress-induced phase transformations Eur. J. Mech. A [75] Arghavani J, Auricchio F, Naghdabadi R, Reali A and
17 789–806 Sohrabpour S 2010 A 3D finite strain phenomenological
[59] Souza A C, Mamiya E N and Zouain N 1998 A new constitutive model for shape memory alloys considering
three-dimensional constitutive model for shape memory martensite reorientation Contin. Mech. Thermodyn.
alloys Computational Mechanics, New Trends and 22 345–62
Applications ed S Idelsohn, E Onate and E Dvorkin [76] Evangelista V, Marfia S and Sacco E 2009 Phenomenological
(Barcelona: CIMNE) pp 1–19 3D and 1D consistent models for shape-memory alloy
[60] Ashrafi M J 2020 Transformation and plasticity of shape materials Comput. Mech. 44 405–21
memory alloy structures: constitutive modeling and finite [77] Choi E, Nam T and Chung Y-S 2010 Variation of mechanical
element implementation J. Mater. Eng. Perform. properties of shape memory alloy bars in tension under
29 5515–24 cyclic loadings Mater. Sci. Eng. A 527 4412–7
[61] Arghavani J, Auricchio F, Naghdabadi R, Reali A and [78] Wang J, Moumni Z and Zhang W 2017 A thermomechanically
Sohrabpour S 2010 A 3D phenomenological constitutive coupled finite-strain constitutive model for cyclic
model for shape memory alloys under multiaxial loadings pseudoelasticity of polycrystalline shape memory alloys Int.
Int. J. Plast. 26 976–91 J. Plast. 97 194–221

18
Smart Mater. Struct. 32 (2023) 075016 C Yang et al

[79] Arghavani J, Auricchio F and Naghdabadi R 2011 A finite crystal plasticity-based model for pseudoelastic NiTi shape
strain kinematic hardening constitutive model based on memory alloy Int. J. Plast. 148 103146
Hencky strain: general framework, solution algorithm and [81] Dornelas V M, Oliveira S A, Savi M A, Pacheco P M C L and
application to shape memory alloys Int. J. Plast. 27 940–61 de Souza L F G 2021 Fatigue on shape memory alloys:
[80] Ju X, Moumni Z, Zhang Y, Zhang F, Zhu J, Chen Z and experimental observations and constitutive modeling Int. J.
Zhang W 2022 A multi-physics, multi-scale and finite strain Solids Struct. 213 1–24

19