Temperature-dependent crystal-plasticity model for magnesium: A bottom-up approach

A crystal-plasticity model is developed to account for temperature-dependent mechanical behaviour of magnesium in this paper. The constitutive description of plastic deformation accounts for crystalline slip and twining as well as their interactions. The temperature dependence is incorporated into the constitu- tive equations for both slip and twin modes based on experimental observations. A bottom-up computational modelling framework is proposed to validate the developed constitutive model. First, the crystal- plasticity model is calibrated with experimental results for plane compression at micro-scale. At meso-scale, a three-dimensional representative element volume was adopted to represent the microstructure of polycrystalline magnesium. In the combination with the proposed constitutive theory, the effects of temperature on mechanical response and evolution of twins and texture in polycrystalline magnesium were predicted. Comprehensive experimental validations at meso-scale were performed to consolidate further the developed crystal-plasticity model incorporating temperature dependence in terms of stress-strain curves, the Hall-Petch relationship and texture evolution. This work provides a useful modelling tool for understanding temperature-dependent behaviour of magnesium, which could be used to improve the formability of this family of materials.


Introduction
With increasing demands for improved fuel efficiency in transportation, there is a strong drive to reduce the weight of vehicles without compromising their structural resilience. Therefore, magnesium (Mg) and its alloys have attracted significant attention in recent years thanks for their high specific strength ( Wei et al., 2015;Zhou et al., 2016 ). However, the widespread structural applications of Mg have been substantially restricted by material's poor ductility and formability. This drawback is primarily due to the underlying hexagonal close-packed (HCP) structure of Mg, which provides a limited number of slip systems for plastic deformation at room temperature ( Mirzadeh, 2014 ). Additional slip systems may be activated at elevated temperatures; consequently, hot processing is advised to overcome the poor formability of Mg ( Figueiredo et al., 2016;Mirzadeh, 2014 ). This needs to be performed with caution, as high temperatures may alter the material's microstructure with a concomitant change in an in-service mechanical response of * Corresponding author.
E-mail address: A.Roy3@lboro.ac.uk (A. Roy). a component ( Yuan et al., 2016 ). Thus, a thorough understanding of hot deformation behaviour of Mg is necessary.
In Mg and its alloys, micro-scale deformation mechanisms include both crystalline slip and deformation twinning. In addition, a significant transition of the dominant deformation mode is observed with a temperature variation. In recent years, crystalplasticity-based approaches that can explicitly implement different deformation modes have been widely used for fundamental investigations on deformation mechanism of various metallic materials, including Mg, further providing guidelines for design of novel materials ( Zhang et al., 2016 ) and formability improvements ( Liu et al., 2016a, c ).
Based on the modelling philosophy adopted in crystal plasticity, the respective modelling technique can be categorised into topdown and bottom-up approaches ( Zhang and Joshi, 2012 ). The topdown approach is well-suited to model polycrystalline behaviour at macro-scale, from which single-crystal parameters are inferred. The commonly used top-down approaches can be classified into isostrain ( Taylor, 1938 ), isostress ( Sachs, 1928 ) and self-consistent ( Eshelby, 1957 ) schemes. There are some successful examples for these different top-down approaches, for instance, the isotraintype ( Ardeljan et al., 2016;Knezevic et al., 2009 ), isostress-type ( Toth et al., 1990 ) ( Beyerlein and Tomé, 2008;Kabirian et al., 2015 ). The top-down approach is relatively easy to implement numerically and can reduce computational cost by incorporating a coarse finite-element (FE) mesh. However, this approach suffers from a high computational time cost at element-level simulation, which impedes parallel computation in FE analysis. Thus, the top-down approach may essentially decrease computational efficiency. Its another drawback is the fact that the choice of a homogenization scheme affects significantly the estimation of single-crystal parameters ( Ardeljan et al., 2016 ). By contrast, the bottom-up approach involves the use of a calibrated (based on experiments) small-scale single-crystalplasticity (SCP) model, incorporating typical deformation modes such as slip, twinning or both. To predict a response of polycrystalline component, individual crystal grains and orientations are represented via the SCP model, which is then employed to assess a stress-strain response and texture evolution during the deformation process. Here, the intra-grain interaction is modelled in a physically representative manner (in contrast to the use of homogenization). Successful implementations of the bottom-up approach were developed for a variety of crystalline materials, FCC ( Cyr et al., 2015 ), BCC ( Lim et al., 2015 ) and HCP metals ( Abdolvand and Daymond, 2013;Cheng and Ghosh, 2015 ) and references there in. Compared to the top-down approach, the bottom-up approach has a higher demand on computational resources, but such an approach is amenable for parallelisation in an FE solver.
In Mg and its alloys, much of the modelling effort involves the use of top-down approach. For example, VPSC models were employed to provide an insightful understanding activity of slip and twin mode at different temperatures in AZ31 alloy ( Kabirian et al., 2015;Zhang et al., 2016 ). In the works, the so-called predominant twin reorientation (PTR) scheme proposed by Tomé et al. (1991) was widely adopted to determine the twin-phase formation. In the PTR scheme, only one twin phase with a highest contribution to total volume fraction was activated in a grain. In these VPSC-based approaches, a critical resolved shear stress (CRSS) needs to be calibrated at different temperatures, which is their primary drawback. Recently, Ardeljan et al. (2016) proposed a Taylor-type modelling scheme, in which the temperature dependence was incorporated into constitutive laws, thus addressing its effect with introduction of appropriate parameters.
For bottom-up approaches, several SCP-based models were developed recently, with their parameters identified through singlecrystal experiments ( Becker and Lloyd, 2016;Gan et al., 2016;Zhang and Joshi, 2012 ). Additionally, these models were also employed to characterise polycrystals at meso-scale ( Chang and Kochmann, 2015;Zhang and Joshi, 2012 ). However, these studies were limited to investigations at room temperature. To date, only some limited attempts were made to capture temperature dependence with different sets of model parameters were used for different temperature conditions ( Hidalgo-Manrique et al., 2015 ). Thus, it is imperative to incorporate temperature dependence into constitutive laws for bottom-up approaches, which will allow for modelling across a wider temperature range exploiting a broader design space.
The aim of this paper is to develop a SCP model to account for the temperature dependence of Mg, henceforth, referred to as T-SCP (temperature-dependent single-crystal plasticity) model. This model was incorporated into a bottom-up modelling framework to investigate the effects of temperature on the mechanical response and texture evolution of single-crystal and polycrystalline Mg. This paper is organized as follows: in Section 2 , a selfcontained description of the governing relations of the proposed T-SCP model was presented. Section 3 presents a modelling strategy of the bottom-up approach based on a commercial FE software package ABAQUS. In Sections 4 and 5 , simulation results and experimental validations are presented and discussed for single-crystal and polycrystalline case at meso-scale, respectively. We end with some concluding remarks in Section 6 .

Constitutive formulas
In this section, a phenomenological T-SCP model is presented to account for the temperature-dependence of single Mg crystals (or grains). In the T-SCP model, four slip and two twin systems were considered for the Mg crystal as listed in Table 1 . Here, four slip planes are considered: basal, prismatic, pyramidal a and pyramidal c + a and two twin planes: tensile twin (TT) and compressive twin (CT) (see Fig. 1 ). Standard notation is adopted here: scalars are in italics, vectors and tensors are indicated with lower-case and upper-case bold letters.

Kinematics
Following a classical crystal plasticity (CP) theory, the deformation gradient F can be decomposed into the elastic and plastic parts, as, where the subscripts 'e' and 'p' denote the elastic and plastic parameters, respectively. The velocity gradient L is introduced following its definition L = ˙ F F −1 , as, For Mg crystal, the plastic deformation is assumed to arise from both crystalline slip and twinning due to its HCP structure with a large aspect ratio. Consequently, the plastic velocity gradient, L p , incorporates contributions from the slip and twin modes as Here L sl p , L tw p and L sl−tw p represent the plastic velocity gradient induced by the slip in the untwined region (or parent phase), deformation twinning in the untwinned region and secondary slip in the twinned region (or child phases), respectively ( Kalidindi, 1998 ). In this paper, the assumption of pseudo slip is adopted for twinning, and its effectiveness has been demonstrated in prior work ( Ardeljan et al., 2016;Gan et al., 2016;Kalidindi, 1998 ). For the sake of clarity, the superscript α is used to represent the slip system in the parent phase, β for the twin system in the parent phase and ˜ α for the secondary slip system in the child phase. The three terms in Eq. (3) can be further expressed as where ˙ γ α is the shear slip rate on the slip system α, f β is the volume fraction of child phase β, and ˙ γ β is the shear strain rate arising from deformation twinning. N s and N tw are the total numbers of slip and twin systems, respectively. The unit vector s represents the direction of slip/twin and m is the unit vector normal to the corresponding slip/twin plane. Furthermore, the velocity gradient can be expressed in terms of a symmetric rate of stretching D and an antisymmetric rate of spin W : From Eqs.

Constitutive laws
Average Cauchy stress σ is evaluated at each material point, by accounting the contributions from both parent and child phases: where σ m and σ tw ( β) denote the Cauchy stress in the parent and child ( β) phases. In each phase, following the work of Huang (1991) , the constitutive law is expressed as the relationship between the elastic part of the symmetric rate of stretching, D e , and the Jaumann rate of Cauchy stress, where, I is the second-order unit tensor, C is the fourth-order, possibly anisotropic, elastic stiffness tensor. The Jaumann stress rate is expressed as ( Liu et al., 2016b ) The temperature dependence of elastic tensor C can be expressed as ( Olsson, 2015 ): where T is temperature (in Kelvin), s ijkl and t ijkl are material constants with the same symmetry as the elastic tensor, C ijkl , and C 0 i jkl corresponds to the elastic tensor in the limit of zero temperature. With the framework of the classical CP theory, the shear strain rate on each slip system, ˙ γ α (or ˙ γ˜ α ) is related to the resolved shear stress τ α (or τ˜ α ) via a well-known power law proposed by Hutchinson ( Hutchinson, 1976 ): Here, ˙ γ 0 is the reference shear rate, τ α c is the slip resistance and n is the rate-sensitivity parameter. On each slip system, the resolved shear stress, τ α , is expressed by the Schmid law τ α = s α m α : σ.
Here, we take pause and observe that determining yield surface in single crystals especially HCP metals is not a trivial matter. We refer the reader to the seminal work of Tomé and Kocks ( Tomé and Kocks, 1985 ) and Ritz and co-workers ( Ritz et al., 2010 ).
Next, assuming twinning to be essentially pseudo slip, the Schmid law is also adopted for twin systems, with a power-law relation employed to describe evolution of the twin volume fraction, where, ˙ f 0 is the reference rate of the twin volume fraction, τ β is the corresponding resolved shear stress, τ β c is the resistance to deformation twinning, m is the rate-sensitivity parameter, and f cr represents the critical twin volume fraction at which lattice reorientation is invoked at the corresponding material point ( Kalidindi, 1998;Zhang and Joshi, 2012 ). That is assumed to be 0.9. The shear strain rate on the twin system, ˙ γ β , may be related to the rate of the twin volume fraction, ˙ f β , via constant shear strain rate, γ tw , as Here γ tw is determined by the aspect ratio of crystal lattice, where χ = 1 . 624 for Mg crystal.
Lattice reorientation is introduced when the accumulated twin volume fraction over all twin systems ( f ) exceeds the critical one ( f cr ) at a material point ( Eq. (13) ). In literature, the PTR scheme is widely used to capture grain reorientation due to twinning. The PTR scheme allows one child phase to evolve in a crystal or grain ( Gan et al., 2016;Zhang and Joshi, 2012 ). In contrast, multiple child phases are allowed in our approach (similar to the work of Ardeljan et al., (2016) ). Here, a child phase is introduced when the volume fraction of the corresponding twin system reaches a threshold value (we assume this to be 1%). Therefore, there can be a maximum of 12 child phases (6 phases resulting from TT and 6 from CT). Though, in practice activating all twin systems simultaneously is difficult. When a new child phase is generated from the twin β, its new orientation can be obtained by The elastic tensor (in indicial notation) of the new child phase, C tw (β ) i jkl , is defined by rotation from that of the parent phase, C ori pqrs , which can be formulated as:

Hardening model
Here, we present hardening laws that capture the evolution of slip resistance (τ α c ) and twinning resistance (τ β c ) introduced in Eqs. (11) and ( 13 ), respectively, accounting also for the temperature dependence. To date, there were limited experimental studies with regard to temperature dependence of different slip and twin modes in single-crystal Mg ( Chapuis and Driver, 2011;Wonsiewicz, 1966 ). It was demonstrated that the basal slip system and TT were temperature-independent, while other slip modes and CT exhibited significant temperature dependence.

Hardening law of slip
Slip resistance is defined as, where τ α 0 , τ α HP and τ α f represent the initial lattice resistance, resistance from the barrier imposed by grain or twin boundaries and forest dislocation interaction, respectively. The first term is obtained from ( Ardeljan et al., 2016;Beyerlein and Tomé, 2008 ): where s α 0 is the reference initial lattice resistance, and T α is an empirical parameter.
The second term, τ α HP , is considered for two cases: with and without child phases in the crystal ( Beyerlein and Tomé, 2008 ), and is presented as a unified expression similar to the classical Hall-Petch effect: Here, μ α and b α are the shear modulus and Burgers vector of slip system α, respectively. H α I is the material parameter depending on the slip mode, with the subscript I indicating three possible conditions: I = 0 -no child phases in the crystal; I = 1a predominant child phase resulting from TT; I = 2 -a predominant child phase resulting from CT. Depending on the choice of I , the parameter d represents grain size ( d g ) when I = 0 , or the mean free path between adjacent child phases ( d α m f p ) when I = 1 or I = 2 ( Beyerlein and Tomé, 2008 ). Clearly, the presence of child phases in a crystal (twin boundaries) introduces an additional barrier for dislocation motion on top of the presence of grain boundary, which manifests in the classical Hall-Petch effect. The meanfree-path, d α m f p , depends on the orientation between the predominant child phase and the slip plane ( Ardeljan et al., 2016 ), which is expressed as where f PTS is the volume fraction of the predominant twin system (PTS) in the crystal, λ represents the ratio of the twin spacing and grain size ( λ= 0 . 2 in this paper), and θ is the angle between the plane of PTS and the slip plane.
where τ α f,sl↔ sl and τ α f,tw → sl represent the slip resistance due to the slip-slip and twin-slip interactions, respectively. As an example, τ α f,α represents the interaction between slip systems α and α , and τ α f,β represents the barrier of twin system β acting on slip system α. The slip-slip system interaction incorporates the effect of active or self-slip-slip interactions (i.e. when α = α ) and latent slip-slip interactions (i.e. when α = α ). We assume that these interactions are related by where q α α is the latent interaction coefficient that generally ranges between 1 and 2. Here, it is assumed that q α α = 1 . As there are some differences between the influences of TT and CT on slip, it is important to distinguish these interactions ( Gan et al., 2016;Zhang and Joshi, 2012 ), as: where N TT and N CT are the number of TT and CT systems, respectively. Similar to slip-slip interactions, the TT-slip interaction is formulated as: Here τ β f,T T is the self-resistance to TT of the system β, and q TT → α is the interaction factor between TT and slip systems. The CT-slip interaction is assumed to follow a Taylor-hardening type form depending on the total volume fraction accumulated over all the N CT CT systems ( Zhang and Joshi, 2012 ), as: where H CT → sl represents the initial hardening parameter of the CTslip interaction. According to Eqs. (22) -( 26 ), the dependence of forest dislocation interaction (i.e. τ α f in Eq. (18) ) on temperature may be expressed by τ α f,α , τ β f,T T (TT) and H CT → sl (CT). Since basal slip and TT are temperature-independent ( Chapuis and Driver, 2011;Wonsiewicz, 1966 ), we can write the following expressions where s α basal and s β T T represent the reference self-resistance for basal slip and TT, respectively. For non-basal slip and the CT-slip interaction, the temperature dependence is defined based on the work of Kocks (1976) , as: where k is the Boltzmann constant, ˙ ε 0 is the reference strain rate ( ˙ ε 0 = 10 7 in this paper), ξ α and ξ CT are the non-dimensional coefficients, s α f is the reference self-resistance of non-basal slip, and the constant H 0 CT → sl is the reference hardening parameter of CT-slip interaction, respectively.
As the growing forest dislocation interaction further increases resistance to slip or results in an increase of the twin volume (e.g. the CT-slip interaction in Eq. (26) ), a hardening law is also required for slip-slip or TT-slip interaction. The hardening rate of slip resistance due to the slip-slip interaction ( τ α f,α ) can be expressed in terms of s α basal or s α f . The basal slip system is observed to follow a linear hardening response based on experiments ( Kelley, 1967 ) For a non-basal slip system, the hardening rate of s α f is formulated in the form proposed by Asaro (1983) : where h α 0 is the initial hardening modulus and s α s is the saturation stress. Similarly, the hardening rate of slip resistance due to the TT-slip interaction ( τ α f,T T ) can be obtained from s β T T (reference selfresistance of TT system) according to Eq. (27) . The hardening law of s β T T is discussed in the Section 2.3.2 .

Hardening law of deformation twinning
The resistance to deformation twinning in Eq. (13) can be also expressed as the sum of three terms (similar to Eq. (18) ), as: The first term is expressed as The second term, τ β HP , the effect of barriers due to grain or twin boundaries, is also expressed in the Hall-Petch-like form ( Beyerlein and Tomé, 2008 ): The resistance to the evolution of TT and CT volume fractions may be attributed to different underlying mechanisms. For example, CT dislocations have lower mobility due to their narrow corewidth, which impedes a CT growth. In contrast, the core-width of TT dislocation is about 3 ∼6 times of their CT counterpart; hence, TT dislocations are easy to nucleate and propagate ( Zhang and Joshi, 2012 ). Consequently, the twin-twin interactions are herein discussed by dividing them into two cases: TT-TT and CT-CT interactions: The resistance due to the TT-TT interaction is expressed as: where q TT is the interaction factor between TT systems, and it is taken as q T T = 1 . 0 in this paper. As discussed before, the TT-TT interaction is temperature-independent with τ β Based on the work of Zhang and Joshi (2012) , to characterise the sluggish kinetics of the CT growth at early stages, the resistance to CT evolution, τ β f,CT , is expressed as: where H CT ↔ CT and η are the empirical parameters controlling the hardening rate of CT-CT interaction. In general, the two parameters should satisfy the conditions of H CT ↔ CT ∼ GPaand η << 1 in order to capture a characteristic of the CT growth. The temperature depen- where H 0 C T ↔ C T is the reference hardening parameter for CT-CT interaction.
The effect of slip on deformation twinning (i.e. slip-twin interaction, τ β f,sl→ tw ) requires further studies. Capolungo et al. (2009b) studied the effect of slip on TT for Mg and demonstrated it to be insignificant. Although there was no direct evidence of the effect of slip on CT in Mg, Capolungo et al. (2009a) concluded that the onset of CT in Zr was insensitive to slip through series of mechanical test and multi-scale modelling. Based on these fundamental studies, we assume the slip-twin interaction in Eq. (34) may be neglected (i.e. τ β f,α = 0 ).

Scheme of bottom-up approach
The T-SCP model proposed in Section 2 was implemented in the commercial FE code ABAQUS/Explicit by employing the user subroutine VUMAT. It is necessary to point out that the stress update algorithm was based on the Green-Naghdi stress rate in ABAQUS/Explicit environment. Therefore, a conversion algorithm was required in order to evaluate a stress update in ABAQUS/Explicit based on the Jaumann stress rate defined in the constitutive law (e.g. Eq. (9) ); one can find more details in our previous work ( Liu et al., 2016c ). First, the model was used to simulate the effect of temperature on the stress-strain response and evolution of microstructure for single-crystal Mg. The model parameters were calibrated using the experimental data reported by Kelley (1967) and Wonsiewicz (1966) . Next, the T-SCP model was employed to predict the overall mechanical properties and microstructure evolution of polycrystalline Mg using a three-dimensional (3D) representative volume element (RVE) modelling approach.  The studies for single-crystal were carried out with the FE modelling strategy illustrated in Fig. 2 . At room temperature, basal slip is the easiest one to be activated among the four slip and two twin systems considered in this paper. A channel-die experimental test may be performed to study an individual slip or twin system ( Kelley, 1967;Wonsiewicz, 1966 ). In these experiments, plane compression loading is imposed on the single crystal in pre-decided orientations, that is, homogeneous compression loading is imposed on a chosen surface of a parallelepiped sample, while one orthogonal surface is held rigid and the third is left free. Seven loading cases were modelled using T-SCP as shown in Fig. 2 and compared with experimental studies. In these test cases, pyramidal < c + a > slip was primarily activated for cases A and B, prismatic slip for cases C and D, tensile twinning for cases E and F, and basal slip for case G. It is obvious that such channel-die tests significantly simplify model parameter calibration, and the involved details are presented in Section 4.1 .
Next, based on the calibrated single-crystal model, a 3D RVE model was used to represent the initial microstructure of polycrystalline Mg as shown in Fig. 3 . Similar RVE models with such idealized grains and meshes were also adopted in the work of Lim et al. ( Lim et al., 2015( Lim et al., , 2011. 200 (8 × 5 × 5) cubic grains were considered in the RVE model, in which each grain was meshed with 64 C3D8 elements and assigned a random initial orientation. A unidirectional tension test case (the loading direction was along the X axis marked in Fig. 3 ) was modelled using this 3D RVE model. The effects of temperature on the stressstrain response and microstructure (i.e. the evolution of twinning and texture) were predicted. In particular, the dependence of yield stress on grain size (i.e. Hall-Petch relationship) was estimated at different tem peratures, and, hence, the effect of temperature on the Hall-Petch relationship could be evaluated. Finally, the presented predictions were compared with the available experimental data ( Ono et al., 2004 ), which are used to validate effectiveness of the T-SCP model at polycrystal meso-scale.

Plane compression and unidirectional tension at room temperature
First, single-crystal deformation mechanism at room temperature was characterised using the T-SCP model. The stress-strain curves, obtained in FE simulations and experiments, are shown in Fig. 4 for the seven loading cases as illustrated in Fig. 2 . Fig. 4 (a) shows the stress-strain curves corresponding to the loading cases with crystal slip dominated. Pyramidal c + a slip was predominant in loading cases A and B with partial activation  of CT. A similar observation was made in the study of Gan and coworkers ( Gan et al., 2016 ). The stress level in loading case B was higher than that of A. This may be understood by analysing Schmid factor, which affects the magnitude of resolved shear stress. Under plane compression, contributions to the Schmid factor come from both the loading direction as well the specific direction in which the constraints are imposed. In loading cases A and B, the loading direction was identical; however, the effects due to constraints lead to a Schmid factor of 0 for loading case A and −0.11 for loading case B ( Zhang and Joshi, 2012 ). This implies that a higher level of compressive stress is required to activate slip in loading case B. Prismatic slip was the easiest to activate in loading cases C and D similar to the conclusions drawn in literature ( Gan et al., 2016;Zhang and Joshi, 2012 ). In contrast to cases A and B, the stressstrain curves are nearly the same for cases C and D due to the identical boundary constraint conditions. For loading case G, plastic deformation was mainly accommodated by basal slip. As shown in Fig. 4 (a), compared to other types of slip systems, basal slip had the lowest slip resistance and the hardening rate is also relatively low. Fig. 4 (b) shows the stress-strain curves of loading cases E and F where TT dominated at the initial stage. For both cases, it is clear that the initial yield stress was low due to the low resistance of TT. At the initial stage, the stress level in case F was higher than that of case E due to a difference in constraint boundary conditions. However, a noticeably large difference occurred after TT-induced reorientation. In loading case E, stress increased rapidly after reorientation, and the strain-hardening phenomenon was significant during the whole loading process. In contrast, for loading case F, stress value saturated after a rapid increase due to TT-induced reorientation.
Reorientation of single-crystal Mg after macro-plane compression is depicted via pole figures in Fig. 5 . Here, the critical twin volume fraction was as assumed to be f cr = 0 . 9 ; thus, orientation with relative low intensity represents the parent phase (or initial orientation) in Fig. 5 . As shown in Fig. 5 (a), the c-axis is rotated   by nearly 90 °through TT-induced reorientation in loading case E, thus, aligning the c-axis parallel to the loading direction. Consequently, subsequent plastic deformation was primarily governed by pyramidal c + a slip and partially by CT, similar to that for loading case A. For loading case F, as shown in Fig. 5 (b), the c-axis was rotated by about 60 °under macro-plane compression. Therefore, prismatic and pyramidal a slip dominated plastic deformation after TT-induced reorientation occurs. In summary, our studies clearly demonstrate different mechanisms involved in TT-induced reorientation, which resulted in significantly different stress-strain responses for loading cases E and F. Next, the stress-strain response and microstructural evolution under unidirectional tension were studied using T-SCP model for Mg crystal, and the predicted results are shown in Fig. 6 . The loading schematic is demonstrated in Fig. 6 (a), with two representative loading directions considered: tension case 1 : loading parallel to c-axis, i.e. along < 0 0 01 > and tension case 2 : loading perpendicular to c-axis, i.e. along < 10 1 0 > . Similar to loading cases E and F in plane compression, TT was initially activated in tension case 1 , followed by a rapid stress increase after the TT volume fraction reached the critical value (see Fig. 6 (b)). For tension case 2 , prismatic slip dominated during the whole loading history; as a result, the corresponding stress-strain curve is similar to that of cases C and D in plane compression. Again, for the two cases under unidirectional tension, the evolution of crystal Table 3 Key slip and twin constitutive equations and related parameters (length in μm, stress in MPa and unavailable parameters with '-').

Constitutive equations Model parameters
570  microstructure was investigated with pole figures as shown in Fig.  6 (c). It is clear that a 90 °rotation of c-axis was observed for tension case 1 , while there is no evidence of TT-induced reorientation in tension case 2 . For tension case 1 , as prismatic slip also dominated plastic deformation after reorientation, its stress-strain curve in the second stage is similar to that of tension case 2 .

Effect of temperature on stress-strain response
As mentioned above, the dependence of basal slip and TT on temperature was not significant according to experimental data ( Chapuis and Driver, 2011;Wonsiewicz, 1966 ). Consequently, in this subsection, only temperature-sensitive slip and CT that are discussed.
The effect of temperature on stress-strain curves, corresponding to slip-dominant plastic deformation, was predicted using the T-SCP model ( Fig. 7 ). The experimental data of Wonsiewicz (1966) was employed to calibrate the related model parameters. Fig. 7 (a) show the stress-strain response at different temperatures for case A under plane comrpession, which is primarily determined by pyramidal c + a slip. With an increase in temperature, initial yield stress, hardening rate and saturated flow stress decreased in loading case A. A similar effect of the temperature increase was also observed for loading case C, which was dominated by prismatic slip. As indicated in Fig. 4 , the stress-strain repsonse dominated by pyramidal a slip (the second stage of loading case F) was nearly the same as that by prismatic slip (loading cases C and D). Therefore, the model parameters of pyramidal a slip were assumed to be same as those of prismatic slip in this paper. A similar strategy was also adopted in the work of Zhang and Joshi (2012) and Gan et al. (2016) . As shown in Fig. 4 , the predicted stress-strain curves are consistent with the experimental data, which indicates the present T-SCP model is capable of capturing the effect of temperature on slip-dominant defomation mechanism.  The experimental data of Chapuis and Driver (2011) was employed to calibrate the model parameters related to the temperature dependence of CT. As there was no systematic experimental data in the form of stress-strain curves, calibration of the model parameters was based on the resistance to CT at different temperatures reported by Chapuis and Driver (2011) , which was not illustrated in this paper. Finally, the model parameters in the T-SCP theory were calibrated for a Mg crystal according to the experimental data at room and higher temperatures. As basal slip and TT were considered to be temperature independent, there were no temperature-dependent parameters involved for basal slip and TT. The parameters related to the temperature dependence of elastic constants are listed Table 2 , and the key constitutive equations and model parameters related to slip and twin modes in the T-SCP theory are summarized in Table 3 .

Simulation of polycrystalline Mg and discussions
In this section, the model parameters calibrated with the experimental data for the single-crystal cases were employed to predict the effective mechanical properties and evolution of microstructure for polycrystalline Mg at meso-scale. An emphasis was placed on the effect of temperature on stress-strain response, the Hall-Petch relationship, relative activity of slip and twin modes, and texture evolution.

Effect of temperature on overall mechanical properties
The stress-strain curves of polycrystalline Mg at different temperatures were predicted by combining the T-SCP theory with the RVE FE model. The simulation results, corresponding to the grain size of 47 μm, are presented in Fig. 8 . Similar to the phenomena observed for single-crystal Mg, both the yield stress and the hardening rate decreased with the increase of temperature for polycrystalline Mg at meso-scale. For example, significant work-hardening was observed with increasing strain at room temperature, while this can be neglected at 523 K. Here, the experimental data reported by Ono et al. (2004) was used for comparisons as shown in Fig. 8 . It is obvious that the present T-SCP model captures the dependence of stress-strain response on temperature for polycrystalline Mg adequately.
Furthermore, the yield stress of polycrystalline Mg was obtained by using the off-set method of 0.2% strain based on the predicted stress-strain curves, which is partially shown in Fig. 8 . Next, the effect of varying the average grain size d from ∼10 μm to ∼100 μm was studied. This allowed us to investigate the relationship between the yield stress and the average grain size (i.e. Hall-Petch relationship), and, in particular, the effect of tem perature on the Hall-Petch relationship. As shown in Fig. 9 , the variations of yield stress with the average grain size were demonstrated for various temperatures. Five different grain sizes, i.e. 35 μm, 47 μm, 70 μm, 10 0 μm and 20 0 μm, were considered in FE simulations, and the corresponding results are illustrated by empty markers in Fig. 9 . A best-fit curve at each temperature shows the linear relation between the yield stress and the inverse square-root of the grain size. With an increase in temperature the slope of the correlation is observed to reduce. Consequently, the Hall-Petch relationship becomes less significant with the growing temperature. As before, the related experimental data ( Ono et al., 2004 ) were used for comparison with the simulation results, showing a good agreement between the prediction and experimental data at each temperature.

Effect of temperature on evolution of twining and texture
In addition to the overall mechanical properties, it was also interesting to study the effect of temperature on the underlying de- formation mechanisms of polycrystalline Mg, which determine texture evolution.
First, the FE model was employed to investigate the activation of child phases during unidirectional tension for polycrystalline Mg. Fig. 10 show the number of activated phases at room temperature and T = 523 K under imposed macroscopic tensile strain of 10%. In Fig. 10 , no child phases being activated is depicted by number of phases equal to 1, which implies that only the parent phase was active. Consequently, the number of phases indicated by 7 implies six activated child phases. As shown in Fig. 10 , the number of activated child phases varies spatially inside different grains as the results of different initial orientations. However, the activation of child phases significantly depends on temperature: child phases were clearly observed in most grains at room temperature but existed only in a few grains at 523 K. Consequently, it was concluded that activation of child phases became less significant as temperature increased.
Based on the FE simulation, the relative activity of each slip and twin system was calculated during the entire loading history in order to better understand the contribution of slip and deformation twinning to macroscopic-plastic deformation. For the slip and twin modes, the relative activity in an activated phase was defined as: where r α and r β represent the relative activity of a slip mode and twin mode, respectively. To obtain the relative activity of a slip or twin mode in the studied sample of polycrystalline Mg, a volume fraction averages of r α or r β were calculated over all phases of all grains in the RVE FE model presented in Section 3 . Evaluation of the relative activity of different slip and twin modes in polycrystalline Mg at room temperature and T = 523 K with strain is shown in Fig. 11 . Here, zero indicates that the particular slip or twin mode was not activated, while the value of one represents that the specific slip or twin mode is the sole contributor to plastic deformation of the component. At room temperature, it was found that the plastic deformation of polycrystalline Mg was predominately accommodated by basal slip due to its low slip resistance, with other slip and twin modes contributing considerably less. Especially, the contribution from pyramidal c + a could be neglected. The prismatic and pyramidal a had nearly equal contribution, which may be due to the same model parameters being used for the two slip systems in this paper. Compared with the slip modes, the twin modes had less contributions to plastic deformation of polycrystalline Mg; however, the relative activity of TT was significant, comparable with those of prismatic and pyramidal slips. As a result, a large number of child phases are observed in Fig. 10 (a). During the loading history, the relative activity of basal slip and TT decreased with tensile strain while other modes exhibited an increasing tendency. The relative activity of TT was much higher than that of CT, though the curves began converging at higher strains. This indicates that the activation of child phases was mainly governed by TT, which explains the high number of activated child phases in Fig 10 (a). The relative activity of slip and twin modes at 523 K is shown in Fig. 11 (b). Here, basal slip played the secondary role in the plastic deformation, while pyramidal a slip became the dominant one, which is consistent with the experimental results of Ono et al. (2004) . That is because basal slip is temperature-independent but slip resistance of the pyramidal a mode decreased significantly with increasing temperature (see Fig. 7 ). Pyramidal c + a slip also contributes substantially at high temperature when compared to room temperature due to a similar reason. It is interesting to note that there was no significant change of the relative activity of prismatic slip with the temperature increase. Since TT is temperature independent, it was difficult to activate at higher temperatures, leading to fewer child phases observed in Fig. 10 (b). This conclusion is in qualitative agreement with that drawn from the experimental results of Figueiredo et al. (2016) .
To better understand the effect of temperature on the microstructural evolution in polycrystalline Mg, pole figure were used to represent texture evolution at different tem peratures in unidirectional tension. The related results are shown in Fig. 12 (a) and (b) for room temperature and T = 523 K, respectively. The initial orientations of grains in the simulations were fully random. As a great number of child phases were induced with accumulation of the TT volume fraction at room temperature (see Fig. 10 (a) and Fig. 11 (a)), a significant difference was observed between deformed and initial textures as shown in Fig. 12 (a). More grains were found to orient perpendicular to the loading direction (i.e. X axis) in the deformed texture. This may be explained with observations from the unidirectional tension studies of single-crystal Mg (results discussed in Fig. 6 ), which indicate that tension along the c-axis resulted in 90 r otation of Mg crystal. Consequently, in polycrystalline Mg, grains with c-axis oriented nearly parallel to the loading direction experienced a similar reorientation during unidirectional tension as well. However, the TT-induced reorientation became less significant with the increase of temperature as TT activation was restricted at higher temperature. As a result, the deformed texture was relatively close to initial texture at T = 523 K as shown in Fig.  12 (b).

Concluding remarks
A temperature-dependent crystal plasticity (T-SCP) model was proposed and implemented to investigate the dependence of mechanical response and texture of Mg on temperature. Different temperature-dependent deformation modes, i.e. crystalline slip and twinning were incorporated into the constitutive model. With the development of this T-SCP model, a bottom-up modelling approach was proposed to study deformation mechanisms in Mg from single crystals to a polycrystalline aggregate at meso-scale. Herein, the relevant model parameters were calibrated by comparing the simulation results against experimental data for single-crystal Mg; the calibrated parameters were subsequently used to predict the effects of tem perature on both mechanical response and texture evolution of polycrystalline Mg by employing a suitable RVE model.
The study indicates that different deformation modes exhibit a distinct dependence on temperature. A physically relevant understanding of polycrystalline deformation can be gleaned from the bottom-up modelling approach. It was found that in polycrystalline Mg at meso-scale, non-basal slip plays a significant role in its plastic deformation, with less twin/child phases at higher temperatures, which is an experimentally verifiable fact. Interestingly, the modelling approach was capable to capture accurately the temperature-dependent Hall-Petch effect. According to the present studies, ductility of Mg indeed improves at higher temperature due to an increase in plasticity from non-basal slip. Using this approach, a representative study on the material formability may be conducted as a function of initial component texture. The proposed T-SCP model is capable of fundamental investigations in temperature-dependent behaviour and intelligent material design for Mg alloys and other HCP metals by employing a similar bottom-up approach as outlined in this paper. ence and Technology (India), project MAST, is gratefully acknowledged.