Process Optimization of the Hot Stamping of AZ31 Magnesium Alloy Sheets Based on Response Surface Methodology

Hot stamping is an important manufacturing process for sheet metal parts. However, it is easy to produce defects such as thinning and cracking in the drawing area during the stamping process. In this paper, the finite element solver ABAQUS/Explicit was used to establish the numerical model of the magnesium alloy hot-stamping process. The stamping speed (2~10 mm/s), the blank-holder force (3~7 kN), and the friction coefficient (0.12~0.18) were selected as the influencing factors. Taking the maximum thinning rate obtained through simulation as the optimization objective, the response surface methodology (RSM) was applied to optimize the influencing factors in sheet hot stamping at a forming temperature of 200 °C. The results showed that the maximum thinning rate of sheet metal was most influenced by the blank-holder force, and the interaction between the stamping speed and the blank-holder force/friction coefficient had a great influence on the maximum thinning rate. The optimal value of the maximum thinning rate of the hot-stamped sheet was 7.37%. Through the experimental verification for the hot-stamping process scheme, the maximum relative error between the simulation and the experimental results was 8.72%. This proves the accuracy of the established finite element model and the response surface model. This research provides a feasible optimization scheme for the analysis of the hot-stamping process of magnesium alloys.


Introduction
With the increasingly severe energy crisis and emission limits, lightweight manufacturing has become the inevitable trend in the future development of automobiles. Magnesium alloys show good application prospects in the automobile field for their excellent properties such as high strength, low density, excellent heat dissipation, and electromagnetic shielding [1,2]. However, as a result of magnesium alloy's hexagonal close-packed (hcp) crystal structure, its plastic deformation ability at room temperature is poor, and tensile cracking and wrinkling defects easily occur [3]. Previous research [4] has shown that when a magnesium alloy sheet is heated to above 200 • C, the first-order conical surface of the crystal structure and {1011} and {1021} slip systems are activated, and the plastic deformation capacity greatly improves. Thus, the hot-stamping technology of magnesium alloy sheets has attracted much attention [5].
Improvements in the hot-stamping process are often based on repeated testing, which is an extremely expensive and time-consuming process [6]. It is therefore necessary to reduce manufacturing-process costs by minimizing design time and physical tests. In recent decades, the combination of finite element simulation [7], experimental techniques [8,9], and optimization methods [10,11] has provided a promising alternative to the optimal hot-stamping process design. For instance, Xiao et al. [12] proposed a multi-objective stochastic method to determine the optimal parameters of the hot-stamping process and Materials 2023, 16, 1867 2 of 13 further studied the effect of these parameters on processing quality in a limited range. Based on the RSM and non-dominated sorting genetic algorithm (NSGA-I), Bao et al. [13] performed a multi-objective optimization procedure for the partition temperature of a hot-stamped steel sheet. Lei et al. [14] developed a novel finite element model coupled with thermoelastic-plastic behavior for the hot stamping of blank patchwork to predict the stamping results and obtain the optimum process parameters. A constitutive equation developed by Namklang et al. [15] was combined with finite element results to assess the local deformability of hot-stamping products. Gao et al. [16] proposed a novel multiobjective optimization method for improving the forming quality and energy consumption in the hot-stamping process. Cui et al. [17] used a three-dimensional forming limit diagram (FLD), which considered phase transformations, to evaluate the forming property of the wrinkling of high-strength steel during hot stamping. A finite element model coupled with a thermal-mechanical phase was established by Quan et al. [18] to examine the effect of hot-stamping process variables on the evolution of the phase field. Hu et al. [19] derived a thermal-mechanical constitutive equation combined with material damage to describe the hot-stamping behavior of high-strength steel and further studied the effects of the blank-holder force and contact relationship on the punch force, crack initiation, and formability.
In the above research, numerous optimization studies mainly focus on the hotstamping process for aluminum alloy/steel sheets, and there is not enough research on a hot-stamping process for magnesium alloy sheets. In addition, there are many factors affecting formability, mainly including the blank-holder force, strain rate, lubrication effects, and die clearance during the hot-stamping process [20]. Due to the mutual restriction of the process variables, it is crucial to obtain the combination of technical parameters that can guarantee the forming quality and performance of the forming parts. Therefore, the aim of this work is to better understand the effects of process parameters on magnesium alloy sheet thinning and optimize the forming process to minimize the sheet thinning of hot-stamped sheets.
The present paper is organized as follows: In Section 2, the constitutive model suitable for magnesium alloy sheet forming is introduced, and a finite element model of the hotstamping process of a magnesium alloy sheet is established and verified. In Section 3, taking the maximum thinning rate of the hot-stamped sheet as the optimization objective, the RSM is applied to investigate the impact of three process parameters on the maximum thinning rate of the sheet. In Section 4, the accuracy of the numerical model and the response surface model is verified using hot-stamping experiments. Finally, the conclusions are summarized in Section 5.

Constitutive Model
Barlat et al. [21] proposed a yield criterion applicable to the plane stress state, which is called the YId2000 yield criterion, to better describe the plastic anisotropy behavior of magnesium, aluminum, and other metal sheets. This criterion represents the anisotropy properties of materials through two linear transformations of the Cauchy stress tensor. The in-plane anisotropy behavior of metal materials can be better described compared with the well-known Hill yield criterion, and its yield surface function is expressed as where where S and S are the linear transformations of stress tensors, and S 1 , S 2 , S 1 , and S 2 are the principal values of the stress tensor.
where α k parameters are the independent coefficients of material (for k from 1 to 8).
The ductile Lemaitre damage model [22], which is based on a thermodynamic framework, is adopted to accurately describe the damaging behavior of the material. Damage can be accurately quantified using the internal variable D, which lies in the range of 0 ≤ D < 1 when represented in its scalar form. This variable (D) represents the ratio of the damaged area of a unit surface (S D ) to the total surface (S): D = S D /S. Based on the hypothesis of strain equivalence, the equivalent stress tensor with damage can be defined as where σ is the stress tensor of the material without considering the damage behavior. The damage strain energy release rate ψ is associated with the damage variable D, and it can be expressed as follows: where E is the Young modulus, ν is the Poisson ratio, q is the equivalent stress, σ H is the hydrostatic stress, and σ H /q represents the stress triaxiality of the material. Accordingly, damage evolution can be expressed as where . λ is the plastic multiplier, F is the plasticity dissipation potential function, and . ε p represents the equivalent plastic strain rate. As soon as equivalent plastic strain exceeds a strain threshold ε D , damage increases by the following formula: The macroscopic fracture is then accounted for by a critical damage value D c . Once the value of D reaches D c , the damage variable D is assigned a value of 1, which signifies the failure of the material.
Since the YId2000 yield criterion can better describe the plastic anisotropy behavior of a light alloy sheet than the traditional yield criterion [23,24], the plastic deformation and damage behavior of the AZ31 magnesium alloy sheet was analyzed with the Lemaitre damage model modified using the YId2000 yield criterion. The constitutive model parameters applicable to describing the deformation of the AZ31 magnesium alloy sheet are discussed in the authors' previous research [25]. The anisotropic yield function of Yld2000 was embedded into the material properties with the help of user subroutine VUMAT, and the numerical calculation for the hot-stamping forming process of the magnesium alloy sheet was realized based on the ABAQUS platform.

Establishment of Finite Element Model
The experimental material was a commercial rolled AZ31 magnesium alloy sheet with a thickness of 0.6 mm, and its chemical composition is listed in Table 1. A rectangular plate with a size of 110 mm × 60 mm was cut off using a wire-cutting machine for the hot-stamping test. The mechanical properties of the AZ31 magnesium alloy sheet are given in the authors' previous study [26]. The hot-stamping process of the AZ31 magnesium alloy sheet was simulated based on the nonlinear finite element solver ABAQUS/Explicit. The three-dimensional coupled thermomechanical finite element model of hot stamping is shown in Figure 1. The whole simulation process of hot-stamping formation was divided into two steps: The holder moved down and made contact with the sheet (step 1), and the punch moved down a fixed distance for hot stamping (step 2). Because changes in temperature and displacement are involved in the forming process, the type of "dynamic, explicit, temperature-displacement" was selected in both analysis steps. The blank magnesium alloy sheet was set as a deformable body, while the die, holder, and punch were set as the rigid body. The Coulomb friction model was adopted for describing the friction between the sheet and dies, and the friction coefficient was 0.15. The sheet was meshed using a four-node thermally coupled doubly curved thin-shell element with reduced integration and hourglass control (S4RT), and five integral points were set in the thickness direction.  Table 2 lists the processing parameters of the AZ31 magnesium alloy sheet during the hot-stamping experiment. Finite element simulation parameters during hot stamping were set in accordance with the actual experimental parameters as follows: The temperature of the environment and the blank were set to 20 °C and 200 °C, respectively; the heat transfer coefficient between the blank and dies was 2000 W/(m 2 °C); the inelastic heat fraction was 0.9; and the film coefficient was 0.025 W/(m 2 °C). In addition, the clearance distance between the punch and die was 0.8 mm.

Validation of Finite Element Model
The hot-stamping test of the magnesium alloy sheet was carried out to verify the  Table 2 lists the processing parameters of the AZ31 magnesium alloy sheet during the hot-stamping experiment. Finite element simulation parameters during hot stamping were set in accordance with the actual experimental parameters as follows: The temperature of the environment and the blank were set to 20 • C and 200 • C, respectively; the heat transfer coefficient between the blank and dies was 2000 W/(m 2 • C); the inelastic heat fraction was 0.9; and the film coefficient was 0.025 W/(m 2 • C). In addition, the clearance distance between the punch and die was 0.8 mm.

Validation of Finite Element Model
The hot-stamping test of the magnesium alloy sheet was carried out to verify the validity and rationality of the established finite element model of hot stamping. The overall dimensions of the hot-stamping dies are shown in Figure 2. The heated magnesium alloy sheet was placed in the dies, and the hot-stamping experiment was performed after the sheet was located.
tance between the punch and die was 0.8 mm.

Validation of Finite Element Model
The hot-stamping test of the magnesium alloy sheet was carried out to verify the validity and rationality of the established finite element model of hot stamping. The overall dimensions of the hot-stamping dies are shown in Figure 2. The heated magnesium alloy sheet was placed in the dies, and the hot-stamping experiment was performed after the sheet was located. Comparing the experimental results with the simulation results, it can be observed that a crack appeared in the sheet metal near the corner region of the die, as shown in Figure 3. The fracture location obtained through simulation matched that of the test piece quite well, which indicates the accuracy of the numerical model. Comparing the experimental results with the simulation results, it can be observed that a crack appeared in the sheet metal near the corner region of the die, as shown in Figure 3. The fracture location obtained through simulation matched that of the test piece quite well, which indicates the accuracy of the numerical model.  The equivalent plastic strain distribution of the hot-stamped sheet obtained through the finite element method is shown in Figure 4. After hot-stamping formation, the fillet region (B) near the punch had a low degree of work hardening, and the material was subjected to greater tensile and compressive stresses in the thickness direction than the material in region A, which was a dangerous area prone to breakage. The straight wall region (C), as the main force transfer region, was affected by radial tensile stress and compressive strain in the thickness direction. In addition, the material in flange region (A) and bottom region (E) was difficult to be replenished in time during hot stamping due to the blankholder force and the friction of the punch. The sheet near the corner region (D) of the die was prone to thinning due to its large strain. Therefore, the maximum thinning rate of the sheet in region D was selected for subsequent analysis. The equivalent plastic strain distribution of the hot-stamped sheet obtained through the finite element method is shown in Figure 4. After hot-stamping formation, the fillet region (B) near the punch had a low degree of work hardening, and the material was subjected to greater tensile and compressive stresses in the thickness direction than the material in region A, which was a dangerous area prone to breakage. The straight wall region (C), as the main force transfer region, was affected by radial tensile stress and compressive strain in the thickness direction. In addition, the material in flange region (A) and bottom region (E) was difficult to be replenished in time during hot stamping due to the blank-holder force and the friction of the punch. The sheet near the corner region (D) of the die was prone to thinning due to its large strain. Therefore, the maximum thinning rate of the sheet in region D was selected for subsequent analysis.
(C), as the main force transfer region, was affected by radial tensile stress and compressive strain in the thickness direction. In addition, the material in flange region (A) and bottom region (E) was difficult to be replenished in time during hot stamping due to the blankholder force and the friction of the punch. The sheet near the corner region (D) of the die was prone to thinning due to its large strain. Therefore, the maximum thinning rate of the sheet in region D was selected for subsequent analysis.

Establishment of Response Surface Model
The response surface methodology (RSM) is an effective method to solve multi-variable and objective optimization problems. Test groups are designed using reasonable experimental design methods, and certain reliable data can be obtained through experiments or simulations. For fitting the functional relationship between various factors and response values, multiple quadratic regression equations are used. By analyzing regression equations, the best parameter combination is determined. This method is widely used

Establishment of Response Surface Model
The response surface methodology (RSM) is an effective method to solve multi-variable and objective optimization problems. Test groups are designed using reasonable experimental design methods, and certain reliable data can be obtained through experiments or simulations. For fitting the functional relationship between various factors and response values, multiple quadratic regression equations are used. By analyzing regression equations, the best parameter combination is determined. This method is widely used in the metal-forming field because of its characteristics of small computation, the interaction between the analyzed factors, and high convergence.
The hot stamping of magnesium alloy sheets is a nonlinear and large deformation process. An excessive number of optimization parameters and constraints can not only make it difficult to test a large number of samples but can also negatively influence the accuracy of the response surface model and disrupt the optimization process. Therefore, three key forming parameters (stamping speed, blank-holder force, and friction coefficient) were selected as the response variables, and the maximum thinning rate (Y) at the corner of the hot-stamped sheet was taken as the response value. In addition, considering the formability of the magnesium alloy sheet, a punching depth of 6 mm was set for simulation calculation and subsequent experimental verification. The factors and levels for the hotstamping simulation are listed in Table 3. Experimental schemes were designed using the Box-Benhnken design (BBD) method in the Design-Expert 12 software after determining the response variables and levels. The hot-stamping process of the magnesium alloy sheet was simulated with the aid of the finite element solver ABAQUS/Explicit. Table 4 shows the finally obtained maximum thinning rate of the stamping sheet under the combined parameters' scheme.

Variance Analysis of Response Surface Regression Model
The simulation results in Table 4 were fitted using multiple quadratic regression. The final regression model equation of the maximum thinning rate in terms of the actual response factors of the magnesium alloy sheet can be obtained as follows: where X 1 is the stamping speed; X 2 is the blank-holder force; X 3 is the friction coefficient; and Y represents the maximum thinning rate at the corner of the hot-stamped sheet.
To further verify whether the established response surface model can accurately express the statistical rule between the variables and the optimization objective, variance analysis was carried out for the above equation, as summarized in Table 5. The confidence interval p-value of the model < 0.0001 for the equation of the maximum thinning rate of the magnesium alloy sheet was considered extremely significant, indicating that this model could be used in this experiment. The lack of fit > 0.05 was not significant, indicating that the model was reasonable in the regression region and could be used to predict the maximum thinning rate of the magnesium alloy sheet. The multivariate correlation coefficient R 2 was 0.9912, and the correction coefficient R 2 (Adj) was 0.9799, which indicated that the model fitted well, and the prediction was close to the actual value.
It can also be concluded from Table 5 that the influence order of the three response factors on the maximum thinning rate was as follows: blank-holder force (X 2 ) > friction coefficient (X 3 ) > stamping speed (X 1 ). In addition, the model also characterized the interaction relationship between the three response factors. The interaction terms X 1×2 and X 1×3 had a significant impact on the maximum thinning rate, indicating that the stamping speed and the blank-holder force/friction coefficient had an obvious interaction. The interaction term X 2×3 had no significant effect on the maximum thinning rate, indicating that there was no significant interaction between the two groups of factors. Figure 5 shows the normal distribution probability of the residual for the maximum thinning rate. The normal distribution probability distribution of the residual is close to a straight line which is fitted well, indicating that the equation between the independent variable and the response variable is highly reliable as a result of regression analysis.   Figure 6 shows the distribution result of the actual and predicted va imum thinning rate of the hot-stamped sheets. The predicted and actual v distributed on the same line, indicating that this model has a certain fore    Figure 6 shows the distribution result of the actual and predicted va imum thinning rate of the hot-stamped sheets. The predicted and actual v distributed on the same line, indicating that this model has a certain fore  Figures 7-9 show the contour plots and corresponding three-dimen surface graphs of the interactions among response factors when the max rate of the magnesium alloy sheet was taken as the response variable. show that the contours on the Y plane are dense, and the contours corre changes in X1 and X2/X3 reveal elliptic characteristics. These results show significant interaction between the stamping speed and the blank-holder between the stamping speed and the friction coefficient.  Figures 7-9 show the contour plots and corresponding three-dimensional response surface graphs of the interactions among response factors when the maximum thinning rate of the magnesium alloy sheet was taken as the response variable. Figures 7 and 8 show that the contours on the Y plane are dense, and the contours corresponding to the changes in X 1 and X 2 /X 3 reveal elliptic characteristics. These results show that there was a significant interaction between the stamping speed and the blank-holder force, as well as between the stamping speed and the friction coefficient.     Figure 9 shows that the contour lines on the Y plane are sparse and appear as regular circles. The slope of the response surface is steeper when one of the response factors changes, and the other is gentler when the response factor changes, which indicates that the interaction between the blank-holder force and the friction coefficient was not significant in the simulation results of the maximum thinning rate of the hot-stamped sheet.    Figure 9 shows that the contour lines on the Y plane are sparse and appear as regu circles. The slope of the response surface is steeper when one of the response fac changes, and the other is gentler when the response factor changes, which indicates t the interaction between the blank-holder force and the friction coefficient was not sign cant in the simulation results of the maximum thinning rate of the hot-stamped sheet.   Figure 9 shows that the contour lines on the Y plane are sparse and appear as regular circles. The slope of the response surface is steeper when one of the response factors changes, and the other is gentler when the response factor changes, which indicates that the interaction between the blank-holder force and the friction coefficient was not significant in the simulation results of the maximum thinning rate of the hot-stamped sheet.

Response Surface Analysis of Regression Models
According to the response surface optimization analysis, the optimal configuration of process parameters was as follows: The stamping speed was 7.067 mm/s, the blank-holder force was 5.533 kN, and the friction coefficient was 0.142. Under these parameters, the optimal maximum thinning rate of the hot-stamped sheet was 7.37%.

Experimental Verification of Hot-Stamping Process
Further experiments on the hot stamping of magnesium alloy sheets were conducted to verify the forming parameters of the sheet reported in Table 2. The equipment used in the stamping forming test was a YT32-200C four-column hydraulic press, with a nominal force of 20 kN and a maximum slide stroke of 710 mm. Before the test, the magnesium alloy sheet was put into the external resistance furnace, heated up and kept for 10 min, and then quickly transferred to the die. Three groups of parameters were randomly selected (i.e., Run 4, Run 7, and Run 13). The hydraulic press stroke was set so that the punch was pressed down to the specified position for the stamping test. Then, the thickness of the sheet after hot stamping (region D in Figure 4) was measured using a Doppler ultrasonic detector.
A comparison of the maximum thinning rates obtained through simulation and testing under the same process parameters was carried out, and the comparative results are presented in Figure 10. The maximum relative error between the simulation and test values for the maximum thinning rate was 8.72%. This may be due to the absence of material damage parameters in the numerical model or the slip error between the fixture and the sheet in the test, which led to some differences between them. Additionally, this also showed that the combined method of the finite element model and response surface analysis had a high accuracy in the study of the magnesium alloy sheet's formability. lected (i.e., Run 4, Run 7, and Run 13). The hydraulic press stroke was set so that the punch was pressed down to the specified position for the stamping test. Then, the thickness of the sheet after hot stamping (region D in Figure 4) was measured using a Doppler ultrasonic detector.
A comparison of the maximum thinning rates obtained through simulation and testing under the same process parameters was carried out, and the comparative results are presented in Figure 10. The maximum relative error between the simulation and test values for the maximum thinning rate was 8.72%. This may be due to the absence of material damage parameters in the numerical model or the slip error between the fixture and the sheet in the test, which led to some differences between them. Additionally, this also showed that the combined method of the finite element model and response surface analysis had a high accuracy in the study of the magnesium alloy sheet's formability.
This work attempts to provide a reliable methodology for the optimal design of the forming process in order to produce the required parts through hot stamping in a limited number of experiments. It was proved that the response surface method combined with the finite element method can be applied to optimize the process parameters of hot stamping. By improving the precision of the numerical model and optimizing the process parameters, the thinning of hot-stamped sheets can be minimized.

Conclusions
In this paper, a constitutive model suitable for describing magnesium alloy sheet forming was introduced. The hot-stamping process of an AZ31 magnesium alloy sheet was numerically analyzed based on the finite element solver ABAQUS/Explicit. Then, the This work attempts to provide a reliable methodology for the optimal design of the forming process in order to produce the required parts through hot stamping in a limited number of experiments. It was proved that the response surface method combined with the finite element method can be applied to optimize the process parameters of hot stamping. By improving the precision of the numerical model and optimizing the process parameters, the thinning of hot-stamped sheets can be minimized.

Conclusions
In this paper, a constitutive model suitable for describing magnesium alloy sheet forming was introduced. The hot-stamping process of an AZ31 magnesium alloy sheet was numerically analyzed based on the finite element solver ABAQUS/Explicit. Then, the response surface methodology was used to study the influence of the key process parameters on the maximum sheet thinning rate, and the optimal combination of the hot-stamping process parameters was obtained. Finally, the accuracy of the numerical model and the response surface model was verified by hot-stamping experiments. The main conclusions are summarized as follows: (1) Numerical simulation results showed that the sheet near the corner region of the die was prone to thinning due to its large strain. With the established model, the numerical simulation of the hot-stamping process of the AZ31 alloy sheet could be achieved with a high level of accuracy. (2) Through analyzing the variance results of the maximum thinning rate obtained using the response surface methodology, it can be concluded that each of the process parameters affected the maximum thinning rate in the order of blank-holder force (X 2 ) > friction coefficient (X 3 ) > stamping speed (X 1 ). Hot-stamped sheets' maximum thinning rate largely depended on the interaction between X 1 and X 2 /X 3 . (3) According to an optimization analysis using the response surface methodology, the optimal process parameters for the hot-stamped AZ31 magnesium alloy sheet were as follows: The stamping speed was 7.067 mm/s, the blank-holder force was 5.533 kN, and the friction coefficient was 0.142. Additionally, the maximum relative error (8.72%) was within a reasonable range after comparing the hot-stamping experiments and simulations.
This study can provide effective guidance for the process improvement and parameter optimization of the hot stamping of magnesium alloy sheets.
Author Contributions: P.Z., conceptualization, methodology, writing-original draft preparation, and funding acquisition; Q.W., formal analysis, and writing-review and editing; Y.-L.Y., resources and project administration; Z.C., supervision and validation. All authors have read and agreed to the published version of the manuscript.