Abdelouahed Chbihi
email: abdelouahed.chbihi@mines-paristech.fr

Pierre-Olivier Bouchard

Marc Bernacki

Daniel Pino Muñoz

P.-O Bouchard

Influence of Lode angle on modelling of void closure in hot metal forming processes

come   L'archive ouverte pluridisciplinaire

with a recent pre-existing mean eld model, named Cicaporo1 model hereafter (Saby et al., 2014b). In comparison with Cicaporo1, the new model, named Cicaporo2 hereafter, uses less constants and is more accurate.

Introduction

Microvoids are frequently detected during large ingots production, which decrease the material quality thereafter. In the industry, hot forming processes are generally used to close these voids and end up with sound products. To insure that the total closure is achieved, the use of 3D X-ray analyses is not possible on such large ingots. Ultrasonic analysis is an alternative solution to detect defects on nal products, but its accuracy is quite low and it does not help in nding solutions to close them.

Therefore, mathematical models are useful to predict void volume evolution and provide an estimation of void closure for given applied plastic strain under given stress states.

To deal with this industrial issue, numerous models were dened in the literature using two principal approaches: macroscopical approach and micro-analytical approach. The rst approach considers a full-eld explicit description of an entire workpiece containing explicit voids. Using adequate numerical features, the description may be accurate for studying void closure in real processes according to a given case study. In micro-analytical analyses, a single void in an innite incompressible matrix is considered. The matrix is usually considered viscoplastic and is dened by a power law type. The evolution of the initial void is studied with respect to various mechanical parameters, for dierent types of material. A large number of studies consider this approach (Zhang et al., 2009), (Saby et al., 2014a) and [START_REF] Feng | A 3-d model for void evolution in viscous materials under large compressive deformation[END_REF]. Voids are typically assumed to be spherical or cylindrical. Most studies do not consider any change of shape during strain and the resulting equations are often inappropriate for large deformation [START_REF] Budiansky | Void growth and collapse in viscous solids[END_REF], [START_REF] Lee | Studies of the growth and collapse of voids in viscous solids[END_REF] and (Zhang et al., 2009).

However, both approaches have signicant limitations in industrial applications.

For instance, macroscopic approaches undermine diculties to take into account, in the same simulations, void-scale characteristics and process-scale parameters. This can lead to the use of heavy meshes, case-dependent description/results and high numerical cost. Regarding micro-analytical approaches, assumptions made for void geometries, loading conditions and material laws are usually too strong for real industrial processes modelling.

It is worth mentionning that the elimination of internal voids is achieved in two principal stages: closure of void volume until the contact between faces is reached and the subsequent bonding of the contact surfaces under suciently high temperature and compressive stresses [START_REF] Park | A study of void crushing in large forgings i: Bonding mechanism and estimation model for bonding eciency[END_REF]. The present work concerns only the void closure stage. [START_REF] Ståhlberg | A study of void closure during plastic deformation[END_REF] presented a micro-analytical criterion for void closure.

They considered square and circular voids in a rigid perfectly plastic material. The proposed model denes the necessary reduction R c to completely close the void to be proportional to the square root of its initial volume V 0 . Experiments were carried out in order to validate the theoretical model and an agreement was found. [START_REF] Tanaka | Factors contributing to crushing of voids during forging[END_REF] suggested that the integral Q of stress triaxiality ratio T x over equivalent strain (Q = T x d ) was a better way to predict void closure. They developed an empirical void closure parameter based on Q and using linear regression to get model constants. Based on the work of Lasne, [START_REF] Saby | Understanding and modeling of void closure mechanisms in hot metal forming processes[END_REF] introduced the Stress Triaxiality Based (STB) model which is implemented in the nite element software FORGE R . The STB model states that the void volume prediction could be computed by the integration of the stress triaxiality T x multiplied by a proportionality coecient, identied by simulations of compression on a spherical void, over the equivalent strain. [START_REF] Budiansky | Void growth and collapse in viscous solids[END_REF] proposed an analytical solution which takes into account the inuence of the stress triaxiality T x and introduces a material parameter m and considered voids as spherical during the whole strain. Based on this solution (Zhang et al., 2009) proposed a semi-empirical solution to predict void closure adding the inuence of the void shape evolution. [START_REF] Lee | Studies of the growth and collapse of voids in viscous solids[END_REF] stated that the non-spherical form of void in the beginning of deformation can signicantly inuence the evolution of void closure. In fact, real voids exhibit complex forms with signicant tortuosity, which makes a sphere far from real voids shape. Therefore, by using ellipsoidal voids, which are closer to reality, the orientation and elongation of voids could also be studied. Recently, [START_REF] Feng | A 3-d model for void evolution in viscous materials under large compressive deformation[END_REF] studied the evolution of dilute ellipsoidal void in power law viscous material under triaxial loading conditions. A semi-analytical expression was dened to predict void closure for large compressive strain. The model accounts for the void's shape to characterize the pores and assumes that the void's orientation is parallel to the principle stress, which is not in agreement with experimental observations. [START_REF] Saby | A geometry-dependent model for void closure in hot metal forming processes[END_REF] used a mesoscale method to develop a model based on a large campaign of simulations on Representative Volume Elements (RVE). This is the only model, to the authors' knowledge, that accounts for the geometry of voids (orientation and form) in addition to mechanical loading parameters (stress triaxiality ratio T x and equivalent strain ). Several congurations of voids were tested in a large range of mechanical loading conditions corresponding to industrial applications. This model (Cicaporo1) presented a great accuracy compared to explicit simulations of real industrial processes (Saby et al., 2014a). Nevertheless, one of its limitations lies on the assumption of axisymmetric loading. Indeed, all void closure prediction models existing in the literature consider this hypothesis, for which loading conditions are dened using only stress triaxiality. However, this parameter does not allow to dene a unique stress state. A detailed review of existing models can be found in (Saby et al., 2014b).

Therefore, the aim of the present work is to propose an enhanced model based on an improved and more accurate description of the stress state in the material. This is achieved by involving the three stress invariants.

Within the context of material damage and failure analysis, many authors have introduced the Lode angle θ, which is a function of the third invariant of deviatoric stress tensor. [START_REF] Danas | Inuence of the lode parameter and the stress triaxiality on the failure of elasto-plastic porous materials[END_REF] modelled the damage of porous elastoplastic materials by two approaches according to the state of the stress triaxiality.

They proved that for triaxial tensile cases, the Lode angle appears to have no eect on the developed model. However, it seems to be an inuencing factor for triaxial compressive cases. [START_REF] Cao | Characterization of ductile damage for a high carbon steel using 3d x-ray micro-tomography and mechanical testsapplication to the identication of a shear modied gtn model[END_REF] studied the calibration of ductile damage models for high strength steels using microtomography. They conrmed that the use of the Lode angle is a good way to represent the stress state for low or negative stress triaxiality conditions. Other works realized by [START_REF] Keshavarz | A new experimental failure model based on triaxiality factor and lode angle for x-100 pipeline steel[END_REF], (Xue and Wierzbicki, 2008) and [START_REF] Mirone | A local viewpoint for evaluating the inuence of stress triaxiality and lode angle on ductile failure and hardening[END_REF] conrmed that the combination of the stress triaxiality ratio T x and Lode angle θ was necessary to dene the stress state accurately.

In this paper, a new void closure model is proposed to improve the prediction accuracy by eliminating the axisymmetric loading hypothesis. The calibration of this model is achieved by using a database of explicit RVE simulations. In the following section, the mesoscale method used to perform explicit simulations is presented.

Several congurations of voids and mechanical loading conditions are considered in order to study their inuence on void closure. The way of dening boundary conditions on the RVE for given stress triaxiality and Lode angle is detailed. The third section concerns the inuence of specic parameters on void closure: void orientation and form as geometry parameters and applied stress triaxiality and Lode angle as mechanical loading parameters. The fourth section is dedicated to the development of the calibration methodology. This methodology is based on an optimization algorithm that enables the identication of the void closure model parameters according to a large number of RVE explicit simulations. This optimization approach also enables to change the void closure analytical function and to proceed to a new calibration in a reduced time. In the fth section, results of the new model are presented and compared to RVE explicit simulations and Cicaporo1 model predictions.

Last section is dedicated to conclusions and perspectives 2 Mesoscale approach and boundary conditions In this section, a short description of the mesoscale approach used to perform explicit simulations is presented. For more details, the reader can refer to [START_REF] Saby | Three-dimensional analysis of real void closure at the meso-scale during hot metal forming processes[END_REF].

Boundary conditions are also detailed in order to take into account the eect of Lode angle θ in addition to stress triaxiality ratio T x .

Representative Volume Element (RVE)

In this work, a mesoscale method is used to perform full-eld simulations. In fact, void dimensions are assumed very small with respect to the workpiece dimensions.

Several works considerd a small void/billet ratio: 0.01 to 0.06 in [START_REF] Kakimoto | Development of forging process design to close internal voids[END_REF]), 0.05 in (Wallerö, 1985) and 0.06 in (Wang et al., 1996). From this condition, two main assumptions can be made: a) the presence of voids has no eect on macroscopic deformation and b) the thermo-mechanical elds that are obtained from the macroscopic scale are considered as locally homogeneous and can be used as remote boundary conditions at RVE scale [START_REF] Saby | Understanding and modeling of void closure mechanisms in hot metal forming processes[END_REF] and [START_REF] Saby | Three-dimensional analysis of real void closure at the meso-scale during hot metal forming processes[END_REF].

The objective is then to perform multiple explicit simulations with stress states representative of industrial conditions (stress triaxiality ratio T x and Lode angle θ). The nite element (FE) simulations were performed using the software Forge R NxT1.0. The void were represented by an ellipsoid of main dimensions (r 1 , r 2 , r 3 ) located at the center of the RVE. A convergence study has been carried out in order to determine the dimensions of the RVE and the optimal mesh size. The RVE dimensions (d 1 , d 2 , d 3 ) are set as four times the ellipsoid initial radii. The maximum mesh size in the RVE is 0.5mm while r 1 , r 2 , r 3 vary between 1mm and 4mm. A mixed velocity-pressure formulation with P + 1 /P 1 tetrahedral elements and local renement are used for an accurate description of void closure.

The viscoplastic material is dened by a Hansel-Spittel law. The equivalent stress is then expressed as:

σ = A( , T ) * ˙ m , (1)
with

A( , T ) = K 0 T m 1 ( + 0 ) n e m 4 /( + 0 ) , (2)
where is the equivalent strain, T is the temperature, m is the strain rate sensitivity, K 0 is the material consistency, m 1 the temperature sensitivity, (n, m 4 ) the hardening and softening coecients, respectively, and 0 is a regularization term.

Boundary conditions

Axisymmetric loading hypothesis

Most of void volume prediction models in the literature [START_REF] Tanaka | Factors contributing to crushing of voids during forging[END_REF]), (Zhang et al., 2009), (Saby et al., 2014a) only take into account the equivalent strain and the stress triaxiality ratio T x , and consider axisymmetric loading conditions. At the RVE scale, the industrial boundary conditions are represented by a velocity applied on the top of the RVE (to obtain a given strain rate) and two identical stresses are applied on the two lateral sides of the RVE (σ 1 = σ 2 ). Three symmetry planes (S 1 , S 2 , S 3 ) are set in the remaining sides (with green, red and blue colors), as illustrated in Figure 1. The velocity V z applied on the top surface is dened so as to apply a given constant strain rate value ˙ zz on the RVE. V z (t) is then dened as -˙ zz D z (t) with D z (t) the current RVE height. It is worth mentioning that the equivalent strain rate has a second order inuence on void closure for our materials and thermomechanical conditions of interest [START_REF] Saby | A geometry-dependent model for void closure in hot metal forming processes[END_REF], therefore, it will not be considered as a parameter of the new model. So a value ˙ zz = 1s -1 is used for all the explicit simulations on RVE.

In order to reproduce real industrial conditions, it is important to address more general loading conditions and then a new approach is proposed in the following.

For such general stress state, it can be shown that the stress triaxiality ratio is not sucient to describe the stress state in a unique way.

Indeed, dierent stress elds can lead to the same state of triaxiality T x and equivalent strain rate ˙ . Figure 2 illustrates two stress states, which can be applied to the RVE as boundary conditions, as presented in Figure 1, and which lead to the values ˙ zz = -1s -1 and T x = -1. Note that the equivalent stress σ is the same in both cases. In the rst case (Figure 2a) the equivalent stress σ is plotted with black solid line while σ 1 and σ 2 are plotted with red triangles and green stars respectively. The evolutions of σ 1 and σ 2 are identical which means that it corresponds to axisymmetric loading. In the second case (Figure 2b), the stresses σ 1 and σ 2 are dierent even if the stress triaxiality ratio is the same (T x = -1), this case does not correspond to axisymmetric loading. Therefore, a third parameter should be added in order to represent the stress state in a unique way.

Non-axisymmetric loading conditions and Lode angle θ

In this work, the axisymmetric hypothesis is removed. The Lode angle θ (0 ≤ θ ≤ π 3 ), which is linked to the third invariant of deviatoric stress tensor (see Eq. 5) by Eq. 7, is added. The stress state is characterized by the 3 eigenvalues of the stress tensor (σ 1 , σ 2 , σ 3 ) and can also be expressed by the rst invariant of the stress tensor (I 1 ) and by the second and third invariants of the deviatoric stress tensor (J 2 and J 3 ), respectively:

σ h = 1 3 (σ 1 + σ 2 + σ 3 ) = I 1 3 , (3)
σ = 1 2 [(σ 1 -σ 2 ) 2 + (σ 2 -σ 3 ) 2 + (σ 3 -σ 1 ) 2 ] = 3J 2 , (4) r = 3 27 2 (σ 1 -σ h )(σ 2 -σ h )(σ 3 -σ h ) = 3 ( 27 2 J 3 ). (5)
Then, the stress triaxiality is obtained by the ratio between the hydrostatic stress σ h and the equivalent stress σ (Eq. 6), while the Lode angle θ is expressed by Eq.

(7):

T x = σ h σ = σ 1 + σ 2 + σ 3 3σ = I 1 3 √ 3J 2 , ( 6
)
cos 3θ = ( r σ ) 3 . (7)
For a better understanding of the meaning of Lode angle θ, reader can consult [START_REF] Cao | Characterization of ductile damage for a high carbon steel using 3d x-ray micro-tomography and mechanical testsapplication to the identication of a shear modied gtn model[END_REF] where more details are presented with a graphical scheme.

A simplied expression of the Lode angle (Eq. 7) was proposed by [START_REF] Lode | Versuche über den einuÿ der mittleren hauptspannung auf das ieÿen der metalle eisen, kupfer und nickel[END_REF] introducing the Lode parameter µ in Eq. ( 9) (-1 ≤ µ ≤ 1). This Lode parameter µ is related to the Lode angle θ by (Eq. 8). It is worth mentionning that the use of the Lode parameter requests that σ 1 ≥ σ 2 ≥ σ 3 , with σ 1 = σ 3 :

µ = √ 3.tan θ - π 6 . (8) µ = 2σ 2 -σ 1 -σ 3 σ 1 -σ 3 . ( 9
)
The stress state is then dened by (T x , σ, µ).

Industrial loading conditions

In order to perform the explicit simulations, the range of mechanical loading conditions observed in industrial applications needs to be determined. Therefore, simulations of industrial processes are performed in order to analyze the variation equivalent strain rate ˙ (red line), stress triaxiality ratio T x (green line) and Lode parameter µ (blue line); which are plotted versus equivalent strain time on gure 3. Rolling and forging processes are modelled and a sensor is located at the ingot's center, which is the place where voids are more likely to exist after casting. This sensor enables to record the mechanical loading conditions during the process. The variation of mechanical loading parameters is plotted in the range of interest, which corresponds to signicant strain rate ˙ zz ≥ 0 or 0.17s ≤ T ime ≤ 0.37s. In this range, T x varies between -0.4 and -0.1 while µ ranges from 0.1 and 1 (Figure 3a). The same procedure is done for the case of forging (Figure 3b). The strain rate is non zero during the whole process. The Lode parameter varies between 0.5 and 1. Therefore, the axisymmetric loading hypothesis (µ = 1) is not consistent with mechanical loading conditions representative of industrial processes.

In fact, as it can be seen from Figure 3, it is dicult to reproduce exactly the stress state variation corresponding to a given forming process. Consequently, the methodology chosen is the following: boundary conditions are iteratively modied to prescribe constant values of strain rate ˙ , stress triaxiality ratio T x and Lode parameter µ during RVE simulations up to complete void closure, in such a way that the loading ranges are covered. This method allows to generate a strain and stress states which could represent a range of mechanical loading at the RVE scale always found in industrial processes as it is presented in Figure 3. In this gure only one case of rolling and forging are presented, however more industrial processes are studied and the same conculsions were reached. The eciency of this method to control the evolution of loading conditions is shown in [START_REF] Saby | Three-dimensional analysis of real void closure at the meso-scale during hot metal forming processes[END_REF].

In the case of axisymmetric loading hypothesis, the applied stresses on the RVE lateral sides are easy to compute. It needs the resolution of a system of two linear equations (Eq. 6) and (Eq. 4). In addition, the equivalent strain rate ˙ is equal to the strain rate in the deformation direction ˙ zz . However, by removing this hypothesis, the resolution becomes more complicated as the equations are in their general forms. The equivalent strain rate ˙ is now expressed by Eq.( 10), so that the system of equations is nonlinear. The idea is to nd the tractions (σ xx , σ yy ) and the velocity V z that should be applied to the RVE in such a way that the desired Triaxiality T x , Lode's Angle µ and strain rate ˙ zz are obtained. Since the material behavior is not linear, then this is done by solving the set of equations ( 4), ( 6) and ( 9). This equations set could be resolved by analytical development, optimization methods or other opproaches. In this work, an optimization algorithm has been developed to compute the boundary conditions to set on the RVE with respect to the constitutive law (Eq. 1), the equivalent strain rate (Eq. 10) and volume conservation (Eq. 11).

˙ = 2 3 ( ˙ 2 1 + ˙ 2 2 + ˙ 2 3 ), (10)
˙ 1 + ˙ 2 + ˙ 3 = 0. (11)
Note that this optimization concerns only the computing of the boundary conditions to set to RVE. Another optimization handling the identication of the mean-eld model's parameters will be discussed later.

Void Morphology

To dene the initial topology of the void during RVE explicit simulations, three orientations (α 1 , α 2 , α 3 ) regarding the canonical basis ( e 1 , e 2 , e 3 )=(x,y,z) and corresponding dimensions (r 1 , r 2 , r 3 ) are dened as shown in Figure 4. However, in order to study the inuence of the void's topology on void closure, dimensionless parameters are prefered. The orientation factors (p 1 , p 2 , p 3 ) are dening the ellipsoid orientation. p i are dened by the vectorial product between the main direction e i of the void and the loading direction e L as expressed in Eq. ( 12). The loading direction e L is computed using a combination of the principal stresses and the three main directions of the canonical basis e i by the expression in Eq. ( 13).

p i = e L ⊗ e i , (12)
e L = σ 1 e 1 + σ 2 e 2 + σ 3 e 3 σ 1 e 1 + σ 2 e 2 + σ 3 e 3 , (13)
Finally, shape factors (γ 1 , γ 2 , γ 3 ) are dened thanks to Eq. ( 14) with V 0 the initial volume of the ellipsoid. Note that the geometrical description of the void is given bi α i and r i , while p i and γ i are used as input for the void closure prediction model.

γ i = 3 √ V 0 r i . ( 14
)
3 Inuence of Lode parameter on void closure

After describing the methodology used to dene the RVE boundary conditions and the range of stress states observed in industrial congurations, the inuence of Lode parameter µ on void closure at RVE scale is studied.

As it was observed in the variation of mechanical loadings in industrial processes (Figure 3), the range of variation of Lode parameter is [0, 1]. However, this study (15)

In Figure 5a, where T x = -1, which is a favorable case for closure [START_REF] Saby | A geometry-dependent model for void closure in hot metal forming processes[END_REF], the axisymmetric loading condition µ = 1 is also favorable for void closure.

The strain dierence is really small (∆ = 0.05) which means almost no inuence of the Lode parameter on void closure. Then, looking at Figure 5b, where T x = -0.66, the gap becomes more important with ∆ = 0.1. The same behavior, concerning the inuence of µ on void closure, is noted again. So that, the gap is still increasing for the case T x = -0.33 (Figure 5c) with ∆ = 0.24. In the last case, T x = 0, an unfavorable case for closure, the gap is even more important (Figure 5d). In fact, ∆ cannot be dened at the threshold value of closure because the void is not completely closed due to its orientation and the mechanical loading condition that represents a generalized shear. However, at the last increment of closure (V /V 0 = 0.4), ∆ = 0.48.

Considering the trend observed in Figure 5, two results could be noted: on the one hand, Lode parameter µ has a signicant inuence on void closure, especially when T x is close to 0. On the other hand, the higher is µ, the faster is the void closure, regardless of the type of geometry neither void geometry. A comparison between dierent void geometry and orientation is given in Figure 6. This bevavior is respected even for a sphere with r = 2 (Figure 6a), an oblate ellipsoid with r i = (2, 2, 1) and α i = (0, 0, 0) (Figure 6b)or a prolate ellipsoidal void with r i =

(2, 1, 1) and α i = (45, 45, 0) (Figure 6c). Consequently, the axisymmetric hypothesis used in the literature represents the most favorable case for the void closure study.

Therefore, it seems quite obvious that the use of this new parameter (µ) must be taken into account in order to improve the accuracy of prediction models used in this eld of research.

Therefore, a large campaign of explicit simulations on RVE needs to be carried out in order to cover the prescribed ranges for all parameters:

T x ∈ [-1, -0.33],
µ ∈ [-0.5, 1], α i ∈ [0, π] and r i ∈ [1, 3]. However, if only 5 values are considered for each parameter, nearly 16000 cases should be carried out, which is too expensive in term of computation time. Thus, a reduced number of cases is used here as following:

• 4 values of µ: -0.5, 0, 0.5, 1

• 3 values of T x : -1, -0.66, -0.33

• 4 shapes of ellipsoids (r 1 , r 2 , r 3 ): (2, 2, 1), (2, 1, 1), (2, 1.5, 1), (2.36, 1.62, 1.34)

• 5 orientations per axis: 0 4 Optimization method and model construction

In the previous section, the inuence of the Lode parameter µ was presented. The new prediction model (Cicaporo2), which includes the Lode parameter µ, is now detailed. In addition, the optimization method required for the model calibration is also introduced. It is worth mentioning that this model takes into account equivalent strain , stress triaxiality ratio T x , Lode parameter µ, orientation factors p i (Eq. ( 12)) and void's shape factors γ i (Eq. ( 14)). This work will not include material parameters eect, as it was shown by [START_REF] Saby | A geometry-dependent model for void closure in hot metal forming processes[END_REF] they have a second order inuence on void closure for our material and thermomechanical conditions of interest.

The choice of the analytical form

To develop a mean eld prediction model, [START_REF] Saby | A geometry-dependent model for void closure in hot metal forming processes[END_REF] and (Zhang et al., 2009) used nonlinear regression to t models constants based on RVE FE simulations. This is the commonly used method due to its simplicity regarding parameters identication. This tting is equivalent to a step by step minimization problem.

Each step handles one parameter while the others are xed. However, if the number of parameters increases, the minimization problem becomes more complex. Additionally, in case of enrichment of one of the studied parameters, the whole tting have to be updated. Thus, this method becomes less ecient and numerically costly.

Therefore, in this work, an optimization method was implemented to simplify the process of model parameters' regression. First, a database of explicit simulations on RVE was built to cover the space of dimensions (mechanical , T x , µ and geometrical p i , γ i ). This database can contain a very large number of simulations, all information coming from the Forge R software outputs (V /V 0 , T x , µ, σ, ˙ , ...) are stored and can be easily enriched if new simulations are necessary. An analytical function is dened by the user to represent the model. This analytical function contains constants (A i ) which must be determined using inverse analysis methods. All types of functions could be used. For example, the following equation, presenting a polynomial form containing all parameters, can be used:

V V 0 = A 0 +(A 1 .T x +A 2 .µ+A 3 .p 1 +A 4 .γ 1 +A 5 .p 2 +A 6 .γ 2 +A 7 .p 3 +A 8 .γ 3 ). +(A 9 .T 2 x + A 10 .µ 2 + A 11 .p 2 1 + A 12 .γ 2 1 + A 13 .p 2 2 + A 14 .γ 2 2 + A 15 .p 2 3 + A 16 .γ 2 3 ). 2 + ... ( 16
)
An oriented object library was coded to manage user-dened functions. In fact, this optimization library ("optimizer") is coupled with the simulations database.

This coupling allows the application of inverse analysis techniques to nd the constants associated with the analytical function dened by the user. These inverse analysis techniques are used to minimize the L 2 -norm of the error between the evolution of the volume of the RVE explicit simulation V V 0 num and predicted volume given by the meaneld model V V 0 mod , as expressed in Eq. ( 17).

f = min V V 0 num - V V 0 mod L 2 . ( 17
)
The advantage of this method lies in its exibility. Indeed, one can easily test dierent models (analytical forms) and thus assess their relevance. 6 & 7). It was also mentionned by [START_REF] Saby | A geometry-dependent model for void closure in hot metal forming processes[END_REF] that this issue is due to the contact between void surfaces and depends on voids tortuosity of the void. In case of equivalent geometries, the importance of this change of slope depends on the shape and orientation of the void (Figure 7). For real industrial voids, this change is more severe and depends a lot on voids tortuosity [START_REF] Saby | Understanding and modeling of void closure mechanisms in hot metal forming processes[END_REF]. Consequently, only the interval V /V 0 ∈ [0.2, 1] will be considered for identication. Thereafter, the model will be compared with explicit simulations until the critical value V /V 0 = 0.2. The optimization algorithm is composed of three main steps as illustrated in Figure 8. The rst step consists in manually providing a mathematical formula that denes the model taking into account all the parameters that inuence the evolution of closure, namely: , T x , µ, γ i and p i . Due to the parabolic form observed in the evolution of void closure for dierent orientations and forms (Figure 7) and for different mechanical loading conditions (Figure 6), a second order polynomial function is proposed, see Eq. ( 18). According to the form of RVE explicit simulation evolution of void closure, a quadratic dependance on equivalent strain is chosen. The analytical function presents also a quadratic dependance to the other parameters as detailed in Eq. ( 18). As it was presented in the general analytical equation (see in Eq. ( 16)), six constants correspond to void's dimension and six constants correspond to orientation properties, which is too weighty regarding the minimization problem. Consequently, in order to reduce the number of model constants, as the Lode parameter µ is also studied here, only two constants are assigned to S and S c which represents a combination of void topology parameters as it is mentionned in Eq. ( 18).

V V 0 = 1 + (A 0 + A 1 T x + A 2 µ + A 3 S) + (A 4 + A 5 T 2 x + A 6 µ 2 + A 7 S c ) 2 , (18)
with: S = 3 i=1 γ i .p i and S c = 3 i=1 (γ i .p i ) 2 . Where A i are the tting constants, p i are the orientation factors and γ i are the shape factors.

The model calibration procedure

The model calibration procedure is illustrated in Figure 8. The BroydenFletcher GoldfarbShanno algorithm (BFGS), which is an iterative method commonly used for solving nonlinear optimization problems, is used within the optimizer to solve the minimization problem [START_REF] Nazareth | A relationship between the bfgs and conjugate gradient algorithms and its implications for new algorithms[END_REF]. Thus, a starting point is required. In order to obtain fast convergence, the initial set of constants should be as close as possible to the optimal set of constants. For this reason, a rst optimization step over a single RVE simulation is carried out. Using only one simulation, randomly chosen, this step enables to nd a starting point very quickly. The result of this initial tting step is plotted in Figure 9. These constants are then used as initial constants for the second stage of optimization.

The second step consists in a loop which contains three optimization iterations (Figure 8). In these three iterations, the various constants of the model are identied sequentially, the red parameters are xed while green ones change. Table 1 presents the dierent sets used for these optimization. Note that for each set, green parameters will vary in the identied range, while the xed ones are randomly chosen in the database. The rst iteration concerns the variation of stress triaxiality ratio T x (set 1). The second iteration is for Lode parameter µ (set 2). The third one is run for the variation of void's morphology and orientation (γ i &p i ) (set 3). Note that in For instance, set 1 of simulations, where T x is varied while µ, p i and γ i are xed, is used to t the constants A 0 , A 1 , A 4 and A 5 while A 2 , A 3 , A 6 and A 7 are considered as xed parameters. These parameters are assigned thanks to the previous step as explained in (Figure . 10). Each iteration step uses a dierent set of simulations (set 1, set 2 or set 3) and stops when the optimizer reaches a very low error, dened by the user, with respect to void volume decrease coming from the database. The set of simulations is chosen so as to give more weight to the modied parameters: T x for set 1, µ for set 2 and (p i & γ i ) for set 3.

On entering the loop, the three optimization iterations will be repeated until the convergence condition is satised. Convergence is reached when the distance between two solutions k and k -1 (as dened in Eq. ( 19)) is lower than a threshold specied by the user.

cvg = A k i -A k-1 i L 2 A k-1 i L 2 . ( 19
)
to show how the model can predict the void volume evolution for dierent void's orientations and shapes and for dierent mechanical loading conditions. A total of ninety cases are compared (Table 2), four of these cases are represented in Figure 11. The void volume evolution is plotted with respect to equivalent strain. The new mean eld model is plotted with red dots while the explicit simulation is plotted with a solid black line. The straight line y = 0.2 represents the critical value of void closure. The value E r (in the legend) represents the mean error (Eq. 20), between the new model ( V V 0 ) mod and explicit RVE simulations ( V V 0 ) num over the whole range of strain. The title above each subgure contains informations about the simulation, r i represents the initial ellipsoid's dimensions, α i represents the initial orientations of the ellipsoid in the RVE. Thus, the void's shape factors γ i and orientations factors p i are calculated as detailed in Eq. ( 12) and Eq. ( 14) using α i and r i , respectively.

E r = 1 N N i V i V 0 num - V i V 0 mod . ( 20
)
Where N is the number of points in the range of strain. As it is shown in Figure 11, all presented cases are well predicted, with small discrepancies. In the interval [0, 0.2], a change of slope is observed in the curve of the explicit simulation. This phenomenon causes more serious discrepancies as the model is dened by a quadratic form, and therefore cannot reproduce this change of slope. However, the results are promising and the prediction is very accurate in the range ∈ [0.2, 1], especially for axisymmetric loading (µ = 1).

Nb of congs

The mean error (Eq. 20) is calculated for all cases and does not exceed 14%.

The maximal error corresponds to a case of orientation which could be considered as an extreme case because the ellipsoid's largest dimension is aligned with the deformation direction which is not desirebale for closure. In addition, the error average is 3.29% while the standard deviation is equal to 2.49% which indicates that error points are clustered closely around the average. As it is seen in Figure 12, the error of 85.56% cases is less than 5.78% which corresponds to the sum of the average E r and the standard deviation σ (E r + σ). In statistics, the Poisson distribution is usually used to study the rare events. In order to apply the Poisson distribution, some conditions have to be satised: (i) the occurence of events is independent, (ii) events occur with constant density λ, and (iii) the event is something that can be counted in integers [START_REF] Grais | Methodes statistiques[END_REF] and (Whittle, 1970). Let us consider the error vector for the ninety simulations used for comparison as a random variable.

The errors will be truncated in order to have integers. Consequently, the error distribution could be approximated by a Poisson distribution as mentionned with red solid line in Figure 12 with λ = 3 which is the average and the variance (truncated) of the error distribution.

Comparison with the Cicaporo1 model

After validating the model with RVE explicit simulations, this paragraph compares the results with the model Cicaporo1 which was dened for axisymmetric loading conditions (Saby et al., 2014b) expressed as:

V V 0 = 1 + B + C 2 , (21)
with

B = 3 i=1 2 j=0 1 k=0 b jk (T x ) k (γ i ) j p i , C = 3 i=1 2 j=0 2 k=0 c jk (T x ) k (γ i ) j p i ,
where the numerical values of six constants b jk and the nine constants c jk are tted using linear regression.

For axisymmetric loading conditions

Dierent congurations are studied for dierent states of triaxiality and for µ = 1, because the Cicaporo1 was calibrated for axisymmetric loading conditions only. In Figure 13, explicit RVE simulations are plotted with solid lines, the Cicaporo2 model and CicaPoro1 model are plotted with dots and stars, respectively.

In Figures 13, both models present good accuracy. The error Er 0.2 (see Eq. 22), which is the relative dierence between model and explicit RVE simulation at the critical value of closure V /V 0 = 0.2, is lower than 5.11%.

with sim and mod are the equivalent strains of the explicit RVE simulation and the predicted model respectively for V V 0 = 0.2.

For general stress states

In order to show the improvement of the new Cicaporo2 model compared to the Cicaporo1 formulation, comparisons are now made with more general stress states (Lode parameter equal to 0, 0.5 and 1). The error Er 0.2 (Eq. 22) is considered for this comparison. In Table 4, values of the error Er 0.2 are presented for both Cicaporo2 and Cicaporo1 models for µ = 1, µ = 0.5 and µ = 0. Note that for Cicaporo1 model, the value mod is the same for the three cases of µ since the model is calibrated for µ = 1. simulations for axisymmetric loading conditions using the expression Eq. ( 22).

In Figure 14, the comparison between four cases is presented. Explicit RVE simulations are plotted with solid lines. Cicaporo1 and Cicaporo2 are plotted with stars and dots respectively. The cases µ = 1, µ = 0.5 and µ = 0 are plotted with red, blue and green, respectively. Note that Cicaporo1 is plotted only for µ = 1, the condition for which it was calibrated. However, it is compared with explicit simulations for other cases to evaluate its discrepancy for other stress states. In the rst case (Figure 14a), the prediction is acceptable for the Cicaporo2 model, especially in the range [0.2, 1] of void closure. The error e 0.2 is reduced from 7.59% for Cicaporo1 to 3.32% for Cicaporo2 model for the case µ = 1. For non-axisymmetric loading conditions, the Cicaporo2 model errors are equal to 8.33% and 10.31% for the cases µ = 0.5 and µ = 0 respectively while the Cicaporo1 model errors are 13.3%

and 22.68% respectively. Even though the model is meant to be used for µ = 1, if it were compared to other loading, the error would be much higher. In table 4, errors of all presented cases are detailed.

For the presented cases, as for the all tested congurations, the case µ = 0 represents the biggest discrepancy, while the case µ = 1 always leads the best prediction.

This could be due to the chosen simulations for the calibration process. Indeed, all explicit RVE simulations used for constants identication when Lode parameter µ is xed, µ is equal to 1. In order to investigate this point and capitalize on the optimization procedure exibility, a database enrichment was made and results are presented in the next section.

Inuence of the set of explicit RVE simulations chosen for the calibration stage

A new set of explicit RVE simulations was used to replace the set 2 (see Table 1) which handle the variation of T x . Lode parameter is now xed to 0. The optimization was run and an improvment of the mean error (Eq. 20) is stated. In fact, the maximal error between the explicit RVE simulation and the prediction of Cicaporo2

Figure 1 :

 1 Figure 1: The boundary conditions set to the RVE in the canonical basis ( e 1 , e 2 , e 3 ).

Figure 2 :

 2 Figure 2: Dierent stress states for the same values ˙ zz = -1s -1 and T x = -1.

Figure 3 :

 3 Figure 3: Variation of stress state and strain rate during industrial processes.

Figure 4 :

 4 Figure 4: Topological denition of ellipsoids with Euler Angles: r i = (2.36, 1.62, 1.34), α i = (45 • , 45 • , 0 • ).

Figure 5 :

 5 Figure 5: Inuence of Lode parameter µ on void closure for 4 dierent stress triaxiality.

Figure 6 :

 6 Figure 6: Inuence of Lode parameter µ on void closure for dierent void geometry.

  Testing a new model just requires the modication of the selected analytical function and a new otimization stage based on the datababase which already exists. It is worth mentioning that the optimization is performed until the critical value V /V 0 = 0.2. The reason for this is the change of slope observed in the last stages of closure (Figures

Figure 7 :

 7 Figure 7: The inuence of void's shape and orientation on void closure, for T x = -0.33, µ = 1.

Figure 8 :

 8 Figure 8: Main steps of the optimization algorithm used to develop the mean eld model.

Figure 9 :

 9 Figure 9: Optimization results to get initial conditions for the loop.

  Dimensions r i orientations α i

Figure 11 :

 11 Figure 11: Comparison between the new model and explicit RVE simulation for dierent void congurations (shape & orientations) and mechanical loadings (T x & µ).

Figure 12 :

 12 Figure 12: Histogram of the probability of the mean error of the new model prediction with the Poisson distribution (λ = 3).

Figure 13 :

 13 Figure 13: Comparison between Cicaporo2 model, explicit RVE simulation and Cicaporo1 model for four dierent void congurations (shape & orientations) and mechanical loadings.

Figure 14 :

 14 Figure 14: Comparison between the Cicaporo2 model and Cicaporo1 model with explicit RVE simulation for dierent void congurations (shape & orientations) and mechanical loadings.

Figure 15 :

 15 Figure 15: Histogram of the error and Poisson law distribution with the new optimization.

  • , 45 • , 90 • , 135 • , 180 • (All similar cases are elimi-

	nated).

A database is used to store all these congurations (695), and it can be easily enhanced if new congurations are needed to enrich the model. More details are presented in section 4.

Table 2 :

 2 Dierent congurations used for comparison with explicit simulations.

Table 3

 3

	presents errors

Table 3 :

 3 Comparison between Cicaporo2 and Cicaporo1 models with explicit RVE

	µ	Congurations	T x	Er 0.2 %	Er 0.2 %
				Cicaporo2 Cicaporo1
			-0.33	5.06	2.12
		r i = (2, 2, 1), α i = (45, 45, 0)
			-0.66	5.06	2.12
	1
			-0.33	2.15	2.89
		r i = (2, 1, 1), α i = (0, 0, 0)
			-0.66	0.82	1.16

Acknowledgement

The authors would like to gratefully thank ArcelorMittal Gandrange & Industeel, Areva TN, Asco Industries, Aubert&Duval, Constellium CRV and Timet Savoie for funding this research.

Varied parameters

Set 1

[-1,0] 1 (2.5, 1.8, 1.1) (45 Variation of µ -0.33 1

(2.5, 1.8, 1.1)

Set 3

(2, 1, 1)

(2, 2, 1) model decreases to 10.76% instead of 13.39% for µ = 1 in set 2. The average of the error remains almost the same with 3.47% instead of 3.29%. This distribution can also be approximated by a Poisson law with λ = 3 (see Figure 15.

Consequetly, one can conrm that the model is a case-dependent as long as the optimization result changes within the change of the sets of explicit RVE simulation used for model calibration. However, the solution given by the second optimization is close to the rst one:

= 0.20, with A µ=1 i and A µ=0 i are the model constants given by the rst (µ = 1 in set 2) and second (µ = 0 in set 2) optimizations, respectively.

Conclusions

The Lode angle was investigated in order to study void closure for more general stress states to avoid the hypothesis of axisymmetric loading always used to study void closure in the literature. Its inuence on void volume evolution was demonstrated and two conclusions were raisen: 1-the axisymmetric loading (µ = 1) represents the most favourable case of closure, 2-the higher is T x (in our interval of interest), the more important is the use of µ in void closure models.

dependance.

Further work will be dedicated to a tensorial form of the model which could predict void shape evolution in addition to void volume closure. This is necessary since void closure is achieved through multi-stages forming operations in the industrial processes where the workpiece is usually rotated between each forming stage. As the shape and orientation of voids are rst order parameters regarding void closure, it is essential to keep track of these parameters during all multistages void closure processes.