JMEPEG ÓASM International

DOI: 10.1007/s11665-016-2149-6 1059-9495/$19.00

Damage and Failure Analysis of AZ31 Alloy Sheet

in Warm Stamping Processes
P.J. Zhao, Z.H. Chen, and C.F. Dong

(Submitted March 17, 2016; in revised form May 3, 2016)

In this study, a combined experimental-numerical investigation on the failure of AZ31 Mg alloy sheet in the
warm stamping process was carried out based on modified GTN damage model which integrated Yld2000
anisotropic yield criterion. The constitutive equations of material were implemented into a VUMAT sub-
routine for solver ABAQUS/Explicit and applied to the formability analysis of mobile phone shell. The
morphology near the crack area was observed using SEM, and the anisotropic damage evolution at various
temperatures was simulated. The distributions of plastic strain, damage evolution, thickness, and fracture
initiation obtained from FE simulation were analyzed. The corresponding forming limit diagrams were
worked out, and the comparison with the experimental data showed a good agreement.

It has been well known that metallic ductile failure process

Keywords AZ31 Mg alloy, damage failure, GTN model, warm
stamping, Yld2000 usually consists of three stages, i.e., micro-voids nucleation,
micro-voids growth, and micro-voids coalescence. McClintock
(Ref 7) took the lead in proposing a micro-void model with
cylindrical voids in the ductile matrix and deduced an
analytical expression of damage evolution in order to describe
1. Introduction this damage and fracture behavior. Rice and Tracey (Ref 8)
investigated the evolution of the growth of spherical voids and
conducted numerous studies on the voids growth mechanism.
Magnesium (Mg) alloys have many merits, such as good
Then Gurson (Ref 9) proposed a classical voids model
thermal conductivity, high specific strength, and low density
through a micromechanical approach, which describes the
(Ref 1). In recent years, they have received increased attention
growth of spherical cavities and its impacts on the plastic
for their potential as structural materials in a variety of
material formability. Tvergaard and Needleman (Ref 10)
applications, from automotive and aerospace industries to
corrected the original Gurson damage model by adding a
consumer products (Ref 2). Mg alloy with a hexagonal close-
specific formulation for the material ductile fracture with the
packed (HCP) crystal structure generally exhibits very poor
consideration of interactions between the neighboring voids.
formability at room temperature because the critical resolved
Zhang (Ref 11) proposed a complete model based on Gurson
shear stress for a basal slip is far below those for nonbasal
model, which can predict both the homogenous and the
slips (Ref 3). Zhang (Ref 4) investigated the mechanical
localized dilatational deformation phases. Now the GTN
properties of twinned AZ31 alloy plates at lower elevated
model has already become one of the most classical
temperature by means of tensile testing. Although Mg alloy
micromechanical damage models for ductile metal material.
shows poor formability under lower temperature, the plastic
Many improved damage models have been proposed based on
formability of the alloy can be enhanced significantly by
original Gurson model and successfully used to explain the
increasing the forming temperature. It is reported that Mg
damage and failure behavior of metal. Nahshon and Xue (Ref
alloy sheets have good stretchability and drawability at
12) modified original GTN model to study the crack
temperature between 200 and 300 °C (Ref 5). In addition to
development under the shear-dominated loading. Le Maout
the processing parameters, the initial texture also effects the
et al. (Ref 13) studied on the hemming process by applying
formability of Mg alloys. Zhang (Ref 6) conducted the
the modified GTN model. Jackiewicz (Ref 14) used the
investigation of drawing deformation behavior, microstructure
modified GTN model to simulate the damage distribution and
evolution, and basal slip change of AZ31 alloy sheets with
ductile fracture under the different stress triaxiality conditions.
different initial textures during deep drawing. Nowadays, an
It can be found that the GTN model has been widely applied
increasing number of Mg alloys have been applied to
to explain the damage evolution and failure behavior in
manufacture mechanical components, of which AZ31 alloy
ductile metal material. But due to the assumption on the
has been considered as one of the most suitable materials for
isotropic matrix material, it is not precise enough for the Mg
the sheet stamping process.
alloy sheet deformation process.
In this paper, Yld2000 anisotropic yield criterion was
integrated with GTN damage model to describe the plastic
anisotropy formability of the metal material. The anisotropic
P.J. Zhao and Z.H. Chen, School of Mathematics and Physics,
University of Science and Technology Beijing, Beijing 100083, China; plastic behavior of AZ31 metal sheet was considered in the
and C.F. Dong, Institute for Advanced Materials and Technology, original damage model. The modified model was programmed
University of Science and Technology Beijing, Beijing 100083, China. by the user-defined subroutine VUMAT provided by the
Contact e-mails: [email protected] and [email protected].

Journal of Materials Engineering and Performance

commercial nonlinear FE software pack ABAQUS/Explicit to where
predict the damage and failure behavior of AZ31 alloy in the
warm stamping progress at elevated temperature. The com- ( a
bined approaches of explicit dynamic solver and implicit stress U0 ¼ X10  X20
a a ðEq 5Þ
integration algorithm provide a better prediction on the damage U00 ¼ 2X200 þ X100 þ 2X100 þ X200
evolution and the forming process of metal material. The
distributions of plastic strain, damage evolution, thickness, and with the parameter a ¼ 8 for FCC materials. The principal
0 00
fracture initiation obtained from FE simulation were analyzed values X1;2 and X1;2 are dependent on the stress deviator X 0
and compared with experiments. 00
and X , respectively. The linear transformation can be ex-
pressed as follows

2. Constitutive Model X0 ¼ L0  r
ðEq 6Þ
X00 ¼ L00  r
2.1 Anisotropic Yield Function The coefficients of matrix L0 and L00 representing linear
In 1989, Barlat and Lian (Ref 15) proposed a yield function transformations of the stress tensor are
suitable for describing the behavior of anisotropic sheet metals 0 1 0 1
under a plane stress state by taking the planar anisotropy and L011 2= 0 0
B 0 C B 3
plane stress state into account. Hence, this function could be B L12 C B 1 C0 1
B C B  =3 0 0CC a1
used to evaluate the effect of yield surface shapes and the B 0 C B CB C
behavior of materials formability during forming processes for BL C ¼ B 0 1=3 0 C @ a2 A;
B 21 C B C
plane stress problems. A simple yield function can be expressed B 0 C B C
BL C @ 0 2= a7
in the following form: @ 22 A 3 0A
L066 0 0 1 ðEq 7Þ
0 00 1 0 10 1
M M M M L11 2 2 8 2 0 a3
ajK1 þ K2 j þ
U¼ ajK1  K2 j þcj2K2 j 2q ¼ 0 ðEq 1Þ B L00 C B CB C
B 12 C B 1 4 4 4 0 CB a4 C
where K1 and K2 are defined as B 00 C 1 B CB C
B L21 C ¼ B 4 4 4 1 0 CB a5 C
B C 9B CB C
B 00 C B CB C
8 @ L22 A @ 2 8 2 2 0 A@ a6 A
> 1 
>
>
< K1 ¼ rx þ hry L0066 0 0 0 0 9 a8
2
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
 2 ðEq 2Þ The independent coefficients ak (for k : 1  8) are the
>
> rx  hry 
>
: K2 ¼ þp2 r2xy material parameters used to describe the anisotropy of the
2 material. It should be noted that these coefficients all reduce to
a, c, 
and  h , and 
p are the anisotropy coefficients. one for the isotropic material. These coefficients employed in
the new anisotropic yield function and the coefficients of
8 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Barlat89 yield function can be transformed each other (Ref 17),
>
> 
h¼r0 ð1 þ r90 Þ=½ð1 þ r0 Þr90  it may be written as
>
> pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
>
< c ¼ 2 r0 r90 =½ð1 þ r0 Þð1 þ r90 Þ 8  
> a ¼ 2  c
ðEq 3Þ >
> a1 ¼ 2   h ðcÞ1=a
>
> >
>
> >
>  
>
:    1 >
>
> a2 ¼ 2 h  1 ðcÞ1=a
p ¼ rp =rb ¼ ðq=ss Þ 2= 2a þ 2M c M >
>
>
< a3 ¼ a4 ¼  aÞ1=a
hð
with r0 , r45 , r90 as the plastic anisotropy parameters in the ðEq 8Þ
rolling, 45°, and transverse directions, respectively. q denotes >
> aÞ1=a
a5 ¼ a6 ¼ ð 
>
>
effective stress. rb , rp denote Cauchy principal stress under >
>
>
uniaxial tensile state and under biaxial tensile state, respec-
>
>
> a ¼pðcÞ1=a
7
>
>
:
tively. ss is yield shear stress under pure shear state. M is a a8 ¼  aÞ1=a
pð
yield criterion index, M ¼ 6 for face-centered cubic (FCC)
materials.
In 2003, Barlat et al. (Ref 16) proposed a new anisotropy 2.2 The GTN Damage Model
yield criterion which is suitable for the plane stress problem
The yield function of GTN damage material model (Ref 10)
with a higher precision and more convenient application. The
can be written as
expression of yield function used two linear transformations
about the Cauchy stress tensor to introduce the anisotropy. In
 
the case of plane stress state, the yield function is expressed in q2  3q2 p  
the general form U¼ þ 2q 1 f cosh  1 þ q3 f 2 ¼ 0 ðEq 9Þ
r2y 2ry

where ry is the equivalent stress of the matrix material which

U ¼ U0 þ U00 ¼ 2qa ðEq 4Þ
represents the actual microscopic stress state, the hydrostatic

Journal of Materials Engineering and Performance

stress p¼ 13 rkk . Three additional fitting parameters q1 , q2 , q3 tensor, and incremental strain tensor, respectively.
introduced by Tvergaard can better predict void interaction Hta ða ¼ 1; 2Þ denotes a state variable with two parame-
effects by Gurson damage equations. ters f and epl
m . Void volume fraction f and matrix equiva-
The damage parameter f  denotes the total effective void lent plastic strain epl
m are treated as two scalar internal
volume fraction, which is a function of the void volume fraction variables.
f . It accounts for the gradual loss of stress carrying capability of (2) Trail elastic state is obtained by assuming that the strain
the material due to void coalescence. If f  is equal to zero, the increment is purely elastic
GTN yield condition reduces to the original yield condition. If retþDt ¼ rt þ C : DetþDt ðEq 15Þ
f  is equal to unity, the material is completely a failure and has
lost load carrying capacity. The porosity is defined as follows where superscript e indicates trail elastic state. C is the
fourth-rank elastic modulus.Calculate the hydrostatic
8 stress.
>
< f f  fc 1

f ¼ fc þ dðf  fc Þ fc < f  fF ðEq 10Þ petþDt ¼  retþDt : I ðEq 16Þ
> 3
: fF f  fF Calculate the equivalent stress of Yld2000-GTN model.
1

e 1 0 00
a
in which qtþDt ¼ ðU þ U Þ ðEq 17Þ
2
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi (3) Calculate the yield function to determine the current
fF  fc q1 þ q21  q3
d¼ 
; fF ¼ ðEq 11Þ state.
fF  fc q3    
UetþDt ¼ U retþDt ; Hta ¼ U petþDt ; qetþDt ; Hta ðEq 18Þ
where fc is a critical value of the void volume fraction at
which void coalescence first occurs, fF is the fracture value of If UetþDt  0, the current time step is elastic,
void volume fraction at which the material will lose the load rtþDt ¼ retþDt , go to step 5.If UetþDt > 0, the step is plas-
carrying capacity completely and the cracks occur. tic, go to Step 4.
The increased rate of total void volume fraction f_ is partly (4) Plastic calculation (omit the subscript t þ Dt to simplify
due to the growth of existing voids f_growth and partly due to the the functions)
nucleation of new voids f_nucl as
(a) Flow direction:
3
f_ ¼ f_growth þ f_nucl ðEq 12Þ n ¼ e Se ðEq 19Þ
2q
The growth rate of voids f_growth is proportional to the (b) Assume crh and creq are the correction of Derh and
hydrostatic component of the plastic strain rate e_ pkk , as follows Dereq , respectively.
(
A11 crh þ A12 creq ¼ b1
f_growth ¼ ð1  f Þ_epkk : ðEq 13Þ ðEq 20Þ
A21 crh þ A22 creq ¼ b2
The nucleation rate of new voids can be expressed by a
plastic strain-controlled nucleation rule by assuming that voids where Aij is stiffness matrix coefficient.
nucleate at second phase particles and there exists a normal (c) Ensure the following functions were satisfied.
distribution of nucleation strain for the total population of @U kþ1 @U
f1 ¼ kþ1 Dep þ Dep ¼0 ðEq 21Þ
particles (Ref 18) @q @p
 
  f2 ¼ U kþ1 p; kþ1 q; kþ1 H a ¼ 0 ðEq 22Þ
fN 1 ep  eN
f_nucl ¼ pffiffiffiffiffiffi exp  e_ pl
m ðEq 14Þ
sN 2p 2 sN drkþ1
kþ1
rm ¼ m
epl
m ðEq 23Þ
where fN represents the volume fraction of void-nucleating depl
m
particles, eN and sN are the average and standard deviation of
Newton-Raphson iterative method is used to solve the nonlin-
the strains at which particles nucleate voids.
ear equations, Eq 22 and 23. Iteration continues until both
2.3 Numerical Algorithms jf1 j and jf2 j < Tolerance (setting to 1:08 ), which means the
Newton method convergence is achieved, then go to Step 5,
The modified GTN damage model for ductile metal material otherwise end the calculation.
was implemented into the solver software pack ABAQUS/Explicit
by developing the user material subroutine VUMAT, which was in (5) Update p, q, H a
line with the the implicit stress update algorithm introduced by (
Aravas (Ref 19). The details of implementation procedures of the p ¼ pe þ KDeP
ðEq 24Þ
stress update algorithm are illuminated as follows: q ¼ qe  3GDeq

where E and v are YoungÕs modulus and PoissonÕs ratio,

E
(1) Obtain initial values rt , et , Hta , DetþDt ; ðt ¼ 0; . . . ; ti Þ respectively. The shear modulus G ¼ 2ð1þ vÞ, the bulk
where rt , et , and DetþDt denote the stress tensor, strain E
modulus K ¼ 3ð12 vÞ

Journal of Materials Engineering and Performance

 Start of VUMAT

Known σ t , ε t , H tα , Δε t +Δt

Calculate σ t +Δt , pt +Δt , qt +Δt

e e e

Yes
φ ( pte+Δt , qte+Δt , H tα ) ≤ 0
No
Initial value k = 0, ΔH α = 0

Calculate f1 , f 2 , k +1σ m
Newton
iterative f1 ≤ toler Yes Update stress σ
t +Δt
and f 2 ≤ toler
Calculate Δε m ( t +Δt ) , Δf t +Δt
pl
No
Yes k≤N Fig. 2 Three kinds of tensile specimens
No Update state variables
Iteration not converge Δε mpl( t +Δt ) , f t +Δt , f t*+Δt
stop ABAQUS Table 1 Chemical composition of the AZ31alloy (mass%)
End of VUMAT
Element Mg Al Zn Mn Fe Others
Fig. 1 VUMAT subroutine flow chart diagram for modified damage
Content 95 2.9 0.6 0.2 0.005 0.45
model

2
r ¼ pI þ qn ðEq 25Þ
3
8
>
< pDep þ qDeq
DH 1 ¼ Depl
m ¼
ð1  f Þry ðEq 26Þ
>
: DH 2 ¼ Df ¼ ð1  f ÞDe þ ADepl
p m
(
H 1 ¼ epl epl
m ¼ epl
mðtÞ þ Dm
ðEq 27Þ
H 2 ¼ f ¼ ft þ Df
If H 2 ¼ f > fc , go to calculate f  ¼ fc þ jðf  fc Þ
(6) Go to the next time step.
The flow chart diagram of stress integration algorithm
for modified damage model is shown in Fig. 1.

3. Experimental Work

3.1 Mechanical Properties Tests at Elevated Temperature

The tensile test is the one of the most fundamental
experiments to study the formability for sheet metal material.
In the present paper, tensile tests at different temperatures were
carried out. The tensile specimens cut from of the AZ31 alloy
sheet with 0.6 mm thickness were made to determine the stress-
strain relationship. The tensile specimens were prepared from a Fig. 3 MTS810 test machine
cold-rolled AZ31 alloy sheet cutting along rolling direc-
tion(RD), 45° direction, and transverse direction(TD), as The true stress-strain curves of AZ31 Mg alloy at elevated
illustrated in Fig. 2. Table 1 summarizes the chemical compo- temperature under the fixed strain rate of 0.001 s1 are shown
sition of the AZ31 alloy. in Fig. 4. It can be observed that the yield stress decreased
A heating furnace equipment was mounted on the MTS810 significantly with the increasing of temperature, and the AZ31
test machine in order to conduct tensile tests at various alloy exhibits the better formability at higher temperature. It is
temperatures, as shown in Fig. 3. Before the tensile tests were also clearly seen that the harden effect of material becomes
performed, Mg alloy specimens were heated to 100, 200 °C, more insignificant with the increasing of test temperature. The
and finally 300 °C due to the poor formability for Mg alloy at mechanical properties of AZ31 alloy sheet at a variety of
room temperature. The temperature of AZ31 sheet kept temperatures under the fixed strain rate of 0.001 s1 are shown
stable until the specimens were stretched to crack. in Table 2. The following Table 3 lists the parameters of the

Journal of Materials Engineering and Performance

anisotropic yield criterion Yld2000 of AZ31 Mg alloy obtained tion of the model parameters, the modified GTN model
by the algebraic manipulation in section 2.1. implemented in the commercial finite element software solver
GTN damage parameters of AZ31 alloy are listed in Table 4 ABAQUS/Explicit was used to predict the formability and
according to the reference (Ref 20). To verify the parameters, failure behavior of mobile phone shell. The schematic illustra-
the whole tensile test process for AZ31 Mg sheet was simulated tion of warm stamping test is shown in Fig. 6. The geometry
using the modified damage model. Fig. 5 shows the compar- size of the AZ31 alloy sheet is 40 9 140 mm2 with the
ison of the experimental and numerical results on the displace- thickness of 0.6 mm. Table 5 lists the processing parameters of
ment-load curves of the AZ31 alloy sheet at room temperature. sheet warm stamping. The other simulation parameters are
It can be observed that the modified GTN model with these blank holder force of 2.5 kN, stamping velocity of 3 mm/s, and
parameters are effective for simulating metal tensile process. die clearance of 0.7 mm on each side. In order to reduce CPU
time, only one half of the specimen is modeled due to the
3.2 Warm Stamping Forming of Mobile Phone Shell symmetry of the specimen geometry and loading conditions.
The reduced integration 4-node shell elements S4R was
Warm stamping is a widely applied processing technology
employed to model the sheet. There were five integration
especially in the field of metal sheet forming. After identifica-
points through the thickness. The blank is set as deformable
body and the other tools including the punch, the die, the holder
were taken as rigid bodies. Movements of the three tools are
represented by the corresponding reference points.
Table 3 Parameters of the anisotropic yield criterion
Yld2000 of AZ31 alloy
Temperature, °C a1 a2 a3 a4 a5 a6 a7 a8

100 1.154 0.879 0.831 0.831 0.909 0.909 1.009 0.864

200 1.123 0.926 0.862 0.862 0.919 0.919 1.004 0.873
300 0.955 1.099 1.045 1.045 0.997 0.997 0.952 0.947

Table 4 Damage model parameters of AZ31 alloy sheet

q1 q2 q3 fN fC fF eN sN f0

1.5 1 2.25 0.04 0.0732 0.1279 0.2 0.1 0.001

Fig. 4 True stress-strain curves at elevated temperature

Table 2 Mechanical properties of AZ31 alloy sheet

Temperature, °C Elastic modulus, GPa PoissonÕs ratio Initial yield stress, MPa r0 r45 r90

100 43 0.35 150 1.91 2.40 3.66

200 41 0.35 105 1.90 2.11 2.92
300 38 0.35 50 1.13 0.99 0.94

Fig. 5 Displacement-load curves of AZ31 alloy sheet and tensile test

Journal of Materials Engineering and Performance

4. Results and Discussion its critical value, material would lose the load-carrying capacity
and the crack occurs.
Warm stamping experiments at various forming tempera- In order to study the development of damage for AZ31 Mg
tures was performed. The experimental results indicate that the alloy sheet during the warm stamping process, three groups of
sheet surface becomes severely oxidized when the initial simulations at different temperatures were performed. Figure 9
temperature exceeds 350 °C, which leads to a decline in shows the distributions of equivalent plastic strain and void
performance and surface quality of the plate. Given this, three volume fraction at elevated temperature. It is found that at
representative forming temperatures were selected and ana- 100 °C (Fig. 9a) damage and fracture occur at the fillet corner
lyzed, as shown in Fig. 7. From Fig. 7(a), a wide range of crack which makes contact with the punch, as void volume fraction of
area was formed at 100 °C, where the blank in the fillet corner AZ31 alloy sheet is beyond the fracture value fF. The failure
cracks extensively. As is shown in Fig. 7(b), the crack appears region is relatively larger, and crack band propagates along the
in the fillet region as the forming temperature increases to sheet corner wall with the increase of punch displacement. At
200 °C. Figure 7(c) reveals that the qualified box-shaped pieces 200 °C (Fig. 9b), the value of void volume fraction near the
could be achieved as the temperature rises to 300 °C. It can be corner exceeds its critical value extensively. Moreover, fracture
concluded that the most suitable stamping temperature of the bands are found at the workpiece corners. But it is obvious that
AZ31 alloy sheet is 300 °C. the failure region at 200 °C is smaller than that at 100 °C. At
Figure 8 shows the morphology near the area of crack 300 °C (Fig. 9c), the critical value of void volume fraction is
observed by SEM at 100 °C. Before the stamping process, the not achieved and therefore material failure does not occur
matrix of material contains a small amount of voids. It is noted during the whole stamping process. It is concluded that the
that there are a series of micro-voids coalescence when the distribution of the damage and failure region obtained by
crack initiation occurs, which indicates that the damage in this numerical simulation at three kinds of temperatures accurately
region accumulates rapidly. If the accumulated damage reached agrees with the experimental results, as shown in Fig. 7.
Figure 10(a) and (b) illustrates the curves of equivalent
plastic strain and damage value versus curvilinear distance in
rolling and transverse directions at the forming temperature of
100 °C. It can be seen that the maximum plastic strain and
damage value in transverse direction predicted by modified
GTN model are lightly higher than those in transverse direction.
The damage could be insignificant at the bottom of the metal
shell body. Damage growth and plastic deformation occur
rapidly adjacent to the punch corner and then exceed the critical
value near the region. Fig. 10(c) shows the thickness distribu-
tion obtained by the experiment and numerical simulation. The
tendency of the percentage changes in thickness predicted by
simulation matches with experimental prediction. Nevertheless,
Fig. 6 Schematic illustration of sheet stamping process minimum thickness appeared at the corner walls since the
plastic flow of material was more restricted compared to that at
straight edges. The maximum thinning of 27% in experiment
Table 5 Processing parameters of warm stamping was observed as compared to 21% predicted by simulation.
Friction coefficient between dies and blank 0.15 The region of fracture obtained from experiment was
Temperature of environment, °C 20 compared with FE simulation at forming temperature of
Temperature of punch, °C 50 100 °C, as shown in Fig. 11. It is shown that a good agreement
Initial temperature of dies and blank, °C 100 200 300 can be found between the experimental and numerical results.
Heat transfer coefficient (blank-dies), W/(m2°C) 2000 The potential failure areas can be identified according to the
Heat transfer coefficient (blank-environment), W/(m2°C) 20 predicted damage distribution in the blank.
Inelastic heat fraction 0.9
In the warm stamping test, the risk of fracture can be
Film coefficient, W/(m2°C) 0.025
predicted by the forming limit curves (FLC). In order to

Fig. 7 Warm stamping results at initial temperature of (a) 100 °C; (b) 200 °C; (c) 300 °C

Journal of Materials Engineering and Performance

Fig. 8 SEM morphology near the crack area. (a) The crack initiation state of material; (b) The original state of material

Fig. 9 Equivalent plastic strain and void volume fraction at initial temperature of (a) 100 °C; (b) 200 °C; (c) 300 °C. (a) Equivalent plastic
strain; void volume fraction. (b) Equivalent plastic strain; void volume fraction. (c) Equivalent plastic strain; void volume fraction

Journal of Materials Engineering and Performance

Fig. 10 (a) Curves of equivalent plastic strain of Mg sheet at 100 °C; (b) curves of void volume fraction of Mg sheet at 100 °C; (c) thickness
distribution obtained from FE simulation and experiment at 300 °C

Fig. 11 Fracture of AZ31 alloy sheet at 100 °C: comparison be-

tween simulation and experiment

evaluate the probability of crack occurrence in the risky zones, Fig. 12 Forming limit curves at three temperatures
the forming limit diagrams (FLD) of the AZ31 alloy sheet at
three kinds of temperatures were measured by Li (Ref 21) as fracture region at 100 °C. In another words, the material failure
seen in Fig. 12. The numerical points obtained from the crack occurs in this region. At 200 °C, formability of AZ31 alloy
band were drawn. The points located above the FLC indicate sheet is significantly improved and a small amount of points are
that the material has failed. Figure 12 demonstrates that most of still in the fracture region. But as the forming temperature
strain values in the fillet zone of AZ31 Mg alloy sheet are in the increases to 300 °C, the whole zone including fillet corner is

Journal of Materials Engineering and Performance

qualified in that all points are located in the safe region of FLD. 4. H. Zhang, Y. Liu, J. Fan, H.J. Roven, W. Cheng, B. Xu, and H. Dong,
The comparison between experimental data and the FLD Microstructure Evolution and Mechanical Properties of Twinned AZ31
Alloy Plates at Lower Elevated Temperature, J. Alloys Compd., 2014,
findings shows a good agreement. 615, p 687–692
5. F.-K. Chen and T.-B. Huang, Formability of Stamping Magnesium-
Alloy AZ31 Sheets, J. Mater. Process. Technol., 2003, 142(3), p 643–
647
5. Conclusions 6. H. Zhang, G. Huang, J. Fan, H.J. Roven, B. Xu, and H. Dong, Deep
Drawability and Drawing Behaviour of AZ31 Alloy Sheets with
In this research, the modified GTN model combined with Different Initial Texture, J. Alloys Compd., 2014, 615, p 302–310
7. F.A. McClintock, A Criterion for Ductile Fracture by the Growth of
Yld2000 anisotropic plastic yield criterion was developed via Holes, J. Appl. Mech., 1968, 35(2), p 363–371
the user material subroutine VUMAT provided by the com- 8. J.R. Rice and D.M. Tracey, On the Ductile Enlargement of Voids in
mercial finite element solver ABAQUS/Explicit. During the Triaxial Stress Fields, J. Mech. Phys. Solids, 1969, 17(3), p 201–217
warm stamping process of the mobile phone shell, the fillet 9. A.L. Gurson, Continuum Theory of Ductile Rupture by Void
corner region of AZ31 alloy sheet was most prone to crack, due Nucleation and Growth: Part I—Yield Criteria and Flow Rules for
to the fact that the damage evolution was accompanied with the Porous Ductile Media, J. Eng. Mater. Technol., 1977, 99(1), p 2–15
10. A. Needleman and V. Tvergaard, An Analysis of Dynamic, Ductile
micro-voids growth and coalescence which led to the fracture Crack Growth in a Double Edge Cracked Specimen, Int. J. Fract.,
in specimen. The formability of AZ31 Mg alloy sheet could 1991, 49(1), p 41–67
reach to the optimal level as the forming temperature of the 11. Z.L. Zhang, C. Thaulow, and J. Ødegård, A Complete Gurson Model
specimen was elevated to 300 °C. The distributions of plastic Approach for Ductile Fracture, Eng. Fract. Mech., 2000, 67(2), p 155–168
strain, damage evolution, thickness, and fracture initiation 12. K. Nahshon and Z. Xue, A modified Gurson Model and Its Application
obtained from FE simulation were analyzed. The comparison to Punch-Out Experiments, Eng. Fract. Mech., 2009, 76(8), p 997–1009
13. N. Le Maoût, S. Thuillier, and P.-Y. Manach, Aluminum Alloy Damage
between the experimental findings and the data through FLC
Evolution for Different Strain Paths–Application to Hemming Process,
shows a good agreement. Thus, the modified GTN damage Eng. Fract. Mech., 2009, 76(9), p 1202–1214
model provides a promising tool to predict the formability and 14. J. Jackiewicz, Use of a Modified Gurson Model Approach for the
failure behavior of AZ31 alloy sheet at various temperatures by Simulation of Ductile Fracture by Growth and Coalescence of
describing both anisotropic plastic deforming and damage Microvoids Under Low, Medium and High Stress Triaxiality Loadings,
Eng. Fract. Mech., 2011, 78(3), p 487–502
evolution during the metal sheet forming process.
15. F. Barlat and K. Lian, Plastic Behavior and Stretchability of Sheet
Metals. Part I: A Yield Function for Orthotropic Sheets Under Plane
Stress Conditions, Int. J. Plast., 1989, 5(1), p 51–66
Acknowledgment 16. F. Barlat, J.C. Brem, J.W. Yoon, K. Chung, R.E. Dick, D.J. Lege, F.
Pourboghrat, S.H. Choi, and E. Chu, Plane Stress Yield Function for
The present study was supported by the National Natural Aluminum Alloy Sheets—Part 1: Theory, Int. J. Plast., 2003, 19(9), p
Science Foundation of China under Grant No. 51222106. 1297–1319
17. F. Barlat, J.W. Yoon, and O. Cazacu, On Linear Transformations of
Stress Tensors for the Description of Plastic Anisotropy, Int. J. Plast.,
2007, 23(5), p 876–896
References 18. C. Chu and A. Needleman, Void Nucleation Effects in Biaxially
Stretched Sheets, J. Eng. Mater. Technol., 1980, 102(3), p 249–256
1. H. Furuya, N. Kogiso, S. Matunaga, and K. Senda, Applications of 19. N. Aravas, On the Numerical Integration of a Class of Pres-
Magnesium Alloys for Aerospace Structure Systems, Mater. Sci. sure—Dependent Plasticity Models, Int. J. Numer. Methods Eng.,
Forum, 2000, 350–351, p 341–348 1987, 24(7), p 1395–1416
2. R. Beals, C. Tissington, X. Zhang, K. Kainer, J. Petrillo, M. Verbrugge, 20. R.-Z. Wang, Z.-H. Chen, Y.-J. Li, and C.-F. Dong, Failure Analysis of
and M. Pekguleryuz, Magnesium Global Development: Outcomes from AZ31 Magnesium Alloy Sheets Based on the Extended GTN Damage
the TMS 2007 Annual Meeting, JOM, 2007, 59(8), p 39–42 Model, Int. J. Min. Met. Mater., 2013, 20(12), p 1198–1207
3. H. Zhang, G. Huang, J. Fan, H. Jørgen Roven, F. Pan, and B. Xu, Deep 21. W. Li, G. Zhao, X. Ma, and J. Gao, Study on Forming Limit Diagrams
Drawability and Deformation Behavior of AZ31 Magnesium Alloy of AZ31B Alloy Sheet at Different Temperatures, Mater. Manuf.
Sheets at 473 K, Mater. Sci. Eng. A, 2014, 608, p 234–241 Process., 2013, 28(3), p 306–311

Journal of Materials Engineering and Performance