**3. Hysteresis modeling: basics**

A one-dimensional numerical approach to simulate the stress-strain and strain-temperature hysteresis behavior of SMAs is described in this section. In particular, in a purely phenomenological way, the hysteresis loop is modeled by a Prandtl-Ishlinskii hysteresis operator *H* [26]; the basic idea of this approach consists in modeling the non-linear hysteretic behavior by a weighted superposition of many elementary hysteresis operators, such as the backlash operators *Hr*, as schematically illustrated in Figure 14.

$$H = \{w\}^T \{H\_r\} \tag{3}$$

where {*Hr*} is the vector of backlash operators and {*w*} is the corresponding vector of weights. As shown in Figure 14.a, each backlash operator *Hri* is characterized by its dead band width *dwi*, while the corresponding weight *wi* represents the slope of the oblique lines of the operator. As illustrated in Figure 14.b, which represent a generic *x*(*t*) − *y*(*t*) hysteretic

10 Will-be-set-by-IN-TECH

It is worth noting that Figure 11.b illustrates the nominal stress-strain behavior of a pseudoelastic SMA, while NiTi alloys exhibit a marked evolution of the stress-strain hysteretic behavior in the first mechanical cycles, up to a stable response, due to the formation of stabilized martensite [32], which causes a reduction of the pseudoelastic recovery of the SMA. In particular, Figure 12 illustrates the evolution of the stress-strain hysteretic behavior of the material in the first 25 mechanical cycles for a fixed value of total strain *εtot* = 3.5%. These data were obtained from experimental testing of a commercial pseudoelastic NiTi alloy [30]; the figure clearly shows a marked reduction of the pseudoelastic recovery, from 3.5% to about 3%, but the stress-strain loops becomes stable after 20 cycles. Furthermore, as reported in section

**Figure 12.** Evolution of the stress-strain response of a commercial pseudoelastic NiTi alloy during the

2.1.2, the critical stresses for direct and inverse phase transformation are significantly affected by the temperature, according to the Clausius-Clapeyron constant (Equation 1), as illustrated in the experimentally measured curves in Figure 13 [30]. These curves have been obtained from isothermal loading unloading cycles, carried out at increasing values of the testing temperature for *T* > *Af* (303*K* − 328*K*), by using a SMA with a stable pseudoelastic response. The analysis of the data in Figure 13 allows to obtain the value of the Clausius-Clapeyron

A one-dimensional numerical approach to simulate the stress-strain and strain-temperature hysteresis behavior of SMAs is described in this section. In particular, in a purely phenomenological way, the hysteresis loop is modeled by a Prandtl-Ishlinskii hysteresis operator *H* [26]; the basic idea of this approach consists in modeling the non-linear hysteretic behavior by a weighted superposition of many elementary hysteresis operators, such as the

where {*Hr*} is the vector of backlash operators and {*w*} is the corresponding vector of weights. As shown in Figure 14.a, each backlash operator *Hri* is characterized by its dead band width *dwi*, while the corresponding weight *wi* represents the slope of the oblique lines of the operator. As illustrated in Figure 14.b, which represent a generic *x*(*t*) − *y*(*t*) hysteretic

*<sup>H</sup>* <sup>=</sup> {*w*}*T*{*Hr*} (3)

first 25 mechanical cycles for a maximum applied deformation *εtot* = 3.5% [30]

backlash operators *Hr*, as schematically illustrated in Figure 14.

constant (*CM* = *CA* = 8.7*MPaK*−1).

**3. Hysteresis modeling: basics**

**Figure 13.** Isothermal stress-strain hysteresis loops of a commercial pseudoelastic NiTi alloy as a function of the testing temperature [30].

**Figure 14.** Schematic depiction of the hysteresis operators: a) elementary backlash operator *Hr* and b) Prandtl-Ishlinskii hysteresis operator *H* given by weighted superposition of elementary hysteresis operators.

behavior, the proposed approach consists of modeling the hysteretic loop by a linear piecewise discretization. The accuracy of the model can be improved by increasing the total number of linear pieces, which represent the total number of the backlash operators. The problem of modeling the hysteretic behavior, starting from the experimental measurements, is now reduced to the determination of the deadband width vector {*dw*} of the backlash operators and the associated gain vector {*w*}. In particular, the parameters of the model can be easily identified by the outer loop of the hysteretic region by using the following simple relation:

$$y\_k = \sum\_{i=1}^k (dw\_{k+1} - dw\_i) w\_i \tag{4}$$

where *yk* is the output value of the lower branch of the loop in the generic point of discontinuity *k*, as shown in Figure 14.b. The vector {*dw*} is a user defined discretization of

#### 12 Will-be-set-by-IN-TECH 132 Smart Actuation and Sensing Systems – Recent Advances and Future Challenges

the total amplitude of the input signal. Equation 4 can be rewritten in matrix form as follows:

$$\{y\} = [A]\{w\} \tag{5}$$

where the matrix [*A*] is constructed, for a given {*dw*} vector, by using equation 4; the unknown vector {*w*} can be found by solving a system of *N* linear equations, where *N* is the total number of backlash operators, as follows:

$$\{w\} = [A]^{-1} \{y\} \tag{6}$$

The main drawback of the Prandtl-Ishlinskii approach consists in the fact that only loops with an odd symmetry to the relative center can be modeled; in fact the symmetry of the elementary hysteresis operator, with respect to the center of the loop, persists under linear superposition. However, this drawback can be overcome by using a modified Prandtl-Ishlinksii operator, as described in details in [26], in which a weighted superposition of saturation operators is combined with the hysteresis operator. The parameters of this sub model, such as saturation limits {*S*} and associated gains {*ws*}, can be identified by using a procedure similar to that described above.

#### **3.1. Modeling of pseudoelastic effect**

The numerical model described in this section is able to simulate the pseudoelastic effect of a shape memory alloy [30], *i*.*e*. the the stress-strain (*σ* − *ε*) hysteretic behavior, based on the Prandtl-Ishlinksii operator and on the assumptions reported in the following.

#### *3.1.1. Basic assumptions*

Figure 15 shows the stabilized stress-strain hysteretic behavior of a commercial NiTi alloy, *i*.*e*. the response of the material after the first training cycles (see Figure 12), for different values of the applied deformation. The figure also illustrates the Young's moduli of austenite and detwinned martensitic structures, *EA* and *EM*, together with the generic young's modulus of the alloy, *E*(*ξM*), corresponding to an incomplete stress induced martensitic transformation, *i*.*e*. as a function of the martensite fraction *ξ<sup>M</sup>* (0 < *ξ<sup>M</sup>* < 1). In particular, *EA* represents the slope of the early stage of the stress-strain loading curve, *EM* is measured from the unloading curve of a complete martensitic transformation (*ξ<sup>M</sup>* = 1), while *E*(*ξM*) is obtained from the unloading path of an incomplete phase transformation.

The total strain *ε* can be decomposed in elastic and a transformation strain components, *εel* and *εtr*, respectively:

$$
\varepsilon = \varepsilon\_{el} + \varepsilon\_{tr} \tag{7}
$$

where the elastic strain can be expressed as a function of the applied stress, *σ*, and of the Young's modulus, *E*(*ξM*), of the material:

$$
\varepsilon\_{el} = \frac{\sigma}{E(\mathfrak{f}\_M)} \tag{8}
$$

As schematically shown in Figure 15 the Young's modulus changes during stress-induced phase transformation between austenite and martensite, *i*.*e*. it decreases from *EA* to *EM*, and

12 Will-be-set-by-IN-TECH

the total amplitude of the input signal. Equation 4 can be rewritten in matrix form as follows:

where the matrix [*A*] is constructed, for a given {*dw*} vector, by using equation 4; the unknown vector {*w*} can be found by solving a system of *N* linear equations, where *N* is

The main drawback of the Prandtl-Ishlinskii approach consists in the fact that only loops with an odd symmetry to the relative center can be modeled; in fact the symmetry of the elementary hysteresis operator, with respect to the center of the loop, persists under linear superposition. However, this drawback can be overcome by using a modified Prandtl-Ishlinksii operator, as described in details in [26], in which a weighted superposition of saturation operators is combined with the hysteresis operator. The parameters of this sub model, such as saturation limits {*S*} and associated gains {*ws*}, can be identified by using a procedure similar to that

The numerical model described in this section is able to simulate the pseudoelastic effect of a shape memory alloy [30], *i*.*e*. the the stress-strain (*σ* − *ε*) hysteretic behavior, based on the

Figure 15 shows the stabilized stress-strain hysteretic behavior of a commercial NiTi alloy, *i*.*e*. the response of the material after the first training cycles (see Figure 12), for different values of the applied deformation. The figure also illustrates the Young's moduli of austenite and detwinned martensitic structures, *EA* and *EM*, together with the generic young's modulus of the alloy, *E*(*ξM*), corresponding to an incomplete stress induced martensitic transformation, *i*.*e*. as a function of the martensite fraction *ξ<sup>M</sup>* (0 < *ξ<sup>M</sup>* < 1). In particular, *EA* represents the slope of the early stage of the stress-strain loading curve, *EM* is measured from the unloading curve of a complete martensitic transformation (*ξ<sup>M</sup>* = 1), while *E*(*ξM*) is obtained from the

The total strain *ε* can be decomposed in elastic and a transformation strain components, *εel*

where the elastic strain can be expressed as a function of the applied stress, *σ*, and of the

*<sup>ε</sup>el* <sup>=</sup> *<sup>σ</sup>*

As schematically shown in Figure 15 the Young's modulus changes during stress-induced phase transformation between austenite and martensite, *i*.*e*. it decreases from *EA* to *EM*, and

Prandtl-Ishlinksii operator and on the assumptions reported in the following.

{*w*} = [*A*]

the total number of backlash operators, as follows:

**3.1. Modeling of pseudoelastic effect**

unloading path of an incomplete phase transformation.

Young's modulus, *E*(*ξM*), of the material:

described above.

*3.1.1. Basic assumptions*

and *εtr*, respectively:

{*y*} = [*A*]{*w*} (5)

*ε* = *εel* + *εtr* (7)

*<sup>E</sup>*(*ξM*) (8)

<sup>−</sup>1{*y*} (6)

**Figure 15.** Stress-strain hysteretic behavior of a commercial NiTi alloy together with the Young's moduli of austenite and detwinned martensite [30].

it is assumed to be dependent on the volume fraction of martensite *ξ<sup>M</sup>* according to the Reuss formula [20]:

$$\frac{1}{E(\mathfrak{f}\_M)} = \frac{\mathfrak{f}\_M}{E\_M} + \frac{1 - \mathfrak{f}\_M}{E\_A} \tag{9}$$

The evolution of martensite is assumed to be a linear function of the stress in the stress-strain transformation curves, *i*.*e*. in the range (*σAM <sup>s</sup>* , *σAM <sup>f</sup>* ) in the loading stage and (*σMA <sup>s</sup>* , *σMA <sup>f</sup>* ) in the unloading stage. In particular, the evolution rule for a complete transformation, can be expressed as follows:

$$
\tilde{\zeta}\_M = \begin{cases}
\frac{\sigma - \sigma\_s^{AM}}{\sigma\_f^{AM} - \sigma\_s^{AM}} & \text{Loading path} \\
1 - \frac{\sigma - \sigma\_s^{AM}}{\sigma\_f^{MA} - \sigma\_s^{AM}} & \text{Unloading path}
\end{cases}
\tag{10}
$$

The assumptions given by equations (9) and (10) have been validated by experimental measurements of the Young's modulus, as reported in Figure16. In particular, the figure shows the measured values of the Young's modulus as a function of the applied stress, on the direct stress-strain transformation plateau, together with the predictions of the Reuss formula, and a satisfactory agreement is observed.

However, it's worth noting that the evolution of martensite is characterized by a hysteretic behavior, *i*.*e*. it is stress path dependent, and equations (9) and (10) do not correctly predict the effects of incomplete transformations. These effects are simulated mathematically, by using a backlash operator, as schematically illustrated in Figure 17; in particular, Figure 17.a shows the stress *vs* time path, while Figure 17.b illustrates the evolution of *ξ<sup>M</sup>* and 1/*E*(*ξM*) *vs* the applied stress. Furthermore, the continuous lines in Figure 17.b are relative to a complete martensitic transformation, *i*.*e*. *ξ<sup>M</sup>* increases from 0 to 1, while the dashed lines show the effects of an incomplete transformation. Due to the modification in Young's modulus during the phase transformation between austenite and martensite, as shown in Figure 17.b, both elastic and transformation strain components are represented by a hysteretic behavior and it can be calculated using equations (7-10); as an example in Figure 18 a typical *σ* − *ε* hysteresis loop for a complete phase transformation is compared with the corresponding computed elastic strain (*σ* − *εel*) and transformation strain (*σ* − *εtr*) hysteresis loops.

**Figure 16.** Evolution of the Young's modulus as a function of the applied stress: experimental measurements *vs* simulations [30].

**Figure 17.** Evolution rule of the martensite fraction *ξ<sup>M</sup>* in the tension cycle: a) example of stress-time path and b) simple hysteresis model to predict the Young's modulus *E*(*ξM*) [30].

**Figure 18.** Stress strain (*σ* − *ε*) hysteresis loops: elastic strain *εel* and transformation strain *εtr* [30].

#### *3.1.2. Numerical flowchart*

14 Will-be-set-by-IN-TECH

**Figure 16.** Evolution of the Young's modulus as a function of the applied stress: experimental

**Figure 17.** Evolution rule of the martensite fraction *ξ<sup>M</sup>* in the tension cycle: a) example of stress-time

**Figure 18.** Stress strain (*σ* − *ε*) hysteresis loops: elastic strain *εel* and transformation strain *εtr* [30].

path and b) simple hysteresis model to predict the Young's modulus *E*(*ξM*) [30].

measurements *vs* simulations [30].

Based on the Prandtl-Ishlinskii hysteresis modeling approach together with the model assumptions described in the previous section, a one-dimensional numerical model can be easily developed by using commercial software packages. In particular, Figure 19 illustrates an implementation of the model in the *MatlabSimulinkTM* platform, in which the two strain components, *εel* and *εtr*, are treated separately. As clearly illustrated in Figure 19, both stress

**Figure 19.** Flowchart of a *MatlabSimulinkTM* model to simulate the pseudoelastic effect in SMAs [30].

and temperature are taken as input variables, and the critical stresses are calculated, based on the current value of the temperature, by using the Clausius-Clapeyron relation (equation 1). In the model the saturation operators are used to identify the stress range where the material exhibits the hysteretic behavior (*σMA <sup>f</sup>* <sup>&</sup>lt; *<sup>σ</sup>* <sup>&</sup>lt; *<sup>σ</sup>AM <sup>f</sup>* ), and the dead zones of transformation where only elastic deformation of austenite (*σ* < *σMA <sup>f</sup>* ) or martensite (*<sup>σ</sup>* <sup>&</sup>gt; *<sup>σ</sup>AM <sup>f</sup>* ) occurs. Three different sub-models are highlighted in Figure 19 which simulate the *σ* − *εtr* and *σ* − *εel* hysteresis loops, and the linear elastic response of the material in austenitic and martensitic conditions. In particular, a Prandtl-Ishlinskii hysteresis operator, was used to model the *σ* − *εtr* hysteretic behavior, a single backlash operator was adopted to model the *σ* − *εel* loop (see Figure 17.b), and a gain operator was used to describe the linear elastic response for *σ* < *σMA f* and for *σ* > *σAM <sup>f</sup>* . The parameters of the Prandtl-Ishlinskii operator describing the *σ* − *εtr* hysteretic behavior, *i*.*e*. the deadband width vector {*dw*} and the associated gain vector {*w*} are determined from an experimentally measured stress-strain hysteresis loop by using the procedure described in section 3. In particular, as shown in Figure 20, the generic input and output variables (*x* and *y*) can be regarded as the stress and strain values (*σ* and *ε*), respectively. As a consequence, the vector {*dw*} represents a user defined discretization of

#### 16 Will-be-set-by-IN-TECH 136 Smart Actuation and Sensing Systems – Recent Advances and Future Challenges

the stress amplitude while the vector {*w*} can be determined by using equation (6), where the vector {*y*} is the vector of strain values {*ε*}. The computational procedure to obtain the aforementioned model parameters from an experimentally measured hysteretic loop can be easily implemented, which generates the vector {*dw*}, by a partition of the input stress amplitude, and calculates the unknown vector of weight {*w*} by solving a system of *N* linear equations, where *N* is the number of backlash operators.

**Figure 20.** Loading branch of the stress-strain hysteresis loop and linear piecewise discretization obtained by a weighted superposition of three backlash operators [30].

#### *3.1.3. Numerical results vs experiments*

The accuracy of the numerical method is illustrated by comparisons with experimentally measured hysteresis loops, by using he thermo-mechanical parameters reported in [30] (*EA*=39 GPa, *EM*=20 GPa, *σAM <sup>s</sup>* =260 MPa, *σAM <sup>f</sup>* =385 MPa, *<sup>σ</sup>MA <sup>s</sup>* =250 MPa, *σMA <sup>f</sup>* =125 MPa, *CA* = *CM*=10.3 MPa/K). Figure 21 shows the stress-strain hysteretic behavior of the SMA for a stress path which involves several incomplete stress-induced martensitic transformations (*A* → *M*), by repeated isothermal tension cycles (*T* = 303*K*) carried out between *σmin* = 0 and decreasing values of *σmax* < *σAM <sup>f</sup>* . In particular, the comparison between numerical simulations and experimental measurements, illustrated in Figure 21.a, clearly shows a satisfactory accuracy of the model in predicting the non-linear stress-strain hysteretic behavior of the material; the figure also shows that the model is able to capture the modification of Young's modulus in the stress-strain transformation curve, as it correctly predicts the change in the slopes of the unloading curves. Furthermore, Figure 21.b shows the evolutions of the transformation strain, *εtr*, and elastic strain, *εel*, for the same input stress path of Figure 21.a. As illustrated in section 3.1.1, both *εtr* and *εel* are characterized by hysteretic behaviors, which are due to the mismatch between the critical stresses in the stress induced transformations and the modification of Young's modulus in the stress-strain transformation curve, respectively.

Figures 22 show comparisons between numerical predictions and experimental measurements for two input isothermal stress paths (*T* = 303*K*) which involve incomplete *M* → *A* transformations (Figure 22.a) and both incomplete *A* → *M* and *M* → *A* transformations (Figure 22.b). In particular, Figure 22.a shows the hysteretic response of the material for a stress path which is composed of some subsequent tension cycles between increasing values of *σmin* < *σMA <sup>f</sup>* and *σmax* = *const*, while Figure 22.b is relative to a stress

**Figure 21.** Numerical simulation for isothermal stress cycles with incomplete *A* → *M* transformations: a) comparison with experimentally measured loops and b) evolution of transformation strain, *εtr*, and elastic strain, *εel* [30].

path which involves different subsequent tension cycles carried out between increasing values of *σmin* > *σMA <sup>f</sup>* and decreasing values of *<sup>σ</sup>max* <sup>&</sup>lt; *<sup>σ</sup>AM <sup>f</sup>* ; both figures shows good agreements between experiments and numerical simulations.

**Figure 22.** Comparison between numerical simulations and experimentally measured hysteresis loops for isothermal stress cycles with: a) incomplete *M* → *A* transformations and b) incomplete *A* → *M* and *M* → *A* transformations [30].

#### **3.2. Modeling of two-way shape memory effect**

16 Will-be-set-by-IN-TECH

the stress amplitude while the vector {*w*} can be determined by using equation (6), where the vector {*y*} is the vector of strain values {*ε*}. The computational procedure to obtain the aforementioned model parameters from an experimentally measured hysteretic loop can be easily implemented, which generates the vector {*dw*}, by a partition of the input stress amplitude, and calculates the unknown vector of weight {*w*} by solving a system of *N* linear

**Figure 20.** Loading branch of the stress-strain hysteresis loop and linear piecewise discretization

*<sup>s</sup>* =260 MPa, *σAM*

The accuracy of the numerical method is illustrated by comparisons with experimentally measured hysteresis loops, by using he thermo-mechanical parameters reported in [30]

*CA* = *CM*=10.3 MPa/K). Figure 21 shows the stress-strain hysteretic behavior of the SMA for a stress path which involves several incomplete stress-induced martensitic transformations (*A* → *M*), by repeated isothermal tension cycles (*T* = 303*K*) carried out between *σmin* = 0

simulations and experimental measurements, illustrated in Figure 21.a, clearly shows a satisfactory accuracy of the model in predicting the non-linear stress-strain hysteretic behavior of the material; the figure also shows that the model is able to capture the modification of Young's modulus in the stress-strain transformation curve, as it correctly predicts the change in the slopes of the unloading curves. Furthermore, Figure 21.b shows the evolutions of the transformation strain, *εtr*, and elastic strain, *εel*, for the same input stress path of Figure 21.a. As illustrated in section 3.1.1, both *εtr* and *εel* are characterized by hysteretic behaviors, which are due to the mismatch between the critical stresses in the stress induced transformations and the modification of Young's modulus in the stress-strain transformation curve, respectively. Figures 22 show comparisons between numerical predictions and experimental measurements for two input isothermal stress paths (*T* = 303*K*) which involve incomplete *M* → *A* transformations (Figure 22.a) and both incomplete *A* → *M* and *M* → *A* transformations (Figure 22.b). In particular, Figure 22.a shows the hysteretic response of the material for a stress path which is composed of some subsequent tension cycles between

*<sup>f</sup>* =385 MPa, *<sup>σ</sup>MA*

*<sup>s</sup>* =250 MPa, *σMA*

*<sup>f</sup>* . In particular, the comparison between numerical

*<sup>f</sup>* and *σmax* = *const*, while Figure 22.b is relative to a stress

*<sup>f</sup>* =125 MPa,

equations, where *N* is the number of backlash operators.

obtained by a weighted superposition of three backlash operators [30].

*3.1.3. Numerical results vs experiments*

and decreasing values of *σmax* < *σAM*

(*EA*=39 GPa, *EM*=20 GPa, *σAM*

increasing values of *σmin* < *σMA*

The numerical model described in this section is able to simulate the two-way shape memory effect of a trained NiTi based shape memory alloy, *i*.*e*. the the strain-temperature (*ε* − *T*) hysteretic behavior [27–29]. Furthermore, the model is able to capture the effects of applied stress on the thermal recovery of the alloy.

#### 18 Will-be-set-by-IN-TECH 138 Smart Actuation and Sensing Systems – Recent Advances and Future Challenges

#### *3.2.1. Basic assumptions*

In a pure phenomenological way, the variation of the two-way shape memory strain, *εtw*, with increasing the applied stress (see Figure 10.b) can be attributed to two different mechanisms: *i*) the variation of Young's modulus in the thermal hysteresis behavior between martensite and austenite, *ii*) the increased volume fraction of favorably oriented martensite variants with increasing external stress.

$$
\Delta \varepsilon\_{tw} = \varepsilon\_{tw} - \varepsilon\_{tw0} = \left(\frac{1}{E\_M} - \frac{1}{E\_A}\right)\sigma + c\sigma \tag{11}
$$

where subscript 0 indicate the stress-free condition, *EM* and *EA* represent the Young's moduli of martensite and austenite, respectively, while *c* can be assumed as a material constant. In particular, the first term at the right end side of equation 11, indicated as Δ*εmech*, describe the first effect while the last terms, namely Δ*εmem*, take into account the second effect. Young's moduli *EM* and *EA* can be measured by isothermal tensile tests carried out at temperatures *T* < *Mf* and *T* > *Af* , respectively, while the parameter *c* can be obtained from experimental measurements of the two-way shape memory strain carried out at different value of the applied stress. In the following the values of the parameters reported in [28] have been used (*EM* <sup>=</sup> <sup>36</sup> <sup>∗</sup> 103*MPa*, *EA* <sup>=</sup> <sup>67</sup> <sup>∗</sup> <sup>10</sup>3*MPa*, *<sup>c</sup>* <sup>=</sup> 8.5 <sup>∗</sup> <sup>10</sup>−5*MPa*−1). Two simple mathematical functions, obtained by a numerical fitting of the experimental data, are used to describe the phase transformation kinetics. In particular, as reported in the model by Tanaka and Nagaki [13], the heating and cooling branches of the hysteresis loop can be represented by two exponential curves:

$$\varepsilon = \begin{cases} \varepsilon\_{tw} e^{a\_A(A\_f - T + \delta T\_M)} & \mathbf{M} \to \mathbf{A} \\ \varepsilon\_{tw} [1 - e^{a\_M(M\_f - T + \delta T\_A)}] & \mathbf{A} \to \mathbf{M} \end{cases} \tag{12}$$

where *εtw*, *Ms*, and *As* are functions of the applied stress, as reported in Equations 11 and 1; *aM*, *aA*, *δTM*, and *δTA*, which define the shape of the heating and cooling branches of the loop, can be identified by a numerical fitting of the experimental data. Figure 23.a shows a comparison between experimental measurement and exponential curves; in the figure, points *P*<sup>1</sup> and *P*<sup>2</sup> represent the range where the numerical fitting is executed to identify the parameters of the heating branch of the loop, while points *P*<sup>3</sup> and *P*<sup>4</sup> are relative to the cooling branch. In Figure 23.b, a linear fitting between points *P*<sup>1</sup> and *P*<sup>2</sup> of the experimental data in the *T* − *logε* plane is shown, where the slope of the line defines the parameter *aA* and the intersection with the *logε* axis allows us to obtain the parameter *δTM*. If the loop is characterized by an odd symmetry with respect to its center, as is quite well observed in the investigated material, the same values can be assumed for the constants *aA* and *aM* and *δTA*, and *δTM*. The two exponential curves describe the two branches of the hysteresis loop in a parametric way for a generic value of the applied stress by using Equations 1 and 11. Starting from the curves *ε* − *T*, the numerical method based on the Prandtl-Ishlinksii operator can be developed, which is able to predict the output response for a generic temperature path, as decribed in the following section.

#### *3.2.2. Numerical flowchart*

The numerical method described above can be easily implemented in a *MatlabSimulinkTM* model, as shown in Figure 24, by a modified Prandtl-Hishlinkii hysteresis operator. The

18 Will-be-set-by-IN-TECH

In a pure phenomenological way, the variation of the two-way shape memory strain, *εtw*, with increasing the applied stress (see Figure 10.b) can be attributed to two different mechanisms: *i*) the variation of Young's modulus in the thermal hysteresis behavior between martensite and austenite, *ii*) the increased volume fraction of favorably oriented martensite variants with

> 1 *EM*

*<sup>ε</sup>tweaA*(*Af* <sup>−</sup>*T*+*δTM* ) <sup>M</sup> <sup>→</sup> <sup>A</sup>

where *εtw*, *Ms*, and *As* are functions of the applied stress, as reported in Equations 11 and 1; *aM*, *aA*, *δTM*, and *δTA*, which define the shape of the heating and cooling branches of the loop, can be identified by a numerical fitting of the experimental data. Figure 23.a shows a comparison between experimental measurement and exponential curves; in the figure, points *P*<sup>1</sup> and *P*<sup>2</sup> represent the range where the numerical fitting is executed to identify the parameters of the heating branch of the loop, while points *P*<sup>3</sup> and *P*<sup>4</sup> are relative to the cooling branch. In Figure 23.b, a linear fitting between points *P*<sup>1</sup> and *P*<sup>2</sup> of the experimental data in the *T* − *logε* plane is shown, where the slope of the line defines the parameter *aA* and the intersection with the *logε* axis allows us to obtain the parameter *δTM*. If the loop is characterized by an odd symmetry with respect to its center, as is quite well observed in the investigated material, the same values can be assumed for the constants *aA* and *aM* and *δTA*, and *δTM*. The two exponential curves describe the two branches of the hysteresis loop in a parametric way for a generic value of the applied stress by using Equations 1 and 11. Starting from the curves *ε* − *T*, the numerical method based on the Prandtl-Ishlinksii operator can be developed, which is able to predict the output response for a generic temperature path, as

The numerical method described above can be easily implemented in a *MatlabSimulinkTM* model, as shown in Figure 24, by a modified Prandtl-Hishlinkii hysteresis operator. The

*<sup>ε</sup>tw*[<sup>1</sup> <sup>−</sup> *<sup>e</sup>aM*(*Mf* <sup>−</sup>*T*+*δTA* ) <sup>A</sup> <sup>→</sup> <sup>M</sup> (12)

where subscript 0 indicate the stress-free condition, *EM* and *EA* represent the Young's moduli of martensite and austenite, respectively, while *c* can be assumed as a material constant. In particular, the first term at the right end side of equation 11, indicated as Δ*εmech*, describe the first effect while the last terms, namely Δ*εmem*, take into account the second effect. Young's moduli *EM* and *EA* can be measured by isothermal tensile tests carried out at temperatures *T* < *Mf* and *T* > *Af* , respectively, while the parameter *c* can be obtained from experimental measurements of the two-way shape memory strain carried out at different value of the applied stress. In the following the values of the parameters reported in [28] have been used (*EM* <sup>=</sup> <sup>36</sup> <sup>∗</sup> 103*MPa*, *EA* <sup>=</sup> <sup>67</sup> <sup>∗</sup> <sup>10</sup>3*MPa*, *<sup>c</sup>* <sup>=</sup> 8.5 <sup>∗</sup> <sup>10</sup>−5*MPa*−1). Two simple mathematical functions, obtained by a numerical fitting of the experimental data, are used to describe the phase transformation kinetics. In particular, as reported in the model by Tanaka and Nagaki [13], the heating and cooling branches of the hysteresis loop can be represented by

<sup>−</sup> <sup>1</sup> *EA* 

*σ* + *cσ* (11)

Δ*εtw* = *εtw* − *εtw*<sup>0</sup> =

*ε* = 

*3.2.1. Basic assumptions*

increasing external stress.

two exponential curves:

decribed in the following section.

*3.2.2. Numerical flowchart*

**Figure 23.** Numerical fitting of the experimental data: a) comparison between exponential curves and experimental measurements; b) numerical fitting in the *T* − *logε* plane to identify the parameters of the exponential curve in the heating branch of the hysteresis loop [28]

**Figure 24.** Flowchart of a *MatlabSimulinkTM* model to simulate the two way shape memory effect in SMAs.

Prandtl-Hishlinkii operator is implemented in the submodel #1 of Figure 24, by a weighted superposition of several backlash operators, and the corresponding parameters, *i*.*e*. the deadband width vector {*dw*} and the associated gain vector {*w*}, are determined from the exponential curves of equation 12, which in turn are obtained from a fitting of experimental data (*T* − *ε*) as illustrated in Figure 23. In particular, the generic input and output variables (*x* and *y*) can be regarded as the temperature and strain values (*T* and *ε*), respectively. As a consequence, the vector {*dw*} represents a user defined discretization of the temperature amplitude while the vector {*w*} can be determined by using equation (6), where the vector {*y*} is the vector of strain values {*ε*} obtained from equation 12 (*ε<sup>i</sup>* = *ε*(*Ti*)). The saturation operator in the submodel #1 is used to simulate the dead zones of transformation, *i*.*e*. the

#### 20 Will-be-set-by-IN-TECH 140 Smart Actuation and Sensing Systems – Recent Advances and Future Challenges

material behavior when the temperature is above *Af* during heating and below *Mf* in the cooling branch of the hysteresis loop. In particular, this operator imposes upper and lower bounds on the temperature, which are *Af* and *Mf* , respectively, so that when the temperature is outside these bounds the signal is clipped to the upper or lower bounds. In Figure 25.a, the response of the submodel #1, in terms of strain-temperature loop, is illustrated for a thermal cycle between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* ; the figure clearly shows that the hysteretic behavior of the material is properly described in the range of temperatures between *Mf* and *Af* , and the dead zone of transformation, when the temperature exceeds *Af* or falls below *Mf* , are also simulated.

**Figure 25.** Numerically simulated loop for a thermal cycle between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* obtained by: a) Prandtl-Ishlinskii model and b) modified Prandtl-Ishlinskii model [28].

Unfortunately, when comparing the experimental results with the numerically simulated loops, high errors are observed in the extremity of the hysteretic region, *i*.*e*. when the temperature is below *Mf* during cooling and above *Af* during heating. To overcome this limitation, a modified *Simulink* model can be implemented by including the submodel #2 of Figure 24. This latter uses two subsystems, for the heating and cooling branches of the loop, which modifies the output response of the system when the temperature is near *Mf* and *Af* . In particular, each subsystem implements a weighted superposition of several dead band operators, which is executed by a series of a dead band block and a gain block, while the saturation block assures that the correction is carried out only in a limited range of temperatures near *Mf* and *Af* . Figure 25.b shows a simulated hysteresis loop, obtained by the modified model, between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* ; the figure clearly shows that the model allows a better simulation of the extremity of the hysteretic region with respect to Figure 25.a.

#### *3.2.3. Numerical results vs experiments*

In this section, the accuracy and efficiency of the 1D numerical model are illustrated by comparing some experimentally measured hysteresis loops with the corresponding numerical predictions. The simulations have been carried out by using a model with 20 backlash operators and 5 dead zone operators to modify the loops in the extremity of the hysteretic region. Figure 26.a shows a comparison between the experimentally measured hysteresis loop for a stress-free martensitic transformation, between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* , and the numerically simulated one; the figure clearly shows a good accuracy of the numerical model with very small errors. In Figure 26.b, a comparison between numerical predictions and experimental results when the material is subjected to a tensile stress *σ* = 100*MPa* is shown. Also, in this case a satisfactory agreement is observed, but the errors increase with respect to the stress free condition. However, it is important to point out that the model parameters were identified by using the measured hysteresis loop under stress-free conditions, and by applying Equations 1 and 11 to modify both the PTTs and *εtw*.

20 Will-be-set-by-IN-TECH

material behavior when the temperature is above *Af* during heating and below *Mf* in the cooling branch of the hysteresis loop. In particular, this operator imposes upper and lower bounds on the temperature, which are *Af* and *Mf* , respectively, so that when the temperature is outside these bounds the signal is clipped to the upper or lower bounds. In Figure 25.a, the response of the submodel #1, in terms of strain-temperature loop, is illustrated for a thermal cycle between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* ; the figure clearly shows that the hysteretic behavior of the material is properly described in the range of temperatures between *Mf* and *Af* , and the dead zone of transformation, when the temperature exceeds *Af* or falls

**Figure 25.** Numerically simulated loop for a thermal cycle between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* obtained by: a) Prandtl-Ishlinskii model and b) modified Prandtl-Ishlinskii model [28].

Unfortunately, when comparing the experimental results with the numerically simulated loops, high errors are observed in the extremity of the hysteretic region, *i*.*e*. when the temperature is below *Mf* during cooling and above *Af* during heating. To overcome this limitation, a modified *Simulink* model can be implemented by including the submodel #2 of Figure 24. This latter uses two subsystems, for the heating and cooling branches of the loop, which modifies the output response of the system when the temperature is near *Mf* and *Af* . In particular, each subsystem implements a weighted superposition of several dead band operators, which is executed by a series of a dead band block and a gain block, while the saturation block assures that the correction is carried out only in a limited range of temperatures near *Mf* and *Af* . Figure 25.b shows a simulated hysteresis loop, obtained by the modified model, between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* ; the figure clearly shows that the model allows a better simulation of the extremity of the hysteretic region with

In this section, the accuracy and efficiency of the 1D numerical model are illustrated by comparing some experimentally measured hysteresis loops with the corresponding numerical predictions. The simulations have been carried out by using a model with 20 backlash operators and 5 dead zone operators to modify the loops in the extremity of the hysteretic region. Figure 26.a shows a comparison between the experimentally measured hysteresis loop

below *Mf* , are also simulated.

respect to Figure 25.a.

*3.2.3. Numerical results vs experiments*

**Figure 26.** Comparison between experimental measurements and numerical predictions for a thermal cycle between the temperatures *T*<sup>0</sup> < *Mf* and *T*<sup>1</sup> > *Af* under a) stress-free condition and b) tensile strees of 100 MPa [28].

The accuracy of the model was also analyzed when the material is subjected to partial thermal cycles, *i*.*e*. to incomplete martensitic transformations. Figures 27.a and 27.b show the hysteretic behavior of the material for two different temperature-time paths in stress-free conditions; in particular, Figure 27.a shows incomplete *A* → *M* transformations, while Figure

**Figure 27.** Comparison between experimental measurements and numerical predictions for two different temperature-time paths: a) incomplete *A* → *M* transformations and b) incomplete *M* → *A* transformations [28].

#### 22 Will-be-set-by-IN-TECH 142 Smart Actuation and Sensing Systems – Recent Advances and Future Challenges

27.b illustrates incomplete *M* → *A* transformations. As shown in the figures, the comparison between experimental measurements and numerical predictions show a good accordance in both cases. It is worth noting that the same model parameters of the first example (Figure 26.a) were used in these numerical simulations.
