A new inelastic mechanism, detwinning -induced plasticity (DIP), is proposed to model the response of NiTi SMAs to cyclic loading , based on thermodynamic considerations . DIP is incorporated into a constitutive framework for NiTi SMA. The constitutive framewo rk also includes well -established inelastic mechanisms of phase transformation, transformation -induced plasticity, residual martensite, and detwinning . The model is constructed at the single crystal scale using the framework of thermodynamics and a crystal plasticity formulation. An explicit scale -transition rule is adopted for the simulation of polycrystalline materials, allowing direct comparison of the model predictions with published experimental test data. Thermodynamic considerations result in a stron g contribution of DIP for cyclic loading regimes where compressive stress occurs during part of the loading cycle. However, the contribution of DIP is negligible for cyclic loading regimes that result exclusively in tensile stress. This predicted dependenc e of DIP on compression, but not on tension, is strongly supported by experimental cyclic loading results. Inclusion of DIP results in improved prediction of experimentally observed NiTi SMAs behavior, including strain -controlled cyclic compression -unloadi ng and cyclic tension -unloading tests and stress -controlled cyclic tension -compression and tension -unloading tests. During the first loading cycle the contribution of DIP is not significant in any loading regimes. However, in cases where compressive stress occurs during part of the loading cycle, DIP contributes strongly to the material response from the second cycle onwards . In strain -controlled cyclic compression -unloading tests DIP leads to a less negative peak stress and a more negative residual strain following several loading cycles. In stress -controlled tension -compression cyclic loading DIP leads to a reduction of peak and residual strains. Keywords: NiTi; Cyclic loading; Detwining -induced plasticity; Constitutive modeling; Crystal plasticity 1 Currently visiting researcher at National University of Ireland, Galway *Corresponding author: arghavani@sharif.edu 2 1. Introduction Shape memory alloy s (SMAs) have unique characteristics such as super -elasticity , shape memory effect, and high damping capacity . Such characteristics make them suitable for diverse applications such as biomedical devices, sensors , and actuators in micro -electro -mechanical systems, aerospace, automotive , and robotics structures (Hartl and Lagoudas, 2007 ; Otsuka and Ren, 1999 ). Of the different categories of SMAs, NiTi alloys are the most common due to high energy density, corros ion resistance , and biocompatibility (Hartl and Lagoudas, 2007 ). The behavior of NiTi SMA is completely dependent on the testing temperature. In the absence of stress, a reduction in temperature will transform austenite to twinned martensite. At temperatures lower than the martensite finish -temperatu re the alloy consists of the original tw inned martensite phase. Application of sufficient stress to twinned martensite will cause the production of inelastic strain which is not recovered by removal of the stress, but is recovered during subsequent heating process. This behavior is known as shap e memory effect in NiTi SMAs. At temperatures higher than the austenit e finish -temperature, the alloy consists of austenite phase. The stress -induced phase transformation leads to strain generation during loading, and subsequent strain recovery upon unload ing. This behavior is known as the super -elasticity effect in NiTi shape memory alloys. Most applications of super -elastic NiTi entail cyclic loading , which can induce changes in the microstructure and can cause macroscopically observable changes in the material behavior . Many experiments have investigate d the effects of cyclic loading on the behavior of NiTi SMAs. The start-stress of forward transformation , and the dissipation energy, are found to decrease with increasing loading cycles, whereas t ransformation hardening is found to increase (Eggeler et al., 2004 ; Kang et al., 2009 ). The peak and residual strains are found to increase up to a saturated level with increasing loading cycles (Kang et al., 2009 ; Nemat -Nasser and Guo, 2006 ; Saleeb et al., 2013 ; Song et al., 2014 ). Microstructural changes of NiTi SMA under cyclic loadings have also been extensively investigated (Ammar et al., 2017 ; Brinson et al., 2004 ; Chen et al., 2019 ; Delville et al., 2011 ; Liu and Xie, 2007 ; Liu et al., 1998 ; Xie et al., 1998 ). Under cyclic loading, when the strain is positive similar mechanisms are involved in the generation of compressive and tensile stresses, as discussed by Liu and Xie (2007). However, when the strain is negative, different mechanisms are involved in the generation of tensile and compressive stresses (Liu, 2002 ; Liu and Xie, 2007 ; Liu et al., 1998 ). Detwinning is defined as the movement of interfaces/twin boundaries between lattice 3 correspondent variants ( Thamburaja, 2005; Thamburaja et al., 2005) . Each martensite habit pl ane variant (or simply variant) is composed of two lattice correspondent variants. For simplicity, lattice correspondent variants are referred to as sub -variants (Yu et al., 2014 ). The microscopic observations show a high density of dislocations at the interfaces of sub -variants during compressive loading while no significant plastic deformation has been detected in these areas under tension (Liu and Xie, 2 007; Liu et al., 1998 ; Xie et al., 1998 ). Since one of the main inelastic deformation mechanisms for super -elastic NiTi is detwinning , and this kind of plasticity occurs at the interface of sub -variants, this phenomenon is referred to as detwinning -induced plasticity (DIP) hereafter. In an effort to develop an improved model for NiTi SMA s, a constitutive law for DIP is proposed and implem ented in the current study . Phenomena such as detwinning -induced plasticity and transformation -induced plasticity generate irrecoverable strains which will accumulate during cyclic loading. In contrast, austenite to martensite transformation results in strains that are recoverable, either by unloading (super -elasticity effect) or by temperature increase (shape memory effect) (Lagoudas, 2008) . Based on experimental observations , many constitu tive models have been derived to describe the behavior of SMAs . These models can be classified into two main groups: macroscopic models and micromechanical models . The macroscopic models are suitable for finite element implementation and are readily applic able to structural analysis . However, the underlying microscale physical mechanisms of SMA behavior are not reflected in macroscale phenomenological constitutive model s. Representative macroscopic SMA models include Auricchio et al. (2007 ), Lagoudas et al. (2012 ), Panico and Brinson (2007 ), Zaki and Moumni (2007 ), Arghavani et al. (2010 ), and Xu et al. (2019 ). Many micromechanical constitutive models have been constructed to incorporate the microscopic physical nature of deformation of NiTi SMA s. Micromechanical models based on crystal plasticity formulation are popular since martensite variants with different morphological feature s and their evolution laws during deformation can be implemented in such a framework . These models are constructed at the scale of a single crystal and then transformed to a polycrystalline scale using either the finite element method , or scale -transition rules. Typical constitutive models include Thamburaja and Anand (2001 ), Anand and Gurtin (2003 ), Thamburaja (2005 ), Thamburaja et al. (2005 ), Wang et al. (2008 ), Yu et al. (2012 ; 2013 ; 2014 , 2015 ), and Dhala et al. (2019 ). Thamburaja and Anand have devel oped a constitutive model for austenite to martensite phase transformation in textured super -elastic NiTi SMA s (Thamburaja and Anand, 2001 ). The thermo - 4 mechanically coupled behavior of super -elastic SMA s consi dering phase transformation ha ve been investigated by Anand and Gurtin (2003 ). In order to consider the effects of detwinning and reorientation, a constitutive model was proposed by Thamburaja (2005 ). This model was modified to include austenite to martensite phase transformation (Thamburaja et al., 2005 ). The interaction s between phase transformation and plasticity in austenite and martensite phases have also been investigated (Dhala et al., 2019 ; Manchiraju and Anderson, 2010 ; Wang et al., 2008 ; Yu et al., 2014 ). In order to simulate the behavior of SMAs under cyclic loadings, Yu et al. (2013 , 2015 ) proposed a constitutive law that incorporates a number of inelastic mechanisms , including : austenite to martensite phase transformation, reorientation, transformation -induced plasticity, reorientation -induced plasticity , and accumulation of residual martensite . While previous studies have considered effects such as phase transformation, transfor mation -induced plasticity, detwinning , and accumulation of residual martensite (Dhala et al., 2019 ; Thamburaja and Anand, 2001 ; Thamburaja et al., 2005 ; Wang et al., 2008 ; Xiao et al., 2018 ; Yu et al., 2014 ; Yu et al., 2013 ; Yu et al., 2015 ), the effect of DIP has not been previously considered . The goal of the current study is to incorporate the effect of detwinning -induced plasticity into a constitutive model for NiTi SMAs . To this end , a constitutive model is proposed based on crystal plasticity and, for the first time, the inelastic mechanism of detwinning -induced plasticity is incorporated . The model is developed at the single crystal scale and also includes phase transformation, transformation -induced plasticity, detwinning, and accumulation of residual martensite. Adopting a new Helmholtz free energy and using the framework of thermodynamics, the driving forces for inelastic mechanisms are derived. An explicit scale -transition r ule (Yu et al., 2014 ; Yu et al., 2013 ; Yu et al., 2012 , 2015 ) is used to obtain the polycrystalline response from that of single crystal. The improved predictive capability of the proposed DIP constitutive law is demonstrated by comparing model results with published cyclic stress -controll ed and strain - controlled ex perimental data . 2. Micro mechanical constitutive model 2.1. Single crystal constitutive model in the framework of thermodynamics 5 Considering small deformation, the total strain tensor in a representative volume element (RVE) of a single crystal is decomposed into five components : elastic strai n tensor eε , transformation strain tensor trε , transformation -induced plastic strain tensor TRIPε , detwinning strain tensor detε , and detwinning -induced plastic strain tensor DIPε . The total stain tensor , ε, and the inelastic strain tensor, inε , are given as: e in=+ε ε ε (1.a) in tr TRIP det DIP= + + +ε ε ε ε ε (1.b) 24 variants are present in the transformed martensite of NiTi shape memory alloy. The transformation strain of each martensite variant can be calculated from the normal and transformation direction of its habit plane . Therefore, the total transformation strain can be calculated considering the contribution of each variant (Patoor et al., 2006 ): 24 1tr trg  ==εP (2.a) ( )1 2    =  + P m n n m (2.b) where trε is the total transformation strain caused by austenite to martensite transformation,  is the volume fraction of the -th martensite variant, and trg is the magnitude of shearing deformation caused by martensite transformation. P is the orientation tensor of the -th martensite variant which is related to transformation direction m , and habit plane normal n , of the -th martensite variant . The magnitudes of habit plane normal n and transformation direction m of 24 variants are listed in the Appendix. As experimentally observed (Brinson et al., 2004 ; Kang et al., 2009 ), during cyclic loading some part of martensite will not transform back to austenite even when the load is completely removed. This part of martensite , which is pinned by the increase of defects , is called residual martensite. Therefore , the volume fraction of each martensite variant can be d ivided into two parts, i.e. the reversible martensite volume fraction rev and the irreversible ( residual ) martensite volume fraction ire . The total volume fraction of martensite is obtained by summing the reversible and irreversible parts of each variant. rev ire    =+ (3.a) 6 24 24 24 1 1 1, , rev rev ire ire           = = == = =   (3.b) Considering reversible and residual parts for each martensite variant, the transformation strain can be re-written as: 24 24 11tr tr tr tr tr rev ire rev ire gg     === + = + ε ε ε P P (4) Due to the misfit at some interfaces, a local stress field is induced at these areas. Plastic deformation will occur to relax this local stress. As discussed by Kang et al. (2009), transformation -induced plasticity (TRIP) is caused by irreversible friction slip at austenite - martensite interfaces. Assuming the same orientation for transformation and transformation - induced plasticity, the strain caused by TRIP can be d ivided into 24 components (Yu et al., 2013 ): 24 1TRIP TRIP  ==εP (5) where TRIPε is the transformation -induced plastic strain tensor and TRIP is the amount of slipping caused by transformation -induced plasticity at the interface of austenite and the -th martensite variant. As discussed by Thamburaja (2005 ) and Thamburaja et al. (2005 ), austenite to martensite phase transformation will occur in the absence of stress if temperature is decrea sed to a value lower than the martensite finish -temperature. This transformation will cause martensitic plates or habit plane variants (HPVs) to form. The martensite within the HPVs comprises of two lattice correspondent variants (LCVs), which are called s ub-variants in this paper. The HPVs and LCVs are separated by interfaces. The movement of interfaces/twin boundaries between LCVs, which is motivated by the application of stress, is referred to as detwinning. As detwinning occurs inside each variant , this contribution is divided into 24 components: ()24 0 1det det       ==−εP (6.a) ( )1 2det    =  + P a w w a (6.b) where the initial volume fraction s of two sub -variants in the -th martensite variant are 0 and 01- .  and 1- are current volume fractions of two sub -variants , and ()0− represents 7 the amount transformed from one sub -variant to the other. detP is the orientation tensor for the -th detwinning system which is related to shear direction for detwinning w , and unit vector normal to detwinning plane a , of the -th martensite variant. The crystallographic data for different detwinning systems are listed in the Appendix . Experimental observations of SMAs undergoing cyclic loading show that some part of the developed strain is not recovered after unloading and accumulates with every transformation cycle. Such plasticity can be activated even when the applied stress is much lower than the yielding stress of austenite and martensite phases (Lagoudas, 2008 ). DIP is related to the part of the irrecoverable strain induced to relax the local stress field at the interfaces of sub -variants . Using the same orientation tensor for detwinning and DIP, a DIP strain tensor is proposed as follows: 24 1DIP DIP det  ==εP (7) where DIP is the amount of shearing deformation caused by detwinning -induced plasticity in the -th martensite variant. The Helmholtz free energy of RVE is considered to be composed of five parts : int e ch tr det rev      = + + + + (8) where e represents the elastic energy of RVE and is dependent on the elastic strain tensor , such that: ()1::2e e e e=ε ε C ε (9.a) where C is the fourth -order elasticity tensor , which is considered to be isotropic for simplicity . ch is the chemical energy of RVE which reflects the effect of ambient temperature on the transformation : ()() ()00 0, lnch TT c T T T T TT    = − − + −    (10) where 0T is ref erence temperature, c is specific heat capacity and  is a constant. 8 int is the part of free energy which reflects the effect of internal stress and is adopted similar to the work done by Yu et al. (2013 ). The role of int is given as int:in=−Bε (11) where B is the internal stress which is introduced to describe the decrease in the start-stress of martensite transformation during cyclic loading . The internal stress facilitates the start of transformation and therefore decreases the start -stress of phase transformation. As the internal stress is caused by martensite variants, it can be d ivided into 24 components. 24 1 ==BB (12) where B is the internal stress caused by the -th martensite variant. The orientation of B is assumed to be the same as the corresponding martensite variant transformation direction and is expressed as:   =B B P (13) where B is the norm of B . For simplicity , it is assumed that the internal stress is a function of accumulated martensite volume fraction c (Yu et al., 2013 ): sat 1B 1 expcb   −=− B (14.a) c rev= (14.b) 1b is the parameter which controls the saturation rate of B . tr rev and det are the energetic contributions due to hardening and softening effects during cyclic loading for transformati on and detwinning, respectively. The rates of tr rev and det are given as: 24 1X tr rev rev  == (15.a) 24 det 1Xdet   == (15.b) where X and detX are the transformation and detwinning resistance s of the -th martensite variant, respectively. 9 As transformation hardening occurs during cyclic deformation, X is considered as a function of the reversible martensite volume fraction. XrevH  = (16.a) ()00 21 expc sat H H H Hb   = + − − −  (16.b) where H is a hardening modulus and a function of accumulated martensite volume fraction c . 0H and satH are initial and saturation values of H , respectively. The parameter 2b controls the evolution rate of H during cyclic loading. The h ardening equation for detwinning resistance is based on (Yu et al., 2014 ) and is expressed as: det detX h= (17) where deth is the hardening modulus for detwinning. Using the 2nd law of thermodynamics and the Helmholtz free energy, the Clausius dissipation inequality can be written as: .:0TTT= − − − qσε (18) By the assumption of quasi -static loading, the temperature in a single crystal can be considered to be homogeneous. Therefore, neglecting the effect of transformation latent heat and the heat generated from internal dissipation, the Clausius dissipation in equality is expressed as: :0 T= − − σε (19) Using equations (1), (9-11), and (15) yield s: () ()() ( )24 24 det 11: : : : : ,,: : X X 0 e tr TRIP det DIP ee ch ch e e tr TRIP det DIP revTTTT T              === + + + +  − − −   + + + + − − −  σ ε σ ε σ ε σ ε σ ε ε εε Bε ε ε ε (20) 10 Considering the special case when 0rev ire TRIP DIP        = = = = = and using the methodology of Coleman and Noll (1974 ), since eε and T can assume either positive or negative values, Eq. (20) can only be satisfied thro ugh the requirement that: (): ee e e ==ε σ C εε (21.a) () 0,lnchT TcTT  =− = −   (21.b) The dissipation inequality can be re -written as: () ()() ()24 24 24 1 1 1 24 24 24 0 1 1 1 24 24 det 11: : + , X X 0 tr tr rev ire TRIP det det DIP det ch revgg T                                    = = = = = = === + + + + + − + − − −       σ B P P P σ B P P P (22) Since reversible and residual martensite transformation, TRIP, detwinning , and DIP are all dissipative phenomena, the strong form requires that each part is non-negative such that: ()() ( )()24 ,0 1,: X 0ch tr tr rev det revTg         =  = + + − − −  σ B P P (23.a) ()() ( )()24 ,0 1,: 0ch tr tr ire det ireTg        =  = + + − − σ B P P (23.b) ()24 1: 0TRIP TRIP  = = + σ B P (23.c) ()24 det 1: X 0det det        = = + − σ B P (23.d) ()24 1: 0DIP det DIP  = = + σ B P (23.e) It is assumed that the dissipation inequality is non -negative for all variants . Therefore, each dissipation part in Eq. (23) can be divided into 24 components, i.e., 11 ()() ( ) ()0 ,: 0 =1,2,...,24,Xtr det chtr rev revg T       + + − = −−σ B P P (24.a) ()() ( )() ,0,: 0 =1,2,...,24ch tr tr ire det ireTg           = + + − − σ B P P (24.b) (): 0 =1,2,...,24TRIP TRIP    = + σ B P (24.c) ()det : X 0 =1,2,...,24det det          = + − σ B P (24.d) (): 0 =1,2,...,24DIP det DIP    = + σ B P (24.e) In order to satisfy the constraint of thermodynamics, it is assumed that: ()( ) sign( )=sign :det +σ B P (25) Thus, the thermodynamics drivin g forces for internal variables can be expressed by: ()() ( )(), 0 0 : X =1,2,...,24tr tr rev det g T T         = + + − − − −σ B P P (26.a ) ()() ( )(), 0 0 : =1,2,...,24 tr tr ire det g T T        = + + − − −σ B P P (26.b) (): =1,2,...,24TRIP=+σ B P (26.c) ()det : X =1,2,...,24det det     = + −σ B P (26.d) (): =1,2,...,24 DIP det=+σ B P (26.e) The work conjugate quantities associated with the thermodynamic driving forces are listed below : ,tr rev rev (27.a) ,tr ire ire (27.b) TRIP TRIP (27.c) det (27.d) DIP DIP (27.e) 2.2. Evolution equations for internal variables 12 2.2.1. Evolution equations for reversible phase transformation Considering the Clausius inequality for reversible martensite volume fraction, the product of rev and ,tr rev should always be non -negative. During forward tr ansformation 0rev , and for Eq. (24.a) to be satisfied, ,tr rev should be positive. If reverse transformation occurs, 0rev and ,tr rev should be negative to satisfy the 2nd law of thermodynamic s. 0rev= during elastic loading and unloading without any martensite transformation , so that Eq. (24.a) is automatically satisfied. Using the concept of transformation functions (Lagoudas and Entchev, 2004 ; Patoor et al., 2006 ), the evolution equation s for the reversible part of martensite volume fraction are proposed here in the form : ( ) ( ), , , , if 0, 1, 1 if 0, , 0 other conditionstr trn tr rev rev tr rev n tr rev rev tr rev ire ire tr revY K                    =      =−     =   (28.a) (28.b) (28.c) where trK is constant , and Y is a positive variable which evolves during cyclic loading. The rate sensitivity of martensite transformation is controlled by trn ; a high value of trn enforces a rate- independent response. The parameter Y controls the width of hysteresis loop , and is dependent on martensite volume fraction : 0 d rev Y Y Y =+ (29) where 0Y is the initial value of Y , and dY describes the dependence of Y on the current martensite volume fraction. Based on (Lagoudas and Entchev, 2004 ), the balance temperature 0T can be calculated from: 00 sYTM=+ (30) where 0sM is the martensite start temperature under the stress -free condition. 13 2.2.2. Evolution equations for transformation -induced plasticity Since transformation -induced plasticity is an irreversible dissipation process, TRIP is always non - negative, i.e., 0TRIP . Considering Eq. (24.c) , the constraint for TRIP can be described as: 0 0 0 no constraintTRIP TRIP TRIP  = (31.a) (31.b) Using E q. (31), an evolution equation for TRIP (Yu et al., 2013 ) can be expressed as: ( )()()() TRIP rev TRIP TRIP TRIP rev TRIP rev TRIP TRIPGT HD            − = + − (32.a) ()()n TRIP rev revG  = () TRIP TRIP TRIP TRIPTd  = (32.b) (32.c) where x and () Hx are McCauley bracket and step function, respectively: when 0x , xx= and ()1 Hx= ; when 0x , ()0 x H x== . () TRIP revG is the generalized force and a function of reversible martensite volume fraction. () TRIP TRIPT is the resistance of the -th friction system. The parameter TRIPD is the reference resistance. The parameter TRIP describes the difference between evolution features of TRIP in forward and reverse transformation. As shown in Eqs. (32. b) and (32. c), the parameters () TRIP revG and () TRIP TRIPT depend on reversible volume fraction of the -th martensite variant rev and TRIP , respectively. TRIPd is a material parameter. () TRIP H is used to ensure that TRIP is always positive and the 2nd law of thermodynamics is satisfied. 2.2.3. Evolution equations for irreversible phase transformation Since the formation of residual martensite is an irreversible dissipative process, the condition for irreversible martensite volume fraction to satisfy the 2nd law of thermodynamic s is: 14 , 0 0 0 no constraintire tr ire ire  = (33.a) (33.b) Using Eqs. (33.a) and (33.b) , an evolution law for residual martensite is expressed as (Yu et al., 2013 ): ( )()()(),ire rev ire ire ire rev ire rev tr ire ireGT HD            − = + − (34.a) ()()n ire rev revG  = (34.b) () ire ire ire ireTd  = (34.c) where () ire revG is the generalized force and () ire ireT is the resistance of residual martensite. ireD is a reference resistance and ire is a material constant. As shown in Eqs. (34.b) and (34.c), () ire revG and () ire ireT are functions of reversible and irreversible parts of martensite variants, respectively. ired is a material parameter. 2.2.4. Evolution equations for detwinning Consider ing Eq s. (26.d) and (27.d) and using a power -law form , an evolution equation for detwinning of the -th martensite variant is expressed as: ()( ) sign : if 0 0 other conditionsdetn det det detK      = +   =σ B P (35.a) (35.b) where detK is the detwinning resistance and detn represents the rate sensitivity of detwinning. Although the rate dependent formulation is used, adopting a large value for detn , a rate independent response can be obtained. It should be noted that , depending on the sign of the term ():det+σ B P , the volume fraction of sub -variants can decrease or increase during loading . 2.2.5. Evolution equations for detwinning -induced plasticity 15 Since DIP is an irreversible dissipative process, DIP will always be non -negative, i.e., 0DIP . Therefore, the condition for DIP to satisfy the 2nd law of thermodynamics is: 0 if 0 0 no constraintDIP DIP DIP  = (36.a) (36.b) According to the experimental observations (Song et al., 2014 ), the residual and peak strains accumulate during cyclic loading of NiTi and will be saturated after some cycles. Using this fact and Eqs. (24.e), (36.a) , and (36.b) , an evolution equation for DIP is proposed here as : () 3sat DIP DIP DIP rev DIP bH       =− (37) where 3b is a parameter controlling the saturation rate of detwinning -induced plasticity , sat DIP is the saturat ion value of DIP . Consid ering the dissipation inequality for DIP , it is obvious that DIP and DIP should have the same sign. () DIP H ensures that the 2nd law of thermodynamics is satisfie d, so that , DIP >0 if 0DIP , and 0DIP= if 0DIP . According to Eq. (26.e), as most of the components of detP are negative, a positive value of DIP is more probable in compressive loading , whereas a negative value is more probable in tensile loading. Therefore, using () DIP H in Eq. (37) causes DIP , and subsequently the DIP strain DIPε , to be zero in tension. A summary of the proposed constitutive model at the single crystal scale is presented in Table 1. Table 1 : Summary of the proposed constitutive model e in=+ε ε ε in tr TRIP det DIP= + + +ε ε ε ε ε 24 24 11tr tr tr rev iregg     ===+ε P P 24 1TRIP TRIP  ==εP ()24 0 1det det       ==−εP 24 1DIP DIP det  ==εP Strain decomposition ()() ( )(), 0 0 :Xtr tr rev det g T T        = + + − − − −σ B P P Thermodynamic deriving forces 16 ()() ( )(), 0 0 :tr tr ire det g T T       = + + − − −σ B P P ():TRIP=+σ B P ()det :Xdet det   = + −σ B P ():DIP det=+σ B P ( ) ( ), , , , if 0, 1, 1 if 0, , 0 other conditionstr trn tr rev rev tr rev n tr rev rev tr rev ire ire tr revY K                   =    =−    =         Evolution equations for reversible phase transformation ( )()()(),ire rev ire ire ire rev rev tr ire ireireGT H D            − = + − Evolution equation for residual martensite ( )()()() TRIP rev TRIP TRIP TRIP rev rev TRIP TRIPTRIPGT H D            − = + − Evolution equation for transformation -induced plasticity ()() sign : if 0 0 other conditionsdetn det det detK     = +  =    σ B P Evolution equations for detwinning ()3sat DIP DIP DIP rev DIP bH       =− Evolution equation for detwinning -induced plasticity 2.3. Explicit scale -transition rule To transform the single crystal response into polycrystalline version, an explicit scale -transition rule (Yu et al., 2013 ) is used in this work. In comparison to implicit self -consistent modelling approach proposed by Hill (1965 ), the explicit scale -transition approach proposed by Kröner (1961 ) leads to a convenient “a priori” interaction law that facilitates easy application and comparison with experimental results. In Kröner’s model, each crystal of the polycrystalline material is tr eated as an inclusion inside a matrix which is made from other crystals of the material. This approach is based on the assumption of elastic interaction between the matrix and the inclusion. An extension of this model by Berveil ler and Zaoui (1978 ) allows for elasto -plastic intergranular accommodation and plastic internal stress relaxation. The method of Berveiller and Zaoui has previously been used to describe ratchetting of polycrystalline materials (Cailletaud and Sai, 2008 ; Kang et al., 2011 ). It has also been adopted to describe the thermo -mechanical properties 17 of NiTi shape memory a lloys (SMAs) (Yu et al., 2014 ; Yu et al., 2013 ; Yu et al., 2012 , 2015 ). With the assumption of uniform isotropic elasticity for the polycrystalline aggregates, the relation between stress tensor in a single crystal σ and the applied uniform macroscopic stress tensor Σ can be ex pressed by: ()in inD= + −σ Σ E ε (38) where in in=Eε is the inelastic strain of polycrystalline aggregates obtained by the volume average of each single crystal strain tensor . The symbol []• denotes the volume average over the entire grains in a polycrystalline aggregate . D is a material parameter reflecting the heter ogeneity of stress field at the inter -granular scale. The value of D can be set to a value between zero and () ()4 4 5 15 1G − − where G is the shear modulus of the material and  is the Poisson’s ratio . If 0D= the heterogeneity at the intergranular scale is neglected and a uniform stress field is obtained in the entire polycrystalline aggregates. As the value of D is increased , the heterogeneity of the polycrystalline material is increase d. 3. Material parameter s identification Since the response of NiTi shape memory alloy is rate -independent and viscosity is not significant , the parameters trn and detn are set to be sufficiently large so that material behavior is effect ively rate independent. The parameter  controls the start -stress of forward transformation at a specific temperature and can be calibrated by fitting the curves related to the start -stress of transformation at different temperatures. Since the elastic constants for monoclinic single crystal martensite in NiTi SMA are not available, we assume the elastic modulus to be isotropic. As we are focusing on inelastic response of shape memory alloys, the elastic modulus is considered the same for austenite and martensite phases. The Poisson’s ratio is based on the study of Thamburaja (2005 ). Two different parameters, Y and trK are adopted to describe the evolution law of reversible martensite volume fraction during forward and reverse transformation. Both parameters control the width of hysteresis loop. Parameters 0Y and dY control the width of hysteresis loop during forward transformation in the first and steady -state cycles, respectively. The role of trK is to control the width of hysteresis loop during reverse transformation . The parameters deth and detK 18 can be obtained by fitting the curve s related to detwinning modulus and start -stress of detwinning, respectively. As commented by Liu et al. (1998 ), the effect of DIP is not obvious during tensile loading. On the other hand, Kang et al. (2009 ) commented that the accumulated peak strain during uniaxial cyclic deformation is mainly caused by transformation -induced plasticity while the accumulated residual strain comes from both the transformation -induced plasticity and residual martensite. Thus, the parameters n , DTRIP , dTRIP , and TRIP can be obtained from the data of accumulated peak strain in tensile cyclic loading. The parameters dire , Dire , and ire can be measured by fitting to the data related to the difference of peak and residual strains in tensile cyclic loading. The p arameter 3b controls the saturation rate of detwinning -induced plasticity and sat DIP is the saturated value of detwinning -induced plasticity . These parameters are calibrated by fitting the evolution curves of residual strain under compressive load ing. The parameters 0H and satH control transformation hardening and can be calibrated by measuring transformation modulus at the first and steady -state cycles, respectively. The fact that the start-stress of martensite transformation decreases with increasing number of cycles is taken into account by using satB . This p arameter can be obtained by measuring the start -stress of martensite transformation in the first and steady -state cycles. Parameters 1b and 2b control the saturation rates of B and H , respectively. They can be calibrated by fitting the evolution curves of decreased start -stress of martensite transformation and hardening modulus, respectively. The parameter trg controls the maximum transformation strain and can be calibrated from the measured value of transformation strain. The value of the parameter D in the scale -transition rule is set to be the same as that used by Yu et al. (2015 ). It should be noted that since all parameters are related to a single crystal, they should ideally be obtained from stress -strain curves of a single crystal. As the experimental data for single crystal is not available, the parameters are calibrated by a tr ial-and-error method from polycrystalline experimental data. In order to explore the predictive capability of the proposed constitutive model, two sets of experimental test data are considered (Kang et al., 2009 ; Yu et al., 2013 ). Both experimental tests entail cyclic loading of super -elastic NiTi SMA with 50.2 at.% Ni. The finish -temperature of austenite for this alloy is 268 K. The original phase of the material is austenite, given that both tests are performed at room tempera ture. In the experiments of Kang et al. (2009) stress -controlled 19 cyclic loading regimes are applied to the SMA, including tension -compression tests and tension - unloading tests. In the experiments of Yu et al. (2013) strain -controlled cyclic loading regimes are applied to the SMA, including compression -unloading tests and tension -unloading tests. The stress -controlled cyclic tests, up to 160 cycles, are related to tension -compression tests with two mean stresses of 100 MPa and 0 MPa, and stress amplitudes of ±450 MPa, as well as tension - unloading tests with peak stress of 550 MPa. The strain -controlled cyclic tests, up to 50 cycles, are related to compression -unloading tests with two minimum strains of -2.5% and -3.5%, and cyclic tension -unloading tests with peak strain of 3.5%. The calibrated parameters used in this study are listed in Table 2. Table 2: Material parameters used in the prediction of cyclic loading considering DIP Parameter Symbol Value Unit Young's modulus E 72 GPa A co nstant controlling the start -stress of forward transformation  0.733 MJ K Poisson’s ratio  0.3 - The martensite start -temperature 0sM 265 K Constants controlling the rate sensitivity of martensite transformation and detwinning , respectively trn, detn 30 , 30 - Saturation value of internal stress satB 50 MPa The initial value of Y, a constant controlling the dependence of Y on the cur rent martensite volume fractio n 0Y, dY 4 , 5 MPa A constant controlling the evolution law of martensite volume fraction during reverse transformation trK 25 MPa Initial and saturation values of hardening modulus 0H, satH 80 , 400 MPa A constant controlling the dependence of TRIP and residual martensite on the applied stress n 3 - Reference resistance for TRIP , a constant controlling the resistance of the -th friction system DTRIP, dTRIP 3.5e-3 , 5e-3 - Reference resistance for residual martensite, a constant controlling the resistance of residual martensite Dire, dire 1e-3 , 1e-3 - Constant s describing the difference between evolution features in forward and reverse transformation s TRIP, ire 0.2 - 20 Constants controlling saturation rate of internal stress, evolution rate of hardening modulus, and saturation rate of DIP, respectively 1b, 2b , 3b 2 , 2, 1 - Hardening modulus for detwinning, detwinning resistance deth, detK 1000, 20 MPa The saturation value of shearing deformation caused by DIP sat DIP 0.05 - A constant reflecting the heterogeneity of stress field in the inter-granular scale D 2 GPa 4. Results and discussion In this section we first establish the role of texture on tension -compression asymmetry in our modeling framework. We then investigate the role of DIP in cyclic tension and compression loading conditions. 4.1. Influence of texture on asymmetry in a singl e loading cycle Depending on manufacturing process, polycrystalline shape memory alloys can have a strong texture (Thamburaja and Anand, 2001 ). Since the samples used in the experiments are manufactured from extruded solid bars, the texture formed in the process of extrusion will affect the polycrystalline response of the material. Based on the method used for forming NiTi shape memory alloys, three different textures with different effects on the overall response of the allo y can be formed : a 100 texture with stiffer mechanical response in tension than in compression, a 111 texture with softer response in tension than in compression , and a 110 texture with almost no effect on the mechanical responses i n tension or compression (Gall et al., 2001 ). In order to show the effect of t exture on asymmetrical response of shape memory alloys, 40 single crystals in the 3 aforementioned textures are considered. Employing rotation matrix in three dimensions, Euler angles for each of 40 crystals are selected such that a specific direction of a ll grains is aligned with the loading direction. Predicted difference s between stress -strain curves in tension and compression for different textures are shown in Fig. 1. For randomly oriented grains, grains with 100 texture, and grains with 110 texture, more hardening occurs in tension than in 21 compression , although the difference is not pronounced for the 110 and random textures. In contrast, the 111 texture results in more hardening in compression than in tension . Fig. 1. Stress -strain curves of NiTi shape memory alloy subjected to strain -controlled tension -unloading for different textures: (a) [100] texture, (b) [111] texture, (c) [110] texture, (d) random texture. 4.2. Influence of DIP in cyclic loading 22 Based on the asymmetrical tension /compression response s observed in the experimental tests of Kang et al. (2009 ) and Yu et al. (2013 ), 40 single crystals with 111 texture are adopted to model the behavior of NiTi SMA in this work . To model the 111 texture, the rotation matrix between local coordinate system at a single crystal and global coordinate system at a polycrystal line material is calculated such that the 111 direction of each grain is aligned with the loading directio n. Before presenting comparison s of model predictions with reported experimental data , we first present a sensitivity analysis of the effect of parameters 3b (controlling the saturation rate of detwinning -induced plasticity) and sat DIP (the saturation value of shearing deformation caused by detwinning -induced plasticity) on the computed stress -strain response , as shown in Fig. 2 . Computed s tress-strain loops for the 20th cycle of stress -controlled tension -compression loading of a single crystal NiTi are shown. As most of the components of detP have negative values, a decrease of either 3b or sat DIP results in less negative values of DIP. Therefore, the computed strain throughout the loading cycle increase s. Fig. 2. Effect of parameters (a) 3b (controlling the saturation rate of DIP ) and (b) sat DIP (the saturation value of shearing deformation caused by DIP ) on stress -strain curves of a single crystal in the 20th cycle for a loading case with a mean stress of 100 MPa and a stress amplitude of 450 MPa. 23 Comparisons between m odel predictions and the stress -controlled cyclic loading experiments of Kang et al. (2009 ) are next presented. Model parameters are listed in Table 2. In Fig . 3 we consider the case of applied mean stresses of 100 MPa, with stress amplitudes of 450 MPa . In the first cycle, DIP affects the part of the cycle with negative stress very slightly, but this effect becomes more pronounced in subsequent cycles . 24 Fig. 3. Experimental and predicted stress -strain curves for stress -controlled cyclic tension -compression tests with a mean stress of 100 MPa and a stress amplitude of 450 MPa ( 100 450 MPa): (a) 1st cycle, (b) 5th cycle, (c) 10th cycle, (d) 20th cycle, (e) 50th cycle, (f) 160th cycle. The corresponding curves of residual strain, peak strain , and dissipation energy as a function of number of cycles are shown in Fig. 4. In stress -controlled cyclic loadings DIP causes a reduction in both peak strain and residual strain. Fig. 4. Experimental and predicted curves for stress -controlled cyclic tension -compression tests with a mean stress of 100 MPa and a stress amplitu de of 450 MPa: (a) peak strain vs. number of cycles, (b) residual strain vs. number of cycles. 25 Predicted and experimental stress -strain curves for the loading case with zero mean stress and a stress amplitude of 450 MPa are shown in Fig. 5. As the density of dislocations increases during cyclic loading, the effect of DIP becomes more pronounced with increasing cycles . Fig. 6 shows a reduction in both residual and peak strains due to DIP. 26 Fig. 5. Experimental and predicted stress -strain curves for stress -controlled cyclic tension -compression tests with zero mean stress and a stress amplitude of 450 MPa ( 0 450 MPa): (a) 1st cycle, (b) 5th cycle, (c) 10th cycle, (d) 5 0th cycle, (e) 160th cycle. 27 Fig. 6. Experimental and predicted curves for stress -controlled cyclic tension -compression tests with zero mean stress and a stress amplitude of 450 MPa: (a) peak strain vs. number of cycles, (b) residual strain vs. number of cycles. Fig. 7 shows the predicted and experimental curves for stress -controlled tensile -unloading cycles, where the applied stress cycles between a maximum value of 550 MPa and a minimum value of 0 MPa . As shown in Fig. 7, DIP does not contribute if compressive stress is not generated during part, or all, of the loading cycle . 28 Fig. 7. Experimental and predicted curves for stress -controlled cyclic tension -unloading tests, where the applied stress cycles between a minimum value of 0 MPa and a maximum value of 550 MPa : (a) simulated stress -strain curves, where dashed and solid lines represent simulation without and with consid ering DIP, respectively, (b) experimental stress -strain curves, (c) peak strain vs. number of cycles, (d) residual strain vs. number of cycles, (e) dissipation energy vs. number of cycles. Model predictions are next presented for the strain -controlled exp eriments of Yu et al. (2013). The set of model parameters listed in Table 2 is once again used, as the same poly crystalline NiTi SMA alloy used in the strain -controlled experiments of Yu et al. (2013) and the stress -controlled experiments of Kang et al. (2 009) presented above. Fig. 8 shows stress -strain curves for cyclic compression -unloading with a maximum applied strain of 0% and a minimum applied strain of - 29 2.5%. The effect of DIP is not significant in the 1st cycle but it becomes significant from 2nd cycle onwards. The incorporation of DIP provides improved predictions of the experimental results of Yu et al. (2013) , with lower peak compressive stress and higher residual compressive strain . 30 Fig. 8. Experimental and predicted stress -strain curves for strain -controlled cyclic compression -unloading tests where the applied strain cycles between a minimum value of -2.5% and a maximum value of 0% : (a) 1st cycle, (b) 2nd cycle, (c) 5th cycle, (d) 10th cycle, (e) 20th cycle, (f) 50th cycle. Dissipated energy over 50 loading cycles is shown in Fig. 9(a). Again, the inclusion of DIP in the model provides an improved prediction of the experimental data. DIP also provides an improved predict ion of the evolution of peak compressive stress over 50 loading cycles, as shown in Fig. 9(b). 31 Fig. 9. Experimental and predicted curves for strain -controlled cyclic compression -unloading tests where the applied strain cycles between a minimum value of -2.5% and a maximum value of 0%: (a) dissipation energy vs. number of cycles, (b) peak stress vs. numbe r of cycles . Fig. 10 shows stress -strain curves for cyclic compression -unloading with a maximum applied strain of 0% and a minimum applied strain of -3.5% . Again, the effect of DIP is not significant in the 1st cycle but it becomes significant from 2nd cycle onwards. As shown in Fig . 11, similar to the case of a minimum applied strain of -2.5%, DIP also provides improved predictions of experimentally measured dissipated energy and peak stress for the higher appl ied compressive strain of -3.5%. 32 Fig. 10. Experimental and predicted stress -strain curves for strain -controlled cyclic compression -unloading tests where the applied strain cycles between a minimum value of -3.5% and a maximum value of 0%: (a) 1st cycle, (b) 2nd cycle, (c) 5th cycle, (d) 10th cycle, (e) 20th cycle, (f) 50th cycle. 33 Fig. 11. Experimental and predicted curves for strain -controlled cyclic compression -unloa ding tests where the applied strain cycles between a minimum value of -3.5% and a maximum value of 0%: (a) dissipation energy vs. number of cycles, (b) peak stress vs. number of cycles . Fig. 12 shows predicted and experimentally measured stress -strain curves for cyclic tension - unloading with a maximum applied strain of 3.5% and a minimum applied strain of 0%. DIP is not predicted to contribute to the material response, due to the absence of compressive stress throug h the loading history. Predict ed stress -strain curves and predicted evolution of dissipated energy, residual strain and peak tensile stress are shown to follow a similar trend to the experimental measurements of Yu et al. (2013). 34 35 Fig. 12. Experimental and predicted stress -strain curves for strain -controlled cyclic tension -unloading tests where the applied strain cycles between a minimum value of 0% and a maximum value of 3.5% : (a) simulated stress -strain curves, where dashed and solid lin es represent simulation without and with considering DIP, respectively, (b) experimental stress -strain curves , (c) peak stress vs. number of cycles, (d) residual strain vs. number of cycles, (e) dissipation energy vs. number of cycles. Based on the experimental work done by Liu et al. (1998 ), the effect of DIP is significant under compressive loading while in tension no significant difference is observed with and without considering DIP . The model results presented in this section demonstrate that the inclusion of DIP provides improved predictions of experimentally observed r esponse of NiTi SMAs to cyclic stress - and strain - controlled loading where compressive stress is generated for part, or all, of the loading history. 5. Conclusion s Based on thermodynamic considerations a new inelastic mechanism, detwinning -induced plasticity (DIP), is proposed for modelling the response of NiTi SMAs to cyclic loading . DIP is incorporated into a constitutive framework for NiTi SMA. The constitutive framework also includes well -established inelastic mechanisms of ph ase transformation, transformation -induced plasticity, residual martensite, and detwinning . The model is constructed at the single crystal scale using the framework of thermodynamics and a crystal plasticity formulation. An explicit scale - 36 transition rule is used to extend the framework for the analysis of polycrystalline materials, thus facilitating simulation of published experimental test data. Thermodynamic considerations are used to motivate the contribution of DIP for cyclic loading regimes whe re compressive stress occurs during part of the loading cycle. However, the contribution of DIP is negligible for cyclic loading regimes that result exclusively in tensile stress . This thermodynamically motivated dependence of DIP on compression, but not o n tension, is strongly supported by experimental cyclic loading results of Liu et al. (1998 ). Inclusion of DIP in a modelling framework results in improved prediction of experimentally observed NiTi SMAs behavior , including strain -controlled cyclic compression -unloading and cyclic tension -unloading tests by Yu et al. (2013 ) and stress -controlled cyclic tension -compression and tension -unloading tests by Kang et al. (2009 ). During the first loading cycle the contribution of DIP is not significant in any loading regime. However, in cases where compressive stress occurs during part of the loading cycle, DIP contributes strongly to the material response from the second cycle onwards . In strain -contro lled cyclic compression -unloading tests DIP leads to a less negative peak stress and a more negative residual strain following several loading cycles . In stress -controlled tension - compression loading , DIP leads to a reduction of peak and residual strains following loading cycles . Future applications of the model should consider non -proportional loading and the effects of reorientation and reorientation -induced plasticity (Yu et al., 2015). Ongoing developments of this modelling framework entail the implementation of DIP in a finite element crystal -plastici ty modelling framework (McCarthy et al., 2014 ; McGarry et al., 2007 ), facilitating the simulation of the mechanical behavior of NiTi SMA stents under long -term cyclic phy siological loading regimes. 37 Appendix Table A.1: Crystallographic data for the 24 martensite variants in NiTi shape memory alloy (Wang et al., 2008 )  1n 2n 3n 1m 2m 3m 1 -0.8889 -0.4044 0.2152 0.4114 -0.4981 0.7633 2 -0.4044 -0.8889 -0.2152 -0.4981 0.4114 -0.7633 3 0.8889 0.4044 0.2152 -0.4114 0.4981 0.7633 4 0.4044 0.8889 -0.2152 0.4981 -0.4114 -0.7633 5 -0.8889 0.4044 -0.2152 0.4114 0.4981 -0.7633 6 0.4044 -0.8889 0.2152 0.4981 0.4114 0.7633 7 0.8889 -0.4044 -0.2152 -0.4114 -0.4981 -0.7633 8 -0.4044 0.8889 0.2152 -0.4981 -0.4114 0.7633 9 0.2152 0.8889 0.4044 0.7633 -0.4114 0.4981 10 0.2152 -0.8889 -0.4044 0.7633 0.4114 -0.4981 11 -0.2152 -0.4044 -0.8889 -0.7633 -0.4981 0.4114 12 -0.2152 0.4044 0.8889 -0.7633 0.4981 -0.4114 13 -0.2152 0.8889 -0.4044 -0.7633 -0.4114 -0.4981 14 -0.2152 -0.8889 0.4044 -0.7633 0.4114 0.4981 15 0.2152 0.4044 -0.8889 0.7633 0.4981 0.4114 16 0.2152 -0.4044 0.8889 0.7633 -0.4981 -0.4114 17 0.8889 -0.2152 0.4044 -0.4114 -0.7633 0.4981 18 -0.8889 -0.2152 -0.4044 0.4114 -0.7633 -0.4981 19 0.4044 0.2152 0.8889 0.4981 0.7633 -0.4114 20 -0.4044 0.2152 -0.8889 -0.4981 0.7633 0.4114 21 0.8889 0.2152 -0.4044 -0.4114 0.7633 -0.4981 22 -0.8889 0.2152 0.4044 0.4114 0.7633 0.4981 23 -0.4044 -0.2152 0.8889 -0.4981 -0.7633 -0.4114 24 0.4044 -0.2152 -0.8889 0.4981 -0.7633 0.4114 38 Table A.2: Crystallographic data for martensite detwinning systems in NiTi shape memory alloy (Thamburaja et al., 2005 )  1a 2a 3a 1w 2w 3w 0 1 -0.5846 0 0.8133 0 0.2804 0 0.271 2 0 -0.5846 -0.8133 0.2804 0 0 0.271 3 -0.5846 0 -0.8133 0 0.2804 0 0.271 4 0 -0.5846 0.8133 0.2804 0 0 0.271 5 -0.5846 0 -0.8133 0 0.2804 0 0.729 6 0 -0.5846 0.8133 0.2804 0 0 0.729 7 -0.5846 0 0.8133 0 0.2804 0 0.729 8 0 -0.5846 -0.8133 0.2804 0 0 0.729 9 -0.8133 -0.5846 0 0 0 0.2804 0.271 10 -0.8133 0.5846 0 0 0 0.2804 0.729 11 -0.8133 0 -0.5846 0 0.2804 0 0.271 12 -0.8133 0 0.5846 0 0.2804 0 0.729 13 -0.8133 0.5846 0 0 0 0.2804 0.271 14 -0.8133 -0.5846 0 0 0 0.2804 0.729 15 -0.8133 0 0.5846 0 0.2804 0 0.271 16 -0.8133 0 -0.5846 0 0.2804 0 0.729 17 0.5846 -0.8133 0 0 0 0.2804 0.729 18 -0.5846 -0.8133 0 0 0 0.2804 0.271 19 0 -0.8133 -0.5846 0.2804 0 0 0.271 20 0 -0.8133 0.5846 0.2804 0 0 0.729 21 -0.5846 -0.8133 0 0 0 0.2804 0.729 22 0.5846 -0.8133 0 0 0 0.2804 0.271 23 0 -0.8133 0.5846 0.2804 0 0 0.271 24 0 -0.8133 -0.5846 0.2804 0 0 0.729 39