Mechanics of Shape Memory Alloy Materials – Constitutive Modeling and Numerical Implications

1. Introduction

Shape Memory Alloys (SMAs) are a unique class of metal alloys which can be deformed severely and afterwards recover their original shape after a thermomechanical cycle (shape memory effect), or a stress cycle within some appropriate temperature regimes (pseudoelasticity, also sometimes called in the literature superelasticity, not to be confused with hyperelasticity). The mechanisms of this recovery are either a diffusionless transformation between the austenite phase (which is a highly ordered phase and is also called the parent phase) and the martensite phase (which is a less ordered one) or the reorientation (detwinning) of martensite variants. Detailed exposures to the physics of the subject may be found in Wayman (1964), Smallman and Bishop (2000) and Bhattacharya (2003). As is shown in these studies the thermomechanical response of SMAs is extremely complex, a fact that in conjunction with the continuously increasing use of SMAs in several innovating applications in many engineering fields results in a greater need for a better understanding of these materials. For the past decades several constitutive models have appeared within the literature (e.g., Raniecki et al., 1992; Abeyaratne and Knowles, 1993; Ivshin and Pence, 1994; Boyd and Lagoudas, 1996; Lubliner and Auricchio, 1996; Panoskaltsis et al., 2004), which within the context of a geometrical linear theory can capture several aspects of the experimentally observed response. Nevertheless, the physics of the problem (e.g., see Smallman and Bishop, 2000), together with some basic results of the crystallographic theory of martensitic phase transformations (e.g., Ball and James, 1987; James and Hane, 2000; Abeyaratne et al. 2001), suggest that a geometrically non – linear approach is more appropriate. Levitas and Preston, (2005) discuss the drawbacks of the infinitesimal models and they report that finite rotations of the crystal lattice can occur at small transformation strains (small strains and finite rotations) and can crucially affect the phase transformation conditions. Rather recently, several researchers have started to develop constitutive models for SMAs within the finite deformation regime. The approaches used for the description of the behavior of these materials are many and almost encompass all branches of mathematics, physics, material science and continuum mechanics. The existing models may be roughly classified in the following categories: (a) Constitutive models based on phase field or Landau – Ginzburg theory, (b) models based on irreversible thermodynamics and (c) models based on plastic flow theories.

The basic idea of the phase field theory is that out of all complexities of statistical mechanics one can reduce the behavior of a system undergoing a phase transformation to that of a few order parameters (i.e., parameters that give a measure of the transformation development), governed by a free energy function, which depends on stress (or deformation), temperature and those parameters. A characteristic example of modeling phase transformations by Landau – Ginsburg theory is provided by Levitas and Preston, 2005.

Also, in the realm of the so – called non equilibrium (or irreversible) thermodynamics several models have been proposed which are based on the use of a set of thermomechanical equations describing the kinetics of the martensitic transformations. The constitutive equations are developed in a non – linear manner on the basis of the laws of thermodynamics. Depending on whether they utilize the full microscopic deformation or the phenomenological one, the thermodynamical models may be classified further as microscopic (e.g., Levitas and Ozsoy, 2009) or macroscopic (e.g., Müller and Bruhns, 2006).

Another approach, which besides being thermodynamically consistent may also furnish a concrete micromechanical justification, is through the employment of plastic flow theories. Recall that the martensite transformation is a diffusionless one during which there is no interchange on the position of neighboring atoms but atom movements resulting in changes in the crystal structure (e.g., see Smallman and Bishop, 2000, pp. 278 – 279). Based on this observation the martensite formation has been explained by a shear mechanism or by a sequence of two shear mechanisms. The shear mechanism can take place either by twinning or by sliding, depending on the composition and on the thermodynamical conditions (Smallman and Bishop, 2000, p. 280). Although in the book of Smallman and Bishop mainly martensitic transformation in steel is described, the authors discuss efforts for the development of a general theory of the crystallography of martensitic transformations. The crystallographic mechanisms of martensite in nickel titanium (NiTi, also known as Nitinol) are similar, i.e., slip or twinning, as in the alloys described in the book of Smallman and Bishop. As a result it can be considered that the role played by the different transformation systems in the martensitic transformations may be suitably parallelized by the role played by the slip systems in crystal plasticity. Models based on this idea have been proposed among others by Diani and Parks (1998), Thamburaja and Anand (2000) and Anand and Gurtin (2003). It should be emphasized that these models are also computationally attractive because a lot of work has been put recently in the algorithms of crystal plasticity, both in their purely algorithmic as well as in their mathematical aspects, resulting in the development of robust algorithms well suited for finite element applications. Accordingly, complex constitutive representations may be considered, since their numerical implementation is no longer intractable, no matter how complex they may be.

An alternative approach, within the context of plastic flow theories, has been proposed by Lubliner and Auricchio (1996) and Panoskaltsis et al. (2004), who developed three – dimensional thermomechanical constitutive models based on generalized plasticity theory in the small deformation regime, and by Panoskaltsis et al. (2011a, 2011b) within finite strains and rotations.

Generalized plasticity is a general theory of rate – independent inelastic behavior which is physically motivated by loading – unloading irreversibility and it may be mathematically founded on set theory and topology (Lubliner 1974, 1984, 1987). Its particular structure provides the theory with the ability to address “non – standard” cases such as non – connected elastic domains.

The objective of this work is twofold: First, to extend the previous works of SMAs modeling based on generalized plasticity, providing a general geometrical framework. This general framework will in turn constitute a basis for the derivation of constitutive models for materials undergoing phase transformations and for arbitrary deformations. Second, as an application, to develop a finite strain model, which can simulate several patterns of the extremely complex response of SMAs under isothermal and non – isothermal loadings.

This chapter is organized as follows: In section 2, a general multi – surface formulation of non – isothermal generalized plasticity, capable of describing the multiple and interacting loading mechanisms which occur during phase transformations (see Panoskaltsis et al., 2011a, 2011b)) is presented within the context of tensor analysis in Euclidean spaces. The derivation of the thermomechanical state equations on the basis of the invariance properties of the local form of the balance of energy equation under some groups of transformations, is attained in section 3; this is a purely geometrical approach. In particular, the fundamental theorem of the covariant constitutive theory of non – linear elasticity (see Marsden and Hughes, 1994, pp. 202 – 203) is revisited and is used in place of the second law of thermodynamics, as a basic constitutive hypothesis for the subsequent derivation of the SMAs thermomechanical state equations. Rate constitutive equations are derived as well. Finally, as an application a specific model is derived within a fully thermomechanical framework in section 4. Computational aspects and numerical simulations are presented in section 5.

2. Generalized plasticity for phase transformations

2.1. Formulation of the governing equations in the reference configuration

Generalized plasticity is a local internal variable theory of rate – independent behavior, which is based primarily on loading – unloading irreversibility. As in all internal – variable type theories, it is assumed that the local thermomechanical state in a body is determined uniquely by the couple (G, Q) where G – belonging to a space $$G$$ – stands for the vector of the controllable state variables and Q – belonging to a space $$Q$$ – stands for the vector of the internal variables, which are related to phase transformations. Following the ideas presented in the review paper of Naghdi (1990) we follow a material (referential) approach within a strain – space formulation. Accordingly, G may be identified by (E, T) where E is the referential (Green – St. Venant) strain tensor and T is the (absolute) temperature. Depending on the nature of the (material) internal variable vector Q, the theory may, in principle, be formulated equivalently with respect to the macro –, meso –, or micro – scale structure of the material.

The central concept of generalized plasticity is that of the elastic range, which is defined at any material state as the region in the strain – temperature space comprising the strains which can be attained elastically (i.e., with no change in the internal variables) from the current strain – temperature point. It is assumed that the elastic range is a regular set in the sense that it is the closure of an open set. The boundary of this set is defined as a loading surface at Q, (see Eisenberg and Phillips, 1971; Lubliner, 1987). In turn, a material state may be defined as elastic if it is an interior point of its elastic range and inelastic if it is a boundary point of its elastic range; in the latter case the material state lies on a loading surface. It should be added that the notion of process is introduced implicitly here. By assuming that the loading surface is smooth at the current strain - temperature point and by invoking some basic axioms and results from set theory and topology, Lubliner (1987) showed that the rate equations for the evolution of the internal variable vector may be written in the form

$$\dot{Q}=HL(G,Q)⟨N\text{:}\dot{G}⟩,$$E1

where < > stands for the Macauley bracket which is defined as:

$$⟨x⟩=\{\begin{array}{cc}x& \text{if }x>0\\ 0& \text{if }x\le 0,\end{array}$$E2

and $$H$$ stands for a scalar function of the state variables. Accordingly, the value of $$H$$ must be positive at any inelastic state and zero at any elastic one. Finally, L stands for a non - vanishing (tensorial) function of the state variables associated with the properties of the phase transformation, N is the outward normal to the loading surface at the current state, while the colon between two tensors denotes their double contraction operation. Furthermore, the set of the material states defined as $$H=H(G,Q)=0,$$ which comprises all the elastic states is called the elastic domain and its projection on the set defined by Q = const. is defined at the elastic domain at Q. In general, the elastic domain at Q is a subset of the elastic range (Lubliner, 1987). The particular case in which the two sets coincide corresponds to classical plasticity and the boundary of the elastic domain, that is the initial loading surface, constitutes the yield surface (see Eisenberg and Phillips, 1971; Lubliner, 1987; Panoskaltsis et al., 2008a, 2008b, 2011c).

It is emphasized that Eq. (1) has been derived under the assumption of a smooth loading surface at the current strain – temperature point, which implies that only one loading mechanism can be considered. On the other hand, phase transformations include multiple and sometimes interacting loading mechanisms, which may result in the appearance of a vertex or a corner at the current strain – temperature point. This fact calls for an appropriate modification of the rate Equation (1).

In order to accomplish this goal we assume that the loading surfaces are defined in the state space by a number – say n – of smooth surfaces, which are defined by expressions of the form

$${\Phi }_{\text{i}}(G,Q)=0,\text{ i=1, 2,}\mathrm{...}\text{, n}\text{.}$$E3

These surfaces can be either disjoint, or intersect in a possibly non – smooth fashion. Each of these surfaces is associated with a particular transformation mechanism which may be active at the current strain – temperature point. Then, by assuming that each equation $${\Phi }_{\text{i}}(G,Q)=0$$ defines an independent (non – redundant) active surface at the current stress temperature point, and in view of Eq. (1), we can state the rate equations for the evolution of the internal variables in the following general form

$$\dot{Q}={\displaystyle \sum _{\text{i=1}}^{\text{n}}{H}_{\text{i}}{L}_{\text{i}}(G,Q)⟨N_{\text{i}}\text{:}\dot{G}⟩},$$E4

where $$H_{\text{i}},L_{\text{i}}$$ and $$N_{\text{i}}$$ are functions of the state variables defined as in Eq. (1) and each set of them – defined by the index i – refers to the specific transformation associated with the part of the loading surface defined by $${\Phi }_{\text{i}}(G,Q)=0.$$ From Eq. (3) one can deduce directly the loading – unloading criteria for the proposed formulation as follows: Let us denote by $${\text{n}}_{\text{adm}}\le \text{n}$$ the number of loading surfaces that may be active at an inelastic state i.e.$$H_{\text{i}}>0$$, and let us denote by $${\text{J}}_{\text{adm}}$$the set of $${\text{n}}_{\text{adm}}$$indices associated with those surfaces, i.e.

$${\text{J}}_{\text{adm}}=\text{{}\alpha \in \text{{1,2,}\mathrm{...}\text{,n}/}{H}_{\alpha }>0\text{}}\text{.}$$E5

Then Eq. (3) implies the following loading – unloading conditions:

$$\begin{array}{l}{\text{If J}}_{\text{adm}}=\varnothing ,\text{ then }\dot{Q}\text{=}0\text{.}\\ {\text{If J}}_{\text{adm}}\ne \varnothing ,\text{ then:}\\ \text{i If }{N}_{\alpha }:\dot{G}\le 0\text{ for all }\alpha \in {\text{J}}_{\text{adm}}\text{ then }\dot{Q}\text{=}0\text{,}\\ \text{ii If }{N}_{\alpha }:\dot{G}>0\text{ for at least one }\alpha \in {\text{J}}_{\text{adm}}\text{ then }\dot{Q}\ne 0\text{.}\end{array}$$E6

Hence, if we denote further by $${\text{n}}_{\text{act}}\le {\text{n}}_{\text{adm}}$$the number of parts for which (ii) holds, and we set:

$${\text{J}}_{\text{act}}=\text{{}\alpha \in {\text{J}}_{\text{adm}}\text{/}{N}_{\alpha }:\dot{G}>0\text{},}$$E7

the loading criteria in terms of the sets $${\text{J}}_{\text{adm}}$$and $${\text{J}}_{\text{act}}$$ may be stated as:

$$\{\begin{array}{l}{\text{If J}}_{\text{adm}}=\varnothing :\text{ elastic state}\text{.}\\ {\text{If J}}_{\text{adm}}\ne \varnothing {\text{ and J}}_{\text{act}}=\varnothing :\\ \text{ i}\text{. If }{N}_{\alpha }:\dot{G}<0\text{ for all }\alpha \in {\text{J}}_{\text{adm}}\text{: elastic unloading,}\\ \text{ ii}\text{. If }{N}_{\alpha }:\dot{G}=0\text{ for at least one }\alpha \in {\text{J}}_{\text{adm}}\text{: neutral loading,}\\ {\text{If J}}_{\text{adm}}\ne \varnothing {\text{ and J}}_{\text{act}}\ne \varnothing :\text{ inelastic loading}\text{. }\end{array}$$E8

2.2. Equivalent spatial formulation

The equivalent assessment of the governing equations in the spatial configuration can be done on the basis of a push – forward operation (e.g., see Marsden and Hughes, 1994, pp. 67 – 68; Stumpf and Hoppe, 1997) to the basic equations. For instance, by performing a push – forward operation onto Eq. (3) the latter is written in the form

$${\text{L}}_{V}q={\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(g,q,F)⟨r_{\text{i}}⟩},$$E9

where F stands for the deformation gradient and g stands for the vector of the controllable variables in the spatial configuration and is composed by the Almansi strain tensor e – defined as the push – forward of the Green – St. Venant strain tensor – and the (scalar invariant) temperature T. Moreover in Eq. (5), q stands for the push forward of the internal variable vector, and $${\text{L}}_V(\cdot )$$ stands for the Lie derivative (e.g., see Marsden and Hughes, 1994, pp. 93 – 104; Schutz, 1999, pp. 73-79; Stumpf and Hoppe, 1997), defined as the convected derivative relative to the spatial configuration. Finally, $$h_{\text{i}}$$ stands for the expression of the scalar invariant functions $$H_{\text{i}}$$ in terms of the spatial variables (e, T, q) and the deformation gradient F, $$l_{\text{i}}$$ stands for the push – forward of the tensorial functions $$L_{\text{i}}$$ and $$r_{\text{i}}$$ denotes the (scalar invariant) loading rates which are written in the form

$$r_{\text{i}}=\frac{\partial {\phi }_{\text{i}}}{\partial e}{\text{:L}}_{V}e+\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>,}$$E10

where $${\phi }_{\text{i}}$$ is the expression for the loading surface associated with the index i, in terms of the spatial variables. The (spatial) loading – unloading criteria flow naturally from Eq. (5) as:

$$\{\begin{array}{l}{\text{If j}}_{\text{adm}}=\varnothing :\text{ elastic state}\text{.}\\ {\text{If j}}_{\text{adm}}\ne \varnothing {\text{ and j}}_{\text{act}}=\varnothing :\\ \text{ i}{\text{. If r}}_{\alpha }<0\text{ for all }\alpha \in {\text{j}}_{\text{adm}}\text{: elastic unloading,}\\ \text{ ii}{\text{. If r}}_{\alpha }=0\text{ for at least one }\alpha \in {\text{j}}_{\text{adm}}\text{: neutral loading,}\\ {\text{If j}}_{\text{adm}}\ne \varnothing {\text{ and j}}_{\text{act}}\ne \varnothing :\text{ inelastic loading, }\end{array}$$E11

where the sets $${\text{j}}_{\text{adm}}$$ and $${\text{j}}_{\text{act}}$$ are now defined in terms of the spatial variables as,$${\text{j}}_{\text{adm}}=\text{{}\alpha \in \text{{1,2,}\mathrm{...}\text{,}\text{n}/}{h}_{\alpha }>0\text{}}$$and$${\text{j}}_{\text{act}}=\text{{}\alpha \in {\text{J}}_{\text{adm}}{\text{/r}}_{\alpha }>0\text{}}\text{.}$$

2.3. Description of rate effects

Rather recent experimental results (see, Nemat – Nasser et al., 2005a, 2005b) on a NiTi shape memory alloy, show that some of the phase transformations depend on the rate of loading. Such a behavior can be accommodated by the (geometrical) framework developed here, by noting that generalized plasticity can be combined consistently with a rate – dependent (viscoplastic) theory. In this case the rate equations for the internal variables may be written in the form

$$\dot{Q}={\displaystyle \sum _{\text{i=1}}^{\text{n}}{H}_{\text{i}}[L_{\text{i}}(G,Q)⟨N_{\text{i}}\text{:}\dot{G}⟩}+M_{\text{i}}(G,Q)],$$E12

where the $$M_{\text{i}}\text{'s}$$ stand for additional functions of the state variables enforcing the rate – dependent properties of the transformation defined by the part of the loading surface associated with the index i. The crucial advantage of this approach lies on the compatibility of the two theories, in the sense that neither viscoplasticity, nor generalized plasticity employs the concept of the yield surface as its basic ingredient.

2.4. Transformation induced plasticity

From a further study of the experimental results of Nemat – Nasser et al. (2005a, 2005b) (see also Delville et al., 2011) it is observed that after a stress cycle within the appropriate limits for pseudo-elastic behavior permanent deformations appear, a fact which implies that a yielding behavior appears within the martensitic transformations.

Such a response can be described within our framework by introducing m additional (plastic) loading surfaces, which control the yielding characteristics of the material. These are assumed to be given by expressions of the form

$$G_{\text{i}}(G,Q,P)=0,\text{ i=1, 2,}\mathrm{...}\text{, m,}$$E13

where P is an additional internal variable vector, which stands for the description of plastic phenomena within the material. In turn, the rate equations for the evolution of the plastic variables within the generalized plasticity context – which includes classical plasticity as a special case (see Lubliner, 1987; Panoskaltsis et al., 2008a, 2011c) – may be stated as

$$\dot{P}={\displaystyle \sum _{\text{i=1}}^{\text{m}}{K}_{\text{i}}{T}_{\text{i}}^{}(G,Q,P)⟨R_{\text{i}}\text{:}\dot{G}⟩},$$E14

where the functions $$K_{\text{i}},T_{\text{i}}^{}$$and $$R_{\text{i}}$$ have an identical meaning with the functions $$H_{\text{i}},L_{\text{i}}$$ and $$N_{\text{i}}$$ which appear in Eq. (3). The constitutive modeling of plasticity phenomena within the martensitic transformations is nowadays a very active area of research. Recent contributions include the phenomenological models by Hallberg et al. (2007, 2010) and Christ and Reese (2009).

A further observation of Eqs. (9) and (10) and their comparison with the basic Eqs. (2) and (3) reveal that both sets of equations show exactly the same qualitative characteristics. Accordingly, it is concluded that from a geometrical standpoint the phase transformation loading surfaces are indistinguishable from the plastic loading surfaces, which means that the internal variable vector P may be absorbed in Q so that the basic equations can simulate both phase transformation and plasticity phenomena in a unified format. This implies that plastic yielding can be understood within the proposed framework as a phase transition.

3. The invariant energy balance equation and the thermomechanical state equations

The concept of invariance plays a fundamental role in several branches of mechanics and physics. In particular, within the context of continuum mechanics the invariance properties of the balance of energy equation, under some groups of transformations, may be systematically used in order to derive the conservation laws, the balance laws and/or to determine some restrictions imposed on the equations describing the material constitutive response (e.g., Ericksen, 1961; Green and Rivlin, 1964; Marsden and Hughes, 1994, pp. 163 – 167, 200 – 203; Yavari et al., 2006; Panoskaltsis et al., 2011c). For instance, Marsden and Hughes (1994, pp. 202 – 203) by studying the invariant properties of the local form of the material balance of energy equation, under the action of arbitrary spatial diffeomorphisms, determined the thermomechanical state equations for a non – linear elastic material. The basic objective of this section is to revisit the approach given in Marsden and Hughes (1994, pp. 202 – 203), within the context of the Euclidean space used herein and to show how this can used as a basic constitutive hypothesis in place of the second law of thermodynamics for the derivation of the constitutive response of the SMA material in question.

3.1. Revisiting Marsden and Hughes’ theorem

Unlike the original approach of Marsden and Hughes where manifold spaces are used and the invariance of the local form of the material balance of energy equation is examined under the action of arbitrary spatial diffeomorphisms, which include also a temperature rescaling, we examine the invariance properties of the local form of the spatial balance of energy equation under the action of the same kind of transformations, within the context of a Euclidean space. Within this framework the basic axioms of Marsden and Hughes (1994, pp. 202 – 203) may be modified as follows:

Axiom 1 (Local energy balance or first law of thermodynamics): For a spatial point with coordinates $x_1,x_2,x_3,$ and a given elastic thermomechanical process (e, T) the balance of energy holds

$$\rho \dot{e}+divh=\sigma :{\text{L}}_{V}e+\rho r,$$E15

where $$e$$ is the energy density, $$\rho$$ is the spatial mass density, $$\sigma$$ is the Cauchy stress tensor, h is the heat flux vector per unit of surface of the spatial configuration, $$r$$ is the heat supply per unit mass and a superimposed dot indicates material time derivative. By introducing the Helmholtz free energy function $$\psi ,$$ obtained by the following Legendre transformation

$$\psi =e-\eta \text{T,}$$E16

where $$\eta$$ is the specific entropy, the local form of the energy balance can be written in the form

$$\rho (\dot{\psi }+\dot{\eta }\text{T}+\eta \text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>)}+divh=\sigma :{\text{L}}_{V}e+\rho r.$$E17

Axiom 2 (Invariance of Helmholtz free energy): We denote by $$S$$ the ambient space, by $$\varphi$$ the deformation mapping, by $$g$$ the space of the control variables and by $$f$$ the set of the $${\text{C}}^{\text{S}}$$ scalar fields all expressed in the spatial configuration. Next we assume the existence of a map $$\widehat{\widehat{\psi }}\text{ : }(S,g,R^{+})\to f$$such that for any diffeomorphism which includes also a temperature rescaling, that is: $(\xi ,\delta )\text{: }(S,R^{+})\to (S,R^{+}),$ the following expression holds:

$$\widehat{\psi }(\varphi \text{,}e,\text{T)=}\widehat{\widehat{\psi }}(\xi \circ \varphi \text{,}{\xi }_{\ast }e,\delta \text{T,}\delta \text{),}$$E18

where $${\text{(}\cdot \text{)}}_{\ast }$$ denotes the push – forward operation.

Axiom 3 (Invariance of the energy balance): For curves $${\xi }_t\text{:}S\to S$$ and $${\delta }_t(x)\in R^{+},$$$${\varphi }_t{}^{\prime }={\xi }_t\circ {\varphi }_t$$, $${\text{T}}_t{}^{\prime }\text{=}{\delta }_t\text{T}{}_t$$ and by assuming that $$\rho ,\psi$$ and $$\eta$$ are transformed as scalars, the heat flux vector is transformed as $$h_t{}^{\prime }={\delta }_t{\xi }_t{}_{*}{h}_t$$ and the “apparent heat supply” due to entropy production, $$r_t{}^{\prime }-{\text{T}}_{\text{t}}{}^{\prime }{\dot{\eta }}_t{}^{\prime },$$ is transformed as $${r^{\prime }}_t-{\text{<msup>T <mo>′</mo></msup>}}_t{\dot{\eta }}_t{}^{\prime }={\delta }_t(r_t-{\text{T}}_t{\dot{\eta }}_t)$$ the balance of energy holds, that is

$${\rho }^{\prime }(\dot{{\psi }^{\prime }}+{\dot{\eta }}^{\prime }\text{<msup>T <mo>′</mo></msup>}+{\eta }^{\prime }\text{<mover accent="true"><msup>T <mo>′</mo></msup><mo>˙</mo></mover>)}+div{h}^{\prime }={\sigma }^{\prime }:{\text{L}}_{V^{\prime }}{e}^{\prime }+{\rho }^{\prime }{r}^{\prime }.$$E19

Then the basic theorem of Marsden and Hughes (Theorem 3.6 p. 203), takes in our case the form:

Theorem 1: Under axioms 1, 2 and 3 the thermomechanical state equations for the Cauchy stress tensor and the entropy density are given as:

$$\begin{array}{l}\sigma =\rho \frac{\partial \psi }{\partial e},\eta =-\frac{\partial \psi }{\partial \text{T}}.\\ \end{array}$$E20

Proof: The proof relies crucially, as in the case examined by Marsden and Hughes, on the evaluation of Eq. (13) at time $$t_0,$$ when $${\xi |}_{t=t_o}=1$$(identity),$$w={\frac{\partial \xi }{\partial t}|}_{t=t_0}\text{ ,}$$ and $${\delta |}_{t=t_0}=1,$$$$\text{u}={\frac{\partial \delta }{\partial t}|}_{t=t_0}$$where u is the velocity of $$\delta$$ at $$t_0\text{.}$$ Then, since at $$t=t_0|\to v^{\prime }\text{=}w\text{+}v\text{, <msup><mover accent="true"><mi>T</mi><mo>˙</mo></mover><mo>′</mo></msup>=uT+<mover accent="true"><mi>T</mi><mo>˙</mo></mover>,}$$ the time derivative of the Helmholtz free energy yields:

$${{\dot{\psi }}^{\prime }|}_{t=t_0}=\dot{\psi }+\frac{\partial \psi }{\partial e}:L_{w}e+\frac{\partial \psi }{\partial \text{T}}(\text{uT)},$$E21

in which $$L_W(\cdot )$$ stands for the autonomous Lie derivative (e.g., see Marsden and Hughes, 1994, pp. 96 – 98; Yavari et al., 2006). Furthermore it holds that

$${{\sigma }^{\prime }:{\text{L}}_{v^{\prime }}{e}^{\prime }|}_{t=t_o}=\sigma :{\text{L}}_{v}e+\sigma :L_{w}e,$$E22

since (see Marsden and Hughes, 1994, p. 98)

$${\text{L}}_{V+w}e=\dot{e}+L_{V+w}e=\dot{e}+L_{V}e+L_{w}e={\text{L}}_{V}e+L_{w}e.$$E23

In light of Eqs. (15) and (16) and the transformation formulae for the heat flux vector and the “apparent heat supply” due to entropy production, Eq. (13) at $$t=t_0$$ can be written as:

$$\begin{array}{l}\rho \dot{\psi }+\rho \frac{\partial \psi }{\partial e}:L_{w}e+\rho \frac{\partial \psi }{\partial \text{T}}(\text{uT)+}\\ \text{+}\rho \eta \text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>+}\rho \eta (\text{uT)+div}h\text{=}\sigma :({\text{L}}_{V}e+L_{w}e)+\rho (r-\text{T}\dot{\eta }\text{)}\text{.}\end{array}$$E24

Subtracting Eq. (12) from Eq. (17) gives the identity

$$\rho \frac{\partial \psi }{\partial e}:L_{w}e+\rho \frac{\partial \psi }{\partial \text{T}}(\text{uT)+}\rho \eta (\text{uT)}-\sigma :L_{w}e=0,$$E25

or

$$(\rho \frac{\partial \psi }{\partial e}-\sigma ):L_{w}e+(\rho \frac{\partial \psi }{\partial \text{T}}\text{+}\rho \eta )(\text{uT)}=0,$$E26

from which and by noting that $$L_{w}e$$ and u can be arbitrarily specified, the thermomechanical state Equations (14) follow. By performing a pull – back operation to Eqs. (14) the following relations are derived

$$S={\rho }_{\text{ref}}\frac{\partial \Psi }{\partial E},\eta =-\frac{\partial \Psi }{\partial \text{T}},$$E27

where S stands for the second Piola – Kirchhoff stress tensor, $${\rho }_{\text{ref}}$$ for the material mass density and $$\Psi$$ for the expression of the Helmholtz free energy in the material configuration. It is concluded that Eqs. (20) are identical to the thermomechanical state equations of Marsden and Hughes (1994, p. 203). Thus, we can state the following proposition:

Proposition 1: The invariance of the local form of the balance of energy equation under the superposition of arbitrary spatial diffeomorphisms, which also include a temperature rescaling, gives identical results with respect to the thermomechanical state equations, irrespectively of whether the energy balance equation is considered in its material or its spatial form.

3.2. Thermomechanical state equations for a SMA material

Building on the previous developments we will derive the thermomechanical state equations for a shape memory material with internal variables, which obeys the rate Equations (5). The development relies crucially on establishing a new set of axioms which will incorporate the presence of the internal variables and their evolution in the course of the phase transformations. We proceed as follows:

Axiom 1: Since the internal variables are not involved explicitly in the balance laws, for a process which is either elastic or inelastic ($${\text{j}}_{\text{adm}}\ne \varnothing$$and$${\text{j}}_{\text{act}}\ne \varnothing$$) axiom 1 remains unaltered, that is Eq. (12) holds.

Axiom 2 is modified as follows:

Axiom 2: In addition to the adopted notation, we denote by $$q$$ the space of the internal variables in the spatial configuration and we assume the existence of a map $$\widehat{\widehat{\psi }}\text{ : }(S,g,q,R^{+})\to f$$such that for any diffeomorphism which includes also a temperature rescaling, that is: $(\xi ,\delta )\text{: }(S,R^{+})\to (S,R^{+}),$ the following expression holds:

$$\widehat{\psi }(\varphi \text{,}e,\text{T,}q\text{)=}\widehat{\widehat{\psi }}(\xi \circ \varphi \text{,}\delta \text{T,}{\xi }_{\ast }e,q({\xi }_{\ast }e,\delta \text{T),}\delta \text{)}\text{.}$$E28

Axiom 3: In addition to the energy invariance axiom it is assumed that under the application of the diffeomorphism $(\xi ,\delta )\text{: }(S,R^{+})\to (S,R^{+}),$ loading surfaces are transformed as scalars, that is $${\phi }_i(\varphi \text{,}e,\text{T,}q\text{)=}{{\phi }^{\prime }}_i(\xi \circ \varphi \text{,}\delta \text{T,}{\xi }_{\ast }e,q({\xi }_{\ast }e,\delta \text{T)},\delta \text{)}\text{.}$$

In this case, as in the previous one, the derivation procedure is the following:

We evaluate Eq. (13) at time $$t_0,$$ when $${\xi |}_{t=t_o}=1$$(identity) and $${\delta |}_{t=t_0}=1$$ with w and u being the velocities of$$\xi$$and $$\delta$$ at $$t_0.$$ Then

$${{\dot{\psi }}^{\prime }|}_{t=t_0}=\dot{\psi }+\frac{\partial \psi }{\partial e}:L_{w}e+\frac{\partial \psi }{\partial \text{T}}(\text{uT)}+\frac{\partial \psi }{\partial q}:L_{w}q.$$E29

The critical step is the evaluation of the loading rates at $$t=t_0,$$which yields

$$\begin{array}{l}{{r^{\prime }}_{\text{i}}|}_{t=t_0}=(\frac{\partial {{\phi }^{\prime }}_{\text{i}}}{\partial e^{\prime }}:{\text{L}}_V{e}^{\prime }+\frac{\partial {{\phi }^{\prime }}_{\text{i}}}{\partial \text{<msup>T <mo>′</mo></msup>}}\text{<msup><mover accent="true"><mi>T</mi><mo>˙</mo></mover><mo>′</mo></msup>}{)|}_{t=t_0}\text{=}\\ \text{ =}\frac{\partial {\phi }_{\text{i}}}{\partial e}:{\text{L}}_{V}e+\frac{\partial {\phi }_{\text{i}}}{\partial e}:L_{w}e+\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>+}\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}(\text{uT)}\text{.}\end{array}$$E30

Accordingly, the rate equation for the internal variables evaluated at time$$t=t_0,$$yields

$$\begin{array}{l}{{\text{L}}_{V^{\prime }}{q}^{\prime }|}_{t=t_0}={{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{}^{\prime }{l}_{\text{i}}{}^{\prime }(e^{\prime },\text{<msup>T <mo>′</mo></msup>},q^{\prime },F^{\prime })r_{\text{i}}{}^{\prime }}|}_{t=t_0}=\\ ={\displaystyle \sum _{\text{i=1}}^{\text{n}}[h_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)(\frac{\partial {\phi }_{\text{i}}}{\partial e}:{\text{L}}_{V}e+\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>)]}}+{\displaystyle \sum _{\text{i=1}}^{\text{n}}[h_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)[\frac{\partial {\phi }_{\text{i}}}{\partial e}:L_{w}e+\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}(\text{uT)],}}\end{array}$$E31

which in view of the rate Equations (5), reads

$$L_{W}q={\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)[\frac{\partial {\phi }_{\text{i}}}{\partial e}:L_{w}e+\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}(\text{uT)]}}.$$E32

Upon substitution of Eq. (24), Eq. (21) takes the form

$${{\dot{\psi }}^{\prime }|}_{t=t_0}=\dot{\psi }+\frac{\partial \psi }{\partial e}:L_{w}e+\frac{\partial \psi }{\partial \text{T}}(\text{uT)}+\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)[\frac{\partial {\phi }_{\text{i}}}{\partial e}:L_{w}e+\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}(\text{uT)]}}.$$E33

Now, by working in a similar manner as in the previous (i.e., the elastic) case, in light of Eq. (25) the basic Eq. (13) evaluated at $$t=t_0,$$yields

$$\begin{array}{l}\rho \dot{\psi }+\rho \frac{\partial \psi }{\partial e}:L_{w}e+\rho \frac{\partial \psi }{\partial \text{T}}(\text{uT)+}\rho \frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)[\frac{\partial {\phi }_{\text{i}}}{\partial e}:L_{w}e+\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}(\text{uT)]}}\text{+}\\ \text{+}\rho \eta \text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>+}\rho \eta (\text{uT)+div}h\text{=}\sigma :({\text{L}}_{V}e+L_{w}e)+\rho (r-\text{T}\dot{\eta }\text{);}\end{array}$$E34

from which by subtracting the balance of energy Eq. (12) we can derive the identity

$$\begin{array}{l}\rho \frac{\partial \psi }{\partial e}:L_{w}e+\rho \frac{\partial \psi }{\partial \text{T}}(\text{uT)+}\rho [\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)\frac{\partial {\phi }_{\text{i}}}{\partial e}]:L_{w}e+}\\ +[\rho \frac{\partial \psi }{\partial q}\text{:}{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}](\text{uT)}}\text{+}\rho \eta (\text{uT)-}\sigma :L_{w}e=0,\end{array}$$E35

or equivalently

$$\begin{array}{l}\{\rho [\frac{\partial \psi }{\partial e}+\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)\frac{\partial {\phi }_{\text{i}}}{\partial e}\text{]-}}\sigma \}\text{:}{L}_{w}e+\\ +\rho [\frac{\partial \psi }{\partial \text{T}}+\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}}+\eta ](\text{uT)=}0,\end{array}$$E36

from which and by noting that$$L_{w}e$$ and u can be specified arbitrarily, we arrive at the expressions:

$$\begin{array}{l}\sigma =\rho [\frac{\partial \psi }{\partial e}+\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)\frac{\partial {\phi }_{\text{i}}}{\partial e}\text{]}}\text{,}\\ \eta =-[\frac{\partial \psi }{\partial \text{T}}+\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}}].\end{array}$$E37

Therefore, unlike the classical elastic case, for the SMA material considered, the invariance of the local form of the energy balance under superposed spatial diffeomorphisms does not yield the standard thermomechanical state equations unless a further assumption is made, namely that an unloading process from an inelastic state (i.e., a process with$${\text{j}}_{\text{adm}}\ne \varnothing$$and$${\text{j}}_{\text{act}}=\varnothing$$) is quasi-reversible, which means that in such a process both the mechanical and the thermal dissipations, defined as

$$d_{mech}=-\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)(\frac{\partial {\phi }_{\text{i}}}{\partial e}{\text{:L}}_{V}e)},d_{therm}=-\frac{\partial \psi }{\partial q}:{\displaystyle \sum _{\text{i=1}}^{\text{n}}{h}_{\text{i}}{l}_{\text{i}}(e,\text{T},q,F)\frac{\partial {\phi }_{\text{i}}}{\partial \text{T}}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>},$$E38

vanish. If this is the case, the classical thermomechanical state equations (Eqs. (14)) can be derived, as in the classical elastic case, directly from Eqs. (29). Thus, we can state the following theorem:

Theorem 2: For the rate – independent SMA material with internal variables whose evolution in the course of martensitic transformations is described by the rate equations (5), (or equivalently by Eqs. (3)), the invariance of the spatial local balance of energy equation under superimposed diffeomorphisms, which also include a temperature rescaling, does not yield the standard thermomechanical state equations, unless further assumptions are made.

It is interesting to note that in the classical theory of thermodynamics with internal variables Lubliner (1974, 1987) has arrived at a similar result by working entirely in the reference configuration and on the basis of the second law of thermodynamics expressed in the form of the Clausius – Plank inequality, which is a stronger (i.e., less general) form of the Clausius - Duhem inequality since it ignores dissipation due to heat conduction. In order to obtain the standard thermomechanical state equations, Lubliner modifies further the Clausius – Planck inequality, by assuming that it holds as an equality for elastic unloading and neutral loading. It is remarkable to note that by working with the covariance axiom we do not have to ignore dissipation due to heat conduction. Also, in comparing the two approaches we note that while in the second law of thermodynamics we focus on all processes, in the covariance axiom we focus on all transformations of a given process (Marsden and Hughes, 1994, p. 201).

4. A constitutive model

Up to now, the proposed formulation was presented largely in an abstract manner by leaving the kinematics of the problem and the number and the nature of the internal variables unspecified. The basic objective of this section is the introduction of a material model that will help make the application of the generalized plasticity concept in modeling phase transformations clearer. The model is based on a geometrically linear model proposed earlier within a stress space formulation by Panoskaltsis and co-workers (Panoskaltsis et al., 2004, Ramanathan et al., 2002) and which has been extensively used in several applications of engineering interesting (e.g., see Freed et al., 2008; Videnic et al., 2008; Freed and Aboudi, 2008; Freed and Banks – Sills, 2007).

There are two fundamental assumptions underlying the new model which is developed here. The first consists of the additive decomposition of the material strain tensor $$E$$ into elastic $$E^e$$ and inelastic (transformation induced) $$E^{Tr}$$ parts, i.e.,

$$E=E^e+E^{Tr}.$$E39

Such a decomposition has its origins in the work of Green and Naghdi (1965). The second fundamental assumption is that the response of the material is isotropic. Accordingly, it is assumed that it can be described in terms of a single scalar internal variable Z, which, as it is common within the literature (e.g., Boyd and Lagoudas, 1996; Lubliner and Auricchio, 1996; Panoskaltsis et al., 2004; Müller and Bruhns, 2006; Thamburaja, 2010), is identified by the fraction of a single (favorably oriented) martensite variant. In turn, and in view of Eq. (31), the internal variable vector is assumed to be composed by the transformation strain tensor $$E^{Tr}$$and the martensite fraction Z.

By noting that the martensitic transformations to be considered are accompanied by variations of the elastic properties of the SMA material and in view of the additive decomposition of strain (Eq. (31)), the Helmholtz free energy can be additively decomposed in elastic and inelastic (transformation) parts, as follows

$$\Psi \text{=}{\Psi }_e\text{(}E-E^{Tr}\text{(Z),T,Z)+}{\Psi }_{Tr}\text{(Z,T)}\text{.}$$E40

It is emphasized that this is not the conventional decomposition of the free energy function performed within the classical inelastic theories (e.g., plasticity, viscoelasticity, viscoplasticity), since the elastic part $${\Psi }_e$$ depends on the internal variable Z. In this sense the decomposition (32) resembles the decompositions employed within the thermomechanical treatment of damage (see Panoskaltsis et al., 2004). The elastic part of the Helmholtz free energy is assumed to be given as

$$$$

$${\Psi }_e(E-E^{Tr}\text{(Z),T,Z)=}U\left(E-E^{Tr}\text{(Z)},Z\right)+\Theta \left(\text{T}\right)+M\left(E-E^{Tr}\text{(Z),Z,T}\right),$$E41

where the terms U, $$\Theta$$ and M will be defined next. U is the mechanical part of $${\Psi }_e^{}$$ and is assumed to be given by a similar expression to the stored energy function of a St. Venant – Kirchhoff material (e.g., see Holzapfel, 2000, pp. 250 – 251), that is

$$U\left(E-E^{Tr}\text{(Z)},Z\right)=\frac{\lambda (\text{Z)}}{\text{2}}\{tr{[E-E^{Tr}\text{(Z)}}}^{\text{2}}+\mu (\text{Z)}tr\text{[}{(E-E^{Tr}\text{(Z)]}}^{\text{2}},$$E42

where $$\lambda$$ and $$\mu$$ are $$\text{Lam<msup>e <mo>′</mo></msup>}$$ type of parameters $$(\lambda >0,\mu \text{ > 0)}$$ and tr denotes the trace operator. These parameters are assumed to be dependent on the martensite fraction of the SMA, according to the following (power) law

$$\lambda (\text{Z)}={\lambda }_A+{\text{Z}}^n({\lambda }_M-{\lambda }_A),\mu (\text{Z)}={\mu }_A+{\text{Z}}^m({\mu }_M-{\mu }_A),$$E43

where $${\lambda }_A$$, $${\mu }_A$$ are the $$\text{Lam<msup>e <mo>′</mo></msup>}$$ type of parameters when the material is fully austenite,$${\lambda }_M$$, $${\mu }_M$$ are these when the material is fully martensite and n and m are two additional model parameters. For the particular case n = m = 1 the rule of mixtures, which has been used extensively within the literature (e.g., Anand and Gurtin, 2003; Hallberg et al., 2007) is derived.

For the thermal part of the stored energy function, that is for the functions $$\Theta \left(\text{T}\right)$$ and $$M\left(E-E^{Tr}\text{(Z)}\right)$$ we consider the following expressions:

$$\begin{array}{l}\Theta \left(\text{T}\right)=c[(\text{T}-{\text{T}}_{\text{0}})-\text{Tln}\left(\frac{\text{T}}{{\text{T}}_{\text{0}}}\right)],\\ M\left(E-E^{Tr}\text{(Z),Z,T}\right)=-[3\lambda (\text{Z})+2\mu (\text{Z)]}{\alpha }_t(\text{T}-{\text{T}}_{\text{0}})tr[E-E^{Tr}\text{(Z)],}\end{array}$$E44

where $${\text{T}}_{\text{0}}$$ is the reference temperature, $$c$$ is the specific heat and $${\alpha }_t$$ the linear expansion coefficient, which may be assumed varying within the phase transformations according to expressions analogous to those given in Eq. (35).

Finally, the transformation part of the Helmholtz free energy is given as

$${\Psi }_{Tr}=-\text{T}{\eta }_{Tr}(\text{Z)}+u_{Tr}(\text{Z)},$$E45

where $${\eta }_{Tr}(\text{Z)}$$ and $$u_{Tr}(\text{Z)}$$ stand for the configurational entropy and the configurational internal energy and for which we assume two expressions justified in the work of Müller and Bruhns (2006) (see also the thermomechanical theory of Raniecki et al., 1992; Raniecki and Lexcellent, 1998), namely

$$\begin{array}{l}{\eta }_{Tr}(\text{Z)=}{s}_0^{\ast A}+\text{Z}\Delta s^{*}-\text{Z(1-Z)}{\overline{s}}_0,\\ u_{Tr}(\text{Z)=}{u}_0^{\ast A}+\text{Z}\Delta u^{*}-\text{Z(1-Z)}{\overline{u}}_0,\end{array}$$E46

where $$s_0^{\ast A}$$, $$\Delta s^{*}$$, $${\overline{s}}_0$$, $$u_0^{\ast A}$$, $$\Delta u^{*}$$ and $${\overline{u}}_0$$ are the model thermal parameters.

Then in light of the first of Eqs. (20) the second Piola – Kirchhoff stress tensor, after extensive calculations, is found to be

$$S=\lambda tr(E-E^{Tr})1+2\mu (E-E^{Tr}\text{)}-\text{(3}\lambda \text{+2}\mu \text{)}{\alpha }_t(\text{T}-{\text{T}}_{\text{0}})1,$$E47

where the dependence of the involved quantities on Z has been dropped for convenience.

The loading surfaces are assumed to be given in the stress – space as a two parameter family of von - Mises type of surfaces, that is

$$\text{F(}S\text{,T) = }\left|DEVS\right|-C\text{T}-R\text{ = 0,}$$E48

where $$|.|$$ denotes the Euclidean norm, $$DEV(.)$$ stands for the deviatoric part of the stress tensor in the reference configuration and $$C$$ and $$R$$ are parameters. On substituting from Eq. (39) into Eq. (40) the equivalent expression for the loading surfaces in the strain – space may be derived as

$$\Phi \text{(}E\text{, }{E}^{Tr}\text{,T) =2}\mu \left|DEV(E-E^{Tr})\right|-C\text{T }-R\text{ = 0}\text{.}$$E49

For the evolution of the transformation strain we assume a normality rule in the strain – space which is given as

$$\text{2}\mu {\dot{E}}^{Tr}\text{ =}{\epsilon }_L\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>}\frac{\partial \Phi }{\partial E},$$E50

where $${\epsilon }_L$$ is a material constant, which is defined as the maximum inelastic strain (e.g., Boyd and Lagoudas, 1996; Lubliner and Auricchio, 1996; Panoskaltsis et al., 2004; Ramanathan et al., 2002), which is attained in the case of one – dimensional unloading in simple tension when the material is fully martensite.

The rate equation for the evolution of the martensite fraction Z, is determined on the basis of the geometrical framework described in section 2 as follows:

For the austenite to martensite transformation $$\text{(A}\to \text{M)}$$ we consider the $${\Phi }_M-$$loading surfaces as:

$${\Phi }_M\text{(}E\text{, }{E}^{Tr}\text{,T) =2}\mu \left|DEV(E-E^{Tr})\right|-C_M\text{T}-R_M\text{ = 0,}$$E51

where $$C_M$$ is a material parameter which can be determined by means of the well – known (e.g., see Lubliner and Auricchio, 1996; Panoskaltsis et al., 2004; Ramanathan et al., 2002; Christ and Reese, 2009) critical stress – temperature phase diagram for the SMAs transformation. Moreover we consider

$$\begin{array}{l}{\Phi }_{Mf}\text{(}E\text{, }{E}^{Tr}\text{,T) =2}\mu \left|DEV(E-E^{Tr})\right|\text{ + }{R}_{Mf}\text{ = 0,}\\ {\Phi }_{Ms}\text{(}E\text{, }{E}^{Tr}\text{,T) =2}\mu \left|DEV(E-E^{Tr})\right|\text{ + }{R}_{Ms}\text{ = 0,}\end{array}$$E52

where

$$R_{Mf}\text{ = }-{\text{ σ}}_{Mf}-C_M(\text{T}-{\text{M}}_f),R_{Ms}\text{ = }-{\text{ σ}}_{Ms}-C_M(\text{T}-{\text{M}}_s),$$E53

where the parameters $${\text{M}}_f$$ and $${\text{M}}_s$$ stand for the martensite finish and martensite start temperatures respectively, and $${\text{σ}}_{Mf}$$and$${\text{σ}}_{Ms}$$ are two additional parameters which may be determined from experimental results. Since $${\Phi }_{Mf}$$is related to the finish values and $${\Phi }_{Ms}$$to the starting values of the $$\text{A}\to \text{M}$$ transformation, the loading surfaces $${\Phi }_{Mf}=0$$ and $${\Phi }_{Ms}=0$$ may be considered as the boundaries of the set of all states for which the$$\text{A}\to \text{M}$$transformation can be active. Then the constant $$H_1=H_M$$ (see Eq. (3)) may be defined as

$$H_M=\frac{⟨-{\Phi }_{Mf}{\Phi }_{Ms}⟩}{\left|{\Phi }_{Mf}{\Phi }_{Ms}\right|}.$$E54

For the function $$L_1=L_M$$ several choices are possible (see Panoskaltsis et al., 2004). In this work, we use a linear expression (see Lickachev and Koval, 1992), which within the present strain – space formulation may be written in the form

$$L_M(E,E^{Tr},\text{T, Z)}=-\frac{\text{1}-\text{Z}}{{\Phi }_{Mf}-2\mu {\epsilon }_L(1-\text{Z)}}.$$E55

In view of Eqs. (45) and (46) the rate equation for the evolution of the martensite fraction of the material during the$$\text{A}\to \text{M}$$transformation may be written in the form

$$\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>}=-\frac{⟨-{\Phi }_{Mf}{\Phi }_{Ms}⟩}{\left|{\Phi }_{Mf}{\Phi }_{Ms}\right|}\frac{\text{1}-\text{Z}}{{\Phi }_{Mf}-2\mu {\epsilon }_L(1-\text{Z)}}⟨L_M⟩,$$E56

where $$L_M=L_1$$ stands for the loading rate in the material description, that is

$$L_M=N_M:\dot{G}=\frac{\partial {\Phi }_M}{\partial E}:\dot{E}+\frac{\partial {\Phi }_M}{\partial \text{T}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>}\text{.}$$E57

Similarly, for the inverse $$\text{M}\to \text{A}$$ transformation we define the $${\Phi }_A-$$loading surfaces as follows

$${\Phi }_A\text{(}E\text{, }{E}^{Tr}\text{,T) =2}\mu \left|DEV(E-E^{Tr})\right|-C_A\text{T}-R_A\text{ = 0,}$$E58

$$\begin{array}{l}{\Phi }_{Af}\text{(}E\text{, }{E}^{Tr}\text{,T) =2}\mu \left|DEV(E-E^{Tr})\right|\text{ + }{R}_{Af}\text{ = 0,}\\ {\Phi }_{As}\text{(}E\text{, }{E}^{Tr}\text{,T) =2}\mu \left|DEV(E-E^{Tr})\right|\text{ + }{R}_{As}\text{ = 0,}\end{array}$$E59

where

$$R_{Af}\text{ = }-{\text{σ}}_{Af}-C_A(\text{T}-{\text{A}}_f),R_{As}\text{ = }-{\text{σ}}_{As}-C_A(\text{T}-{\text{A}}_s),$$E60

and the parameters $$C_A$$,$${\text{A}}_f$$,$${\text{A}}_s$$,$${\text{σ}}_{Af}$$ and $${\text{σ}}_{As}$$ are material parameters, all related to the $$\text{M}\to \text{A}$$ transformation. By applying analogous to the $$\text{A}\to \text{M}$$ transformation case arguments, we derive the rate equation for the evolution of Z for the $$\text{M}\to \text{A}$$transformation as

$$\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>}=-\frac{⟨-{\Phi }_{Af}{\Phi }_{As}⟩}{\left|{\Phi }_{Af}{\Phi }_{As}\right|}\frac{\text{Z}}{{\Phi }_{Af}+2\mu {\epsilon }_L\text{Z}}⟨L_A⟩,$$E61

where

$$L_A=N_A:\dot{G}=-(\frac{\partial {\Phi }_A}{\partial E}:\dot{E}+\frac{{\Phi }_A}{\partial \text{T}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>)= }-L_M\text{.}$$E62

As a result, the final form for the rate equation for the evolution of the internal variable Z (see Eq. (3)) takes the form

$$\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>}=-\frac{⟨-{\Phi }_{Mf}{\Phi }_{Ms}⟩}{\left|{\Phi }_{Mf}{\Phi }_{Ms}\right|}\frac{\text{1-Z}}{{\Phi }_{Mf}-2\mu {\epsilon }_L(1-\text{Z)}}⟨L_M⟩-\frac{⟨-{\Phi }_{Af}{\Phi }_{As}⟩}{\left|{\Phi }_{Af}{\Phi }_{As}\right|}\frac{\text{Z}}{{\Phi }_{Af}+2\mu {\epsilon }_L\text{Z}}⟨L_A⟩.$$E63

The thermomechanical coupling phenomena, which occur during the martensitic transformations may be studied on the basis of the energy balance equation. It should be mentioned here that with the aid of the fundamental concept of energy it is possible to relate different physical phenomena to one another, as well as to evaluate their relative significance in a given process in mechanics and more generally in physics (Lubliner, 2008, p. 44). This will be accomplished as follows:

The energy balance Eq. (12) can be written in a material setting as

$${\rho }_{\text{ref}}(\dot{\Psi }+\dot{\eta }\text{T}+\eta \text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>)}+DIVH=S:\dot{E}+{\rho }_{\text{ref}}\text{R,}$$E64

where $$DIV(.)$$ is the divergence operator, H is the heat flux vector and R is the heat supply per unit mass, all expressed in the material description. By taking the time derivative of the Helmholtz free energy function and inserting it in Eq. (54) we obtain

$${\rho }_{\text{ref}}(\frac{\partial \Psi }{\partial E}:\dot{E}+\frac{\partial \Psi }{\partial E^{Tr}}:{\dot{E}}^{Tr}+\frac{\partial \Psi }{\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>+}\frac{\partial \Psi }{\partial \text{T}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>)+}{\rho }_{\text{ref}}\eta \text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>+}{\rho }_{\text{ref}}\dot{\eta }\text{T+}DIVH={\rho }_{\text{ref}}\text{R+}S\text{:}\dot{E}.$$E65

This equation in turn, upon substitution of the thermomechanical state Eqs. (20), yields

$${\rho }_{\text{ref}}(\frac{\partial \Psi }{\partial E^{Tr}}:{\dot{E}}^{Tr}+\frac{\partial \Psi }{\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>)+}{\rho }_{\text{ref}}\dot{\eta }\text{T+}DIVH={\rho }_{\text{ref}}\text{R}\text{.}$$E66

The time derivative of the entropy density is determined by the second of Eqs. (20) as

$$\dot{\eta }=-\frac{{\partial }^2\Psi }{\partial \text{T}\partial E}:\dot{E}-\frac{{\partial }^2\Psi }{\partial \text{T}\partial E^{Tr}{}^{}}:{\dot{E}}^{Tr}-\frac{{\partial }^2\Psi }{\partial \text{T}\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover> }-\frac{{\partial }^2\Psi }{\partial {\text{T}}^{\text{2}}}\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>}\text{.}$$E67

Upon definition of the specific heat $$c$$ as

$$c=-\frac{{\partial }^2\Psi }{\partial {\text{T}}^{\text{2}}}\text{T,}$$E68

and upon substitution of Eqs. (57) and (58) into Eq. (56), the latter yields the temperature evolution equation as

$$c\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>}=-(\frac{\partial \Psi }{\partial E^{Tr}}:{\dot{E}}^{Tr}+\frac{\partial \Psi }{\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>)+(}\frac{{\partial }^2\Psi }{\partial \text{T}\partial E}:\dot{E}+\frac{{\partial }^2\Psi }{\partial \text{T}\partial E^{Tr}}:{\dot{E}}^{Tr}+\frac{{\partial }^2\Psi }{\partial \text{T}\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>)T+(R}-\frac{\text{1}}{{\rho }_{\text{ref}}}DIVH).$$E69

If we now define the heating due to thermoelastic effects as

$$\begin{array}{l}{\dot{Q}}_e=\text{T(}\frac{{\partial }^2\Psi }{\partial \text{T}\partial E}:\dot{E}+\frac{{\partial }^2\Psi }{\partial \text{T}\partial E^{Tr}}:{\dot{E}}^{Tr})\\ \end{array}$$E70

and the inelastic (transformation) contribution to heating as

$${\dot{Q}}_{Tr}=-(\frac{\partial \Psi }{\partial E^{Tr}}:{\dot{E}}^{Tr}+\frac{\partial \Psi }{\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>)+T}\frac{{\partial }^2\Psi }{\partial \text{T}\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>=}{D}_{Tr}\text{+T}\frac{{\partial }^2\Psi }{\partial \text{T}\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>,}$$E71

where $$D_{Tr}=-(\frac{\partial \Psi }{\partial E^{Tr}}:{\dot{E}}^{Tr}+\frac{\partial \Psi }{\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>)}$$ is the inelastic dissipation due to phase transformations, the temperature evolution equation takes the following, remarkably simple, form (see also Rosakis et al., 2000)

$$c\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>}={\dot{Q}}_e+{\dot{Q}}_{Tr}\text{+ (R }-\frac{\text{1}}{{\rho }_{\text{ref}}}DIVH).$$E72

This expression has the obvious advantage of decoupling the elastic and inelastic contributions to material heating and is well suited for computational use.

It is noted that in an adiabatic process, that is in a process with $$\text{R}-\frac{\text{1}}{{\rho }_{\text{ref}}}DIVH=0,$$

Eq. (62) takes the form

$$c\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>}={\dot{Q}}_e+{\dot{Q}}_{Tr}$$E73

,

from which and by assuming that the temperature evolution due to structural heating $$H_{str}$$, defined as $$H_{str}=-\text{(}\frac{{\partial }^2\Psi }{\partial \text{T}\partial E}:\dot{E}+\frac{{\partial }^2\Psi }{\partial \text{T}\partial E^{Tr}}:{\dot{E}}^{Tr}+\frac{{\partial }^2\Psi }{\partial \text{T}\partial \text{Z}}\text{<mover accent="true"><mi>Z</mi><mo>˙</mo></mover>)T}$$ (e.g., see Simo and Miehe, 1992), is negligible in comparison to that due to inelastic dissipation $$D_{Tr},$$ the temperature evolution equation takes the following simple form

$$c\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>}=D_{Tr}.$$E74

Finally, as a constitutive law for the heat flux vector we assume the standard Fourier’s law (e.g., Simo and Miehe, 1992; Müller and Bruhns, 2006):

$$H=-kGRAD\text{T}\text{.}$$E75

5. Computational aspects and numerical simulations

As a final step we examine the ability of our model in simulating qualitatively several patterns of the extremely complex behavior of SMAs under simple states of straining. Isothermal and non – isothermal problems are considered.

5.1. Isothermal problems

Focusing our attention first in the isothermal case we note that when the total strain tensor E is known, the rate equations for the evolution of the internal variables (Eqs. (42) and (53)) and the mechanical state (thermoelastic stress-strain law) equation (Eq. (39)) together with the appropriate initial and boundary conditions form a system of three equations in the three unknowns $$E^{Tr}$$, Z and S. The numerical solution of this system of equations and accordingly the numerical implementation of the proposed model relies crucially on the general loading – unloading criteria (see Eq. (4)), which can be expressed in a remarkably simple form, based on the following observation:

As it has been mentioned the $$\text{A}\to \text{M}$$ transformation is active when $$L_M>0,$$ while the inverse transformation is active when $$L_A>0.$$ Since we always have $$L_M=-L_A,$$ it is clear that only one phase transformation can be active at a given time of interest. Then we can treat the two phase transformations as two different inelastic processes and replace the general loading – unloading criteria by the following decoupled ones:

$$\text{A}\to \text{M}$$Transformation:

$$\{\begin{array}{l}\text{If }{H}_M=\frac{⟨-{\Phi }_{Mf}{\Phi }_{Ms}⟩}{\left|{\Phi }_{Mf}{\Phi }_{Ms}\right|}=0:\text{ elastic state,}\\ \text{ If }{H}_M=\frac{⟨-{\Phi }_{Mf}{\Phi }_{Ms}⟩}{\left|{\Phi }_{Mf}{\Phi }_{Ms}\right|}=1\text{ then }\\ \text{ i}\text{. If }{L}_M\text{<0: elastic unloading,}\\ \text{ ii}\text{. If }{L}_M\text{=0: neutral loading,}\\ \text{ iii}\text{. If }{L}_M\text{>0: inelastic loading}\text{. }\end{array}$$E76

$$\text{M}\to \text{A}$$Transformation:

$$\{\begin{array}{l}\text{If }{H}_A=\frac{⟨-{\Phi }_{Af}{\Phi }_{As}⟩}{\left|{\Phi }_{Af}{\Phi }_{As}\right|}=0:\text{ elastic state,}\\ \text{ If }{H}_M=\frac{⟨-{\Phi }_{Af}{\Phi }_{As}⟩}{\left|{\Phi }_{Af}{\Phi }_{As}\right|}=1\text{ then }\\ \text{ i}\text{. If }{L}_A\text{<0: elastic unloading,}\\ \text{ ii}\text{. If }{L}_A\text{=0: neutral loading,}\\ \text{ iii}\text{. If }{L}_A\text{>0: inelastic loading}\text{. }\end{array}$$E77

Then the governing equations, along with the aforementioned loading – unloading criteria, can be solved by a time discretization scheme based on backward Euler. The resulting system of the discretized equations is solved by means of a three step predictor - corrector algorithm, the steps of which are dictated by the time discrete loading - unloading criteria. Algorithmic details regarding the enforcement of the time discrete loading – unloading criteria and the solution of the system, within the framework of large deformation generalized plasticity in the case of a single loading surface, can be found in Panoskaltsis et al. (2008a, b).

To this end it is emphasized that predictor – corrector algorithms work well in case of domains which are connected. The commonly used predictor – corrector algorithms for elastoplasticity employ an elastic predictor and an inelastic corrector. The most important assumption is that the solution is unique for a particular set of values of the state variables. The predictor step freezes the plastic flow and checks for an elastic solution. The yield criterion then is checked and if it is satisfied the elastic solution is acceptable, otherwise the inelastic corrector is activated. In the cases of elastic – plastic analysis there exists a set of consistency conditions the enforcement of which “returns” the (wrong) elastic solution onto the exact solution point on the evolving yield surface. However, in the case of disconnected elastic zones separated by inelastic zones the predictor – corrector algorithm is very sensitive on the strain step used, while going from an inelastic zone to an elastic one.

This is the case of SMAs, which have a transformation (inelastic) zone separating the fully martensite and fully austenite zones (being treated as elastic zones). During forward or reverse transformation, the predictor strain step is very important as we near the elastic – plastic (i.e. transformation) boundary. If the predicted solution lies within the transformation zone (i.e., outside the elastic range) the corrector step is activated and the resulting set of non –linear equations are solved. However, as we approach the end of the transformation zone and therefore the boundary between the inelastic and the elastic zones, the predictor could predict an elastic solution, which the algorithm accepts as a valid one, but which is within the next elastic zone, achieved without the transformation being fully complete (i.e., achieved while the state is still inelastic) and is therefore an unacceptable solution. This would cause errors in the minimization process and results in jumps in the solution and kinks in the stress strain curve. This problem is resolved here by making the strain step very small and by checking the limits of the transformation.

The first problem we study is a standard problem within the context of finite inelasticity and is that of finite shear, defined as

$$x_1=X_1+\gamma X_2,x_2=X_2\text{, }{x}_3=X_3,$$E78

where $$X_1,X_2,X_3$$ are the material coordinates and γ is the shearing parameter. For this problem the model parameters are set equal to those reported in the work of Boyd and Lagoudas (1994), that is:

$$\begin{array}{l}{\lambda }_{\text{M}}=9,486.95\text{ MPa},{\mu }_{\text{M}}=4,887.22\text{ MPa, }{\lambda }_{\text{A}}=21,892.97\text{ MPa},{\mu }_{\text{A}}=11,278.20\text{ MPa,}\\ {\text{M}}_{\text{f}}=5^o{\text{C, M}}_{\text{s}}={23}^o{\text{C, A}}_{\text{s}}={29}^o{\text{C, A}}_{\text{f}}={51}^o{\text{C, C}}_{\text{M}}=11.3{\text{ MPa/}}^o{\text{C, C}}_{\text{A}}=4.5{\text{ MPa/}}^o\text{C,}\\ {\sigma }_{\text{Mf}}\text{=}{\sigma }_{\text{Ms}}\text{=}{\sigma }_{\text{Af}}\text{=}{\sigma }_{\text{As}}\text{=0 MPa, }{\epsilon }_{\text{L}}\text{=0}\text{.0635}\text{.}\end{array}$$E79

All numerical tests that performed start with the specimen in the parent (austenite) phase, (Z=0).

The first simulation demonstrates the pseudoelastic phenomena within the SMA material. In this case the temperature is held constant at some value above$${\text{A}}_{\text{f}}$$. The purpose is to study a complete stress – induced transformation cycle. The results for this finite shear problem are shown for constant material stiffness $$(\lambda ={\lambda }_{\text{A}}=21,892.97\text{ MPa},\mu ={\mu }_{\text{A}}=11,278.20\text{ MPa),}$$ as well for a linear (n = m = 1) and a power (n = m = 5) type of stiffness variation in Figures 1, 2 and 3. On loading, the material initially remains austenite (elastic region and straight shear stress – strain curve). As loading is continuing and the shear strain attains the value at which the material point crosses the initial loading surface for the $$\text{A}\to \text{M}$$ transformation $$({\Phi }_{Ms}=0)$$, the transformation starts (inelasticity and curvilinear shear stress – strain curve; coexistence of the two phases). If the loading continues and the strain crosses the final loading surface for the $$\text{A}\to \text{M}$$ transformation $$({\Phi }_{Mf}=0),$$ the material is completely transformed into martensite and on further loading since the state of the material is elastic the shear stress – strain diagram is straight. Then, during unloading, the material is fully martensite (elastic region and straight shear stress – strain curve) until the strain crosses the initial loading surface $$({\Phi }_{\text{As}}=0)$$ of the$$\text{M}\to \text{A}$$transformation, which is subsequently activated (phase coexistence, inelasticity and curvilinear shear stress – strain curve). On further unloading and when the strain meets the last boundary surface for the $$\text{M}\to \text{A}$$transformation $$({\Phi }_{\text{Af}}=0)$$, the material becomes fully austenite and on further unloading the stress – strain curve is straight going back to zero, which means that no permanent deformation exists and the austenite is completely recovered. This is expected as the martensite phase is not stable at a temperature above $${\text{A}}_{\text{f}}$$at zero stress level.

Figure 1.

Finite shear. Isothermal one – dimensional behavior under monotonic loading. Shear stress $${\text{S}}_{\text{12}}$$ vs. shear strain $$\gamma .$$

Figure 2.

Finite Shear. Isothermal one – dimensional behavior under monotonic loading. Normal stress $${\text{S}}_{\text{11}}$$vs. shear strain $$\gamma .$$

Figure 3.

Finite Shear. Isothermal one –dimensional behavior under monotonic loading. Normal stress $${\text{S}}_{\text{22}}$$vs. shear strain $$\gamma .$$

Next, the model is tested under multiple shear stress cycles, by subjecting it to partial unloading (incomplete $$\text{M}\to \text{A}$$ transformation) and partial reloading (incomplete $$\text{A}\to \text{M}$$ transformation). The results for linear stiffness variation are illustrated in Figure 4. A series of loops appears inside the complete loading – unloading cycle. These loops exhibit slight ratcheting which stabilizes in a few cycles. The response of the model is absolutely compatible with that described by other investigators (e.g., see Ivshin and Pence, 1994; Lubliner and Auricchio, 1996). In view of Figure 4 and since the dissipated energy can be estimated by the area of the $${\text{S}}_{\text{12}}-\gamma$$ loop, the dissipated energy in the case of partial unloading and reloading is the area of the loop times the number of the loops. This explains the important property of the high internal damping of SMA materials. (For a discussion of the relation between areas of stress-strain diagrams and dissipated energy see Lubliner and Panoskaltsis, 1992.)

The ability of the model to simulate phase transformations and the corresponding stiffness variations under cyclic loading is demonstrated further by three additional tests. The first one illustrates the case of partial loading with complete unloading, the second the case of partial unloading with complete loading and the third the case of a series of partial loading and partial unloading. The results are shown in Figures 5, 6 and 7 respectively.

Figure 4.

Finite Shear. Pseudoelasticity with partial loading and unloading. Shear stress $${\text{S}}_{\text{12}}$$vs. shear strain $$\gamma .$$

Figure 5.

Finite Shear. Partial loading followed by complete unloading. Shear stress $${\text{S}}_{\text{12}}$$vs. shear strain $$\gamma .$$

Figure 6.

Finite Shear. Partial unloading followed by complete loading. Shear stress $${\text{S}}_{\text{12}}$$vs. shear strain $$\gamma .$$

Figure 7.

Finite Shear. Series of partial loading and partial unloading. Shear stress S₁₂ vs. shear strain γ.

5.2. Non – isothermal problems

In this section we examine the ability of the model in predicting pseudoelastic phenomena under non – isothermal conditions. In general, the numerical treatment of the coupled thermomechanical problems is performed on the basis of a split of the governing equations (equations of motion, constitutive equations, energy balance equation and the appropriate boundary conditions) into their mechanical and thermal parts. Most popular among the several computational schemes which have been proposed within the literature is the isothermal split proposed in the work of Simo and Miehe (1992). However, this rather straight forward numerical scheme has the disadvantage of not being unconditionally stable. An alternative time integration algorithm relies on the so-called adiabatic split (see Armero and Simo, 1993). In this approach, the problem is divided in a mechanical phase during which the entropy is held constant, followed by a thermal phase in which the configuration is held constant, leading to an unconditionally stable algorithm.

Nevertheless, since our objective is to discuss the proposed framework in its simplest setting, we consider two rather simple problems, namely a simple shear and a plane strain problem, where the equations of motion and the (mechanical) boundary conditions are trivially satisfied. Accordingly, within our simulations, a simultaneous solution of the remaining governing equations, namely the constitutive equations and the balance of energy equation, is performed.

First, an adiabatic test in finite simple shear is considered. We assume that due to the dynamic rates resulting in adiabatic response, heat exchanges due to conduction, convection and radiation can be neglected in comparison to the temperature changes induced by inelastic (transformation) dissipation, which leads to thermomechanical processes that can be considered as homogeneous. The elastic constants, the mass density and the thermal parameters used in this simulation are those considered in the work of Müller and Bruhns (2006), that is:

$$\begin{array}{l}\lambda ={\lambda }_{\text{A}}=25,541.80\text{ MPa},\mu ={\mu }_{\text{A}}=13,157.90\text{ MPa,}\\ {\rho }_{\text{ref}}^{}=6.45\times {10}^{-3}{\text{ k/m mm}}^{\text{2}},{\alpha }_{\text{t}}=8.8\times {10}^{-6}1{/}^{\text{o}}\text{K, c=837}\text{.36 J/kg}\text{.K}\\ \Delta u^{\text{*}}=16800.0\text{ J/kg, }\Delta s^{\text{*}}=64.50\text{ J/kgK, }{\overline{u}}_0=4264.5\text{ J/kg, }{\overline{s}}_0=11.5\text{ J/kgK, }\end{array}$$E80

while the other parameters are set equal to those used in the isothermal problems studied before. The shear stress – strain curves predicted by the model, for both adiabatic and isothermal cases, are shown in Figure 8. It is observed that the stress – strain curves have similar qualitative characteristics with the adiabatic and the isothermal curves of a perfect gas in a pressure – volume diagram, with the adiabatic stress curve being above the corresponding isothermal one. This fact has to be attributed to material heating due to inelastic dissipation during the $$\text{A}\to \text{M}$$transformation, which shifts the stress – strain curve upwards. Moreover, due to the higher stress attained during the $$\text{A}\to \text{M}$$ transformation, the initial loading surface for the inverse transformation$$({\Phi }_{\text{As}}=0)$$ is triggered at a higher stress level, a fact which results in a corresponding higher stress – strain unloading curve. The corresponding temperature – shear strain curve for the adiabatic specimen is shown in Figure 9 (for constant stiffness). Consistently with the experimentally observed adiabatic response of a SMA material, the model predicts heating of the material during the forward $$\text{A}\to \text{M}$$ transformation and cooling during the inverse $$\text{M}\to \text{A}$$ transformation.

Figure 8.

Finite Shear. Adiabatic and isothermal one – dimensional behavior under monotonic loading. Shear stress $${\text{S}}_{\text{12}}$$vs. shear strain $$\gamma .$$

Figure 9.

Finite Shear. Adiabatic one – dimensional behavior under monotonic loading. Temperature$$\text{T}$$vs. shear strain $$\gamma .$$

Figure 10.

Plane strain (restrained tension). Monotonic loading at various temperatures. Normal stress$${\text{S}}_{\text{11}}$$vs. axial displacement $$\lambda .$$

Figure 11.

Plane strain (restrained tension). Monotonic loading at various temperatures. Normal stress$${\text{S}}_{\text{33}}$$ vs. axial displacement $$\lambda .$$

Figure 12.

Plane Strain (biaxial extension). Monotonic loading at various temperatures. Normal stress$${\text{S}}_{\text{11}}$$ vs. axial displacement $$\lambda .$$

Figure 13.

Plane strain (biaxial extension). Monotonic loading at various temperatures. Normal stress$${\text{S}}_{\text{33}}$$ vs. axial displacement $$\lambda .$$

Next, we study a plane strain model, that of the biaxial extension of a material block. The straining occurs along $X_1\text{ and }{X}_2$ axes while the block is assumed to be fixed along the $X_3$ direction. This problem is defined as

$$x_1=(1+\lambda )X_1,x_2^{}=(1+\omega )X_2,x_3=X_3^{},$$E81

where λ and ω are the straining parameters.

The isothermal stress – displacement curves for the limiting cases $$\omega =0$$(restrained tension) and $$\omega$$=$$\lambda$$ (biaxial tension) for three different material temperatures ($${\text{T}}_{\text{0}}={60}^o\text{C}>{\text{A}}_{\text{f}},$$ $${\text{A}}_{\text{f}}>{\text{T}}_{\text{0}}={45}^o\text{C}>{\text{A}}_{\text{s}},$$ $${\text{A}}_{\text{s}}>{\text{T}}_{\text{0}}={25}^o\text{C}$$) are shown in Figures 10, 11, 12 and 13. By referring to these figures for $${\text{T}}_{\text{0}}>{\text{A}}_{\text{f}}$$, we easily verify the ability of the model in predicting pseudoelastic phenomena in two dimensions.

The isothermal tests for $${\text{A}}_{\text{f}}>{\text{T}}_{\text{0}}>{\text{A}}_{\text{s}}$$ and $${\text{A}}_{\text{s}}>{\text{T}}_{\text{0}}$$ are conducted in order to show the ability of the model in predicting the shape memory effect. In the first of them, upon loading the $$\text{A}\to \text{M}$$ transformation is activated, but since the temperature is less than the temperature required for the complete reverse transformation at zero stress, upon unloading the two phases coexist and permanent deformations appear. However, these deformations are recovered after increasing the temperature. In the second test the temperature initially is kept constant at a value less than the austenite start temperature at zero stress. As a result, at the end of the stress cycle the material is completely in the martensite phase and large permanent deformation appears. Nevertheless, like in the previous test, this deformation may be eliminated upon heating. For these new non – isothermal (i.e. heating) problems we assume thermal boundary conditions corresponding to convective heat exchange between the specimen and the surrounding medium on the free faces (with area $$A$$) of the specimen. In this case the normal heat flux is given by Newton’s law of cooling (e.g., see Simo and Miehe, 1992) as: $$H_u=hA({\text{T}}_{\infty }-{\text{T}}_{\text{0}}),$$ with $$h$$ being the constant convection coefficient, which is chosen as $$h={17.510}^{-3}{\text{ N/mm}}^{\text{o}}\text{K}$$, and $${\text{T}}_{\infty }$$ is the surrounding medium temperature. By assuming that the size of the tested material is small, the contribution to the material heating due to heat conduction can be neglected, so that the temperature evolution equation (see Eq. (62)) can be written in the form

$$c\text{<mover accent="true"><mi>T</mi><mo>˙</mo></mover>}={\dot{Q}}_e+{\dot{Q}}_{Tr}-\frac{\text{1}}{{\rho }_{\text{ref}}}{H}_u.$$E82

The results of these tests are illustrated in Figures 14 and 15, where the elongation along $$X_1$$ axis is plotted versus the surrounding medium temperature. The slight increase of the elongation of the SMA material due to the (elastic) thermal expansion occurring prior to the activation of the $$\text{M}\to \text{A}$$transformation, for initial temperature $${\text{T}}_{\text{0}}{\text{<A}}_{\text{s}}$$, is noteworthy (Figure 15).

Figure 14.

Plane strain. Shape memory effect $${\text{(A}}_{\text{s}}{\text{<T}}_{\text{0}}{\text{<A}}_{\text{f}}).$$ Axial displacement $$\lambda$$ vs. surrounding medium temperature $${\text{T}}_{\infty }.$$

Figure 15.

Plane strain. Shape memory effect $${\text{(T}}_{\text{0}}{\text{<A}}_{\text{s}}).$$ Axial displacement $$\lambda$$ vs. surrounding medium temperature $${\text{T}}_{\infty }.$$

6. Concluding remarks

In this chapter we developed a geometrical framework for the establishment of constitutive models for materials undergoing phase transformations and in particular for shape memory alloys. The proposed framework has the following characteristics:

1. It is quite general for the derivation of the kinetic equations governing the transformation behavior and it can describe multiple and interacting loading mechanisms.

2. It formulates general loading – unloading criteria, in both their material and spatial settings, that can be systematically employed for the numerical implementation of the derived constitutive models.

3. It can describe rate effects.

4. It can model non-isothermal conditions.

5. It can model transformation induced plasticity by considering it as an additional phase transformation.

6. It employs the invariance of the spatial balance of energy equation under the superposition of arbitrary spatial diffeomorphisms – that is spatial transformations which can change the Euclidean metric – as a basic constitutive hypothesis, in place of the second law of thermodynamics.

As an application a specific three – dimensional thermomechanical constitutive model for SMA materials is derived. The model can simulate several patterns – under isothermal and non-isothermal conditions – of the extremely complex behavior of these materials such as:

1. The pseudoelastic behavior observed under monotonic loading.

2. The pseudoelastic behavior observed under several cyclic loadings.

3. The stiffness variations occurring during phase transformations.

4. The shape memory effect.

Additionally, the basic differences between the classical return mapping algorithms and the one used here for the case of not connected regions, have been outlined.