Developing a Zener – Hollomon Based Forming Limit Surface for Warm Sheet Forming of Magnesium Alloy

The concept of Zener–Hollomon (Z) based Forming Limit Surface (Z-FLS)—which has minor strain, major strain, and ln(Z) as its three axes—was initially proposed in a previous study. In the Z-FLS diagram, strain rate and temperature effects on the major limit strain are reflected by ln(Z). In the current study, the concept of Z-FLS is revisited to provide a practical approach to construct Z-FLS. A Z-FLS then is constructed for magnesium alloy AZ31B sheet material using available experimental forming limit curves. The constructed Z-FLS is used to identify fracture in a non-isothermal warm cup forming process, which was modeled as a coupled thermo-mechanical process. Based on the Z-FLS, the determined limiting draw ratio (LDR) matches well with the published experimental results.


Introduction
As the lightest structural metal, magnesium alloys are attractive candidates for lightweighting automotive components [1,2].However, magnesium alloys have poor formability at room temperature, which is attributed to the fact that basal slip and twinning of the HCP crystal structure of magnesium are the only two active deformation mechanisms at room temperature.At elevated temperatures, the formability of magnesium alloys is improved substantially due to the activation of additional slip systems, e.g., <c + a> slip, and thermally activated dynamic recrystallization [3].This is the reason for developing warm forming processes for manufacturing magnesium alloy parts [4].Among the wrought magnesium alloys, AZ31 is the most frequently considered sheet material for automotive applications.Considerable increase in the formability of the AZ31 alloy was observed at temperatures ranging from 200 to 350 • C in warm forming of conical cups by [5], circular cups by [6][7][8] and square cups by [9].With improved formability, vehicle deck lid inner panels formed from AZ31 sheet alloy using hot gas pressure forming process at 450 • C have been reported by Carter et al. [10].
In warm forming magnesium alloy sheet materials, the stretching failure is generally regarded due to ductile fracture mechanism, in which fracture happens as a result of void nucleation, growth, and coalescence.Figure 1 gives a SEM morphology of the AZ31 sheet material at temperature of 100 • C [11]. Figure 1a shows that small amount of voids can be observed before deformation.Due to plastic deformation, the voids grow and coalesce at the initiation of fracture as shown in Figure 1b.Finite element method based forming simulations are widely used to design and optimize sheet metal forming processes [12].Fracture and necking failure in sheet metal forming are usually predicted using strain-based Forming Limit Curves (FLCs), in which major local strains at the initiation of localized necking and fracture are plotted against minor local strain [13].Although it is sensitive to the strain path, the convenience of its implementation and lack of a better criterion has made the use of FLC a widely accepted method in both academia and industry for the identification of localized necking and fracture in sheet metal forming [14,15].
In various studies on formability of magnesium alloys in warm forming conditions, it is observed that the levels of FLCs increase with increasing temperature, but decrease with increasing strain rate [5,7,11,[16][17][18][19].Due to variation in tooling speed and heat transfer between the sheet and the tool, warm forming is, in general a non-isothermal process, which creates spatially and temporally different distribution of strain rates and temperatures on the sheet during the forming operation [6][7][8][9].Thus, the prediction of fracture in warm forming requires a complete representation of forming limits that vary with strain rates and temperatures.Three-dimensional forming limit diagrams with either temperature or strain rate as the third axis have been proposed [20,21].A three-dimensional forming limit surface that takes into account both temperature and strain rate was proposed in [1].In this concept, the Zener-Hollomon parameter (Z), or the so-called temperaturecompensated strain rate [22], is used to represent the effects of strain rate and temperature on the forming limit and forms the third axis in the Forming Limit Surface diagram.The three-dimensional FLS, which is termed Z-FLS in the current paper, was based on the observation that the major limit strains on the forming limit curves of aluminum alloy 5083-O sheet at different strain rates and temperatures can be described as a function of ln(Z).In the current study, the concept of Z-FLS is revisited from following two aspects: (1) Introduction of a normalized Z parameter in natural logarithmic scale _____ ) ln( Z , which makes it compatible with the strain values and thus provides a practical and user friendly interface for identifying failure.
(2) Introduction of a construction procedure, in which firstly, major limit strain is represented in a space of _____ ) ln(Z and strain ratio ( ρ ), which is a ratio between minor strain ) ( 2 ε and major limit strain ) ( 1 ε of the measure data point at forming limit; then the major limit strain is transformed into a space of minor strain ) ( 2 ε and _____ ) ln(Z .The construction procedure can avoid reduced amount of data points measured at larger minor strain levels, which happens when Z-FLS is constructed in the space of ) ) ln( , ( 2 Z ε directly [1].For example, in the test data of sheet material Al 5083 [1,23], at
Finite element method based forming simulations are widely used to design and optimize sheet metal forming processes [12].Fracture and necking failure in sheet metal forming are usually predicted using strain-based Forming Limit Curves (FLCs), in which major local strains at the initiation of localized necking and fracture are plotted against minor local strain [13].Although it is sensitive to the strain path, the convenience of its implementation and lack of a better criterion has made the use of FLC a widely accepted method in both academia and industry for the identification of localized necking and fracture in sheet metal forming [14,15].
In various studies on formability of magnesium alloys in warm forming conditions, it is observed that the levels of FLCs increase with increasing temperature, but decrease with increasing strain rate [5,7,11,[16][17][18][19].Due to variation in tooling speed and heat transfer between the sheet and the tool, warm forming is, in general a non-isothermal process, which creates spatially and temporally different distribution of strain rates and temperatures on the sheet during the forming operation [6][7][8][9].Thus, the prediction of fracture in warm forming requires a complete representation of forming limits that vary with strain rates and temperatures.Three-dimensional forming limit diagrams with either temperature or strain rate as the third axis have been proposed [20,21].A three-dimensional forming limit surface that takes into account both temperature and strain rate was proposed in [1].In this concept, the Zener-Hollomon parameter (Z), or the so-called temperature-compensated strain rate [22], is used to represent the effects of strain rate and temperature on the forming limit and forms the third axis in the Forming Limit Surface diagram.The three-dimensional FLS, which is termed Z-FLS in the current paper, was based on the observation that the major limit strains on the forming limit curves of aluminum alloy 5083-O sheet at different strain rates and temperatures can be described as a function of ln(Z).In the current study, the concept of Z-FLS is revisited from following two aspects: (1) Introduction of a normalized Z parameter in natural logarithmic scale ln(Z), which makes it compatible with the strain values and thus provides a practical and user friendly interface for identifying failure.(2) Introduction of a construction procedure, in which firstly, major limit strain is represented in a space of ln(Z) and strain ratio (ρ), which is a ratio between minor strain (ε 2 ) and major limit strain (ε 1 ) of the measure data point at forming limit; then the major limit strain is transformed into a space of minor strain (ε 2 ) and ln(Z).The construction procedure can avoid reduced amount of data points measured at larger minor strain levels, which happens when Z-FLS is constructed in the space of (ε 2 , ln(Z)) directly [1].For example, in the test data of sheet material Al 5083 [1,23], at ε 2 = −0.3,only four major limit strains can be measured from tested FLCs at temperatures of 423, 473, and 573 • K and strain rates of 10 −1 /s, 10 −2 /s, and 10 −4 /s while 12 data points (major limit strains) can be measured at ε 2 = −0.1.
Based on those further developments, a Z-FLS is constructed for AZ31B sheet material by using available forming limit curves and used to identify fracture failure by using FEM simulations of a non-isothermal round cup drawing process.

Zener-Hollomon Based Forming Limit Surface (Z-FLS)
The Zener-Hollomon based forming limit surface (Z-FLS) is a three dimensional representation of forming limits of sheet materials that can be expressed in the following equation form.
In this equation, ε 1 and ε 2 are the major and minor local strains, and Z represents the Zener-Hollomon parameter [22] given by Equation (2). where, .
ε is the strain rate, Q is the activation energy, R is the gas constant, and T is the sheet temperature in K.
The dependency of the effective strain at fracture on ln(Z) was reported in hot tensile/torsion tests on bulk aluminum alloys [24] and warm forming of magnesium sheet material [25].Alexandrov et al. [24] used a polynomial function of the form to represent the effective limit strain at uniaxial tension where A, B, and C are material parameters that are determined by regression analysis of effective limit strains from tensile tests at a few selected temperatures and strain rates.Such an empirical relationship can be observed at other strain paths on the FLCs measured from different temperatures and strain rates on Al5083 sheet material [1,23].Based on this, the major limit strain ε * 1 can be expressed as where A, B, and C are functions of strain ratio ρ, which is a ratio between ε 2 and ε * 1 .They are determined by fitting Equation (4) to experimentally determined FLCs at a few selected temperatures and strain rates and certain strain path, which is represented by ρ.The use of ρ is based on an assumption that sheet material deforms in proportional and quasi-proportional strain paths.Since strains in sheet forming range from −2 to 1 and ln(Z) is a much larger number, a normalized value of ln(Z) was used to make it compatible with the strain values.The normalized value, represented by ln(Z), is calculated by using the equation where, ln(Z 0 ) is the value of ln(Z) at the lowest strain rate and the highest temperature in the selected warm forming experiments, i.e., at 300 • C and strain rate of 10 −4 /s in [26].Then, replacing ln(Z) with ln(Z), Equation (4) becomes Equation ( 6) suggests a Z-FLS in a space of (ε 1 , ρ, ln(Z)).Based on test data in [26], Figure 2 gives major limit strain plots as functions of Equations ( 4) and ( 6) under three critical strain paths of Uniaxial Tension (UT), Plane Strain (PS), and Equal Biaxial Tension (EBT).The analysis shows that the same level of correlation coefficient (R 2 ) can be achieved using either ln(Z) or ln(Z) along three different strain paths.This rationalizes the usage of ln(Z).
the same level of correlation coefficient (R 2 ) can be achieved using either ln(Z) or   6), Z-FLS can be transformed to a space of When working with a finite element forming simulation model to predict forming limit, the risk of fracture can be determined by considering the difference of the major limit strain * 1 ε on the surface of the Z-FLS diagram and the calculated major strains 1 ε in various elements on the model.
In equation form, the difference is given by Δ is higher than or equal to zero, fracture will occur.Given a safe margin of 0.05, if 1 ε Δ is between 0 and −0.05, the element is considered to be at risk of fracture; otherwise, the element is considered safe.

Construction of Z-FLS for AZ31 Magnesium Alloy
Forming limit curves (FLCs) reported in the literature at four different temperatures and two different strain rates were used to construct the Z-FLS for an AZ31 magnesium alloy as shown in Figure 3. Three forming limit curves from [5,9], were generated using Nakajima dome test data on 1.2 mm thick rectangular specimens of AZ31 magnesium alloy sheet under isothermal temperature conditions of 100, 200, and 300 °C and a constant punch speed of 0.1 mm/s [5,9].The obtained FLCs were fitted from measured data from three orientations to the rolling direction, i.e., 0°, 45°, and 90°.In constant speed mechanical stretching tests [26], such as the Nakajima dome test, the strain rate increases with punch travel and its value depends on the punch speed.Albakri et al. [26] have shown that, under balanced biaxial test conditions in a Nakajima dome test, the strain rate in AZ31B sheets reached 0.008 s −1 at a punch speed of 0.1 mm/s and punch displacement of 32 mm.Since the strain rate value was not given in [5,9], it is assumed that the strain rate in their test is approximately 0.008 s −1 .The strategy of strain rate determination conforms to the method used by Naka et al. [23], which used the strain rate close to the failure in their tests on determining FLCs.The fourth forming limit curve from [18], was developed using average values from 0° and 90° to the rolling direction  Using the major limit strain ε * 1 (ρ, ln(Z)) in Equation ( 6), Z-FLS can be transformed to a space of (ε 1 , ε 2 , ln(Z)) as When working with a finite element forming simulation model to predict forming limit, the risk of fracture can be determined by considering the difference of the major limit strain ε * 1 on the surface of the Z-FLS diagram and the calculated major strains ε 1 in various elements on the model.In equation form, the difference is given by ∆ε If ∆ε 1 is higher than or equal to zero, fracture will occur.Given a safe margin of 0.05, if ∆ε 1 is between 0 and −0.05, the element is considered to be at risk of fracture; otherwise, the element is considered safe.

Construction of Z-FLS for AZ31 Magnesium Alloy
Forming limit curves (FLCs) reported in the literature at four different temperatures and two different strain rates were used to construct the Z-FLS for an AZ31 magnesium alloy as shown in Figure 3. Three forming limit curves from [5,9], were generated using Nakajima dome test data on 1.2 mm thick rectangular specimens of AZ31 magnesium alloy sheet under isothermal temperature conditions of 100, 200, and 300 • C and a constant punch speed of 0.1 mm/s [5,9].The obtained FLCs were fitted from measured data from three orientations to the rolling direction, i.e., 0 • , 45 • , and 90 • .In constant speed mechanical stretching tests [26], such as the Nakajima dome test, the strain rate increases with punch travel and its value depends on the punch speed.Albakri et al. [26] have shown that, under balanced biaxial test conditions in a Nakajima dome test, the strain rate in AZ31B sheets reached 0.008 s −1 at a punch speed of 0.1 mm/s and punch displacement of 32 mm.Since the strain rate value was not given in [5,9], it is assumed that the strain rate in their test is approximately 0.008 s −1 .The strategy of strain rate determination conforms to the method used by Naka et al. [23], which used the strain rate close to the failure in their tests on determining FLCs.The fourth forming limit curve from [18], was developed using average values from 0 • and 90 • to the rolling direction based on pneumatic stretching tests on 0.98 mm thick rectangular specimens of AZ31B sheet at 400 • C and a constant strain rate of 0.005 s −1 [18].It is to be noted that the sheet thickness in the second set of FLC is different than that in the first set.Antoniswamy et al. [27] have shown that initial sheet thickness did not affect the limit strain based on their tests on AZ31B sheet material at two different thicknesses (1.28 mm and 2.0 mm) at a temperature of 350 • C and different strain rates.
based on pneumatic stretching tests on 0.98 mm thick rectangular specimens of AZ31B sheet at 400 °C and a constant strain rate of 0.005 s −1 [18].It is to be noted that the sheet thickness in the second set of FLC is different than that in the first set.Antoniswamy et al. [27] have shown that initial sheet thickness did not affect the limit strain based on their tests on AZ31B sheet material at two different thicknesses (1.28 mm and 2.0 mm) at a temperature of 350 °C and different strain rates.[5,9,18] and used for constructing Z-FLS shown in Figure 5. (a) FLCs at 100, 200, and 300 °C and strain rate ~0.008/s [5,9]; (b) FLCs at 400 °C and strain rate 0.005/s [18].Table 1 lists the major limit strain values corresponding to strain ratios ranging from −0.25 to 0.6 based on above discussed FLCs.The major limit strain values in white boxes were measured from FLCs in Figure 3.The values of A, B, and C as a function of minor strain were determined using regression analysis of the data given in Table 1 and are plotted in Figure 4. Using Equation ( 6), major limit strains corresponding to several other minor strains and _____ ) ln( Z were then calculated.These values are in the blue boxes in Table 1.

Table 1.
Major limit strains at different strain ratios, temperatures and strain rates.   (The experimental data in white boxes were obtained from published FLCs in [5,9,18]. As shown in Table 1, among the four FLCs selected, the lowest ln(Z) is at 400 °C and strain rate of 0.005 s −1 .It was calculated using R as 8.31 J mol −1 K −1 and an average activation energy Q of 135 kJ mol −1 [12].Its value is 18.83, which was used as ln(Zo) for normalizing ln(Z).The upper limit of _____ ) ln( Z on the Z-FLS diagram corresponds to a combination of 50 °C temperature and 0.01 s −1 strain rate.Its value is 2.425, which corresponds to a ln(Z) value of 45.67.The selection of this upper limit is based on an assumption that the strain rate effect is trivial at room temperature.Using Equation (7), Table 1 lists the major limit strain values corresponding to strain ratios ranging from −0.25 to 0.6 based on above discussed FLCs.The major limit strain values in white boxes were measured from FLCs in Figure 3.The values of A, B, and C as a function of minor strain were determined using regression analysis of the data given in Table 1 and are plotted in Figure 4. Using Equation ( 6), major limit strains corresponding to several other minor strains and ln(Z) were then calculated.These values are in the blue boxes in Table 1.

Table 1.
Major limit strains at different strain ratios, temperatures and strain rates. (1)  .
As shown in Table 1, among the four FLCs selected, the lowest ln(Z) is at 400 • C and strain rate of 0.005 s −1 .It was calculated using R as 8.31 J mol −1 K −1 and an average activation energy Q of 135 kJ mol −1 [12].Its value is 18.83, which was used as ln(Z o ) for normalizing ln(Z).The upper limit of ln(Z) on the Z-FLS diagram corresponds to a combination of 50 • C temperature and 0.01 s −1 strain rate.Its value is 2.425, which corresponds to a ln(Z) value of 45.67.The selection of this upper limit is based on an assumption that the strain rate effect is trivial at room temperature.Using Equation (7), data in Table 1 was transformed into the space of (ε 1 , ε 2 , ln(Z)) and a Z-FLS was constructed as shown in Figure 5.In this figure, the major strain ε f 1 from Table 1 is plotted along the vertical axis and the corresponding minor strain ε 2 , calculated from the strain ratio ρ, is plotted along one of the horizontal axes.The third horizontal axis represents ln(Z), which has its lowest value at 1.It represents the combined effect of temperature and strain rate.
shown in Figure 5.In this figure, the major strain 1 ε from Table 1 is plotted along the vertical axis and the corresponding minor strain 2 ε , calculated from the strain ratio ρ, is plotted along one of the horizontal axes.The third horizontal axis represents   The use of Equations ( 6) and ( 7) to construct of Z-FLS was validated by determining FLCs from the constructed Z-FLS and comparing them with the experimental FLCs reported by Bruni et al.  data in Table 1 was transformed into the space of Z ε ε and a Z-FLS was constructed as shown in Figure 5.In this figure, the major strain f 1 ε from Table 1 is plotted along the vertical axis and the corresponding minor strain 2 ε , calculated from the strain ratio ρ, is plotted along one of the horizontal axes.The third horizontal axis represents   The use of Equations ( 6) and ( 7) to construct of Z-FLS was validated by determining FLCs from the constructed Z-FLS and comparing them with the experimental FLCs reported by Bruni et al.

Fracture Prediction Using Z-FLS
The Z-FLS constructed in Figure 5 is now used to predict fractures in an AZ31B sheet material during a non-isothermal round cup drawing process conducted by Dröder [28] and investigated numerically by Palaniswamy et al. [29].In the experiments, Dröder determined the limiting draw ratio (LDR) of AZ31B at blank holder/die temperatures ranging from 150 to 300 °C and an initial punch temperature of 25 °C.The punch diameter and initial sheet thickness were 100 mm and 1.3 mm, respectively.At 5 mm/s punch speed and blank holder/die temperature of 200 °C, the maximum LDR was 2.52.In the test, the sheet blank was initially heated to a uniform temperature that is equivalent to blank holder/die temperature.During drawing, the cup bottom lost heat to cold punch while the flange area stayed at a high level due to the heating effect of blank holder/die.In order to predict failure, a thermo-mechanically coupled finite element simulation of the cup drawing process was conducted using the implicit and explicit codes of LS-DYNA.The plastic strains, thickness distribution in the sheet, and punch load were calculated using the explicit dynamic algorithm, while the thermal phenomenon was calculated using the implicit algorithm.

Flow Stress Modeling
The flow stress at different temperatures and strain rates are expressed by a set of baseline flow stress curves while the strain rate effect is calculated by the Cowper-Symonds equation [30] where, 0 σ is the baseline flow stress and P eff ε is the effective plastic strain rate.C and P are temperature-dependent material parameters and are determined by minimizing the difference between the calculated flow stress and the flow stress at the target high strain rate.
In order to determine C and P in Equation ( 9), flow stress curves are needed at several temperatures and strain rates that are representative of the warm cup forming process being analyzed.Before localized necking, the strain rate in the cup wall during the cup drawing process with a punch speed of 5 mm/s is estimated to be in a range from 0.01 s −1 to 0.1 s −1 .Thus, flow stress curves at temperatures ranging from 25 °C to 300 °C at strain rate of 0.01 s −1 were chosen for baseline while those at strain rate of 0.1 s −1 were chosen as target for determining C and P in Equation (9).The baseline flow stress curves and target flow stress curves were generated using the flow stress

Fracture Prediction Using Z-FLS
The Z-FLS constructed in Figure 5 is now used to predict fractures in an AZ31B sheet material during a non-isothermal round cup drawing process conducted by Dröder [28] and investigated numerically by Palaniswamy et al. [29].In the experiments, Dröder determined the limiting draw ratio (LDR) of AZ31B at blank holder/die temperatures ranging from 150 to 300 • C and an initial punch temperature of 25 • C. The punch diameter and initial sheet thickness were 100 mm and 1.3 mm, respectively.At 5 mm/s punch speed and blank holder/die temperature of 200 • C, the maximum LDR was 2.52.In the test, the sheet blank was initially heated to a uniform temperature that is equivalent to blank holder/die temperature.During drawing, the cup bottom lost heat to cold punch while the flange area stayed at a high level due to the heating effect of blank holder/die.In order to predict failure, a thermo-mechanically coupled finite element simulation of the cup drawing process was conducted using the implicit and explicit codes of LS-DYNA.The plastic strains, thickness distribution in the sheet, and punch load were calculated using the explicit dynamic algorithm, while the thermal phenomenon was calculated using the implicit algorithm.

Flow Stress Modeling
The flow stress at different temperatures and strain rates are expressed by a set of baseline flow stress curves while the strain rate effect is calculated by the Cowper-Symonds equation [30] where, σ 0 is the baseline flow stress and .
ε P e f f is the effective plastic strain rate.C and P are temperature-dependent material parameters and are determined by minimizing the difference between the calculated flow stress and the flow stress at the target high strain rate.
In order to determine C and P in Equation (9), flow stress curves are needed at several temperatures and strain rates that are representative of the warm cup forming process being analyzed.Before localized necking, the strain rate in the cup wall during the cup drawing process with a punch speed of 5 mm/s is estimated to be in a range from 0.01 s −1 to 0.1 s −1 .Thus, flow stress curves at temperatures ranging from 25 • C to 300 • C at strain rate of 0.01 s −1 were chosen for baseline while those at strain rate of 0.1 s −1 were chosen as target for determining C and P in Equation (9).The baseline flow stress curves and target flow stress curves were generated using the flow stress model [31] given by Equation (10) and the tensile test data from [28,29].Based on this model, the strain rate and temperature effect on flow stress is expressed in terms of log (Z/2) as From the tensile stress-strain diagrams shown in Figure 7a,b, K and n in Equation ( 10) were determined by regression as model [31] given by Equation ( 10) and the tensile test data from [28,29].Based on this model, the strain rate and temperature effect on flow stress is expressed in terms of log (Z/2) as From the tensile stress-strain diagrams shown in Figure 7a (c) calculated at strain rate of 0.01 s −1 using Equation ( 10), (d) calculated at strain rate of 0.1 s −1 using Equation (10).
Then C and P at different temperatures were generated by minimizing the difference between the target flow stress and the baseline flow stress obtained from Equation (10).Table 2 gives the C and p values obtained by using the process described above.

FEM Modeling of Round Cup Draw
Figure 8 shows the FEM model for the round cup draw.Table 3 lists the processing conditions and Table 4 lists the material properties used in the simulation.The blank holder pressures vary with different setups to be at a minimum level to prevent wrinkling.It is assumed that no phase transformation happens during plastic deformation [4].Anisotropic behavior of AZ31 sheet material was studied by several studies [3,17,18,30].Based on polycrystal simulation study, the anisotropy can be explained by the increased activity of non-basal <c + a> and dislocations [3].In their study, (b) at temperature 200 • C and strain rates from 0.002 to 2.0 s −1 of 1.3 mm thick AZ31B sheet material [28,29]; (c) calculated at strain rate of 0.01 s −1 using Equation ( 10); (d) calculated at strain rate of 0.1 s −1 using Equation (10).
Then C and P at different temperatures were generated by minimizing the difference between the target flow stress and the baseline flow stress obtained from Equation (10).Table 2 gives the C and p values obtained by using the process described above.

FEM Modeling of Round Cup Draw
Figure 8 shows the FEM model for the round cup draw.Table 3 lists the processing conditions and Table 4 lists the material properties used in the simulation.The blank holder pressures vary with different setups to be at a minimum level to prevent wrinkling.It is assumed that no phase transformation happens during plastic deformation [4].Anisotropic behavior of AZ31 sheet material was studied by several studies [3,17,18,30].Based on polycrystal simulation study, the anisotropy can be explained by the increased activity of non-basal <c + a> and dislocations [3].In their study, Anew and Duygulu [3] observed that normal anisotropy (r) values at the transverse direction (TD) and 45 • direction drop and become closer to those of the rolling direction (RD), at which its r-values vary around 1 at temperatures ranging from 25 • C to 200 • C. Such a trend has also been observed by Tari and Worswick [30] in their study at the temperatures ranging from 200 • C to 250 • C. Based on the observed tendency and to reduce the complexity of material modeling, material deformation is assumed as isotropy in this study.
The tooling components were modeled as rigid, while the sheet blank was modeled as an elasticviscoplastic material with thermal effects.Material model 106 (MAT 106) in LS-DYNA, which uses von Mises criteria for isotropic yielding, was used to describe the elastic-viscoplastic behavior of the sheet material.The punch was modeled by using tetrahedron solid elements with 10 nodes and an average mesh size of 3 mm and had 8 elements distributed over the punch radius to calculate the heat transfer between the punch and the blank.Since the blank holder and the die temperatures are higher than the punch and blank temperatures, they function as heat source and were modeled by shell elements with an average mesh size of 2.5 mm, which ensured nine elements distributed over the die radius.The sheet blank was meshed with mixed quad and triangular shell elements with an average size of 2 mm.The quad elements were uniformly distributed from the edge of the blank to the start of the bottom flat surface of the punch to avoid high stiffness introduced by the triangular elements.The die displacement was constrained while the punch and the blank holder were allowed to move only in a vertical direction.Anew and Duygulu [3] observed that normal anisotropy (r) values at the transverse direction (TD) and 45° direction drop and become closer to those of the rolling direction (RD), at which its r-values vary around 1 at temperatures ranging from 25 °C to 200 °C.Such a trend has also been observed by Tari and Worswick [30] in their study at the temperatures ranging from 200 °C to 250 °C.Based on the observed tendency and to reduce the complexity of material modeling, material deformation is assumed as isotropy in this study.The tooling components were modeled as rigid, while the sheet blank was modeled as an elastic-viscoplastic material with thermal effects.Material model 106 (MAT 106) in LS-DYNA, which uses von Mises criteria for isotropic yielding, was used to describe the elastic-viscoplastic behavior of the sheet material.The punch was modeled by using tetrahedron solid elements with 10 nodes and an average mesh size of 3 mm and had 8 elements distributed over the punch radius to calculate the heat transfer between the punch and the blank.Since the blank holder and the die temperatures are higher than the punch and blank temperatures, they function as heat source and were modeled by shell elements with an average mesh size of 2.5 mm, which ensured nine elements distributed over the die radius.The sheet blank was meshed with mixed quad and triangular shell elements with an average size of 2 mm.The quad elements were uniformly distributed from the edge of the blank to the start of the bottom flat surface of the punch to avoid high stiffness introduced by the triangular elements.The die displacement was constrained while the punch and the blank holder were allowed to move only in a vertical direction.Thermal conductivity of the sheet material (W/m K) 95.8 Heat capacity of the sheet material (J/Kg K) 1767.5 Thermal conductivity of the tool material (W/m K) 60.5 Heat capacity the tool material (J/Kg K) 3410 Interface contact heat conductance (W/m 2 K) 4500 Factor to convert plastic deformation energy to heat 0.95 Table 3. Drawing conditions and material used in the simulation of round cup drawing [6,29].

Tooling Setup
Punch diameter (mm) 100 Punch and die corner radii (mm  4. Thermal properties of AZ31 sheet material and tool material [29]. Thermal conductivity of the sheet material (W/m K) 95.8 Heat capacity of the sheet material (J/Kg K) 1767.5 Thermal conductivity of the tool material (W/m K) 60.5 Heat capacity the tool material (J/Kg K) 3410 Interface contact heat conductance (W/m 2 K) 4500 Factor to convert plastic deformation energy to heat 0.95 The Coulomb friction coefficient at the interfaces between the tooling components and the sheet material was assumed as 0.1 at 25-150 • C, 0.2 at 200 and 250 • C, and 0.4 at 300 • C. The selection of relatively high values of friction coefficients at higher temperatures is based on friction study on magnesium alloys at different lubrication conditions at elevated temperatures [32][33][34][35].Local increase of friction coefficient at high temperatures is due to the possible failure of lubricant and friction coefficient selection strategy discussed in [36].

Results and Analysis
In this study, cup drawing at four different temperatures of sheet blank/die/blank holder namely 150, 200, 250, and 300 • C-were simulated.The initial punch temperature was assumed to be 25 • C. The finite element models were first validated by comparing the calculated temperature distribution, punch load, and thickness distribution on the formed cups with the experimental measurements.Then, the Z-FLS in Figure 4 was used to predict fracture and determine the LDRs of the numerically formed cups.

Temperature Distribution
During simulations, the blank holder and die temperatures were maintained at the initial temperature, which caused the flange area to remain at a higher temperature than the rest of the cup.As the sheet blank was drawn into a cup, the cup wall was cooled due to heat loss to the cold punch.The temperature distribution on the cup wall is given in Table 5.The difference between the predicted values in this study and the reported values in [29] is small (<3 • C).The Coulomb friction coefficient at the interfaces between the tooling components and the sheet material was assumed as 0.1 at 25-150 °C, 0.2 at 200 and 250 °C, and 0.4 at 300 °C.The selection of relatively high values of friction coefficients at higher temperatures is based on friction study on magnesium alloys at different lubrication conditions at elevated temperatures [32][33][34][35].Local increase of friction coefficient at high temperatures is due to the possible failure of lubricant and friction coefficient selection strategy discussed in [36].

Results and Analysis
In this study, cup drawing at four different temperatures of sheet blank/die/blank holder namely 150, 200, 250, and 300 °C-were simulated.The initial punch temperature was assumed to be 25 °C.The finite element models were first validated by comparing the calculated temperature distribution, punch load, and thickness distribution on the formed cups with the experimental measurements.Then, the Z-FLS in Figure 4 was used to predict fracture and determine the LDRs of the numerically formed cups.

Temperature Distribution
During simulations, the blank holder and die temperatures were maintained at the initial temperature, which caused the flange area to remain at a higher temperature than the rest of the cup.As the sheet blank was drawn into a cup, the cup wall was cooled due to heat loss to the cold punch.The temperature distribution on the cup wall is given in Table 5.The difference between the predicted values in this study and the reported values in [29] is small (<3 °C).The Coulomb friction coefficient at the interfaces between the tooling components and the sheet material was assumed as 0.1 at 25-150 °C, 0.2 at 200 and 250 °C, and 0.4 at 300 °C.The selection of relatively high values of friction coefficients at higher temperatures is based on friction study on magnesium alloys at different lubrication conditions at elevated temperatures [32][33][34][35].Local increase of friction coefficient at high temperatures is due to the possible failure of lubricant and friction coefficient selection strategy discussed in [36].

Results and Analysis
In this study, cup drawing at four different temperatures of sheet blank/die/blank holder namely 150, 200, 250, and 300 °C-were simulated.The initial punch temperature was assumed to be 25 °C.The finite element models were first validated by comparing the calculated temperature distribution, punch load, and thickness distribution on the formed cups with the experimental measurements.Then, the Z-FLS in Figure 4 was used to predict fracture and determine the LDRs of the numerically formed cups.

Temperature Distribution
During simulations, the blank holder and die temperatures were maintained at the initial temperature, which caused the flange area to remain at a higher temperature than the rest of the cup.As the sheet blank was drawn into a cup, the cup wall was cooled due to heat loss to the cold punch.The temperature distribution on the cup wall is given in Table 5.The difference between the predicted values in this study and the reported values in [29] is small (<3 °C).9.The cups drawn at 150 and 200 • C do not exhibit any fracture, while the cups drawn at 300 • C fracture at 33 mm punch displacement.The punch loads reach the peak values at about 40 mm stroke at 150 and 200 • C, then decrease due to a reduction of restraining force from the work piece with reduced flange area clamped between the blank holder and the die.At 300 • C, the peak load is reached at a stroke of about 25 mm and stays at this level till the stroke reaches 33 mm, then drops due to the fracture of the cup wall.These trajectories are similar in trend to the punch load trajectories obtained in round cup draw [28,29].At 150 • C, the predicted punch load is very close to the experimental value; however, at 200 • C and 300 • C, the predicted punch loads are higher at larger punch travels and the peak punch load is higher compared with the experimental values.The possible reasons are: (1) the stress-strain curves were extrapolated after the peak stress in the tensile test data and did not consider the softening stage; and (2) the elements were not deleted when the fracture criterion was met in the case of 300 • C.

Punch Loads
Punch loads predicted by finite element simulations are compared with the measured values in Figure 9.The cups drawn at 150 and 200 °C do not exhibit any fracture, while the cups drawn at 300 °C fracture at 33 mm punch displacement.The punch loads reach the peak values at about 40 mm stroke at 150 and 200 °C, then decrease due to a reduction of restraining force from the work piece with reduced flange area clamped between the blank holder and the die.At 300 °C, the peak load is reached at a stroke of about 25 mm and stays at this level till the stroke reaches 33 mm, then drops due to the fracture of the cup wall.These trajectories are similar in trend to the punch load trajectories obtained in round cup draw [28,29].At 150 °C, the predicted punch load is very close to the experimental value; however, at 200 °C and 300 °C, the predicted punch loads are higher at larger punch travels and the peak punch load is higher compared with the experimental values.The possible reasons are: (1) the stress-strain curves were extrapolated after the peak stress in the tensile test data and did not consider the softening stage; and (2) the elements were not deleted when the fracture criterion was met in the case of 300 °C.

Thickness Distribution
Figure 10 shows a comparison of the predicted thinning of the cup wall with experimental measurements at 200 and 250 °C at a draw ratio (DR) of 2.3 [28,29].The patterns of thinning distribution match with the experimental data.The maximum thinning increases with increases in temperature.Compared with the measurement, the predicted maximum thinning is 3~5% lower while the maximum thickening at the flange is about 8~11% lower in both cases.One possible reason for the difference is our assumption of the constant coefficient of friction at each temperature investigated.The constructed Z-FLS in Figure 5 was used to predict fracture and its location on the numerically formed cup.The maximum blank diameter (Do) that can be drawn in the die without fracture was used to calculate the LDR, which is a ratio of Do and punch diameter Dp.The fracture of numerically formed cups was identified by plotting the data points (ε1, ε2, _____ ) ln(Z ) of concerned elements against the Z-FLS, see Figure 11.The data points of three elements representing maximum thinning 'ET', maximum major strain 'EM', and no risk of fracture 'ES' are plotted at the initiation of fracture, see Table 6.In all cases, the elements with the maximum major strain and maximum thinning fracture at the same depth of draw.At 250 °C, the thinnest element also develops the highest major strain.The fracture location shifts from location close to punch corner at 150 °C to the

Thickness Distribution
Figure 10 shows a comparison of the predicted thinning of the cup wall with experimental measurements at 200 and 250 • C at a draw ratio (DR) of 2.3 [28,29].The patterns of thinning distribution match with the experimental data.The maximum thinning increases with increases in temperature.Compared with the measurement, the predicted maximum thinning is 3~5% lower while the maximum thickening at the flange is about 8~11% lower in both cases.One possible reason for the difference is our assumption of the constant coefficient of friction at each temperature investigated.

Thickness Distribution
Figure 10 shows a comparison of the predicted thinning of the cup wall with experimental measurements at 200 and 250 °C at a draw ratio (DR) of 2.3 [28,29].The patterns of thinning distribution match with the experimental data.The maximum thinning increases with increases in temperature.Compared with the measurement, the predicted maximum thinning is 3~5% lower while the maximum thickening at the flange is about 8~11% lower in both cases.One possible reason for the difference is our assumption of the constant coefficient of friction at each temperature investigated.The constructed Z-FLS in Figure 5 was used to predict fracture and its location on the numerically formed cup.The maximum blank diameter (Do) that can be drawn in the die without fracture was used to calculate the LDR, which is a ratio of Do and punch diameter Dp.The fracture of numerically formed cups was identified by plotting the data points (ε1, ε2, _____ ) ln(Z ) of concerned elements against the Z-FLS, see Figure 11.The data points of three elements representing maximum thinning 'ET', maximum major strain 'EM', and no risk of fracture 'ES' are plotted at the initiation of fracture, see Table 6.In all cases, the elements with the maximum major strain and maximum thinning fracture at the same depth of draw.At 250 °C, the thinnest element also develops the highest major strain.The fracture location shifts from location close to punch corner at 150 °C to the upper portion of cup wall at 200, 250, and 300 °C, which conforms to the observation in [29].This can be explained by increased strength difference between sheet material at the die curvature area and The constructed Z-FLS in Figure 5 was used to predict fracture and its location on the numerically formed cup.The maximum blank diameter (D o ) that can be drawn in the die without fracture was used to calculate the LDR, which is a ratio of D o and punch diameter D p .The fracture of numerically formed cups was identified by plotting the data points (ε 1 , ε 2 , ln(Z)) of concerned elements against the Z-FLS, see Figure 11.The data points of three elements representing maximum thinning 'ET', maximum major strain 'EM', and no risk of fracture 'ES' are plotted at the initiation of fracture, see Table 6.In all cases, the elements with the maximum major strain and maximum thinning fracture at the same depth of draw.At 250 • C, the thinnest element also develops the highest major strain.The fracture location shifts from location close to punch corner at 150 • C to the upper portion of cup wall at 200, 250, and 300 • C, which conforms to the observation in [29].This can be explained by increased strength difference between sheet material at the die curvature area and that at the punch radius area, which is strengthened by the cold punch.The identified depths of cup at the initiation of fracture are close to the observations in the experiment [29].For example, at 300 • C and DR = 2.3, the calculated depth is 33 mm while 32 mm is the measured depth, and at 250 • C and DR = 2.4, the calculated depth is 28.5 mm while 30 mm is the depth reported in the experiment.that at the punch radius area, which is strengthened by the cold punch.The identified depths of cup at the initiation of fracture are close to the observations in the experiment [29].For example, at 300 °C and DR = 2.3, the calculated depth is 33 mm while 32 mm is the measured depth, and at 250 °C and DR = 2.4, the calculated depth is 28.5 mm while 30 mm is the depth reported in the experiment.The evolution of the risk of fracture can be observed by plotting historical traces of the critical elements under the Z-FLS diagram.Figure 12a gives a plot of three critical elements in the case of temperature 200 • C. Element 'ET' and element 'EM' are at about same radial distance from the center of blank and thus are in contact with die and blank holder during most of the punch stroke.Because of this, their temperature remains approximately at 200 • C as shown in Figure 12b.As a result, their major strain traces are close to each other, which reaches the limit strain on the Z-FLS at strain rates of 0.06 s −1 and 0.07 s −1 , respectively, when the punch travel becomes equal to 27 mm as shown in Figure 12a.Element 'ES', on the other hand, makes an early contact with the cold punch, which brings its temperature down at a faster rate (Figure 12b).The decrease in temperature leads to a faster increase in ln(Z) which increases to about 2.03 at a punch travel of 27 mm, see Figure 12d.However, a decrease in temperature also strengthens the local material and major strain has a slower growth ( .ε 1 ≤ 0.02 s −1 ) thereafter, as shown in Figure 12c.The slowly increasing major strain helps to keep the element 'ES' safe under the Z-FLS.
Figure 13 shows that the predicted LDRs at different temperatures match well with those obtained in the experiment.One exception is at temperature of 200 • C, the prediction is higher by 0.03 (which is less than 0.15%) than the experimental value.The discrepancy is quite small.The evolution of the risk of fracture can be observed by plotting historical traces of the critical elements under the Z-FLS diagram.Figure 12a gives a plot of three critical elements in the case of temperature 200 °C.Element 'ET' and element 'EM' are at about same radial distance from the center of blank and thus are in contact with die and blank holder during most of the punch stroke.Because of this, their temperature remains approximately at 200 °C as shown in Figure 12b.As a result, their major strain traces are close to each other, which reaches the limit strain on the Z-FLS at strain rates of 0.06 s −1 and 0.07 s −1 , respectively, when the punch travel becomes equal to 27 mm as shown in Figure 12a.Element 'ES', on the other hand, makes an early contact with the cold punch, which brings its temperature down at a faster rate (Figure 12b).The decrease in temperature leads to a faster increase in _____ ) ln( Z which increases to about 2.03 at a punch travel of 27 mm, see Figure 12d.Figure 13 shows that the predicted LDRs at different temperatures match well with those obtained in the experiment.One exception is at temperature of 200 °C, the prediction is higher by 0.03 (which is less than 0.15%) than the experimental value.The discrepancy is quite small.(a)

Conclusions and Discussion
In this study, the concept of forming limit surface (Z-FLS), which utilizes ln(Z) to represent the combined effect of strain rate and temperature on forming limit strains, was revisited to develop a practical approach constructing Z-FLS.A Z-FLS was then constructed for AZ31 magnesium alloy sheet material.A non-isothermal cup drawing process of AZ31 sheet material was modeled using thermo-mechanically coupled finite element forming simulation.The constructed Z-FLS was used to identify failure by fracture on the numerically formed cups.Results show that the depths at fracture and limiting draw ratios match well with the experimental data.
Since both strain rate and temperature effects are represented by a single parameter Z, the Z-FLS may provide a concise and convenient way to take them into account in predicting formability in warm and hot forming.Its usefulness becomes evident in elevated temperature forming, which involves non-isothermal transient conditions, in which both temperature and strain rate change continuously with punch travel.As discussed, the use of available experimental data from various sources and several approximations made in this study may have introduced some uncertainties in the construction of Z-FLS in this paper, but the validity of Z-FLS has been demonstrated by using it to predict FLCs at different temperatures and LDRs in a non-isothermal cup draw and correlating the results closely with experimental data.In practice, Z-FLS will be constructed using the forming limit data at a few carefully selected temperatures and a wide range of strain rates.The Z-FLS can then be used for representing continuously distributed temperature/ strain rate field in warm forming condition.
A number of assumptions were made in this study to reduce the analytical complexity of the thermo-mechanical phenomenon that occurs during the warm forming process.Although good correlations were obtained with the experimental data, some of these assumptions, such as microstructural differences in the AZ31 alloys used in experimental data used and variation of thermo-physical properties of the blank and tool materials with temperature, may have introduced uncertainties in the results.Planned experiments need to be done to further validate the Z-FLS diagrams so that they can be utilized in applications with less uncertainties.

Conclusions and Discussion
In this study, the concept of forming limit surface (Z-FLS), which utilizes ln(Z) to represent the combined effect of strain rate and temperature on forming limit strains, was revisited to develop a practical approach constructing Z-FLS.A Z-FLS was then constructed for AZ31 magnesium alloy sheet material.A non-isothermal cup drawing process of AZ31 sheet material was modeled using thermo-mechanically coupled finite element forming simulation.The constructed Z-FLS was used to identify failure by fracture on the numerically formed cups.Results show that the depths at fracture and limiting draw ratios match well with the experimental data.
Since both strain rate and temperature effects are represented by a single parameter Z, the Z-FLS may provide a concise and convenient way to take them into account in predicting formability in warm and hot forming.Its usefulness becomes evident in elevated temperature forming, which involves non-isothermal transient conditions, in which both temperature and strain rate change continuously with punch travel.As discussed, the use of available experimental data from various sources and several approximations made in this study may have introduced some uncertainties in the construction of Z-FLS in this paper, but the validity of Z-FLS has been demonstrated by using it to predict FLCs at different temperatures and LDRs in a non-isothermal cup draw and correlating the results closely with experimental data.In practice, Z-FLS will be constructed using the forming limit data at a few carefully selected temperatures and a wide range of strain rates.The Z-FLS can then be used for representing continuously distributed temperature/strain rate field in warm forming condition.
A number of assumptions were made in this study to reduce the analytical complexity of the thermo-mechanical phenomenon that occurs during the warm forming process.Although good correlations were obtained with the experimental data, some of these assumptions, such as microstructural differences in the AZ31 alloys used in experimental data used and variation of thermo-physical properties of the blank and tool materials with temperature, may have introduced uncertainties in the results.Planned experiments need to be done to further validate the Z-FLS diagrams so that they can be utilized in applications with less uncertainties.

Figure 2 .
Figure 2. Plots of major limit strains at strain paths of UT (Square solid dots), PS (Circular solid dots) and EBT (Triangular solid dots) as functions of ln(Z) and _____ ) ln(Z against data points from measurement at strain rates (10 −4 /s to 10 −2 /s) and temperatures (25 °C to 300 °C) in [23].(a) Function of ) ln(Z ; (b) Function of _____ ) ln(Z .

Figure 2 .
Figure 2. Plots of major limit strains at strain paths of UT (Square solid dots), PS (Circular solid dots) and EBT (Triangular solid dots) as functions of ln(Z) and ln(Z) against data points from measurement at strain rates (10 −4 /s to 10 −2 /s) and temperatures (25 • C to 300 • C) in [23].(a) Function of ln(Z); (b) Function of ln(Z).
Z , which has its lowest value at 1.It represents the combined effect of temperature and strain rate.

Figure 4 .
Figure 4. Plots of material constants (A-C) at different strain ratios.
[17]   and Hsu et al.[19] by using Nakajima dome tests.The experiments were conducted at 300 °C and two different strain rates, namely 0.001 s −1 and 0.01 s −1 .The corresponding ) ln(Z are 1.14 and 1.26 as shown in Figure6.The forming limit curves determined by using Z-FLS and shown by the solid lines match well with the distribution of experimental data[17,19].

Figure 4 .
Figure 4. Plots of material constants (A-C) at different strain ratios.
Z , which has its lowest value at 1.It represents the combined effect of temperature and strain rate.

Figure 4 .
Figure 4. Plots of material constants (A-C) at different strain ratios.
[17]   and Hsu et al.[19] by using Nakajima dome tests.The experiments were conducted at 300 °C and two different strain rates, namely 0.001 s −1 and 0.01 s −1 .The corresponding ) ln(Z are 1.14 and 1.26 as shown in Figure6.The forming limit curves determined by using Z-FLS and shown by the solid lines match well with the distribution of experimental data[17,19].

Figure 6 .
Figure 6.FLCs constructed using the Z-FLS in Figure 5 and compared with data reported in [17,19].(a) FLC at ) ln(Z = 1.14 (00 °C and strain rate of 0.001/s): the solid line is from Z-FLS, while green stars are from [19] and circular solid dots are from [17]; (b) FLC at

Figure 6 .
Figure 6.FLCs constructed using the Z-FLS in Figure 5 and compared with data reported in [17,19].(a) FLC at ln(Z) = 1.14 (300 • C and strain rate of 0.001/s): the solid line is from Z-FLS, while green stars are from [19] and circular solid dots are from [17]; (b) FLC at ln(Z) = 1.26 (300 • C and strain rate of 0.01/s): the solid line is from Z-FLS, while green stars are from [19] and the solid dots are from [17].

Figure 8 .
Figure 8. Finite element model of round cup draw.(a) Tooling design and setup; (b) meshed blank.

Table 5 .
Predicted temperature distribution on the cup being drawn at 200 °C.

Table 5 .
Predicted temperature distribution on the cup being drawn at 200 °C.

Figure 9 .
Figure 9. Predicted punch load curves at different temperatures for draw ratio of 2.3 compared with experimental data in [28,29].

Figure 9 .
Figure 9. Predicted punch load curves at different temperatures for draw ratio of 2.3 compared with experimental data in [28,29].

Figure 9 .
Figure 9. Predicted punch load curves at different temperatures for draw ratio of 2.3 compared with experimental data in [28,29].
However, a decrease in temperature also strengthens the local material and major strain has a slower growth ( ) thereafter, as shown in Figure12c.The slowly increasing major strain helps to keep the element 'ES' safe under the Z-FLS.

1 difference 1 major
C and P material parameters at different temperatures in the Cowper and Symonds equation Z Zener-Hollomon parameter Q activation energy R universal gas constant T absolute temperature ε 1 , ε 2major and minor strains on the deformed element in the FEM simulation ∆ε * of the major limit strain on the FLS and the major strain of the concerned element during the forming operation ε f t effective strain at fracture in tensile test ε * limit strain on the FLC

Table 2 .
Material parameters in the Cowper and Symonds equation

Table 2 .
Material parameters in the Cowper and Symonds equation.

Table 5 .
Predicted temperature distribution on the cup being drawn at 200 • C.

) [29]
Punch loads predicted by finite element simulations are compared with the measured values in Figure • C

Table 6 .
Failure status of representative elements in non-isothermal cup draw test.

Table 6 .
Failure status of representative elements in non-isothermal cup draw test