Temperature-dependent yield criterion for high strength steel sheets under warm-forming conditions

In this paper, uniaxial and biaxial tensile tests with cruciform specimens were conducted to investigate the deformation behaviour of dual phase steel sheet with a tensile strength of 590 MPa (DP590) under evaluated warm-forming temperatures (20–190 ◦C). Detailed analyses were then carried out to obtain the corresponding experimental yield loci. For the purpose of describing the temperature-dependent yield behaviour of DP590 appropriately, the Yld2000–2d yield function with temperaturedependent exponent was proposed. The identification procedures of the introduced parameters were then proposed based on Levenberg-Marquardt optimization algorithm. Afterwards, the proposed model was implemented into ABAQUS as user subroutine VUMAT with NICE (Next Increment Corrects Error) explicit integration scheme. The numerical simulations of biaxial tensile tests were then conducted to confirm the validity of the proposed model. It could be concluded that the flexibility and accuracy of the proposed model guarantee the applicability in warm-forming applications.


Introduction
The application of the advanced high strength steel alloys (AHSS) has been increasing in the automobile industry for the weight reduction and collision-safe performance.Among the AHSS, the dual phase (DP) steel has gained great attention since it achieves high strength with relatively good formability [1].However, it was observed that the temperature of the DP steel tends to increase evidently during the plastic deformation [2], which would inevitably affects the plastic flow and formability of the materials [3].This coupled thermal-mechanical behaviour of the DP steel might bring some difficulties in predicting the spring-back.[4] or cracking during the high-speed punching processes, which have been widely accepted in the automobile manufacturing.Therefore, a practical constitutive model is needed to characterize the plastic behaviour of the DP steel sheet under warm-forming conditions.
Recently, some studies have been carried out on the thermal-mechanical behaviours of the sheet metals.Ozturk et al. conducted the uniaxial tensile and V-bending tests under warm-forming temperatures for DP600 [1].The results indicated that the formability and springback of the DP steels decreased with the increase of the temperature.Lee et al. performed the bulge tests for a variety kinds of AHSS and observed the decline of the biaxial stresses as the initial bulge temperature increased [5].Sung et al. proposed a modified Hollomon/Voce model [6], which successfully described the strainrate and temperature sensitivity of the flow stress for the dual phase steels.Merklein et al. performed a series of biaxial tensile test under evaluated temperatures for aluminium alloy AA6016 as well as magnesium alloy AZ31 [7].It was discovered that the yield loci and work hardening behaviour of these two materials were influenced by the temperature significantly.Similar approaches were also carried out by Tetsuo Naka et al. to identify the suitable yield criteria for 5083 aluminium alloy and AZ31 magnesium alloy [8,9].The results revealed that the Logan-Hosford and Yld2000-2d yield criteria could fit the experimental yield loci with better precision.To the best of the author's knowledge, few studies have been carried out on the yield criteria for the AHSS under warm-forming conditions, which would be significant in enhancing the precision of the corresponding FE simulations.In this paper, uniaxial and biaxial tensile tests with cruciform specimens were conducted to investigate the deformation behaviour of DP590 under evaluated warm-forming temperatures (20-190 • C).Detailed analyses were then carried out to obtain the corresponding experimental yield loci.For the purpose of describing the temperature-dependent yield behaviour of DP590, the modified Yld2000-2d yield criterion was proposed.In this model, the exponent of the Yld2000-2d was assumed as the function of the temperature.The identification procedures of the introduced parameters were then proposed based on Levenberg-Marquardt optimization algorithm.Afterwards, the proposed model was implemented into ABAQUS as user subroutine VUMAT.The numerical simulations of biaxial tensile tests were then conducted to confirm the validity of the proposed model.

Material
The materials used in this paper were the dual phase steel sheets with tensile strength of 590 MPa (DP590), which produced by BaoSteel, Inc.The thickness of the sheets was 1.5 mm, and the mechanical properties under evaluated temperatures are listed in Table 1.The 0 , 45 and 90 in Table 1 represent the initial yield stresses (with plastic strain of 0.2%) obtained from the uniaxial tensile tests along the directions of 0 • , 45 • and 90 • from the rolling direction, respectively.

Testing apparatus and procedures
The biaxial tensile testing system used in this study includes a loading test machine, a heating device, and related control unit.It was improved based on the prototype originally designed by Wu et al. [10].The loading test machine, shown in Fig. 1(a), has two pairs of hydraulic cylinders, and the hydraulic pressure of these cylinders is servo-controlled independently.The heating device is based on the resistance heater with closed-loop temperature control.In order to lower the shear stress in the gauge area, the shape and dimension of the cruciform specimen were optimized by orthogonal design with FEM simulations.The biaxial tensile tests with linear loading path were conducted to verify the applicability of the yield criteria for the dual-phase steel under warm conditions.9 fixed ratios were selected as the biaxial loading proportions, e.g.F x :F y = 4 : 0, 4:1, 4:2, 4:3, 4:4, 3:4, 2:4, 1:4, 0:4.For the load ratios of 4:0 and 0:4, the standard uniaxial tensile tests were adopted.

Measurement of the experimental yield loci
In order to investigate the initial and subsequent yield behavior of the materials, the concept of plastic work contour was adopted [11].Plastic work contours are commonly considered as yield loci for simplicity.According to this concept, two different stress states are considered to be on the same yield locus only if their plastic work per unit volume is equal.Therefore, the relationship between the biaxial and uniaxial stress state with the same equivalent plastic strain could be written as where the 1 , p 1 and 2 , p 2 are the stress and plastic strain along the rolling and transverse directions respectively.The ¯ and ¯ are the equivalent stress and equivalent plastic strain, normally referring to the true stress and logarithmic strain obtained by uniaxial tensile test from rolling direction.A detailed description of calculating the plastic work contour was given by Kuwabara et al. [11].

Yield loci under evaluated temperatures
The experimental yield loci of DP590 under different temperatures, shown in Fig. 2(a), (b) and (c), were determined according to Eq. (1).Also depicted in the figures were the theoretical yield loci based on Yld2000-2d yield criterion [12].The exponent m of the Yld2000-2d was chosen as 4, which was reported as the suitable value for DP590 in Ref. [13,14].For the propose of evaluating the theoretical yield criteria, a deviation function was introduced for each yield locus (2) where i 1 and i 2 are the stress components of the experimental yield points measured from the rolling and transverse directions, respectively; d i stands for the distance between the experimental yield point and corresponding theoretical one at the same stress ratio; n stands for the number of the experimental points at each yield locus.The deviation values between the experimental and the theoretical yield loci based on the von Mises, Hill's quadratic (Hill'48) and Yld2000-2d yield criteria were calculated and plotted in Fig. 2(d).Figure 2(d) indicates that, compared with the Mises or Hill'48 yield criterion, the Yld2000-2d with exponent of 4 has the closest agreement with the obtained experimental yield loci under given temperatures.

Analyses on the variation of the yield loci
As indicated in Fig. 2, the outlines of the experimental yield loci in the first quadrant vary evidently with the growth of the temperature, which leads to the decreased accuracy of the Yld2000-2d at the temperatures of 60 • C and 190 • C. Since the shape of the Yld2000-2d yield locus strongly depending on the value of exponent m [14], it could be a feasible idea to determine a suitable exponent for each evaluated temperature to enhance the precision of the model in a wider temperature range.Previous studied carried out by Barlat et al. [12] and Kuwabara [14] also indicated that flexible value of the exponent m could be accepted to characterize the behaviours of variety kinds of sheet metals.
Figure 3(a) presents the deviation's trend of the theoretical yield loci with the increase of exponent m.It confirms that an optimal exponent m opt do exists for each evaluated temperature.For the purpose of obtaining m opt , the following optimization procedure was proposed min ICNFT 2015  where g 1 ( , m) is the objective function of the optimization procedure; P i exp stands for the experimental yield point; P i theory stands for the corresponding theoretical yield point predicted by Yld2000-2d yield function; stands for the set of the anisotropy coefficients determined by the material parameters and the exponent m; T stands for the Celsius temperature; N stands for the total number of the experimental yield points at the current yield locus.The nonlinear equations F ( , m) = 0 and G( , m) = 0 are the constraints of the anisotropy coefficients, which were firstly proposed by Barlat et al. [12] aiming to determine the anisotropy coefficients using eight material parameters.Detailed definition of F ( , m) and G( , m) were adequately elaborated in the Ref. [12].With classic Courant Penalty function method, the constrained minimization problem 3 could be converted into unconstrained one expressed as follows: min where g 2 ( , m) is the objective function of the optimization procedure; p is the penalty parameter and could be set to one in this case.The minimization 4 has been performed using Levenberg-Marquardt algorithm implemented in MATLAB and the results were then plotted in Fig. 3(b).
For the purpose of modeling the variation trend of the optimized exponent, the following function was selected as the evolution model where m 1 , m 2 , k 1 and k 2 are the introduced parameters to describe the evolution of the yield behaviour under warm conditions.Afterwards, the nonlinear fitting was performed on the obtained (T , m opt (T )) to acquire these introduced parameters.The fitting results were then listed in Table 2.

Implementation and verifications of the proposed model
The above presented model has been implemented in a general purpose finite element code via VUMAT subroutine of ABAQUS.The NICE (Next Increment Correct Error) explicit algorithm, firstly proposed by Halilovic et al. [15], was adopted as the integration scheme for the constitutive model.This algorithm was derived from forward-Euler explicit algorithm by replacing the consistency condition with where f represents the yield function of the constitutive model.This amendment not only enhances the accuracy and stability of the algorithm, but also retains the computational efficiency of the explicit scheme, which would be significant to the dynamic explicit simulations.
The FE simulations of biaxial tensile tests were conducted to verify the proposed model for DP590 under warm-forming conditions.Due to the symmetry of the specimen, a quarter FE model with 4-node shell element (ABAQUS S4R) was introduced, shown in Fig. 4(a).The FEM predicted stress-strain curves with different yield criteria were presented in Fig. 4(b).It can be concluded that the simulated results based on the proposed model are in better agreement with the experimental results than the original Yld2000-2d (with exponent of 4).

Conclusions
(1) The experimental yield loci of DP590 under warm-forming condition (20-190 • C) were obtained based on biaxial tensile tests.It was found that the Yld2000-2d yield criterion with exponent of 4 could describe the yield behaviour of the DP590 under evaluated temperatures with relatively acceptable precision.(2) In order to describe the temperature-dependent yield behaviour of DP590, the Yld2000-2d yield criterion with temperature dependent exponent was proposed.The identification procedures of the introduced parameters were then proposed based on optimization algorithm.(3) The proposed model was implemented into ABAQUS as user subroutine VUMAT.The biaxial tensile simulations were then conducted to confirm the validity of the proposed model.

Figure 1 .
Figure 1.(a) The biaxial tensile testing apparatus; (b) geometry of the cruciform specimen (dimensions in mm).

Figure 2 .
Figure 2. (a)-(c) Comparisons between experimental and theoretical yield loci for DP590 under evaluated temperatures; (d) deviations between the experimental and the theoretical yield loci with different yield criteria.

Figure 3 .
Figure 3. (a) Deviation's trend of the theoretical yield loci with the increase of exponent m under evaluated temperatures; (b) optimized exponents of the Yld2000-2d for DP590 under evaluated temperatures.

07008-p. 5 MATECFigure 4 .
Figure 4. (a) FE model for the biaxial tensile simulations; (b) comparison between the experimental and simulated stress-strain curves along the rolling direction of the biaxial tensile tests under 60 • C.

Table 1 .
Mechanical properties of DP590 under different temperatures.

Table 2 .
Nonlinear fitting results on the optimized exponents.