Modeling the ductile fracture and the plastic anisotropy of DC01 steel at room temperature and low strain rates

This paper presents different extensions of the classical GTN damage model implemented in a finite element code. The goal of this study is to assess these extensions for the numerical prediction of failure of a DC01 steel sheet during a single point incremental forming process, after a proper identification of the material parameters. It is shown that the prediction of failure appears too early compared to experimental results. Though, the use of the Thomason criterion permitted to delay the onset of coalescence and consequently the final failure.


Introduction
There are currently several damage models available in the literature in order to describe the material degradation. The most popular micromechanical model for ductile damage is the Gurson-Tvergaard-Needleman (GTN). Gurson [1] first proposed a porous plasticity model with the void volume fraction as internal variable. The model was further improved to account for the loss of load carrying capacity associated with void nucleation and coalescence [2], hence considering the three stages of damage: void nucleation, growth and coalescence. Since then, the GTN model has been extended to cover different applications. For instance, Benzerga and Besson [3] performed the same analysis as Gurson [1] but considering an anisotropic matrix of the Hill'48 type [4]. A fully implicit integration scheme was implemented by Guzmán [5] for this anisotropic GTN model, also including the effect of the kinematic hardening.

Material model
The hardening of the DC01 steel is described by both Swift isotropic and Armstrong-Frederick kinematic hardening models. Swift law is formulated as follows: where is the flow stress, is the equivalent plastic strain and K, n, 0  are material constants identified from a tensile curve in the rolling direction. Armstrong-Frederick law is expressed as: where and sat are the kinematic hardening saturation rate and value, respectively. This backstress X ij is used to define the shifted stress tensor: ̃= − . The equivalent stress of DC01 is defined using Hill'48 criterion [4]: where F, G, H, L, M, N are material parameters. Note that if = = = 1 and = = = 3, Hill'48 yield criterion recovers the isotropic von Mises criterion. The classical GTN yield surface is defined by: where ̃ is the macroscopic mean stress and f is the void volume fraction (also called porosity), defined as the average ratio of the void volume to the total volume of the material. The damage parameters 1 and 2 , originally equal to 1.0 in the initial Gurson model, are usually set to 1.5 and 1.0 allowing the continuum model to be in good agreement with the localization strain for cell analysis [2]. The extended Gurson model used for the simulations has been implemented by Guzmán [5], featuring Hill'48 type anisotropy of the matrix coupled with a mixed hardening law and classical nucleation and coalescence laws proposed in the GTN model. The evolution of voids is additively decomposed in the nucleation n f and growth part g f , hence: with as the equivalent plastic strain in the matrix and , , as material parameters. Coalescence is triggered when the porosity reaches a critical value. Its effect appears as an accelerated evolution of voids through the effective porosity * function: cr cr where the failure porosity F and the critical coalescence porosity cr f are material parameters in the classical GTN model, while the critical porosity is supposed to be reached when the following criterion is no more fulfilled in the Thomason coalescence model [6]:  is the maximum principal stress,  is a material parameter defined as a function of the hardening exponent n,  is equal to 1.24 and  is the void space ratio [6]. Finally, the model was extended to account for the directional growth of voids at low triaxiality (e.g. during shearing), following the Nahshon and Hutchinson shear approach [7].

Parameter identification
The elasto-plastic behavior of the DC01 steel (including isotropic and kinematic hardening parameters) was characterized from experimental tensile, shear and plane strain tests along different directions in the plane of the sheet [5]. The methodology is based on using two optimization techniques: OPTIM software [8] based on inverse modeling and the classical simulated annealing [5]. The obtained parameter sets are given in Table 1. The identification of the GTN model parameters does not follow the same methodology than the one used for the elasto-plastic behavior. The nature of the model, mixing macroscopic and microscopic variables, implies that the identification of the parameters should rely on various experimental measurements [5] (fractography, optical microscopy, strain field measurements on complex geometry samples). The obtained set of parameters is shown in Table 2.

Single Point Incremental Forming (SPIF)
The studied SPIF process consists in the forming of a truncated cone of 30mm in depth from a circular sheet having a diameter of 182mm and a thickness of 1mm [5]. The wall angle is defined by the path of the spherical tool, whose diameter is 10mm. Various wall angles were tested in order to evaluate the formability of the material, quantified as the maximum achievable wall angle without failure and the obtained experimental value is 67° [9].
The finite element (FE) simulations were performed using the updated Lagrangian FE code Lagamine developed by the ArGEnCo Department of the University of Liège. The 8-node 3D Solid-Shell element, called RESS [10], was selected. It uses the enhanced assumed strain method with only one additional deformation mode. It involves a reduced in-plane integration scheme (with only one integration point) and an hourglass control technique, while the number of integration points along the thickness of the sheet can be adapted (a value of 3 was chosen in the present study). The SPIF process was modeled with 1492 RESS element on a 90° angle pie (using rotational boundary conditions) [5].
The goal of this study was to evaluate the accuracy of variants of the GTN damage models to predict the failure of the sheet for different wall angles. Four variants, available in our FE code, were selected:  GTN  Classical GTN model where coalescence is initiated by the parameter cr f .  GTN + Thomason  idem but Thomason criterion is used for coalescence (the critical coalescence porosity cr f is computed from Thomason criterion).  GTN + Shear  GTN model with Nahshon and Hutchinson shear extension [7].  GTN + Shear + Thomason  idem with Thomason coalescence criterion. For each variant of the GTN model, several simulations were performed with increased values of the wall angle (by step of 1°) until failure is predicted. The numerical results for the maximum achievable wall angle are presented in Table 3 and compared to experimental result. For each case, a wall angle of 1° more resulted in the failure of the sheet before the end of the SPIF process. The maximum values on the whole FE mesh of the porosity reached when the coalescence starts (this value is only meaningful for the Thomason coalescence criterion) and the maximum effective porosity reached at the end of the process are also presented in Table 3. It can be noted: (i) The maximum achievable wall angles predicted by the variants of the GTN model are significantly smaller than the experimental value. (ii) The shear extension has a very limited influence on the results. (iii) The Thomason coalescence criterion permitted to increase the maximum achievable wall angle by delaying the onset of coalescence. Indeed, the porosity attained when Thomason criterion is no more fulfilled is way larger than the parameter cr f of the classical GTN model. (iv) The maximum effective porosity exceeded the failure limit F f . However, such values only appear very locally for the simulations in Table 3, which was not considered as a failure of the DC01 sheet in this study. Figure 1 presents numerical results when the failure is predicted (maximum achievable wall angle exceeded by 1°) for the GTN + Shear + Thomason variant. It is shown that the porosity reaches large values only in a zone around the final path of the tool; locally, the porosity can be significantly larger than the failure limit. According to figure 1(b), coalescence appears in a similar zone. Also, the edge effects in figure 1 prove the limitation of the rotational boundary conditions used.