Modeling of anisotropic and asymmetric behaviour of magnesium alloys at elevated temperature coupled with ductile damage

Poor formability of Magnesium alloys at room temperature is due to their Hexagonal Closed Packed (HCP) crystal structure. These materials also have a pronounced Strength Differential (SD) effect. In the present work, an improved constitutive model of thermo-elasto-Dviscoplasticity with mixed nonlinear isotropic and kinematic hardening strongly coupled with isotropic ductile damage is developed. The induced anisotropy as well as tension compression asymmetry are carefully considered including their interaction with thermal effects. The numerical implementation of the developed model into ABAQUS/Explicit FE is made through the user subroutine VUMAT. The proposed model is used to simulate material responses of AZ31 Magnesium alloy during sheet metal forming processes at elevated temperature.


Introduction
Due to the high strength to weight ratio of Hexagonal Closed Packed (HCP) structured metals which are desired in light-weight design, these kinds of materials have become an attractive research focus, meanwhile, modeling the plastic behavior of these materials is becoming a highly challenging task. Plastic deformation of HCP materials can be divided mainly into slip and twinning modes. The activation of these modes is highly depending on the critical stress and loading directions. Twinning is a directional deformation mechanism which means a pronounced strength differential (SD) effect: the compressive strengths are much lower than the tensile strengths [1][2][3][4][5]. Experimental results have shown that mechanical response of magnesium alloys have a strong anisotropy and asymmetry under tension and compression at room temperature [1]. Modeling of the plastic deformation with various phenomenological yield functions has been widely applied to sheet metal forming processes. Various isotropic and anisotropic [2] yield functions have been developed by introducing more coefficients to describe the plastic behavior more accurately. To capture the SD effects in the anisotropic model, Cazacu [3][4][5] has developed two yield functions: (1) by introducing the third stress invariant based on the Drucker's criterion; (2) by introducing a new parameter to control the asymmetry in tension and compression and extend to anisotropy using a liner translation on stress based on Balart's criterion. These two yield functions have been extended by others to describe the anisotropic hardening [6] and pressure sensitive metals [7]. In the meantime, the HCP metals (e.g. Magnesium alloys) often have poor formability at room temperature due to the limited number of active slip systems in their HCP crystal structure. Accordingly, hot sheet metal forming technology is proposed to increase the formability of Magnesium alloys, as a result, the thermomechanical behavior complexity is significantly increased under high temperature. Experimental  temperature affect the anisotropic response and tension compression asymmetry of HCP materials [8,9]. Besides the thermal and SD effect, the induced anisotropy due to the evolution of the texture affects significantly the hardening evolution and yield surface. Few macroscopic approaches have been published to capture this phenomenon, see [10][11][12][13].
In order to describe more accurately the hot sheet metal forming processes of the magnesium alloys, we propose a thermo-elasto-viscoplasticity model with mixed nonlinear isotropic and kinematic hardening coupled with isotropic ductile damage. Distortion of the subsequent yield surface was included for modeling the distortion-induced anisotropy. The formulation of the model is performed in the framework of thermodynamics of irreversible processes with state variables [14] using generalized non-associative theory under finite strains [15]. The coupling with isotropic ductile damage is made in the framework of continuum damage mechanics with the effective variables based on the total energy equivalence assumption.

State potential and associated state relations
In order to ensure thermodynamical admissibility of the formulation, the framework of thermodynamics of irreversible processes is assumed. The Helmholtz free energy, a convex and closed function of strain-like state variables in the effective strain space, is taken as a state potential. Assuming that the plastic strain and hardening do not affect the elastic properties of the material, the state potential can be decomposed into two parts: the thermoelastic part Following the thermodynamics of irreversible processes, the stress-like state variables are obtained: where  is the material density, ET is temperature dependent Young's modulus and  is Poisson ratio, () CT and () QT are the kinematic and isotropic hardening temperature dependent moduli respectively, γ is a material parameter governing the damage effect on the isotropic hardening(see [14] for further details).

Dissipation analysis and evolution equations
In this paper, the dissipative phenomena are described by: (i) a yield criterion or yield surface The usual stress deviator S is replaced by a 'distorted stress'  T is the melting temperature of the material and  is a temperature independent material parameter.

Thermal dissipation analysis
The combination between the first and second principle of thermodynamics supplies the Clausius-Duhem inequality, which leads to the definition of the state relations (equations (3)- (7)) and the residual (or dissipation) inequality defined by equations (32) to (34). It is assumed that the mechanical and the thermal dissipation are uncoupled, so the dissipation analysis can be split into two parts mech  and the  : 0 :: The heat flux vector q can be obtained from Fourier potential using the classical linear heat theory, k is the heat conduction coefficient. q kg  (35) The generalized heat equation can be obtained by using this equation in conjunction with the first law of thermodynamics [14]. It will be used for solving thermal problem. : :

Applications
The developed model has been implemented in ABAQUS/Explicit finite element code using the VUMAT user material subroutine. The numerical integration algorithm developed in this routine is based on elastic prediction and viscoplastic correction with radial-return mapping algorithm and considering adiabatic assumption for the thermal coupling. In the following section this routine is applied within parametric study of the local response using only one element.

Parametric study: effect of the asymmetry and distortion
To examine the effect of asymmetry parameter w of the proposed model without damage full coupling, we assumed the material parameters under isotropic mode, leading to F=G=H=0.5; L=M=N=1.5. The assumed model parameters can be found in Table 1. Figure 1 shows the effects of asymmetry parameter w. When w is greater than zero, the yield stress in tension is greater than that in compression, and with the increase of the value of w, the difference between yield stress in tension and compression increase. Similarly, we can find the same tendency when w is negative. The SD effect is correctly captured with the proposed model. In order to examine the distortional parameters, we use the same conditions as mentioned above, and compare the initial yield surface to the subsequent surface after 10% tension pre-straining on the principal stresses plane. Figure 2a shows the effect of the parameter 1 c l X , it is found that 1 c l X controls the distortional ratio of the yield surface without changing the size of the yield surface. Figure 2b shows that with the value of 2 cross size of the yield surface orthogonal to the loading direction. The details of the parameters study can be found in the work of H.Badreddine [12] . The damage effect is studied by comparing the yield surface evolutions for uncoupled and fully coupled cases in uniaxial tension as shown in Figure 3. The evolution curves of Cauchy stress, kinematic hardening back stress and isotropic hardening stress with plastic strain under coupled and uncoupled conditions are shown in Figure 3a, when the damage effect increases to fracture, all the three stresses decrease. It is remarked that the yield surface reaches a maximum surface at 15% strain after which the surface size decreases and the surface center moves towards the reference origin, for the fully coupled case showed in Figure 3c. When the damage value d = 1 the yield surface will reduce to a single point located at the origin.

Application to AZ31 magnesium alloy
In order to validate the proposed model, the experimental results were taken from the work developed by Khan et al., [8]. These results concern the tests which have been done on Magnesium alloy AZ31 under two different temperatures (25°C, 150°C) at strain rate 0.01s -1 on the rolling direction. Note that, the numerically predicted results of Figure 4 are obtained based on the model using the material parameters given in Table 2. Initial isotropy was assumed, F=G=H=0.5; L=M=N=1.5. For the coefficient  ,  =0.9 for Young's modulus;  =1.08 for damage factor S, for other modules  =1.03.

Conclusions
An advanced constitutive model of thermo-elasto-viscoplasticity fully coupled with nonlinear isotropic and kinematic hardening and isotropic damage is developed. The induced anisotropy by distortion of yield surface is considered and the strength differential effect is carefully taken into account to extend the proposed model for the HCP materials in sheet metal forming processes. The capabilities of the proposed model have been investigated through parametric study and validation with some experimental results. Future applications will be made to hot sheet metal forming with higher temperatures and more complex loading paths.