Linear Shape Memory Alloy Thermomechanical Actuators

Abstract

The thermally activated changes in crystal structure in nickel-titanium shape memory alloy (SMA) material, which produces transformation strains of an order of magnitude higher than and opposite to the thermal strains, have been described using different models. Some of these models are defined by cyclic functions (trigonometric functions) which they become complex to control when they are subjected to a wide range of transformation temperatures. This chapter presents an alternative model to better describe the behavior of SMA for a more general temperature range, which an SMA-powered actuator might be subjected to. The proposed model is then implemented to the analysis of one-dimensional (1D) problem oriented to the two-dimensional (2D) space. The simulation results were then graphically compared to the experiment.

1. Introduction

Shape memory alloys (SMAs) are being described as smart materials due to their ability of memorizing shape. Their ability in memorizing shapes requires some training suitable for user’s application. The trained SMA exhibits two behaviors, that is, superelasticity and shape memory effect depending on the temperature. Those behaviors are triggered by the phenomenon occurring at the microstructural level called phase transformation. Phase transformation can be induced by stress or temperature. When SMA is mechanically loaded at room temperature (or at a temperature below thermal transition), the twinned martensite transforms to detwinned martensite [1]. When the detwinned martensite is subjected to heat, the detwinned martensite transforms to austenite phase. When the fully austenite SMA is cooled down, there will be a temperature range where austenite will transform back to twinned martensite. This whole process is called shape memory effect. When a mechanical load is applied on SMA which is subjected to temperature above austenitic finish temperature, the deformation will take place. The deformation will be reversed completely upon the removal of mechanical load. This phenomenon is called superelasticity. These two physical behaviors are the ones that make SMA a good candidate for actuation applications. There are a number of models reported in the study which are meant to describe the phase transformation kinetics [1]. The first model was developed by Tanaka and Nagaki [2, 3] which was later improved by Boyd and Lagoudas [4]. The firstly developed model was based on exponential functions. Boyd and Lagoudas redefined the material constants so as to suit their application. Later, a model based on cosine functions was developed by Liang and Rogers [5]. Liang and Rogers developed their model with the purpose of applying it to the analysis of acoustic vibrations control [5, 6, 7]. The common factor in these models is the calculation material constants which are calculated based on transformation temperatures. Brinson [8] developed another alternative phase transformation model, where the internal variable measuring the evolution of martensite phase was split into two.

The first internal variable is induced by temperature and the second one is induced by stress. This work proposed that the elastic modulus should be split into austenitic and martensitic modulus. It should be noted that these models developed by the previously listed authors are functions of transformation temperature (martensite start temperature, martensite finish temperature, austenite start temperature and austenite finish temperature). These temperatures are crucial since they are used to calculate certain material constants. More work has been done in improving the already existing models [9, 12, 13, 14, 15]. This chapter proposes a one-dimensional (1D) new model to describe the physical behavior of SMA for all transformation ranges. The proposed model is then used to simulate 1D structural deflection which is orientated in two-dimensional (2D) space. The simulated results are compared to the experimental results.

2. Phase transformation models

This section shows the already available phase transformation models. It should be noted that the law adopted for the proposed model is based on the one adopted by Liang and Rogers [5] with the difference on the phase transformation estimation. The first phase transformation model based on exponential functions was developed by Tanaka et al. [2, 3]. The transformation from martensite to austenite is given by

$$\xi ={\xi }_0\text{exp}\left[a_A\left(T-A_s\right)-b_A\sigma \right]$$E1

And the transformation from austenite to martensite is given by

$$\xi =1-\text{exp}\left[a_M\left(T-M_s\right)-b_M\sigma \right]+{\xi }_0$$E2

where ${\xi }_0$, $a_M$ and $b_M$, $M_s$, $\sigma$, $a_A$, $b_A$ and $A_s$ are internal variables measuring martensite volume fraction evolution, martensitic material constants, martensite start temperature, austenitic material constants and austenite start temperature, respectively. Boyd and Lagoudas [11] developed the calculation of material constants ($a_M$, $b_M$, $a_A$ and $b_A$) as follows:

$$a_M=\frac{2\text{ln}10}{M_s-M_f},b_M=\frac{a_M}{C_M},a_A=\frac{2\text{ln}10}{A_f-A_s},b_A=\frac{a_A}{C_A}$$E3

Later, Liang and Rogers [5] presented another alternative model for martensite volume fraction calculation. Their model was based on trigonometric functions (cyclic functions), that is, cosine functions. The following equation describes the phase transformation during heating (martensite transforms to austenite):

$$\xi =\frac{{\xi }_0}{2}\left\{1+\text{cos}\left[A_A\left(T-A_s-\frac{\sigma }{C_A}\right)\right]\right\}$$E4

This transformation occurs between the stress ranges described as follows:

$$C_A\left(T-A_f\right)<\sigma <C_A\left(T-A_s\right)$$E5

The transformation occurs during cooling (austenite transforms to martensite) and is given by

$$\xi =\frac{1-{\xi }_0}{2}\text{cos}\left[A_M\left(T-M_f-\frac{\sigma }{C_M}\right)\right]+\frac{1+{\xi }_0}{2}$$E6

The stress restriction is described as follows:

$$C_M\left(T-M_s\right)<\sigma <C_M\left(T-M_f\right)$$E7

The material constants involved in Liang and Rogers [5] are the same as those used in Tanaka’s model [2]. The difference comes with the calculation of the two constants which are given as follows:

$$A_M=\frac{\pi }{M_s-M_f};A_A=\frac{\pi }{A_f-A_s}$$E8

It should be noted that the calculation of the two material constants, $A_A$ and $A_M$, involves the difference in transformation temperatures. It was noticed that the previously mentioned models were not compatible with our application due to the size of our temperature range. The new non-cyclic model had to be developed. The proposed phase transformation model is based on the hyperbolic functions. The phase transformation during heating (martensite transforms to austenite) is given by the following equation:

$$\xi ={\xi }_0\left\{1-\text{tanh}\left[A_A\left(T-A_s\right)+{\vartheta }_A\right]\right\}$$E9

The phase transformation during cooling (austenite transforms to martensite) is described by the following equation:

$$\xi ={\xi }_0\left\{1-\text{tanh}\left[A_M\left(T-M_f\right)+{\vartheta }_M\right]\right\}$$E10

The calculation of material constants $A_A\text{and}{A}_M$ is similar to that of Eq. (8). The phase transformation constants ($A_s\text{and}{M}_s$) are those described previously. The constants ${\vartheta }_A\text{and}{\vartheta }_M$ are given by

$${\vartheta }_A=-\frac{A_A}{C_A}\sigma ;{\vartheta }_M=-\frac{A_M}{C_M}\sigma$$E11

The constants $C_A\text{and}{C}_M$ are the same as those used by Liang and Rogers. The developed phase transformation model was then used in performing the 2D analysis on a steel beam deflection explained in the next section. Prior to the performance of simulation analysis, an experiment was performed with the purpose of obtaining simulation parameters and verification of numerical results.

3. Beam deflection test experiment

The experiment was performed using the following equipment:

1. Explorer GLX Data logger

2. Nickel titanium shape memory alloy wire of 0.5-mm diameter

3. ISO-TECH IPS 2303 Laboratory DC Power Supply

4. Linear variable differential transformer (LVDT)

5. Fluke 190C Scope meter

6. Rectangular steel beam

7. Rectangular steel frame

The steel beam used in the performance of this experiment was 300-mm long with a 25-mm breadth and 1.2-mm thickness. Both ends of the steel beam were fixed on top of the vertical rectangular steel frame. The effective length for the steel beam after installation to the vertical rectangular steel frame was 243 mm. Figure 1 shows the complete experimental setup.

Figure 1.

Beam test experiment setup.

A 184-mm long NiTi SMA wire with a 0.5-mm diameter was attached with one end at the center of the steel beam and the other end at the center foot of the steel frame. The attachment of the NiTi SMA wire to the frame and to the beam was isolated using high-melting temperature plastic tubes. The isolation was required since the NiTi SMA wire was heated through the joules heating method. The NiTi SMA wire was attached at the center of the beam so as to get the maximum deflection. During the performance of the experiment, two types of data were being logged, that is, beam deflection and the NiTi SMA wire temperature. The beam deflection was recorded through the LVDT which was powered by the DC power supply (shown in Figure 1). The LVDT data were logged in the form of output voltage and it was translated into millimeters using the method described subsequently. The temperature was recorded using a GLX explorer data logger. A 1.54 A was supplied to the NiTi SMA wire, and the temperature of the wire increased from the room temperature to a maximum temperature. The wire was allowed to cool down after reaching the maximum deflection of the beam, and then the beam could go back to its rest position. This process was repeated several times so as to get stable results.

4. Deflection test results

The following figures show the experimental results for the beam deflection test and these results are used to validate FEA results. Although there were several tests conducted, only one graph is presented to avoid repetition. It should be noted that the results presented are already translated into millimeters using the method explained in the study [10]. Figure 2 shows the time it takes for a beam to finish the full deflection cycle, that is, the deflection from the rest position to the final position or the maximum deflection and from the final position back to the rest position.

Figure 2.

Time-deflection curve.

During heating, the beam deflects from the rest position to the maximum deflection and the reverse during cooling. It is clearly seen from the figure that the beam’s maximum deflection is reached in less than 10 s which is very quick. This quickness is suggested to be related to the wire diameter. The smaller the diameter, the faster the response of the wire. It is noted in the figure that there is a flatness of the graph between 10 and 70 s. This flatness shows the beginning of the cooling cycle. It should be remembered that the driving force behind this graph is the NiTi SMA wire phase transition as discussed in Section 1. The martensite phase transforms to austenite phase between 0 and 10 s upon heating. The austenite phase is the dominant phase between 10 and 70 s, and then the nucleation of martensite starts to occur after 70 s. Martensite start occurs after 70 s and finishes after 140 s. The flatness behavior is seen after 140 s which depicts the domination of martensite phase. The number of data points between 0 and 10 s are few compared to the rest of the graph, and this is caused by the fact that the sampling rate was not easy to control.

The SMA wire force, which produced the maximum deflection of the steel beam, is shown in Figure 3. The maximum force produced by the current SMA wire (184-mm long) was found to be approximately 35.87517 N (~3.5 kg).

Figure 3.

SMA wire force-temperature graph.

Figure 3 shows the deflection of the beam to the negative vertical (y) direction as the temperature increases and the reverse during cooling. It is noted that all the values on the y-axis of each figure are negative. This indicates that the direction of steel beam is to the negative y-direction. It is noted also that the curves (from Figures 2–4) are not smooth, and this is caused by the uncontrolled environmental conditions. The maximum deflection is taken with consideration of the sign of the value since it symbolizes the axis direction.

Figure 4.

Steel beam deflection-temperature graph.

The deflection test results seem to be in agreement with the general behavior of shape memory alloy material. Therefore, these results will be used as the benchmark for simulation results and will be compared with FEA results in the following section.

5. Finite element analysis on a 2D beam setup

This section reflects the application of the proposed NiTi SMA model in predicting the response of the steel beam subjected to mechanical loading from NiTi SMA wire. Prior to the performance of the finite element analysis, an experiment was performed so as to find the simulation parameters like transformation temperature, the beam effective length and breadth, SMA wire effective length, and so on. To perform finite element analysis, a four-noded beam bar structure was constructed as shown in Figure 5. Two elements were steel beam elements and one bar element was NiTi SMA wire.

Figure 5.

Steel beam bar-SMA setup.

6. SMA model implementation

This section deals with the computer implementation of the finite element analysis. The local and the global plane beam bar element systems are being demonstrated schematically and are shown in Figure 6. The system has got three elements where each element has six degrees of freedom (three degrees of freedom per node). The axial degrees of freedom are ignored due to the nature of loading, hence the degrees of freedom per node are reduced from three to two, and therefore each element remains with four degrees of freedom (i.e. vertical displacement and rotation). The element node displacements and forces are given by Eqs. (12)–(14) as follows:

Element 1:

$$\begin{array}{l}\underset{¯}{f^{\left(1\right)}}=\left[\begin{array}{l}{f}_{y1}\\ m_{\mathit{\theta R}1}\\ f_{y2}\\ m_{\mathit{\theta R}2}\end{array}\right]\\ \underset{¯}{u^{\left(1\right)}}=\left[\begin{array}{l}{u}_{y1}\\ {\theta }_{R1}\\ u_{y2}\\ {\theta }_{R2}\end{array}\right]\end{array}$$E12

Element 2:

$$\begin{array}{l}\underset{¯}{f^{\left(2\right)}}=\left[\begin{array}{l}{f}_{y2}\\ m_{\mathit{\theta R}2}\\ f_{y3}\\ m_{\mathit{\theta R}3}\end{array}\right]\\ \underset{¯}{u^{\left(2\right)}}=\left[\begin{array}{l}{u}_{y2}\\ {\theta }_{R2}\\ u_{y3}\\ {\theta }_{R3}\end{array}\right]\end{array}$$E13

Element 3:

$$\begin{array}{l}\underset{¯}{f^{\left(3\right)}}=\left[\begin{array}{l}{f}_{x2}\\ f_{y2}\\ f_{4x}\\ f_{4y}\end{array}\right]\\ \underset{¯}{u^{\left(3\right)}}=\left[\begin{array}{l}{u}_{x2}\\ u_{y2}\\ u_{4x}\\ u_{4x}\end{array}\right]\end{array}$$E14

Figure 6.

Plane bar elements in its local and global system.

For beam stiffness calculations, the only material properties required are the modulus of elasticity E, the second moment of area I and the length L of the beam. For the bar stiffness calculation, the only properties required are the cross-sectional area A, the modulus of elasticity E and the length L of the bar. These properties are taken to be constant throughout the bar and the beam. The solution procedure for our present problem is shown in Appendix A.

The schematic diagram for solution flow is shown in Figure 7. The solution diagram is a general solution flow to be used to solve one-sided temperature and phase coupling of any given 2-D geometry. In our case, we assume that phases are only affected by temperature and not by stress or pressure. The temperature affects the calculation of volume fraction through Eqs. (9) and (10). The linear averaged property (Young’s modulus, E) is then calculated together with other parameters mentioned under phase transformation section.

Figure 7.

Finite element program solution flow process.

The stiffness matrix of beam is formed using material properties of the beam and then that affects the formation of forces due to shape memory effect which is triggered by temperature variation. In our case, there are no external loads. This brings us to the calculation of the linear solution of Eq. (A8) from Appendix A. The process gets repeated until the maximum allowed temperature is reached.

Figure 7 summarizes the solution flow process of the developed simulation tool.

7. Simulation and experimental results

This section reflects the numerical results obtained from 1D finite element analysis. The material parameters used during FEA simulation are tabulated in Table 1. These parameters were derived from experiment using the procedure available in the study [11].

Property	Value	Units
EM	16.8	GPa
EA	31.8	GPa
AM	0.6	^oC^−1
AA	0.06	^oC^−1
${\xi }_M$	0.5	—
As	43	^oC
MF	18.3	^oC
${\epsilon }_L$	0.0485	mm/mm
ASMA	1.9635E−7	m^2

Table 1.

Measured material parameters for 0.5-mm diameter NiTi SMA wire.

Figure 8 shows the deflection-time graph where we wanted to see if the numerical results would match the experimental results. Therefore, based on the results, it can be seen that the simulation curve resembles the experimental curve. Both curves show the maximum deflection of about 4.9 mm in less than 10 s. The two curves follow the same trend and they both do not exceed the maximum deflection of 5 mm.

Figure 8.

Deflection-time curve.

The force-temperature curves are shown in Figure 9. This graph shows the force variation which is effected by the temperature variation. The negative values indicate that the SMA wire was contracting. The maximum deflection (in Figure 8) was produced by a force of approximately 35.0 N (in Figure 9). The numerical results are in agreement with the experimental results with minor negligible difference. It can be seen from the figure that the numerical prediction follows a similar trend with experiment.

Figure 9.

Force-temperature curves.

The maximum deflection for a steel beam as a function of temperature is shown in Figure 10. This figure shows graphical comparison between the numerical and the experimental results. The maximum deflection obtained experimentally agrees to that obtained numerically. The minor difference between the experimental and the numerical results is seen during the cooling path, and this is suggested to be caused by the uncontrolled environmental conditions. It can be concluded that the developed mathematical equations can be used to predict the behavior of smart simple structures.

Figure 10.

Experimental and numerical results for deflection-temperature curve.

8. Conclusion

A finite element analysis based on the proposed SMA model was performed successfully. The finite element analysis was performed on the 1D setup which was oriented on the 2D space. The C++ code was developed in order to perform the 2D numerical analysis. The experiment was performed in order to obtain the parameters to input in the developed code and also to validate the numerical results. The maximum deflection obtained numerically matches that which was measured experimentally. It was verified through the results that the developed SMA model has the ability to capture all the temperature range and not only the intended range.

The element stiffness matrix in global coordinates (X,Y) is explicitly given as follows:

Element 1:

$$\left[K^{\left(1\right)}\right]=\frac{\mathit{EI}}{L^3}\left[\begin{array}{llll}12& 6L& -12& 6L\\ 6L& 4L^2& -6L& 2L^2\\ -12& -6L& 12& -6L\\ 6L& 2L^2& -6L& 4L^2\end{array}\right]\begin{array}{l}{u}_1\\ {\theta }_1\\ u_2\\ {\theta }_2\end{array}$$EA1

Element 2:

$$\left[K^{\left(1\right)}\right]=\frac{\mathit{EI}}{L^3}\left[\begin{array}{llll}12& 6L& -12& 6L\\ 6L& 4L^2& -6L& 2L^2\\ -12& -6L& 12& -6L\\ 6L& 2L^2& -6L& 4L^2\end{array}\right]\begin{array}{l}{u}_2\\ {\theta }_2\\ u_3\\ {\theta }_3\end{array}$$EA2

Element 3:

$$\begin{array}{l}\left[K^{\left(3\right)}\right]=\left[\begin{array}{llll}{\mathit{kc}}^2& \mathit{ksc}& -{\mathit{kc}}^2& -\mathit{ksc}\\ \mathit{ksc}& {\mathit{ks}}^2& -\mathit{ksc}& -{\mathit{ks}}^2\\ -{\mathit{kc}}^2& -\mathit{ksc}& {\mathit{kc}}^2& \mathit{ksc}\\ -\mathit{ksc}& -{\mathit{ks}}^2& \mathit{ksc}& {\mathit{ks}}^2\end{array}\right]=k\left[\begin{array}{llll}{c}^2& \mathit{sc}& -c^2& -\mathit{sc}\\ \mathit{sc}& s^2& -\mathit{sc}& -s^2\\ -c^2& -\mathit{sc}& c^2& \mathit{sc}\\ -\mathit{sc}& -s^2& \mathit{sc}& s^2\end{array}\right]\begin{array}{l}{u}_{2x}\\ u_{2y}\\ u_{4x}\\ u_{4y}\end{array}\\ \theta ={270}^0;c=0;\mathrm{s}=-1;\mathrm{cs}=0\\ \left[K^{\left(3\right)}\right]=\frac{\mathit{AE}}{L}\left[\begin{array}{llll}0& 0& 0& 0\\ 0& 1& 0& -1\\ 0& 0& 0& 0\\ 0& -1& 0& 1\end{array}\right]\begin{array}{l}{u}_{2x}\\ u_{2y}\\ u_{4x}\\ u_{4y}\end{array}\end{array}$$EA3

The next step is the assemblage of the above stiffness matrix to get the global stiffness matrix.

$$\left[K_G\right]=\left[\begin{array}{llllllll}\frac{12\mathit{EI}}{L^3}& \frac{6\mathit{EI}}{L^2}& -\frac{12\mathit{EI}}{L^3}& \frac{6\mathit{EI}}{L^2}& 0& 0& 0& 0\\ \frac{6\mathit{EI}}{L^2}& \frac{4\mathit{EI}}{L}& -\frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& 0& 0& 0& 0\\ -\frac{12\mathit{EI}}{L^3}& -\frac{6\mathit{EI}}{L^2}& \frac{24\mathit{EI}}{L^3}& 0& -\frac{12\mathit{EI}}{L^3}& \frac{6\mathit{EI}}{L^2}& 0& 0\\ \frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& 0& \frac{8\mathit{EI}}{L}+\frac{\mathit{AE}}{L}& -\frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& 0& \frac{\mathit{AE}}{L}\\ 0& 0& -\frac{12\mathit{EI}}{L^3}& -\frac{6\mathit{EI}}{L^2}& \frac{12\mathit{EI}}{L^3}& -\frac{6\mathit{EI}}{L^2}& 0& 0\\ 0& 0& \frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& -\frac{6\mathit{EI}}{L^2}& \frac{4\mathit{EI}}{L}& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& -\frac{\mathit{AE}}{L}& 0& 0& 0& \frac{\mathit{AE}}{L}\end{array}\right]$$EA4

The structural equation would be

$$\left[\begin{array}{l}{F}_{1x}\\ M_1\\ F_{2x}\\ F_{2y}\\ F_{3x}\\ M_3\\ F_{4x}\\ F_{4y}\end{array}\right]=\left[\begin{array}{llllllll}\frac{12\mathit{EI}}{L^3}& \frac{6\mathit{EI}}{L^2}& -\frac{12\mathit{EI}}{L^3}& \frac{6\mathit{EI}}{L^2}& 0& 0& 0& 0\\ \frac{6\mathit{EI}}{L^2}& \frac{4\mathit{EI}}{L}& -\frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& 0& 0& 0& 0\\ -\frac{12\mathit{EI}}{L^3}& -\frac{6\mathit{EI}}{L^2}& \frac{24\mathit{EI}}{L^3}& 0& -\frac{12\mathit{EI}}{L^3}& \frac{6\mathit{EI}}{L^2}& 0& 0\\ \frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& 0& \frac{8\mathit{EI}}{L}+\frac{\mathit{AE}}{L}& -\frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& 0& \frac{\mathit{AE}}{L}\\ 0& 0& -\frac{12\mathit{EI}}{L^3}& -\frac{6\mathit{EI}}{L^2}& \frac{12\mathit{EI}}{L^3}& -\frac{6\mathit{EI}}{L^2}& 0& 0\\ 0& 0& \frac{6\mathit{EI}}{L^2}& \frac{2\mathit{EI}}{L}& -\frac{6\mathit{EI}}{L^2}& \frac{4\mathit{EI}}{L}& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& -\frac{\mathit{AE}}{L}& 0& 0& 0& \frac{\mathit{AE}}{L}\end{array}\right]\left[\begin{array}{l}{u}_1\\ {\theta }_1\\ u_2\\ v_2\\ u_3\\ {\theta }_3\\ u_4\\ v_4\end{array}\right]$$EA5

The boundary conditions for our system are as follows:

$$\begin{array}{l}{F}_{x1}=M_1=F_{x3}=M_3=F_{x2}=F_{x4}=F_{y4}=0\\ u_{x1}={\theta }_1=u_{x3}={\theta }_3=u_{x4}=u_{y4}=0\\ F_{y2}=F_{\mathit{SMA}}\end{array}$$EA6

The structural equation reduces to the following equation:

$$F_{\mathit{SMA}}=\left(\frac{8\mathit{EI}}{L_{\mathit{\text{beam}}}}+\frac{\mathit{AE}}{L_{\mathit{bar}}}\right)U_{2y}$$EA7

The deflection could be calculated as follows:

$$\frac{F_{\mathit{SMA}}}{\left(\frac{8\mathit{EI}}{L_{\mathit{\text{beam}}}}+\frac{\mathit{AE}}{L_{\mathit{bar}}}\right)}=U_{2y}$$EA8

This is the deflection that needs to be solved by the FEA. It should be noted that SMA force is the only variable in this equation, and also this SMA force is calculated inside the FEA code as follows:

$$F_{\mathit{SMA}}=\mathit{EA}{\epsilon }_L\xi$$EA9