A phase-field numerical framework to study ductile damage to fracture transition: An application to material forming processes

. In this work we present a phase-field based method to accurately predict the nucleation and evolution of damage. The damage pattern is then used as a criterion for automatic, remeshing-based, 3D crack insertion and propagation. The proposed framework has been implemented in Forge® finite element software and it offers a robust numerical tool for the modeling of damage to fracture transition in complex industrial processes. A bar shearing simulation will be used to show the robustness and efficiency of the approach. The flexibility of the approach is presented through the use of different damage criteria introduced into the phase-field formulation.


Introduction
Material forming processes take advantage of the ductility of materials, in particular metals, to obtain a given geometry of a part.Depending on the complexity of the part, different stages might be needed.From a mechanical perspective, the material will be submitted to a non-proportional loading path and large deformation.As a result of these processes damage and fracture might appear within the material.Depending on the process, damage, and subsequent fracture, should be either avoided or controlled.In either case the simulation of the kinetics of the evolution of damage is fundamental in order to design new parts and manufacturing processes.The damage to fracture transition becomes particularly necessary in order to study processes where damage should be controlled (e.g., blanking and cutting processes).To this end, the first step consists in being able to predict the damage kinetics.Then the damage prediction should be used to initiate and propagate cracks into the material.
Regarding the prediction of damage, when dealing with brittle materials Griffith theory can be used to explain the evolution of fracture.In fact, this energetic approach can be used to reproduce experiments of fracture of brittle materials [1].On the other hand, when dealing with ductile materials the literature about models that allow to predict the evolution of damage is often restricted to specific material and/or particular loading conditions.It is thus fair to say that there is no unified theory that quantitatively explains the kinetics of ductile damage and fracture phenomena.However, we can say that there are three main stages on the process: 1. void nucleation 2. void growth and distortion 3. void coalescence.
Voids are nucleated due to stress concentrations and the presence of defects in the materials.These voids will grow under plastic deformation and their growth will be significantly affected by the material's stress state.And finally void coalescence is the result of the merging of previously nucleated voids [2,3].The numerical modeling of these ductile damage mechanisms has been done by using different numerical approaches: uncoupled and coupled approaches.Uncoupled approaches can be seen as a post processing tool that allows to compute the local damage state of the material by using a set of thermomechanical variables integrated over time.These approaches are easy to implement and do not affect the convergence of the finite element solver [1].They present an important drawback: uncoupled approaches do not allow to reproduce the drop of stress experimentally observed once damage has been initiated.On the other hand, coupled models link the evolution of damage to the elastoplastic properties of the material.This allows to capture the softening behavior experimentally observed but this brings along some drawbacks.The softening of the materials leads to numerical difficulties that are well documented in literature [9] and the numerical solver must be adapted to tackle these challenges.Extensive literature can be found in the context of ductile fracture applied to material forming as in [5] and [6] (and the references therein).The phase-field community has also proposed remarkable contributions to the field such as [7,8].Such approaches often focus on capturing the drop of loading carrying capacity of the material.However, the degradation of the material is taken into account within a continuum mechanics framework and thus no actual crack is introduced in the computational domain.
Once damage has developed in the material, a discontinuity should be introduced.Early attempts to introduce the discrete cracks based on a remeshing approach in 2D were proposed in [9].Crack propagation direction is given by the maximum damage at different distances from the crack tip, which led to wrong crack directions especially in the case of complex crack patterns (e.g., merging, branching and multiple cracks) and its extension to 3D was complicated.An alternative to this approach based on the marching ridges method was introduced in [10].The idea was to find the local maximum of the damage field and propagate the crack by introducing a segment linking these local maximum points.It is worth mentioning that there exist alternatives to discrete crack insertion methods (such as the XFEM [11] or the GFEM [12]) which suffer from difficulties when it comes to large plastic strain and remeshing.But these methods will not be covered in this work.
In this article, the phase-field based damage modeling strategy proposed in [13] is briefly discussed.The proposed method, based on a phase-field-inspired approach to predict the evolution of damage, will be presented in the next section.Then the phase-field damage is used to drive the crack insertion and propagation following the Crack Insertion and propagation using the Phase Field and Adaptive Remeshing (CIPFAR) algorithm introduced in [14].Finally, some results are presented.

Phase-Field Inspired Damage Formulation
Phase-field approaches to simulate damage and fracture can be seen as an extension of the generalized formulation of Griffith theory introduced by Francfort and Marigo [15].The idea is to use a regularized representation of the crack by introducing a damage variable (d).The evolution of this damage problem is driven by the minimization of an energy functional that links the mechanical behavior of the material to its damage evolution.The energy functional is given by: Where is the elastic strain tensor, is the elastic energy density, is the fracture toughness, is the regularization length and is the degradation function given by: The choice of the degradation function satisfies the following conditions:   (0) = 1,   (1) = 0 and   ′ () = 0.The minimization of the energy functional leads to the following set of equations: . = 0 (Neumann boundary condition) A similar set of equations is found when a non-local continuum damage formulation is used.The main difference is that phase-field evolution includes the term   ′ () which vanishes when the material is completely damaged ( = 1 →   ( = 1) = 0).This vanishing term prevents spurious diffusion of the damage variable [13].
The previous set of equation corresponds to the classic phase-field damage formulation for brittle materials.If we were to extend this formulation to ductile damage modeling while keeping the phase field philosophy, we would need to come up with an energy functional that accounts for the plastic dissipation which should also account for the plasticity-damage coupling.Real material exhibit complex damage behavior often involving a strong coupling damage, stress triaxiality and/or even the Lode parameter [2].Obviously proposing a general energy functional that fulfills all these requirements is not an easy task.
In order to overcome this problem of finding an appropriate energy functional, we simply replace the local history functional by a phenomenological law.The idea is to introduce any of the existing Uncoupled damage criterion into the local functional: is the local damage (damage criterion),   is a material parameter used to control the Where post-peak stress response of the material by controlling the amount of effective energy needed to create the crack surfaces once the crack is initiated and  ℎℎ is a damage threshold value used and  ℎℎ is shown in Fig. 1, which to delay the softening of the material.The effect of corresponds to the 1D homogeneous solution of the phase-field equation.

Fig. 1. The effect of model parameters on the evolution of the phase field in a 1D homogeneous
solution [16].
Once the phase-field is computed, the crack insertion strategy presented in [14] is used in order to introduce the crack discontinuity into the mesh.In the next section, a forming application is presented and discussed.The implementation and validation of the different building blocks of the proposed approach has been done in Forge®.

Results
The proposed strategy has already been validated and tested [13,14].It has been proven that the approach is very robust and flexible as different damage criterion could be used to study a given problem.In fact, as it has been mentioned before, there is not a unified theory that explains the evolution of ductile damage and therefore the choice of the damage criterion is done depending on the material and the loading conditions.Here we will present a simple case of a bar shearing.The boundary conditions of the problem are shown in Fig. 2. Depending on the choice of damage criterion, the failure pattern predicted by the model will be different.The different failure criteria tested are listed in Table 1 and the corresponding properties used in the simulations are shown in Table 2.
Table 1.Damage criteria used for the coupled phase-field simulations [16].

Damage criterion
Driving force (DF) Parameters The different crack patterns are shown in Fig. 3 for each of the models in Table 1.For each model the upper row shows the phase-field damage level ( ) and the lower row represents the cracks that have been inserted.Fig. 3a shows the results obtained when using the Normalized Latham-Cockcroft model (number 1 in Table 1) at different macroscopic loading states.It can be seen that two cracks are initiated at the left and right boundaries of the sheared zone.It is important to highlight that the approach is robust enough to handle the two cracks that propagate independently and that end up merging into a single crack that breaks the rod into two pieces.
Crack patterns obtained when using models 2 and 4 in Table 1 present a similar kinetics: first two cracks near the boundary of the rod are initiated and these cracks merge at the end in order to fully break the rod.Interestingly the crack initiation region from these two models is different from the one observed in model 1.This time, the two branches are initiated in the upper and lower boundaries of the rod.This change on the initiation region can be explained by the driving term of the damage criterion used (Column 3 in Table 1).
The von Mises stress as well as the maximum shear stress are concentrated over the regions that are in contact with the tools in the shearing regions.On the other hand, the normalized maximum principal stress is maximum over the left and right boundaries of the rod, which explains the pattern observed in model 1.Finally, model 3 (triaxiality driven) leads to a completely different crack pattern.In fact, the regions that present high triaxiality are located behind the shearing region.Therefore, the crack nucleates in the upper and lower regions of the rod but inside of the tool.Fig. 2. Geometry and boundary conditions of a bar subjected to shear loading [16].1.

Summary
To conclude, the presented examples confirm the ability of the proposed approach to accurately model the damage to fracture transition under complex loading conditions.In addition, the model provides a general framework for crack initiation and propagation in 3D that can be adopted for different applications.Since there is no a universal model that accurately predicts the ductile damage evolution, the use of custom damage criteria (often developed in the context of uncoupled approaches) combined with the phase-field formulation allows to account for complex damage phenomena experimentally observed.Once the damage has initiated, the CIPFAR algorithm is used to insert cracks into the computational domain while being robust and reliable.Further developments and tests should be carried out in order to study the impact of different damage criteria in the case of multistage non-proportional loading processes.