Modeling dynamic formability of porous ductile sheets subjected to biaxial stretching:

7 This paper investigates the effect of porous microstructure on the necking formability of ductile sheets subjected 8 to dynamic in-plane stretching. We have developed an original approach in which finite element calculations which 9 include actual void distributions obtained from additively manufactured materials are compared with simulations in 10 which the specimen is modeled with the Gurson-Tvergaard continuum plasticity theory (Gurson, 1977; Tvergaard, 11 1982) which considers porosity as an internal state variable. A key point of this work is that in the calculations 12 performed with the continuum model, the initial void volume fraction is spatially varied in the specimen according 13 to the void distributions included in the simulations with the actual porous microstructure. The finite element 14 computations have been carried out for different loading conditions, with biaxial strain ratios ranging from 0 (plane 15 strain) to 0.75 (biaxial tension) and loading rates varying between 10000 s−1 and 60000 s−1. We have shown that 16 for the specific porous microstructures considered, the necking forming limits obtained with the Gurson-Tvergaard 17 continuum model are in qualitative agreement with the results obtained with the calculations which include the 18 actual void distributions, the quantitative differences for the necking strains being generally less than ≈ 25% (the 19 calculations with actual voids systematically predict greater necking strains). In addition, the spatial distribution 20 of necks formed in the sheets at large strains is very similar for the actual porosity and the homogenized porosity 21 models. The obtained results demonstrate that the voids promote plastic localization, acting as preferential sites 22 for the nucleation of fast growing necks. Moreover, the simulations have provided individualized correlations 23 between void volume fraction, maximum void size and necking formability, and highlighted the influence of the 24 heterogeneity of the spatial distribution of porosity on plastic localization. 25


28
Ductile fracture of metals and alloys commonly occurs by plasticity-mediated growth and coalescence of voids, 29 which are generally nucleated from inclusions, second phase particles and micro-cracks formed during materials 30 manufacturing and processing (Tipper, 1949;Goods and Brown, 1979;Marino et al., 1985;Pineau et al., 2016). 31 For instance, solidification defects, such as porosity and hot cracking, commonly observed in metal-based additive 32 manufacturing processes (e.g. Aboulkhair  The mechanical behavior of the material is modeled as elastic-plastic, with yielding defined by the von Mises 130 criterion (Mises, 1928) in the simulations with explicitly resolved porosity, see Section 2.2, and by the  Tvergaard criterion (Gurson, 1977;Tvergaard, 1982) in the calculations with homogenized porosity, see Section 2.3.

132
The Gurson-Tvergaard criterion considers the material to display a periodic porous microstructure approximated 133 by an array of representative volume elements idealized as hollow spheres with a central hole and with matrix 134 described by the von Mises criterion. 135 2.1. General equations for elastic/plastic materials 136 We assume the additive decomposition of the total rate of deformation tensor d into an elastic part d e and a 137 plastic part d p : where the elastic part of the rate of deformation tensor is related to the rate of the stress by the following linear 139 elastic law: where σ =σ+σΩ−Ωσ is the Green-Naghdi objective derivative of the Cauchy stress tensor used by ABAQUS/Explicit 141 (in order to achieve incremental objectivity of the constitutive equations), where Ω =ṘR T with R being the polar 142 rotation tensor. Note that() denotes differentiation with respect to time. Moreover, L is the tensor of isotropic 143 elastic moduli given by: with 1 and I being the unit second-order tensor and the unit deviatoric fourth-order tensor, respectively. Moreover, 145 G and K are the shear modulus and the bulk modulus.

147
Assuming an associated plastic flow rule, the plastic part of the rate of deformation tensor is: whereλ is the rate of plastic multiplier and Φ is the flow potential. The effective plastic strain rate is obtained from the work conjugacy relation: Assuming adiabatic conditions of deformation (no heat flux) and considering that plastic work is the only 165 source of heat, the evolution of the temperature is given by: where ρ 0 is the initial material density, C p the specific heat and β the Taylor-Quinney coefficient. Specialization of the flow potential for the Gurson-Tvergaard yield criterion (Gurson, 1977;Tvergaard, 1982) 169 leads to: where q 1 and q 2 are material parameters, and f is the porosity (void volume fraction). Recall from Section 2.2 171 thatσ and σ h are the effective von Mises stress and the hydrostatic stress, respectively. The flow strength of the 172 matrix material σ Y is defined by equation (8).

174
The effective plastic strain rate in the matrix material is obtained assuming that the rate of macroscopic plastic 175 work is equal to the rate of effective plastic work in the matrix material: Assuming the incompressibility of the matrix material, the evolution of the void volume fraction is defined as: We consider that the evolution of the porosity is due only to the growth of preexisting voids (void nucleation is 179 neglected). If the initial material porosity f 0 is set to zero, the Gurson-Tvergaard model reduces to the von Mises 180 yield criterion. Moreover, the evolution of temperature is computed using equation (10). The values of the material parameters used in the simulations reported in Section 4 are given in Table 1 Material parameter, Eq. (11) 1 β Taylor-Quinney coefficient, Eq. (10) 0.9

Finite element modeling 188
The problem addressed is that of a plate subjected to in-plane biaxial stretching. In order to reduce the withε 0 xx andε 0 yy being the imposed initial strain rates, see Fig. 1. The loading condition is determined by 197 the ratio χ =ε 0 yy /ε 0 xx , which is varied between 0 (plane strain tension) and 0.75 (near equibiaxial stretching) 198 in the calculations presented in this paper (Z-direction is parallel to the thickness of the plate). Hereinafter, 199ε 0 xx andε 0 yy will be referred to as the imposed initial major and minor strain rate, respectively, and χ as the 200 loading path (ε 0 xx will be also called loading rate). We have investigated imposed initial major strain rates varying these strain rates, the dimensionless parameter customarily used to represent the inertial resistance to motion Zhou et al. (2006)  so that modeling a wider plate does not affect the results. Moreover, since for the considered loading conditions 216 and isotropic material behavior the orientation of the necks is perpendicular to the major loading direction (Stören 217 and Rice, 1975), the specimen width does not affect the necking pattern. The specimen includes the actual porous microstructure of three additively manufactured metallic materials: 220 stainless steel 316L, titanium alloy Ti6Al4V and Inconel 718, see Fig. 2. The mechanical behavior of the matrix 221 Figure 1: Schematic of the geometry and boundary conditions of the problem addressed: a strip of initial length L 0 x = 8 mm and square cross-section of initial thickness L 0 y = L 0 z = 0.5 mm, subjected to biaxial stretching. The Lagrangian Cartesian coordinate system associated to the applied velocity field is denoted by (X, Y, Z). The origin of coordinates is located at the rear bottom left corner of the finite element model. material is modeled with the constitutive framework described in Section 2.2, which has been implemented in  Table 1 Table 2. For each of the three porous microstructures, we have generated three realizations of voids 232 size and positions distribution (R1, R2, R3) which fulfill the Log-normal probability function. The goal is to assess 233 the scatter in the finite element results induced by the random spatial distribution of pores and the distribution 234 of void sizes. Moreover, Table 3 shows the number of voids, the initial void volume fraction and the number of 235 voids per mm 3 included in the strips, for the 3 microstructures investigated, and for the 3 realizations that we have 236 generated per specimen. The differences between the void volume fractions and the void densities obtained from 237 the experimental measurements, Table 2, and included in the finite element models, Table 3 where f 0 i and N i are the initial void volume fraction and the number of voids in the cell i, respectively, V c denotes 268 the volume of the cells, and d j is the diameter of a given void. The cell size is taken to be a sub-multiple of the 269 dimensions of the strip, and it has to be larger than the mean spacing between voids, so that the void volume

295
The presentation of results is split up into four sections: the procedure for comparison of the localization 296 patterns obtained from the calculations with actual porosity and homogenized porosity is described in Section 4.1,

297
Section 4.2 investigates the influence of the porous microstructure on the necking formability of the strip, Section  The results presented in this section correspond to the three realizations generated for SS5Z. This microstructure 302 displays the greater void volume fraction and the larger voids, see Table 2, which helps to enlighten the role of  to determine the formation of necks in plates and rings subjected to dynamic stretching. The saturation value of 325 the major strain outside the necks is referred to as the major necking strain ε c xx , and the minor strain measured at 326 the same material point, and at the same loading time, is called the minor necking strain ε c yy (note that ε c yy = 0 for 327 χ = 0, and that major and minor necking strains will be used in Section 4.4 to construct forming limit diagrams).   and Needleman, 1984), which describes the effects of void coalescence (for a coalescence porosity of 0.12). We have 374 found that both models yield the same results in terms of necking strain (results not shown for the sake of brevity).

375
The reason is that, for the microstructures considered in this paper, which all have relatively low values of initial

402
The results for plane strain stretching (χ = 0) are shown in Fig. 9(a). The necking patterns for both 403 modeling approaches are very similar (as in Fig. 7), with four necks developing at different locations of the strip,  The results for χ = 0.4 are shown in Fig. 9(b). The necking time is less for the homogenized porosity approach,   Fig. 9, as for this microstructure the porosity is significantly 436 smaller (≈ 10 times less than for INC1XY, see Table 2) and the loading rate considered in the calculations is 437 greaterε 0 xx = 60000 s −1 .     shows that the relative difference between both modeling approaches increases with the maximum void diameter: Increasing the imposed loading rate to 60000 s −1 raises these figures up to 15% and 20%, respectively, suggesting ≈ 38% for 10000 s −1 , and 39% for the higher strain rate 60000 s −1 .     crease of the necking strain with loading rate is more for plane strain stretching than for χ = 0.4, namely   respectively. These results make apparent that the delay in the necking formation with the increase of the load-618 ing rate depends to some extent on the microstructure. This behavior, however, does not seem to be caused by the dynamic behavior of the voids (microinertia), as it is proven to be negligible before necking localization, see 620 Appendix D. More likely, it is due the interplay between the microstructure (which seems to act as a material 621 imperfection) and macroscale inertia effects.       Fig. 16. The origin of the Cartesian coordinate system (X , Y , Z ) is located at the center of mass of the 648 void, so that the axes X , Y , Z are parallel to the axes X, Y , Z shown in Fig. 1.

654
For the biaxial stretching case, χ = 0.4, Fig. 19 shows the void for t = 0 µs, 10 µs, 20 µs and 30 µs (before 655 necking time, see Fig. 16(b)). The loading times selected are greater than in the case of plane strain, because the

683
The effect of microinertia has been shown to be small prior localized necking, so that the dynamic behavior The loading rate isε 0 xx = 60000 s −1 and the loading path is χ = 0 (plane strain stretching). The origin of the Cartesian coordinate system (X , Y , Z ) is located at the center of mass of the void, so that the axes X , Y , Z are parallel to the axes X, Y , Z shown in Fig. 1. The initial diameter of the void is 43 µm. The void lies within the necked section 4 indicated in red in Fig. 16(a).      x (see Fig. 1) for calculations with explicitly resolved voids. The loading rate is 60000 s −1 and the loading path χ = 0 (see Fig. 6). The results obtained for three different mesh discretizations are compared: Mesh 1 has been used in the simulations reported in Section 4, Mesh 2 increases the number of elements by a factor of 7, and Mesh 3 includes the same number of elements of Mesh 1, but the element type is C3D10 instead of C3D4. note that the initial porosity is relatively low for the three microstructures investigated in this paper, less than 825 0.05%. It is possible that the stabilizing effect of microinertia is more important for materials with larger initial 826 void volume fraction, for which the onset of fracture may precede necking (Zheng et al., 2020).