Forming limit diagram prediction of 6061 aluminum by GTN damage model

. Forming limit diagram (FLD) is one of the formability criteria which is a plot of major strain versus minor strain. In the present study, Gurson-Tvergaard-Needleman (GTN) model is used for FLD prediction of aluminum alloy 6061. Whereas correct selection of GTN parameters ’ is effective in the accuracy of this model, anti-inference method and numerical simulation of the uniaxial tensile test is used for identi ﬁ cation of GTN parameters. Proper parameters of GTN model is imported to the ﬁ nite element analysis of Nakazima test for FLD prediction. Whereas FLD is dependent on forming history and strain path, forming limit stress diagram (FLSD) based on the GTN damage model is also used for forming limit prediction in the numerical method. Numerical results for FLD, FLSD and punch ’ s load-displacement are compared with experimental results. Results show that there is a good agreement between the numerical and experimental results. The main drawback of numerical results for prediction of the right-hand side of FLD which was concluded in other researchers ’ studies was solved in the present study by using GTN damage model.


Introduction
Aluminum alloys are extremely used in different industries such as automotive and aircraft. 6061 alloy is one of these alloys with high strength and corrosion resistance which are used in the automotive, aircraft, highway and marine applications. Therefore, 6061 aluminum sheets are used in a different forming process such as stamping, deep drawing, and hydroforming to convert to the desirable parts. Necking and fracture are the main defects which happen in the sheet metal forming. Therefore, forming limit prediction of sheet metals in the forming processes is important to prevent necking and fracture in the final part. Forming limit diagram (FLD) is used as a criterion for forming limit prediction of sheet metal in the forming processes. It displays in principal strain space (major and minor strains) at the onset of local necking. There are three methods for FLD calculation: experimental, numerical and analytical methods. Whereas, experimental method is time-consuming and expensive, numerical and analytical are useful methods. There are many different analytical models of ductile damage for fracture prediction, but Gurson-Tvergaard-Needleman (GTN) criterion is one of the well-known micromechanical models for ductile fracture which is based on the void growth and coalescence phenomenon. This model was purposed by Gurson [1]. The GTN damage model has been used in many commercial software as ductile damage criterion. This model was used by researchers for fracture prediction in different metal forming processes. There are different parameters in this model which incorrect identification of them has an influence on the fracture prediction. Slimane et al. [2] studied the phenomenon of nucleation by GTN model using the numerical simulation of axisymmetric notched samples. Their results showed that void growth and failure of specimens affected significantly by Tvergaard's parameters. Alegre et al. [3] used GTN model for fracture prediction of 15.5 PH steel in the pre-cracked small punch tests (P-SPT). The effect of damage parameters was analyzed using numerical simulation. He et al. [4] used modified GTN model to identify FLD and forming limit stress diagram (FLSD) of aluminum alloy 5052. There was a good agreement between numerical and experimental results. Two parameters for void growth mechanism and void shear mechanism were introduced into yield function by Jiang et al. [5]. They used modified GTN damage model for crack initiation and growth in the numerical analysis by ABAQUS/Explicit. Numerical results were verified by experimental tests of the tensile bar, flat grooved plate, torsion tube and compression cylinder. Abbasi et al. [6] used GTN damage model for FLD prediction of interstitialfree steel (IF steel). They used finite element method (FEM) and response surface methodology (RSM) for identification of GTN parameters. Chen and Dong [7] developed GTN damage model based on Hill's quadratic yield criterion and isotropic hardening rule. This model was used as a user-defined subroutine for damage prediction of deep drawing of a cylindrical cup. Their results showed that anisotropy of ductile sheet metal has an influence on the damage evolution and deformation behavior of the material. Hu et al. [8] used the simulation of uniaxial tensile test for parameters identification of GTN model by inverse FEM. Numerical results had a good agreement with experimental results for the AA6061 sheet. Influence of GTN damage parameters of f 0 and f n were investigated on the mechanical properties of AA6061. Amaral et al. [9] used three different damage models of Lemaitre [10], Gurson [1] and Johnson and Cook [11] for damage prediction of the AA5182-O sheet. These damage models were imported to the finite element (FE) analysis of ABAQUS for damage prediction. Their results showed that GTN and Lemaitre have a good agreement with experimental results for damage prediction. Yu et al. [12] used the tensile test to study fracture behavior of rolled and aged ultrafine grained AA6061 sheets. Three different damage criteria of Cockcroft and Latham [13], Tresca and GTN were used in the FE simulation of the tensile test. Their results showed that fracture angle decreases with the number of rolling passes increasing. GTN criterion has the best agreement with experimental results for fracture prediction. Safdarian et al. [14] used different numerical criteria for FLD prediction of tailor welded blanks (TWBs) of interstitial-free steel sheets. These numerical methods were Müschenborn-Sonne forming limit diagram, forming limit diagram criterion and ductile fracture criterion (DFCcrt) and numerical method of the second derivative of thinning. Their results showed that numerical methods are successful for prediction of the left-hand side of FLD, but not successful for the right-hand side. Safdarian [15] used numerical criteria of ductile fracture and forming limit diagram for FLD and FLSD prediction of AA6061. Numerical results were compared with experimental ones. Results showed that DFCcrt has good accuracy for FLD and FLSD prediction. Safdarian [16] used modified Marciniak-Kuczynski (M-K) model for necking prediction of IF TWB by adding bending strain was to the M-K model. This criterion was used as a script in the numerical simulation of ABAQUS using Python programing language. Afshar et al. [17] studied the FLD of tubular material (Al 7020-T6) using numerical and experimental methods. Numerical method of the acceleration of plastic strain was applied to compute the hydroforming strain limit diagram in their study. The numerical results were verified with experimental ones which were a good agreement between numerical and experimental results. Boissiere et al. [18] studied the effects of punch geometry and sample size on the FLD of a 2024 aluminium alloy. Four configurations of flat punch (Marciniak test) or hemispheric punch and decimetric vs. centimetric tooling dimensions were selected. Their results showed that hemispherical punch allows measuring higher strains and is less sensitive to size effect than Marciniak test. Parente et al. [19] investigated the effect of weld line orientation on the formability and FLD level of aluminum TWBs. They used friction stir welding to join aluminum TWBs which consist of AA 6061 and AA 5182 with equal thickness as base metals. Their results showed that the level of FLD and formability of aluminum TWBs will be decreased by increasing of weld line orientation. Gatea et al. [20] used GTN damage model with the consideration of shear for fracture prediction in the single point incremental. The results showed that the shear modified GTN model improved the modeling accuracy of fracture over the original GTN model under shear loading conditions. Teng et al. [21] used GTN model for ductile fracture and the FLD prediction of aluminium alloy 5A06 sheet. The results show that the GTN fracture criterion can give a good prediction under high stress triaxiality. Peng et al. [22] used different model for failure prediction in the hydroforming of sheet metals. Their results showed that the M-K model and the GTN-Thomason model are revealed to be able to accurately predict the ultimate pressure and the height at the onset of failure by comparing to the experimental results.
In the present study, the GTN damage model is used for FLSD and FLD prediction of 6061 aluminum sheet. Inverse FEM of the uniaxial tensile test is used for GTN parameters identification. For this purpose, different sets of GTN parameters are used to simulate the uniaxial tensile test in the ABAQUS software. Then, the optimum parameters are selected based on the comparison of experimental and numerical stress-strain curve. The identified GTN parameters are used in the numerical simulation of forming limit tests. FLD of AA 6061 sheet is characterized experimentally using Nakazima FLD test. FLD, FLSD and punch's load-displacement of FEM are compared with experimental results.

Experimental materials and properties
Many industries such as automotive and aircraft industry use aluminum alloy 6061 sheets which have high strength to weight ratio for production of body parts. In this study, FLD and FLSD of AA6061 with the thickness of 1 mm are investigated using the numerical method of GTN. The chemical composition of aluminum alloy 6061 is shown in Table 1. Figure 1 shows engineering stress-strain curve of aluminum alloy 6061. Mechanical properties of this aluminum are shown in Table 2. These mechanical properties are yield stress (YS), ultimate tensile strength, work hardening exponent (n), work hardening coefficient  [24] for FLD calculation in the experimental tests. As Figure 2 shows, these samples have equal length of 175 mm and different width from 25 to 175 mm. Specimens were grid marked with circles of 2.5 mm by an electrochemical etching method to measure major and minor strains after samples deformation.

Experimental set up for FLD
FLD or forming limit curve (FLC) is a criterion for formability prediction in the sheet metal forming. This curve is drawn using minor principle strain and major principal strain. When the strain state of any region of sheet metal in the forming process is under the FLC, that region is safe. One of the standard methods for FLD calculation is Nakazima [25] test. In this study, Nakazima standard method and die was used for FLD calculation in the experimental part of this study. Based on this method a hemispherical punch of 101.6 mm diameter was used on a 200 kN hydraulic press to form the uniform blanks of Figure 2 until fracture. The experimental setup (punch, die, blank holder and data acquisition system) is shown in Figure 3.
Oil was used as lubrication between punch and sheet surface contact in all tests. Punch speed was 20 mm/min. An optimum blank holding force in the range of 60-100 kN was applied on the upper die. The press was equipped with load and displacement sensors and experiments were stopped when forming load decreased suddenly [15].

GTN model
There are many different models for material degradation which called damage. Ductile damage is the failure mode that occurs when a structure is subjected to an increasing uniform loading, wherein the material can endure plastic deformations. Ductile damage starts with nucleation of micro-cavity which mostly happens at weak material points such as grain boundaries. Then, cracks were created with growth and coalesce of voids. The GTN is one of the micromechanical damage models which is used for fracture prediction of a ductile material. At first, Gurson [1] introduced a method for considering the effect of the void in the yield surface. This method was based on the modeling a spherical void in a unit cell and the following yield surface was suggested by Gurson: where s eq is the Von Mises equivalent stress, s m and s y are the mean stress and YS of the material, respectively. f is the void volume fraction (VVF) which is the ratio of the total volume of cavities to the volume of the body. Tvergaard [26] considered the effect of cavities interaction by parameters of q 1 , q 2 and, q 3 for better agreement between experimental results and prediction. Tvergaard and Needleman [27] substituted the VVF f in equation (1) by modified VVF f * to consider the loss of load bearing capacity because of void coalescence. Finally, the following equation was accrued by considering the above modification.
The relationship between f * and f is given as follows: where f c is critical VVF and f f is VVF at total failure. Total VVF increasing is caused by the growth of existing voids and nucleation of new voids. Therefore, rate of VVF can be shown by the following equation: The plastic incompressibility of the circumambient matrix of the voids based on the mass balance in representative volume elements is used to specify the void growth rate.
where e pl kk is the plastic hydrostatic strain. Normal distribution of void nucleation was proposed by Chu and Needleman [28] using the following equation: Length: 175mm Length: 175mm  where _ e pl is the equivalent plastic strain of the material, S N standard deviation, f N the volume fraction of void nucleation, e N mean value of the distribution of plastic strain. According to the equations which were mentioned for the GTN model, nine coefficients require to be identified. These coefficients are: f 0 , f c , f N , f f ,e N , S N , q 1 , q 2 and q 3 . Correct selection of these coefficients is necessary for successful prediction of fracture using GTN model. In the present study design of experiment of RSM is coupled with numerical simulation of tensile test to identify the parameters of GTN model.

Parameters identification of GTN model
There are nine parameters in the GTN model that should be correctly identified. The values of some of these parameters are different for each material and there is not in the literature. Therefore, some values are considered for unknown parameters based on the literature. Then, a design of experiment and anti-inference method is used to identify proper parameters of GTN model.
Based on the literature [12], 12 runs were prepared with different parameters for f 0 , f c , f N , f f as presented in Table 3. Other parameters of GTN are presented in Table 4 as known parameters of GTN damage model which are selected based on the literature [12]. These parameters were constant in the anti-inference analysis. Therefore, 12 models of uniaxial tensile test with different parameters for GTN were prepared and simulated in the ABAQUS software.

Numerical simulation
In the numerical part of this study, two types of simulation were done. In the first type, 12 uniaxial tensile tests were modeled and simulated for identification of GTN parameters under the same condition of the experiment using the ABAQUS/Explicit. Anti-inference method [4] through the combination of FE simulation with the experiment of the uniaxial tensile test was used for GTN parameter identification. VVF criterion was used for identification of the first element which starts to damage. When the VVF exceeds the f f , element start to damage and stress-strain curve of the damaged element was compared with experimental results. Maximum stress and strain at the maximum stress of the first damaged element in the numerical tensile test were compared with experimental results. Then the set of parameters which caused the minimum of difference between numerical and experimental stress-strain curve were selected for the second type of simulation.
In the second type of simulation, seven model FLD test were prepared based on the Nakazima experimental standard test [25]. The Numerical model consisted of a hemispherical punch, blank holder, die and the blank as shown in Figure 4. Punch, die and blank holder was modeled as analytical rigid parts because they have negligible deformation. The blank was modeled as a deformable part using four nodes Kirchhoff thin shell elements (S4R). The circular draw-bead model is obtained by constraint forces applied on a circular partition of the sheet at a distance of 66 mm from the center of the die. The die was fixed and the punch was moved downward with a numerical speed of 1000 mm/s. This speed was selected based on the quasi-static condition of forming process [15]. 3 Results

Parameters identification of GTN
The real stress-strain curve of numerical tensile samples is compared with experimental one in Figure 5. As this figure shows, the variation of GTN parameters influence on the stress-strain curve, maximum stress and strain at the maximum stress of tensile samples. Moreover, GTN parameters influence on the fracture strain and elongation of tensile samples. For a better understanding of GTN parameters' effect on the mechanical properties of AA6061, maximum stress and strain at the maximum stress of numerical samples are compared with experimental results in Table 5. Two last column of this table shows error values of stress and strain comparison of experiment and numerical results. As Figure 5 and error values of Table 5 show, GTN parameters of run 2, 4, 5 and 11 cause elongation decrease of tensile samples and more difference between experimental and numerical stress-strain curve. Although stressstrain curve of run 1, 3, 6, 7, 9, 10 and 12 are near to the experimental results with acceptable error values, the best fracture trend just happens for run 10. Error values for maximum stress and strain at the maximum stress of run 10 are 6.375 and 0.065 percent, respectively. Stress-strain curve of this run is separately compared with experimental results in Figure 6.
Based on the results of anti-inference method and comparison of numerical and experimental stress-strain curve of the uniaxial tensile test, GTN parameters of run 10 are the best parameters for damage prediction of 6061 aluminum alloy. Therefore, these parameters were used in the second part of FE analysis which is FLD prediction of AA6061 by GTN model. Table 6 shows these parameters.

Load-displacement curves
Punch's load-displacement curve of numerical samples of FLD tests was compared with experimental ones in Figure 7. As this figure shows, numerical curves have a good agreement with experimental results. For both numerical and experimental curves, load increase until necking starts to happen in the samples and suddenly decrease at necking point. As this figure shows, for experimental samples force reduction happen suddenly.
For the better understanding of numerical and experimental comparison, punch's maximum force and related displacement of FLD samples have been summa-  Experimental and numerical samples after fracture have been compared in Figure 8. As this figure shows, for all numerical samples when the VVF exceeds the f f (0.6) fracture happen. For all numerical and experimental samples fracture happen near the pole.

Forming limit diagram (FLD)
Numerical FLD of present study which is based on the GTN model has been compared with experimental FLD and reference [15] in Figure 9. FLD of the present study was extracted from numerical samples using VVF criterion of GTN damage model. When the VVF of an element of FLD sample exceed f f , major strain and minor strain of element which has the maximum major strain was extracted. This method was repeated for all FLD samples (seven samples) and extracted data was imported to the FLD. DFCcrt is s fracture criterion in the ABAQUS software which was used by Safdarian [15] for FLD prediction of AA 6061. As Figure 9 shows, this criterion has good accuracy for prediction of the left-hand side of FLD, but for the righthand, it is not useful. Predicted FLD by GTN model which was used in this study has a good accuracy for prediction of both sides of FLD.

Forming limit stress diagram (FLSD)
One of the main drawbacks of FLD criterion for formability prediction is its dependency to strain path.  Table 6. GTN parameters of 6061 aluminum alloy.  This criterion is valid only for cases of proportional loading, where the ratio between the principal stresses remains constant throughout the forming process. Complex parts in the industry are usually manufactured in multi-step and influence of non-proportional strain history can influence on the FLD. Therefore, FLSD is used to overcome this drawback. When strain-based FLCs are converted into stress-based FLCs, the resulting stressbased curves have been shown to be minimally affected by changes to the strain path [29]. Experimental and numerical FLSD of the present study are compared in Figure 10. Experimental FLSD was calculated from FLD by using the equation of yield function and some other plasticity equations which were presented in reference [15]. The Numerical curve of FLSD was obtained by GTN model. When the VVF of GTN model exceed f f , major and minor in-plane stresses of sample's element which has the maximum stress were extracted. This method was repeated for all samples and extracted data were imported to the FLSD. As Figure 10 shows, numerical results are in the safe region and have good agreement with experimental results.

Conclusion
In the present study, GTN damage model was used for FLD and FLSD prediction of 6061 aluminum alloy. Antiinference method and numerical simulation of the uniaxial tensile test by ABAQUS/Explicit were used for identification of GTN parameters. Results of this analysis showed that the best values for GTN model of 6061 aluminum alloy are: f 0 = 0.0001, f N = 0.0005, f f = 0.06 and f c = 0.02. Other parameters of GTN model were selected from the literature. Identified parameters caused the best prediction for fracture stress and strain, maximum stress and, strain at the maximum stress of uniaxial tensile test simulation.
Proper parameters of GTN model was imported to the numerical simulation of FLD tests of 6061 aluminum alloy. Punch's load-displacement curves, FLD and FLSD of these tests were compared with experimental ones. There was a good agreement between the numerical and experimental results. The numerical results using GTN damage criterion has also good accuracy for prediction of fracture position in the FLD samples. Although in other studies [15,30] was indicated that numerical methods have a drawback for prediction of right the hand-side of FLD, the numerical method of present study based on the GTN damage model has a good accuracy for prediction of the right-hand side of FLD. Whereas FLD dependent to the strain path, GTN damage model was used for prediction of FLSD. When VVF parameter of this criterion exceeds f f , major and minor stresses of FLD samples were imported to the FLSD. Results showed that predicted FLSD is in the safe region and near to the experimental results. Whereas experimental FLD tests are expensive and time-consuming, numerical methods such as GTN damage model are a useful alternative method for FLD prediction in the sheet metal forming.