Advanced numerical simulation based on a non-local micromorphic model for metal forming processes

An advanced numerical methodology is developed for metal forming simulation based on thermodynamically-consistent nonlocal constitutive equations accounting for various fully coupled mechanical phenomena under finite strain in the framework of micromorphic continua. The numerical implementation into ABAQUS/Explicit is made for 2D quadrangular elements thanks to the VUEL users’ subroutine. Simple examples with presence of a damaged area are made in order to show the ability of the proposed methodology to describe the independence of the solution from the space discretization.


Introduction
Nowadays, the fully local constitutive equations have been well established to model the induced material softening behaviour due to thermal, damage and other microstructure-dependent phenomena. However, the solutions of the fully local constitutive equations are highly sensitive to the space and time discretization. The mechanics of generalized continua makes possible the straightforward introduction of the characteristic lengths into the constitutive equations of materials with microstructure in order to overcome this drawback [7]. Among others, the most recent and comprehensive review of these generalized continua and their use to solve various problems in mechanics of solids and fluids is found in the recent books and publications ( [22,23], [24], [2], [7]). The micromorphic theory was proposed simultaneously by (Eringen and Suhubi, 1964[25]; Mindlin, 1964[26]). It consists of introducing a general non-compatible full field of micro-deformation as an extra degree of freedom, in addition to the classical displacement field. The micromorphic approach can, in fact, be applied to any macroscopic quantity in order to introduce a character length scale in the original classical continuum model in a systematic way, as presented by ([7], [8]). From the comparison between nonlocal and micromorphic theories, Forest and Aifantis [27] concluded that when the micromorphic variable remains as close as possible to the plastic strain, the micromorphic model reduces to the strain gradient theory.
The main goal of this paper is to show the efficiency of the formulated fully coupled nonlocal micromorphic constitutive equations in getting really mesh independent solutions even with high values of the ductile damage. The concept of internal lengths, that have to be experimentally determined from accurate measurements by using advanced methods, is introduced.
The introduction of the new balance equations in the nonlocal model with the effect of the ductile damage demands the construction of an element based on the assumed strain method combined with the Hu-Washizu variational form. The standard B matrix obtained from the displacement field analytically and including the hourglass control part, is projected to eliminate volumetric and shear locking phenomena. An air bending test of a DP1000 dual phase steel sheet is performed to validate and to show the mesh independent solution predicted by the micromorphic nonlocal model. This kind of nonlocal elastoplastic constitutive equations with damage can be later used for the numerical simulation of various sheet and bulk metal forming processes.

The generalized principle of virtual power
The virtual power of internal forces of the classical local continuum is extended to the micromorphic continuum by using the virtual micromorphic damage ( [2,7]): : The virtual power of the inertia forces are enhanced by the micromorphic phenomena: where, d ] d is a scale factor which maps the local density to the micromorphic level. The generalized virtual power enhanced with the micromorphic damage takes the following form for any given kinematically admissible fields: Substitute Eq. (1) to Eq. (3) into the generalized form of virtual power Eq. (4), and by using the divergence theorem to transform the volume integrals, leads to the two momentum equations: The two partial differential equations are the strong forms with respect to displacement and micromorphic damage with Initial Boundary Value Problem (IBVP).
And the stress-like variables V , Y Y and Y Y Y will be defined by the state relations developed under the framework of thermodynamic irreversible processes.

State relations
In this section, a specific state potential as well as yield function and dissipation potential are constructed for an isothermal isotropic elastoplasticity fully coupled with isotropic local damage, isotropic and kinematic hardening and isotropic micromorphic damage. The damaged state variables are applied to the local state variables according to the assumption of the equivalence of total energy developed in (see [12][13][14][15], [16], [17], [18], [19], [7], [21,22]). And for simplification, we also postulate that only the local damage is coupled with the micromorphic damage. With these assumptions, a quadratic the Helmholtz free energy can be expressed as: 2 1 2 1 2 : : : In Eq. (10) J is the parameter governing the effect of the damage on the isotropic hardening and h f a simple function equal to 1 if the hydrostatic part of the stress is positive, else it is equal to 0.
The energy density release rate given by state relation with respect to damage in Eq. (10) can be decomposed into local part and coupled nonlocal part:

Evolution equations
Since there's no dissipation from micromorphic damage, to account for the damage effects, we limit ourselves to the isotropic plasticity for simplicity. The von Mises yield function and dissipation potentials used in ( [7], [12][13][14], [16], [17], [18], [19]) are used:  where, defines the well-known von Mises equivalent stress and dev V is its deviatoric part of the stress tensor. The parameters D and b reflect the nonlinear property of the kinematic and isotropic hardening respectively and S , s and E characterize the nonlinear evolution of the ductile damage. Clearly, the yield function and plastic potential are identical to the local theory. However, the indirect contribution of the micromorphic state variables is introduced from the state relations. Applying the generalized normality rule to the above defined local yield function and dissipation potential lead to the following evolution equations: is the outside normal to the yield surface in the stress space and O O is the plastic multiplier.

Transformation of the micromorphic balance equations
Considering the generalized balance equation of micromorphic momentum Eq. (5b), it can be easily expressed in strain-like variables space by using micromorphic state relations Eq. (11) and Eq. (12) The Eq. (21) is used as a strong form from which the associated weak form will be deduced in order to solve the fully coupled IBVP.

Internal lengths
The development of an 'advanced' modelling of multiphysic thermomechanical coupling in the framework of the generalized continuum mechanics (micromorphic theory) is proposed in order to introduce the concept of internal lengths that are representative of the materials microstructures while accounting for the various initial and induced anisotropies under large deformations. These internal lengths have to be experimentally determined from accurate measurements of highly localized displacement/strain (or velocity/strain rate) fields by using advanced methods to measure the kinematic fields at the relevant scales. For that, the electronic speckle pattern interferometry (ESPI) method will be used.
We seek to locally measure the displacement and velocity fields in order to access to the local strain/strain rate fields inside the localized zones. An example of a localization pattern obtained with the ESPI technique on a 460 steel sheet is shown in Figure 1 [6]. The tensile direction is horizontal and the interference fringes due to the displacement field are shown in the left where each grey level corresponds to 2.8 nm and the total relative displacement is about 3 µm. The corresponding strain rate field (tensile component) extracted with the geometrical model is shown in the right. The widths B1 and B2 of the localization bands are obtained directly from the experimental measure. These experimental measurements will provide inputs for the advanced micromorphic constitutive equations and mainly to determine the internal lengths related to the targeted micromorphic phenomena [5]. with: x a , the inclination of the band with respect to the transverse direction. x b , the location of the band.
x B , the integral width of the band.
In all cases of the present study, the specimens failed in the same direction, along the band that was oriented almost perpendicular to the rolling direction or along band 1 as illustrated in Fig.1. Thereby to simplify the explications, this band will be called dominant band and the other one disappearing band, and the superscript "dom" and "dis" will be used to indicate them. For the identification of the internal lengths, an effort is under progress. According to the analyticalexperimental study described above, an experimentalnumerical study is evolving in order to determine the values of H H and g H g H and verify the characteristic length for each material and mesh size in total agreement with the experiments and the information extracted by Fig.2 (bandwidths, total average strain, estimated error between numerical and experimental strain rate etc).
Let us note here, that the above study concerns a 460 (X8Cr17) steel specimen.
For this paper and the air bending application presented in section 5, a DP1000 dual phase sheet is used. The internal length for this type of material is determined empirically and fixed at ~0.2mm by following the rule that Mesh size (Ms) d l d .

Numerical aspects
From the numerical point of view, an advanced numerical procedure is presented. The finite element code that is used to implement the micromorphic model is developed in ABAQUS®. The numerical integration in mixed one/four quadrature(s) (or Gauss) point for all the coupled constitutive equations is necessary in order to compute the stress tensor and the overall stress-like internal variables which appear in the weak forms associated to the balance equations. The weak form of the balance equations can be obtained thanks to the weighted residual method [4] in order to get: 2 : : 0 (   ^`^`^`^0 ang et al [1] proposed the construction of a quadrangular element to the field of the assumed strains in a 2D plane strain in order to control the shear and volumetric locking. The isoparametric shape functions for the 4-node quadrilateral can be written in terms of a set of orthogonal base vectors as: