Thermal Effects in a Mechanical Model for Pseudoelastic Behavior of NiTi Wires

A mechanical model for pseudoelastic behavior of NiTi wires is proposed with the aim to predict the behavior of Shape Memory Alloys(SMA) damping wire elements in model structures. We have considered at first a simple linearwise stress-strain relationship to describe the basic isothermal behavior of the SMA members. Then, this basic model is modified in order to include the effect of the strain rate. The model is based on detailed experimental characterization performed on a Ni rich NiTi superelastic wire which included the study of the localized character of the deformation and the local heat generation associated with the stress induced martensitic transformation occurring in these alloys. Heat conduction along the wire and heat interaction with the surroundings was also considered. In that way, the resulting local temperature field around the transformation front is assessed and its effect on the progression of the transformation is evaluated. It is shown how the simple mechanical model reproduces the global mechanical behavior, including the existence of a maximum in the damping capacity with the transformation rate.


Introduction
In certain temperature ranges, shape memory alloys (SMA) as Niquel-Titanium (NiTi) or Copper based (CuZnAl, CuAlNi, CuAlBe, etc) exhibit the pseudoelastic effect (SE).It consists on the elongation of the specimen upon the application of a critical load, accompanied with a total, or almost total, recovery of the strain upon unloading 1 .Maximun recoverable deformation typically ranges from 8 to 10% strain in polycrystalline and single-crystalline specimens, respectively.These recoverable strain regimes are possible because the material, before the yielding point is reached, suffers a martensitic transformation induced by the applied load.This transformation produces a shear type deformation in the grains of the material with a component in the axial direction.In several cases the transformation occurs at approximately constant stress, giving place to the characteristic upper plateu of the transformation during loading and the lower plateau during the unloading.The hysteresis appeared in the pseudoelastic cycles is due to friction associated to interface movement and plastic work associated to the grain accommodation.This implies an energy dissipation that motivates the study the SMA as an interesting alternative for dampers design for seismic or mechanical vibrations control 2 .In this sense, NiTi alloys reached more technological concern than the others SMA.
In comparison with other technologies, a steel device working by tensile-compressive cycles in the plastic range may result more economical than a device with NiTi as the kernel component.However, for the same level of deformation in uniaxial loading, the life of the steel devices in number of cycles is greatly diminished besides the NiTi base devices.Moreover, NiTi allows to design dampers with recentering capacity 3 , being, in this way able to work without an external force source on the damper system.
With the aim to evaluate the mechanical response of the material in the presence of any load event, is common to adopt pseudoelastic behaviors models.There are models of diverse complexity.In a materials science approach, thermodynamics and kinematics relations over elements at the micro or mesoscopic scales (i.e. grain size scales) are formulated, to describe shape memory and pseudoelasticity effects when macroscopic aggregates are analyzed.In general, constitutive laws are derived, using continuum mechanics formulations, where the magnitudes of strains are evaluated from crystallographic observations.Lagoudas et al. 4,5 present a review of modeling techniques in single crystals and polycrystals.This approach allows explaining from first principles the occurrence of the shape memory and SE effects, but the high computational cost of these models limits its use to static loading situations on SMA specimens using homogenizing criteria to treat the high number of martensite variants orientations observed.Moreover, the complexity in the SMA behavior makes that experimental results retrofit continuously the existent models incorporating new particularities observed for each SMA alloy.NiTi crystallographic texture effect coming from heavy cold working was analyzed by Thamburaja et al. 6 .On the other hand Sittner and Novak 7 have investigated the crystallographic orientation effects in the origin of lüders like behavior in NiTi tensile tests, while selfheating effects due the latent heat was studied by Anand et al. 8 .
In structural and control engineering, dynamical applied loads on pseudoelastic members are analyzed, adopting simples models that fit trajectories in the σ−ε space obtained empirically with linear segments 9,10 .It is well known that conditions of strain amplitude, strain rate, temperature of cycling and number of cycles affects the shape of the trajectories.That is why these dependences are incorporated with different criterion, keeping a simple and macroscopic approach.A similar approach constitutes the surface models 11 where a particular device with a superelastic element is rigorously tested in a wide range of temperatures including the finish martensite transformation temperatures.Several isothermal σ−ε curves are evaluated to make up a limit surface for the behavior of the device.Any (σ, ε, T) trajectory can be calculated once the surface was measured, using interpolating numerical tools, and a criterion for the internal trajectories.
In this article the effect of the strain rate on the pseudoelastic behavior for uniaxial cycling of NiTi wires is studied.For this, the thermal problem associated to the latent heat of the transformation is analyzed.This latent can be descomposed in two terms.One is the chemical latent heat and the other is the pseudoelastic work to produce the transformation 12,13 .A localized transformation pattern with moving martensite/austenite interfaces is assumed to occur.Taking a reference isothermal σ−ε cycle, when transformations or retransformations conditions are reached, the interfaces will move giving rise to latent heat exchange with surroundings.This will change the temperature profile along the wire and then, because of the Clasius-Clapeyron relation the stresses needed to proceed with the transformation or retransformation will change 12 .This effect will be increased if the strain rate is high enough to avoid the latent heat to be dissipated.It is considered that σ−ε trajectories for different strain rates could be obtained evaluating the change of temperature at the interfaces region due to latent heat of martensitic transformation exchange, and then considering the change in the stress with respect to an empirical cycle measured at isothermal conditions.
Simulations of strain controlled pseudoelastic cycles were performed for various strain rates, and a comparison with experimental pseudoelastic cycles is made, using the equivalent viscous damping ξ as a magnitude of damping capacity: (1)   where ∆W is the hysteresis area and W is the maximum strain work made during the pseudoelastic loading.

Proposed Model
When considering a pseudoelastic tension cycle on a NiTi wire with very low strain rate and the ambient temperature fixed at a value T, the induced transformation is slow, and the latent heat delivered by the transformation does not change significantly the local temperature of the wire, as it is maintained in thermal equilibrium with the surroundings.These isothermal cycles can be fit with straight lines, being defined by the stress σ and strain ε parameters of start and end of the transformations (Figure 1a).The slope dσ / dε of any generic point "G" at the different regions of the cycle are obtained using the coordinates of the pole P (Figure 1b) with the expression: (2) when no transformations are occurring and: (3) (4) if the material is stress transforming from austenite to martensite (3)  or from martensite to austenite (4), respectively.Therefore, in the present model, the local slope of the σε path is a function of σ, ε and the sign of .ε, and neither the stress nor the strain can be expressed as analytical functions.
A pseudoelastic cycle performed at the same strain rate but at a higher temperature T + ∆T, will show higher levels for both plateau corresponding to the transformation and retransformation stresses, in a quantity: where C is the value of dσ / dT from the Clausius-Clapeyron relation, being in practice constant for each of these alloys.It is assumed that the strains of start and end of the transformations do not change with the temperature, and, in consequence, the stress σ p as the others suffers the same change given by (5).Then an isothermal σ−ε diagram is defined for each temperature.
If during cycling, the temperature around the interface changes due to, for example, the latent heat, a new isothermal diagram must be used to calculate the slope dσ / dε.Depending on the control of the test, the variables of the wire are calculated every time step as follows: Starting with state variables: T, ε, σ a) Definition of isothermal diagram as Figure 2, using ( 5  In this work, strain controlled tests are used to perform simulations.

Temperature evaluation
To evaluate the local temperature of the wire an analysis of the heat generated and dissipated from some control volume must be carried out.This is shown in Figure 3a where the martensite-austenite interface are considered as localized heat sources as is explained at the end of this section.The volume is assumed to be small enough to have a uniform temperature T.
Then, the rate of heat accumulated .Q acum in the volume is given by the rate of heat generated .Q gen and the rate of heat exchanged with surroundings .
The change in temperature is given by the thermal energy rate inside the volume: (7)   where ρ is the density, C e is the specific heat.These properties are considered to be equal for both phases.
The generated heat rate depends on interfaces velocity and has two contributions.One is related with the work produced by the external stress along the transformation branches goes on (i.e. the work for interfaces movement).This heat generation rate is thus, related with the stretching velocity V and is stress dependent.It is positive during the transformation (V and σ positives) and negative during the retransformation (V negative, σ positive).The other contribution is the heat generated due the chemical latent heat of the transformation.It is independent of the applied stress and it is considered to be independent of temperature.Its rate of generation must be evaluated with the interface velocity, or for convenience its absolute value ν.Then the total heat generation is: (8)   Here, L is the chemical latent heat per volume unit.x f is an interface position, and is time dependent.The dirac integral over the volume points that the generations takes place only if there is any interface inside it.The sign of the chemical heat generation is positive for the transformation (exothermic) and negative for the retransformation (endothermic).
For the exchanged heat rate conduction with adjacent volumes and convection with surroundings are considered: (9)   Being k the thermal conductivity identical in both phases, and h the convection coefficient.
Rewritting (6), with ( 7), ( 8) and ( 9) Now, dividing by ρ.A.∆x .C e , and taking limits when ∆x → 0: This equation can be expressed as a differences equation with a suitable time and length discretization: (12)   This can be evaluated using an explicit finite difference method that would be stable with 14 (13) The position and velocity of interfaces can be determined after a transformation pattern during pseudoelastic cycling is defined.For this work is considered that martensite is nucleated simultaneously at both ends, where the clamp system generate a complex stress state.The interfaces go towards the center of the wire during loading, and they go back to the ends.This pattern of localized transformation with two or one interface fronts is a feature observed in polycrystalline samples of NiTi 15 .
For the calculation of position and velocity of each interface the actual martensite fraction Z is considered.It can be calculated with the Figure 1b giving: (14)   where E A and E M are the slopes dσ / dε of austenite and martensite respectively, evaluated directly from isothermal diagrams of Figure 1b: The two interfaces are placed symmetrically in such a way that a Z portion of the wire length be behind them.When solving (11)  the velocity ν is not evaluated directly but the displacement ν∆t is calculated using the difference between the interface positions of an interface after a time step.
For boundary conditions, an effective conductivity k eff was adjusted for agreement between numerical and experimental results:

Size of discretization
The usefullness of any numerical model is limited by its computational cost.The stability condition (13) relates the size of the space (∆X) and time (∆t) domain discretization.Therefore, the sensitivity analysis is restricted to one of these variables.The size of the discretization should be as large as possible.Figure 4 shows the effect of dividing the wire in different nodes numbers.For this work the simulations were made taking a 1000 nodes discretization, but stable results could be obtained with less nodes.

Wideness of the transforming zone
Even though the interface is modeled as punctual heat sources, in the discrete model, it is considered that the latent heat and the irreversible work are produced homogeneously in the discrete volume.Moreover, better results are obtained if the heath source is considered as a wider environment including several consecutives volume, giving a gradual heat generation instead of a sharp one.This means that the martensite fraction will present a gradual profile from 0 to 1 along the wire.The length in which the fraction grows from 0 to 1 is determined by the width of the heat source, as is shown in Figure 5 either for continue or space-discretised interpretation.Figure 5a  instantaneous scheme of the transforming zone."tr" is the value of the transforming fraction on each element.This zone moves along the wire so the total area below this curve is equal one.Figure 5b is a scheme of the martensite fraction "f" on each element at the same instant as in Figure 5a.
This picture is a better representation of the experimental marteniste advancing front.In a polycristal this front is not a planar one, instead its interface is composed by a complex arrangements of martensite plates giving the profile shown in Figure 5b.
Complete pseudoelastic cycles with several values for source width (referred to as ∆) are shown in Figure 6.Basically, the effect of considering wider transformation environments is to soften the transition corners increasing the roundness.Besides it, there is better fitting of the σ-ε trajectories towards the end of the martensitic transformation, i.e. in the zone where transforming zones interact.
On other hand, it is observed that the transition between the end of the retransformation plateau and the elastic section present a sharp stress decreases.This effect is hardly observed under usual experimental conditions and its magnitude can be diminished in the simulations by using wider heat sources.

Experimental Data
Uniaxial tension tests, with strain control were carried out in an universal mechanical test machine INSTRON 5567 using a 1 kN load cell, and a temperature control chamber.The material was a NiTi wire of 0.5 mm of diameter, 50.8% atomic Ni, from the memory metalle TM obtained by conventional conforming methods with a final heat treatment of straight annealing 16 consisting in an annealing under load at 773 K for 1-2 minutes.A sample of 75 mm of gage length was used, and, setting the chamber temperature at 35 ° C and the velocity of extension at 1 mm/minutes a number of 120 pseudoelastic cycles was carried.The extension needed for a complete transformation was of 5.8 mm.After these complete cycles, the hysteresis area diminished and the residual strain increased, until the material reached a stable situation 17 .Then the chamber temperature was set to 30 °C and an isothermal cycle was measured in order to obtain the reference parameters.Then, complete cycles at 5, 20, 35, 50 respectively mm/minutes were carried out, and finally, extension limits were changed to record

Numerical Results
From pseudoelastic complete cycles at deformation velocity of 0.1 mm/minutes the parameters σ sm = 326 MPa; ε sm = 0.0125; σ fm = 357 MPa; ε fm = 0.0544; σ sa = 296 MPa; ε sa = 0.0505; σ fa = 257 MPa; ε fa = 0.0083 were obtained.Thermal properties of NiTi are difficult to assess and the information is rather disperse in the technical literature.For the conductivity k = 20 W.mK -1 was taken, the convection coefficient was set at h = 100 W.m -2 K -1 corresponding to the case of forced convection due to the fan of the chamber, The other parameters were L = 156 MJ.m -3 , ρ = 6450 kg.m -3 and Ce = 837.1 J.K -1 •kg.
Taking ∆x as 0.001 times the gage length and assuming that the boundaries have a reduced conductivity of 0.1 k / ρC e , numerical simulations at different strain rates were performed.
In Figure 7 simulated and experimental complete cycles are shown, while in Figure 8 internal cycles are represented.In the case of internal cycles for 1 and 2 mm/minutes, simulations were performed considering only one interface propagating along the wire.This pattern seems to be more suitable for low strain rates even thought the 2 interface pattern results better for higher velocities as shown in Figure 9 in the equivalent viscous damping comparison.Also in Figure 9, the model results show a velocity for which the damping reaches a maximum of both, internal and complete cycles.
Higher velocities simulations do not have their experimental counterpart in this article, but, they are included in Figures 10 and 11 for completeness.It shows a lower limit tendency of the viscous damping capacity.

Conclusions
The model proposed in the present work reproduces the main features observed in experimental pseudoelastic cycles.It also captures the rounded edges of the hysteresis cycles observed in the experimental behavior and these features could be well explained by transient thermal effects.However, it is observed that the simulated cycles do not reproduce correctly the experimentally determined  stress-strain trajectories at the end of the retransformation.The reason for this discrepancy resides on the exact details of the heat boundary conditions that were overlooked in this simple model.The existence of a velocity showing a maximum damping capacity, reported by Piedbouef 18 and which are present in the experimental results is also found in our numerical simulations.Numerical results of cycles at high velocities show a slope that seems to reach an asymptotic lower limit.In test made on prototypes, Castellano 19 also points that the damping capacity does not change with the strain velocity at high rates.This is consistent with the results obtained here.As a conclusion it can be remarked that a numerical tool to analyze pseudoelastic cycling, taking into account the strain rate effect over the σ−ε trajectories, was developed.It gives results in good agreement with the experiments.

Figure 1 .Figure 2 .
Figure 1.Isothermal cycle definition.a) Linearwise representation of an isothermal cycle; and b) Ubication of Pole "P" for slope calculus of generic points "G".

Figure 3 .
Figure 3. Volume control representation.a) Energy balance; and b) Localized strain deformation.

Figure 5 .
Figure 5. a) Transformating portion "tr" of each element in a piece of the wire; and b) Martensite fraction "f" on the same piece that a).

Figure 6 .
Figure 6.Better results are obtained using width heat sources.

Figure 11 .
Figure 11.Lower limit tendency of ξ at high extension velocities.