In this study the mechanisms of transformation -induced plasticity (TRIP) , detwinning -induced plasticity (DIP) , and accumulation of residual martensite, are incorporated into a finite - deformation crystal plasticity model of NiTi SMA for the first time. The constitutive model is constructed at the single -crystal scale and also includes phase transformation and detwinning mechanisms. Using a proposed Helmholtz free energy , the driving forces for inel astic mechanisms are derived within the framework of thermodynamics . The constitutive model has been implemented in the Abaqus/Explicit finite -element program, using VUMAT subroutine to simulate a polycrystalline material . Considering various orientations for crystals, the effect of texture on tension -compression asymmetry is investigated. It is shown that different textures may cause stiffer , softer, or similar response in compression compared to tension. Due to the incorporation of the effect of residual martensite, the model provides accurate predictions of experimentally measured residual strain. The incorporation of the afo rementioned inelastic deformation mechanisms is shown to accurately capture the k ey features related to cyclic loading . Finally, the effect of DIP in compressive cyclic loading of NiTi SMA is investigated. Conside ring the effect of DIP causes a less negative peak stress and a more negative residual strain to be observed i n the strain -controlled cyclic compression -unloading tests. Keywords: NiTi; Large deformation; Detwining -induced plasticity; Constitutive modeling; Crystal plasticity * Corresponding author s: patrick.mcgarry@nuigalway.ie , naghdabd @sharif.edu 1. Introduction Shape memory alloys (SMAs), due to unique characteristics such as super -elasticity and shape memory effect, are used for a range of applications , including biomedical devices (Morgan, 2004; Petrini and Migliavacca, 2011) , robotics (Coral et al., 2018; Furuya and Shimada, 1991) , aerospace (Hartl and Lagoudas, 2007) , actuators (Benafan et al., 2019) , and automotive structures (Williams and Elahinia, 2008) . NiTi is the most widely used SMA due to its high energy density, corrosion resistance, and biocompatibility (Hartl and Lagoudas, 2007; Morgan, 2004) . The behavior of NiTi SMA is based on the temperature -dependent reversible phase transformation from the high symmetry austenite phase to the low symmetry martensite phase and vice versa (Lagoudas, 2008) . At temperatures higher than the austenite fini sh-temperature, application of stress to austenite will cause strain generation, which is recoverable upon unloading (super -elasticity). Upon cooling, in the absence of stress, the high temperature austenite phase will transform to the low temperature mart ensite phase (austenite to twinned - martensite transformation). Application of enough stress to twinned -martensite will cause strain generation during loading , which is not recovered upon unloading, and requires subsequent heating to be recovered (shape -mem ory effect). Many constitutive models have been developed to predict the behavior of SMAs based on experimental observations. The models can be categorize d into two main groups: macroscopic phenomenological models and microstructure -based models. Typical macroscopic constitutive models include Auricchio et al. (Auricchio et al., 2007) , Lagoudas et al. (Lagoudas et al., 2012) , Panico and Brinson (Panico and Brinson, 2007) , Zaki and Moumni (Zaki and Moumni, 200 7), Arghavani et al. (Arghavani et al., 2010) , and Xu et al. (Xu et al., 2019) . Whil e macroscopic models do not consider the microscale mechanisms underlying SMA behavior, they are readily implementable in finite element codes for structural analysis and design. Microstructure -based constitutive models have been proposed to incorporate the underlying microscale physical mechanisms underlying SMA behavior. The popularity of crystal plasticity -based models is due to considering various variants of martensite with dif ferent morphological features and their evolution during deformation . These models are developed at the single -crystal scale and then, adopting finite element method or scale -transition rules, transformed to the poly -crystal scale. Representative constitut ive models include Thamburaja and Anand (Thamburaja and A nand, 2003, 2002, 2001) , Anand and Gurtin (Anand and Gurtin, 2003) , Thamburaja (Thamburaja, 2005) , Thamburaja et al. (Thamb uraja et al., 2005) , Wang et al. (Wang et al., 2008) , Manchiraju and Anderson (Manchiraju and Anderson, 2010) , Yu et al. (Yu et al., 2015, 2014, 2013, 2012) , Xiao et al. (Xiao et al., 2018) , and Dhala et al. (Dhala et al., 2019) . It should be noted that some of the above -mentioned constitutive models do not consider the effect of large deformations in the modeling framework. In practical applications, the strains and rotations applied to SMA can be quite large such that the theories of small deformation cannot be used to reasonably predict the behavior of these materials. In recent years, the behavior of NiTi SMAs under finite deformations has been studied by some researchers. The effect of phase transformation on the mechani cal and thermomechanical response of NiTi SMA has been investigated by Thamburaja and Anand (Thamburaja and Anand, 2003, 2002, 2001) , and Anand and Gurtin (Anand and Gurtin, 2003) . The effects of detwinning and reorientation have been incorporated into constitutive modeling of NiTi SMA by Thamburaja (Thamburaja, 2005) . This model has been modified to consider austenite to martensite phase transformation, as well (Thamburaja et al., 2005) . Manchiraj u and Anderson (Manchiraju and Anderson, 2010) considered the effects of phase transformation and austenite plasticity, which results from phase transformation at high temperatures. It should be noted that the residual strain observed in this modeling is due to austenite plasticity, occurs only at sufficiently high temperatures, and differs from the mechanism associated with the accumulation of residual martensite, which can be activated at any temperature. The effects of austenite -plasticity and martensite -plasticity on the mechanical behavior of super -elastic NiTi have been investigated by Dhala et al. (Dhala et al., 2019) . During phase transformation, some martensite variants are locked, due to increase of defects, and cannot return to austenite phase. This process, called accumulation of residual martensite, is known to be one of the main mechanisms in the deformation of super -elastic NiTi SMA (Kang et al., 2009; Yu et al., 2013) . It should be noted that the effect of residual martensite has not been considered in the above - mentioned finite deformation constitutive models, while it is evident in experimental observations. On the other hand, some kind of local stress field and plasticity can be activated in SMAs which is different from conve ntional plasticity (Lagoudas, 2 008). The plastic deformation which is created to overcome the local stress field at the interface of asutenite and martensie phases is called transformation -induced plasticity (Yu et al., 2015) , whi ch has significant ef fect on predicting the behavior of NiTi SMA under cyclic loading (Lagoudas, 2008) . In addition , the experimental results of Liu et al. (Liu et al., 1998 ) show a high density of dislocations at the interfaces of lattice correspondent variants during compressive loading , but no significant plastic deformation is detected in such areas under tension. Detwinning - induced plasticity is referred to that part of plastic deformation which is induced to relax the local stress field at the interfaces of lattice correspondent variants (Ebrahimi et al., 2020) . As demonstrated by Ebrahimi et al. (Ebrahimi et al., 2020) , DIP has considerable effect on compressive cyclic loading of NiTi SMAs . It should be noted that the cyclic behavior of NiTi SMA under large deformations , considering the crucial mechanisms of TRIP and DIP, has not been investigated in previous studies . The goal of the current study is to predict the behavior of NiTi under mo notonic and cyclic loadings considering large deformations. Using crystal plasticity formulation, diverse mechanisms of phase transformation, transformation -induced plasticity, detwinning, detwinning -induced plasticity, and accumulation of residual martens ite have been incorporated into finite -deformation modeling of NiTi SMA. The single -crystal constitutive model has been implemented in the Abaqus/Explicit finite -element program, using VUMAT subroutine to obtain the polycrystalline response. The model is s hown to capture various features observed during monotonic and cyclic loading of super -elastic NiTi. 2. Microstructure -based constitutive model 2.1. Definition of various deformation mechanisms in the framework of thermodynamics Using multiplicative decomposition , the deformation gradient tensor F can be decomposed into elastic , eF, and inelastic , pF , parts (Lubarda, 2004; Thamburaja, 2005) : ep F F F (1) where pF is composed of several parts related to different mechanisms of phase transformation (tr), detwinning (det), transformation -induced plasticity (TRIP), and detwinning -induced plasticity (DIP). The velocity gradient tensor, L , can also be divided into elastic, eL , and inelastic, pL , parts (Dhala et al., 2019; Erdle and Böhlke, 2017; Fohrmeister et al., 2019) : 1 1 1 1 11( )( )e p e p p e e e e p p e L FF F F F F F F F F F F F F (2) 11, e e e p p p L F F L F F (3) The plastic velocity gradient tensor, pL , comprises various velocity gradient tensors related to phase transformation trL , transformation -induced plasticity TRIPL , detwinning detL , and detwinning -induced plasticity DIPL . Due to experimentally observed a ccumulation of residual martensite (Kang et al., 2009) , the volume fraction of each martensite variant can be divided into the reversible , rev , and the irreversible , ire, volume fractions . The total volume fraction of martensite is considered as the sum of eversible and irreversible parts (Ebrahimi et al., 2020) . rev ire    =+ (4) 24 24 24 1 1 1, , rev rev ire ire           = = == = =   (5) where  , rev , and ire are the total, reversible, and irreversible martensite volume fractions, respectively.  is the martensite volume fraction for the -th martensite variant. The orientation tensor P of the -th martensite variant can be expressed by (Anand and Gurtin, 2003; Thamburaja and Anand, 2001) : P b m (6) where b and m denote transformation direction and habit plane normal of the -th martensite variant , respectively. The magnitudes of theses parameters , related to different martensite variants, are presented by (Wang et al., 2008) . In the absence of stress , by reducing the temperature to lower than the martensite finish - temperature, austenite will transform to twinned martensite. Consequently, some habit plane variants (HPVs) will be formed. Each HPV is composed of two lattice correspondent variants (LCVs) or sub-variants (Ebrahimi et al., 2020; Yu et al., 2014) . The stress -induced movement of interfaces between sub-variants is called detwinning (Thamburaja, 2005) . The orientation tensor of the -th detwinning system, detP , is such that: det P a w (7) where a and w denote shear direction and unit vector normal to detwinning plane for the -th martensite variant , respectively. The magnitudes of theses parameters, related to different martensite variants are presented by (Thamburaja et al., 2005) . According to crystal plasticity (Dhala et al., 2019; Rengarajan an d Rajagopal, 2001; Thamburaja et al., 2005) , the inelastic velocity gradient tensor is proposed here such that : 1 24 24 24 24 24 det det 1 1 1 1 1rev ireP p p tr tr TRIP DIP gg L F F P P P P P (8) The terms on the right -hand side of Eq. ( 8) are related to reversible martensite phase transformation, residual martensite, transformation -induced plasticity, detwinning, and detwinning -induced plasticity, respectively. trg is the magnitud e of shearing deformation due to phase transformation, TRIP and DIP are the amount s of slip due to transformation -induced plasticity and detwinning -induced plasticity, respectively . For the -th martensite variant, t he volume fractions of two sub -variants are denoted by and 1 , and represents the rate change of sub-variants’ volume fractions . Since TRIP is the plastic deformation occurred to relax the local stress field at the austenite -martensite interfaces, its orientation tensor is assumed to be the same as phase transformation mechanism (Yu et al., 2015, 2013) . The Helmholtz free energy of a representative volume element (RVE) is considered to be composed of different parts related to elastic energy, e , chemical energy , ch , the part related to the influence of internal stress, int , and the part reflecting the effects of hardening or softening during phase transformation , trrev : int =trreve ch (9) The elastic part of free energy can be expressed by (Yu et al., 2015) : 1( , ) .( ( ) )2e e e e E E E (10) where () is the fourth -order stiffness tensor , considered as a function of the volume fraction of martensite according to the rule of mixtures, (1 )AM . eE is the Green - Lagrange elastic strain tensor. The chemical energy is used to account for the effect of temperature on the phase transformation (Thamburaja, 2005) : ( , ) ( )ch T T Th (11) T is the phase equilibrium temperature, at which the free energies of the austenite and martensite phases are equal. Th is the absorbed or released latent heat per unit reference volume at the temperature of T during phase transformation . During cyclic loading, some reduction in the start stress of phase transformation will be observed by increasing number of cycles (Kang et al., 2009) . It is concluded that local stress fields due to different inelastic mechanisms will assist the occurrence of phase transformation, and consequently lower the start stress of forward transformation in cyclic loading. Using this fact and considering the work of Yu et al. (Yu et al., 2015) , the form of int is proposed here such that: int:P BL (12) B is an internal variable, known as internal stress. trrev is used to account for the hardening or softening effects during phase transformation. Hardening or softening will be observed in the response of shape memory alloy by any increase or decrease in the amount of trrev , respectively (Yu et al., 2015) . 24 1tr rev rev X (13) X represents the resistance to phase transformation. Since internal stress is related to local stress field caused by different variants of martensite , it can be divided into 24 components (Yu et al., 2015) : 24 6 1, (1 exp( )) sat revBb B B P B (14) where B is the norm of B (the internal stress caused by the -th martensite variant). 6b is the controlling parameter for the saturation rate of B , satB denotes the saturation value of B . rev is the accumulated reversible martensite volume fraction , which is shown in the rate- form by (Yu et al., 2013) : =rev rev (15) The experimentally observed transformation hardening in cyclic loading of SMA (Kang et al., 2009) , is considered in the model by adopting X as a function of the reversible martensite volume fraction (Yu et al., 2015) . 0 0 4, ( )(1 exp( ))rev rev rev sat rev X H H H H H b (16) where 0H and satH are initial and saturation values of hardening modulus , revH , respectivel y. 4b is the control ling parameter for evolution rate of revH . The Clausius inequality , in finite -deformation regimes, can be expressed by (Anand and Gurtin, 2003; Thamburaja, 2005) : 01. : 0 q S F (17) where S is the first Piola –Kirchhoff stress tensor. Considering quasi -static loading and neglecting the transformation latent heat and the heat from internal dissipation, the third term in Eq. ( 17) can be neglected. Therefore , using the equality **: : ( ) :e e P S F T E C T L , Eq. (17) can be re -written as: **: ( ) : 0e e P T E C T L (18) where 1*(det )Tee T F F TF (Thamburaja, 2005) is the second Piola –Kirchhoff stress tensor associated with the elastic part , and T is the Cauchy stress tensor. Using Eqs. ( 10) and ( 11), e and ch are expressed by: :eeee e E E (19.a) ch chch (19.b) where: ,1 .( ), 2eee e e MA e E E E E (20.a) , ( )ch chTT T TThh (20.b) Using Eqs. ( 12), (13), and ( 19), Eq. ( 18) can be re -written as: ** 24 1: ( ) : : :0eee e P e e ch chP revX T E C T L E E BL (21) Since either positive or negative values can be assigned to the parameters eE and , according to the method proposed by Coleman and Noll (Noll et al., 1974) , the necessary condition for Eq. (21) to be satisfied is that: *=ee e TE E (22) chT Th (23) Using Eqs. (2 1-23), the dissipation inequality can be stated by: 24 24 24 24 24 * det det 1 1 1 1 1 24 24 24 24 2 1 1 1 1 1( ) : ( )rev ire rev ire rev ire rev ire rev iree tr tr TRIP DIP e e ch ch revgg X C T B P P P P P 4 0 (24) Adopt ing a stronger from of Eq. (24), since all the remaining phenomena are dissipative, the energy dissipation of each process can be considered non-negative. 24 * 1(( ) : ) 0e che tr rev rev revgX C T B P (25.a) 24 * 1(( ) : ) 0e che tr ire ire ireg C T B P (25.b) 24 * 1(( ) : ) 0e TRIP C T B P (25.c) 24 * det 1(( ) : ) 0e C T B P (25.d) det24 * 1(( ) : ) 0e DIP C T B P (25.e) Considering non -negative dissipation inequality for each variant , Eq. (25) can be written as: *(( ) : ) 0 =1,2,...,24e ch e tr rev revrevgX C T B P (26.a) *(( ) : ) 0 =1,2,...,24e ch e tr ire ireireg C T B P (26.b) *(( ) : ) 0 =1,2,...,24e TRIP C T B P (26.c) * det(( ) : ) 0 =1,2,...,24e C T B P (26.d) * det(( ) : ) 0 =1,2,...,24e DIP C T B P (26.e) Therefore, the dissipation inequalities, in terms of the thermodynamics driving forces and their conjugate variables, are expressed as follows: 0 =1,2,...,24rev rev (27.a) 0 =1,2,...,24ire ire (27.b) 0 =1,2,...,24TRIP TRIP (27.c) det0 =1,2,...,24 (27.d) 0 =1,2,...,24DIP DIP (27.e) where rev , ire , TRIP , det , and DIP are the thermodynamic driving forces related to reversible phase transformation, irreversible phase transformation (residual martensite), transformation -induced plasticity, detwinning, and detwinning -induced plasticity, respectively. *( ) :e che tr rev rev revgX C T B P (28.a) *( ) :e che tr ire ire ireg C T B P (28.b) *( ) :e TRIP C T B P (28.c) * det det( ) :e C T B P (28.d) det*( ) :e DIP C T B P (28.e) 2.2. Evolution laws 2.2.1. Phase transformation of reversible martensite According to Clausius inequality, Eq. (27.a) , rev and rev should always have the same sign. rev captures positive values during forward, and negative values during reverse transformation. Therefore , rev should be non -negative during forward , and non -positive during reverse transformation. During elastic deformation , rev is equal to zero, and the entropy inequality is automatically satisfied. Using the concept of transformation functions (Lagoudas and Entchev, 2004; Patoor et al., 2006) , the evolution equations for the reversible martensite volume fraction are as follows: 5 if 0, 1, 1 if 0, , 0 other conditionstr trn rev rev n revrevrev ire ireY bY (29) where Y is the resistance to phase transformation. 5b reflects the difference in evolution features of reversible martensite volume fraction during forward and reverse transformations . trn and Y are the controlling parameters for rate sensitivity of martensite transformation and width of hysteresis loop, respectively. 0 drev Y Y Y (30) The magnitude of Y at the start of forward transformation is denoted by 0Y . dY describes the evolution of Y during forward and reverse transformations . 2.2.2. Transformation -induced plasticity Considering Eq. (27.c) and the fact that TRIP is an irreversible process with 0TRIP at all times, the constraint conditions for TRIP are such that : 0 0 0 no constraintTRIP TRIP TRIP  = (31.a) (31.b) Therefore, the evolution equation for TRIP can be explained by (Yu et al., 2015) : 2 ()sat TRIP TRIP TRIP TRIP bH (32) where x and () Hx denote McCauley bracket and unit step function, respectively. 2b is the controlling parameter for saturation rate of TRIP , and sat TRIP is the saturated value of TRIP . () TRIP H insures that the 2nd law of thermodynamics is always satisfied. 2.2.3. Phase transformation of residual martensite Due to irreversible dissipative nature of residual martensite phase transformation, the related parameter cannot assign negative values. Therefore, following conditions are needed for residual martensite volume fraction to satisfy the 2nd law of thermodynamics . 0 0 0 no constraintire ire ire  = (33.a) (33.b) Experimental observations reveal an asymmetry in the stress -strain response of NiTi SMA under compressive and tensile loadings. While a flat stress -plate au is observed during tension, a rapid strain hardening is reported during compression . This difference in material behaviour is addressed to the generation and migration of dislocations during compressive loading (Liu et al., 1998) . On the other hand, it is known that residual martensite is referred to the part of martensite which is locked by defects and can not return to austenite phase (Yu et al., 2015) . Considering these points , it is concluded that the amount of residual martensite will be higher in compression than in tension. This conclusion is supported by the experimental results of Thamburaja and Anand (Thamburaja and Anand, 2001) . To reflect the difference of residual martensite volume fraction in tension and compression, separate saturation values of this parameter are assigned for tension and compression. Since the amount of residual martensite in super -elastic NiTi increases with i ncreasing number of cycles , eventually reach ing a saturated value (Kang et al., 2009) , the evolution equation for residual martensite volume fraction is proposed here as: 1( ) for tension ( ) for compressionsat irev irev sat irev irevrev irev irev rev irev bH H (34) where satirev is the saturated value of irev , and 1b controls the saturation value of irev in compressive loading . 2.2.4. Detwinning Considering the constraint imposed by the 2nd law of thermodynamics , the evolution equation for detwinning of the -th martensite variant is expressed as: det det det det()n signK (35) Since can assume negative , zero, or positive values, det()sign guarantees the satisfaction of the second law of thermodynamics . detK represents resistance to detwinning. detn is a parameter used to determine the degree of dependence of the process on the loading rate and the higher its value, the lower the dependence of the detwinning process on the loading rate. 2.2.5. Detwinning -induced plasticity Considering DIP as an irreversible dissipative process, with 0DIP at all times, the necessary conditions for DIP to satisfy the 2nd law of thermodynamics are such that : 0 if 0 0 no constraintDIP DIP DIP  = (36.a) (36.b) Experimental observations reveal accumulation, and finally saturation, of the residual and peak strains during cyclic loading of NiTi SMAs (Song et al., 2014) . Therefore, considering the dissipative inequ ality, the evolution equation for DIP is such that (Ebrahimi et al., 2020) : 3 ()sat DIP DIP DIP DIP bH (37) where sat DIP and 3b are the saturation value of DIP and the controlling parameter for the saturation rate of detwinning -induced plasticity, respectively . () DIP H is used to guarantee the satisfaction of dissipation inequality , so that , DIP >0 if 0DIP , and 0DIP= if 0DIP . Table 1 represents a summary of the formulation used in the proposed constitutive model. Table 1: Summary of the equations used in the proposed constitutive model ep F F F 24 24 11tr tr tr tr rev ire rev iretrgg     === + = + PP L L L 24 1TRIP TRIP LP 24 det det 1, LP Deformation gradient tensor Velocity gradient tensor 24 det 1DIP DIP LP *( ) :e che tr rev rev revgX C T B P *( ) :e che tr ire ire ireg C T B P *( ) :e TRIP C T B P * det det( ) :e C T B P det*( ) :e DIP C T B P Thermodynamic deriving forces 5 if 0, 1, 1 if 0, , 0 other conditionstr trn rev rev n revrevrev ire ireY bY Evolution laws for reversible martensite phase transformation 1( ) for tension ( ) for compressionsat irev irev sat irev irevrev irev irev rev irev bH H Evolution laws for residual martensite 2 ()sat TRIP TRIP TRIP TRIP bH Evolution law for TRIP det det det det()n signK Evolution law for detwinning 3 ()sat DIP DIP DIP DIP bH Evolution law for DIP 2.3 Calculation of stress In order to model the behavior of NiTi SMAs under large deformations, the constitutive model has been implemented in the Abaqus/Explicit finite -element program, using VUMAT subroutine. The following describes how stress is calculated using the relationship s presented in Section 2. Using Eqs. (3) and (8), the plastic deformation gradient tensor is calculated as follows: 124 24 24 24 det det 1 1 1 1) ( ) )((nnP tr P TRIP DIPI g F P P P P F (38) where 1)(nP F is the plastic deformation gradient tensor in the ( 1) thn+− increment , and )(nPF is the corresponding variable in the thn− increment. Using Eq. (1), With F and pF to be known, the elastic deformation gradient tensor can be determined. 1ep F FF (39) The elastic Green -Lagrange strain tensor, in terms of the elastic deformation gradient tensor, is calculated such that (Anand and Gurtin, 2003) : 1()2Te ee EI FF (40) Using eE and *T (Eq. (22)), the Cauchy stress tensor can be calculated. * 1 detTee T F T FF (41) 3. Material parameters identification For austenite phase, with cubic crystal structure, three independent elastic coefficients are used to define the fourth -order stiffness tensor (Jamal et al., 2014) : 11 12 12 12 11 12 12 12 11 44 44 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 = 44 0 0 0 0         (42) The corresponding independent elastic constants for austenite phase are denoted by 11A , 12A , and 44A , and considered to be 130, 98, and 34 GPa, respectively (Brill et al ., 1991; Thamburaja and Anand, 2001) . The modulus of elasticity of austenite is usually two to three times greater than the modulus of elasticity of martensite, while the Poisson’s ratio s are almost the same . Therefore, the elastic constants of martensite ( 11 12 44,,MMM ) are considered equal to half of the corresponding constants of austenite (Thamburaja and Anand, 2001) . Sufficiently large values are selected for the parameters trn and detn , to reproduce the rate- independent response of NiTi SMAs . The evolution law of reversible martensite volume fraction during forward and reverse transformations is described by the parameter Y , while 5b is responsible for difference in evolution of reversible martensite during loading and unloading . The width of hysteresis loop for the first and steady -state cycles is controlled by the parameters 0Y and dY , respectively. Using the curves related to start-stress of detwinning , detK can be obtained. As experimentally observed (Liu et al., 1998) and computationally predicted (Ebrahimi et al., 2020) , DIP has no considerable effect in tension . Therefore, the p arameters related to DIP cannot be obtained by using the data related to tensile loading of SMA. On the other hand, as commented by Kang et al. (Kang et al., 2009) , TRIP affects both the accumulated peak strain and accumulated residual strain during uniaxial cyclic deformation, while the influence of residual martensite is mainly focused on the accumulated residual strain. Thus, the data of accumulat ed peak strain in tensile cyclic loading can be used to obtain the parameters 2b and sat TRIP . Using the data related to the difference of peak and residual strains in tensile and compressive cyclic loading s, the parameters 1b and sat irev can be obtained . The parameter s 3b and sat DIP , which are related to DIP, can be obtained by using the evolution curves of residual strain in compressive loading. Using transformation modulus at the first and steady -state cycles, the parameters 0H and satH can be calibrated, respectively. satB , which is used to demonstrate the reduction in the start - stress of martensite transformation during cyclic loading, can be calibrated by measuring the start-stress of martensite transformation in the first and steady -state cycles. In order to control the saturation rates of B and H , parameters 6b and 4b are used, respectively. Therefore, the evolution curves of decreased start -stress of martensite transformation and hardening modulus can be used to obtain the values of 6b and 4b , respectively. The maximum transformation strain is controlled by the parameter trg . The parameter Th controls the temperature dependence of the driving force. The phase equilibrium temperature can be calculated by 1()2T ms as  + , where ms and as are martensite and austenite start temperatures at stress -free conditions , respectively (Thamburaja and Anand, 2001) . It should be noted that the above -mentioned parameters are related to a single crystal and should be calibrated using the corresponding single crystal experimental data. In the lack of such experimental results for single crystal, a trial -and-error method is used to calibrate the parameters using the polycrystalline experimental data . The adopted parameters f or predicting the behavior of NiTi SMA by considering the effect of residual martensite, are listed in Table 2. Table 2: Parameters used to predict Nitinol behavio ur in tensile and compressive loads considering the effect of residual martensite accumulation Unit Value Parameter GPa 130, 98, 34 11 12 44, , A A A GPa 65, 49, 17 11 12 44, , MMM MPa 10.2, 0 0, d YY - 50 det, trnn MPa 15 detK MPa 1 satB MPa 80, 1 0, sat HH - 0.19 trg - 2.5, 1, 1, 1, 1.5, 1 1 2 3 4 5 6, , , , ,b b b b b b MPa 140 Th K 298, 256 , T - 0.1, 0, 0 , , sat sat sat irev TRIP DIP   The parameters used to predict the cyclic behavior of NiTi, considering the effect s of DIP and TRIP are presented in Table s 3 and 4, respectively . Table 3: Parameters used to predict the cyclic behavior of NiTi by considering the effect of detwinning -induced plasticity Unit Value Parameter GPa 130, 98, 34 11 12 44, , A A A GPa 65, 49, 17 11 12 44, , MMM MPa 10.2, 10 0, d YY - 50 det, trnn MPa 15 detK MPa 50 satB MPa 600, 0 0, sat HH - 0.1708 trg - 1 1 2 3 4 6, , , ,b b b b b MPa 100 Th K 298, 256 , T - 0, 0, 0,0 01 , , sat sat sat irev TRIP DIP   Table 4: Parameters used to predict the behavio ur of NiTi SMA by considering the effect of transformation - induced plasticity Unit Value Parameter GPa 130, 98, 34 11 12 44, , A A A GPa 65, 49, 17 11 12 44, , MMM MPa 10.2, 10 0, d YY - 50 det, trnn MPa 15 detK MPa 50 satB MPa 200, 800 0, sat HH - 0.1708 trg - 1 1 2 3 4 5 6, , , , ,b b b b b b MPa 100 Th K 298, 256 , T - 0, 0.001, 0 , , sat sat sat irev TRIP DIP   4. Results and discussion In polycrystalline shape memory alloys, manufacturing process can affect the crystallographic texture which is the main reason for tension -compression asymmetry observed in such materials (Thamburaja and Anand, 2001) . As commented by Gall et al. (Gall et al., 2001) , based on the method used to form an alloy, three different textures with different effects on the overall response can be formed: a 100 texture making the mechanical response stiffer in tension than in compression, a 111 texture causing softer response in tension than in compression, and a 110 texture which almost has no effect on the mechanical responses in tension or compression. It should be noted that the parameters related to phase transformation direct ion, normal to the phase transformation plane, and detwinning are expressed in a crystal local coordinate system. Therefore, a transformation of the deformation gradient tensor from poly -crystal to single - crystal coordinate system is performed at the begin ning of each increment. At the end of calculations, the calculated stress tensor should be transformed back to poly -crystal coordinate system, which is done by using the rotation matrix Q . Appropriate ranges of Euler angles have been calculated in order to model different textures . To this end, the single -crystal and poly -crystal coordinates are considered as local and global coordinate systems, respectively. According to the Bunge system, the orientation of each crystal can be identified using 3 Euler angles,  , , and  . In this system, the rotation matrix in 3D is expressed as (Depriester, 2 018): cos cos sin sin cos sin cos cos sin cos si n sin ( , , ) cos sin sin cos cos sin sin cos cos cos cos sin sin sin co s sin Q                              −+ = − − − + − cos    (43) In order to model crystals with [111] texture, the [111] direction of each crystal should be nearly parallel to the loading direction. A set of randomly generated Euler angles without any preferred orientation, is used to simulate a random texture/untextur ed polycrystalline material. Using the proposed model, the behavior of a single -crystal for different textures under tensile and compressive strain -controlled loading -unloading is presented in Fig. 1. Fig. 1. Comparison of the behavior of NiTi SMA under tensile and compressive loads for different textures of crystals: (a) [100] texture, (b) [111] texture, (c) [110] texture, (d) Random texture. Magnitudes of Cauchy stress and logarithmic strain are presented for direct comparison of tension and compression. As shown in Fig. 1, for [111] texture the material behavior is harder in compression than in tension. In contrast, for the two other textures, and also for the randomly oriented crystals, the behavior in compression is softer in compression than in tension. Tension -compress ion differences are less pronounced for [110] and random textures. Based on the asymmetrical tension/compression response observed in the experimental test of Thamburaja and Anand (Thamburaja and Anand, 2001) , the [111] texture is adopted to model the behavior of NiTi SMA in this work. Since the number of crystals used to model the polycrystalline material can affect the predicted behavior, an appropriate number of crystals should be selected. The number of crystals should be as small as possible to reduce the computational cost and on the other hand should be large enough to prevent considerable change in material behavior by increasing number of crystals afterwards . Therefore, a sensitivity analysis is performed to determine the number of crystals in a representative volume element required for a converged stress -strain r esponse . Results for 1, 8, 27, 64, and 125 crystals are shown in Fig. 2. Converged stress -strain behavior is obtained for 27, 64, and 125 crystals. Based on this sensitivity analysis 64 crystals will be used to model the polycrystalline behavior for the re mainder of this study . Fig. 2. Comparison of stress -strain curves of NiTi SMA for different number of crystals with a minimum strain of 0% and a maximum strain of 6%. Von Mises stress distributions for 27, 64, and 125 crystals are shown in Fig. 3. Although different stress contours are obtained for different number of crystals, the stress -strain curves have the same trend in all three cases. Beside using various number of crystals , altering the randomly assigned orientations of the crystals changes the stress contours and locations of stress concentrations. However, the overall stress -strain response is independent of the randomly assigned orientations in all three cases . The von Mises stress distributions for two various arrangements of crystals in a polycrystalline material with 125 crystals are presented in Fig. 4. It should be mentioned that the same sets of Euler angles have been used to capture [111] texture in both cases , but the Euler angles assigned to eac h crystal are rando mly selected from th at set in each case. Fig. 3. Von Mises stress distribution at the end of loading for tension -unloading with a minimum strain of 0% and a maximum strain of 6% in a polycrystalline material with different number of crystals: (a) 27, (b) 64, and (c) 125 crystals. Fig. 4. Various von Mises stress distributions for two different arrangements of crystals in a polycrystalline material with [111] texture composed of 125 crystals . The results are presented at the end of loading in tension - unloading with a minimum strain of 0% and a maximum strain of 6%). In order to evaluate the proposed model , in Fig. 5 computed results for tensile and compressive loadings are directly compared with corresponding NiTi experimental data reported by Thamburaja and Anand (Thamburaja and Anand, 2001) . The parameters used to model Nitinol's behavior are presented in Table 2. Several key features should be noted: (i) the experimental results show residual strain after unloading, in both compression and tension. This is accurately predicted by the proposed model due to the in corporation of the effect of residual martensite . Experimental results of Liu et al. (Liu et al., 1998) show a stress- plateau in the stress -strain response during tensile loading, while the response is quickly strain hardened and no plateau is observed during compressive loading , which is related to the generation of dislocations (Liu et al., 1998) . On the other hand, the increase of defects causes the residual martensite to be accumulated (Yu et al., 2015) . Using these 2 facts and the experimental results of Thamburaja and Anand (Thamburaja and Anand, 2001) , it is concluded that the amount of residual martensite should be higher in compression than in tension. In order to reflect the difference of residual martensite volume fraction in tension and compression, two separate values are assigned to the saturation value of this parameter for 2 loading cases in Eq. (34). (ii) The proposed model accurately predicts the significant experimentally measured asymmetry between tension and compression, both in terms of loading and unloading behavior. As shown in Fig. 1, adopting [111] texture will cause to get a harder response in compression compared to tension. On the other hand, considering hardening/softening effects in the proposed model (Eq. (13)), the experimental results are predicted more precisely. To highlight the improved predictions of the current model, corresponding predictions from the model of Thamburaja and Anand (Thamburaja and Anand, 2001) are reproduced. This earlier model does not include the effect s of residual martensite and hardening/softening effects during phase transformation , and therefore does not predict residual strain or increased hardening under compressive loading. Fig. 5. Comparison of stress -strain curves of NiTi SMA with experimental results and available analytical results for two loading cases: (a) compression with a minimum strain of -5% and unloading to zero stress, (b) tension with a maximum strain of 6% and unload ing to zero stress . The distributions of von Mises stress es and residual martensite volume fractions for the aforementioned compressive and tensile loading cases are shown in Figs. 6 and 7, respectively. Fig. 6. Distribution of von Mises stress at the end of loading for (a) compression to a minimum strain of -5% and unloading to zero stress, (b) tension to a maximum strain of 6% and unloading to zero stress. Fig. 7. Contour plots of residual martensite volume fraction at the end of unloading for (a) compression to a minimum strain of -5% and unloading to zero stress, (b) tension to a maximum strain of 6% and unloading to zero stress. The ability of the proposed model in predicting the cyclic behavior of NiTi SMA is shown in Fig. 8. The stress -strain curves are plotted for the first 5 cycles in tensile loading -unloading with a minimum stress of 0 MPa and a maximum stress of 550 MPa . Fig. 8. NiTi SMA cyclic stress -strain diagram under tensile loading -unloading with a minimum stress of 0 MPa and a maximum stress of 550 MPa . The proposed model captures the main features observed in the experimental tests (Kang et al., 2009; Song et al., 2014) . In summary, with increasing number of cycles , the start stress of austenite to martensite transformation decreases, the width of the hysteresis loop decreases, and the residual strain approaches a steady state, indicating that the accumulation of residual martensite approaches a saturation level . In order to investigate the effect of DIP in large deformation, the stress -strain diagram s of the material under compressive loading to a minimum strain of -5% and unloading to zero stress are shown for the first, fifth and tenth cycles in Fig. 9. Fig. 9. Stress -strain diagram of NiTi SMA under compressive cyclic loading with a minimum strain of -5% and unloading to zero stress , with and without considering the effect of detwinning -induced plasticity : a) cycle 1, (b) cycle 5, and (c) cycle 10 As shown in Fig. 9, the effect of detwinning -induced plasticity in the first cycle is negligible but is more pronounced with increasing number of cycles. Similar to the results prese nted by Ebrahimi et al. (Ebrahimi et al., 2020) , considering the effect of DIP reduces the absolute value of the maximum stress and increases the absolute value of the residual strain in the strain - controlled cyclic loading of NiTi SMA. The contour plots of shearing deformation caused by detwinning -induced plasticity for two variants of martensite are de picted in Fig. 10 Fig. 10. contour plots of shearing deformation , due to detwinning -induced plasticity , in two variants of martensite. 5. Conclusions In this study a microstructure -based model for NiTi SMAs is developed by incorporating transformation -induced plasticity (TRIP), detwinning -induced plasticity (DIP), and accumulation of residual martensite into a finite -deformation crystal plasticity framework. The constitutive model is constructed at the single -crystal scale and also includ es phase transformation and detwinning mechanisms. Using a proposed Helmholtz free energy , the driving forces for inelastic mechanisms are derived within the framework of thermodynamics. The constitutive model has been implemented in the Abaqus/Explicit fi nite-element program, using VUMAT subroutine to simulate a polycrystalline material. Considering various orientations for crystals, the effect of texture on tension -compression asymmetry is investigated. It is shown that different textures may cause harder , softer, or similar response in compression compared to tension, which is in good agreement with the comments of Gall and Sehitoglu (Gall and Sehitoglu, 1999) , Gao et al. (Gao et al., 2000) , and Thamburaja and Anand (Thamburaja and Anand, 2001) . Due to the incorporation of the effect of residual martensite, the model provides accurate prediction s of experimentally measured residual strain (Thamburaja and Anand, 2001) . The incorporation of the inelastic DIP and TRIP deformation mechanisms and using appropriate form of free energy is shown to accurately capture the key features related to cyclic loading , including decrease in start stress of austenite to martensite transformation, decrease in width of the hysteresis loop, and saturation of residual strain (Song et al., 2014) . Finally, the effect of detwinning -induced plasticity in compressive cyclic loading of NiTi SMA is investigated. DIP affects the results of strain -controlled cyclic compression -unloading simulations by causing the peak stress to be less negative and residual strain to be more negative in different cycles (Ebrahimi et al., 2020) . Future applications of the model should include analysis of the inelastic deformation mechanism s such as DIP on the behavior of stents during crimping and deployment, including detailed 3D crystal plasticity modelling (McGar ry et al., 2007) of the role of grain boundaries on residual martensite and DIP concentrations. Additionally, future studies should incorporate DIP and TRIP into microstructure based fatigue indication parameters (McCarthy et al., 2014, 2014 a). Fatigue failure remains a significant clinical problem for NiTi stents (Scheinert et al., 2005) and the analysis presented in the current study and the role of inelastic microstructural mechanisms in cyclic response of SMAs can potentially be used to guide the material and structural design of next -generation stents.