Finite element prediction of ductile fracture in sheet metal forming processes

Journal of Materials Processing Technology 177 (2006) 278–281

Finite element prediction of ductile fracture in sheet metal forming processes Pedro Teixeira b , A.D. Santos a,∗ , F.M. Andrade Pires a , J.M.A. C´esar de S´a a a

FEUP-DEMEGI–Fac. Engineering University of Porto, Dep. Mech. Eng. Ind. Management, R. Dr. Roberto Frias s/n, 4200-465 Porto, Portugal b INEGI–Institute of Mechanical Engineering and Industrial Management, R. Barroco 174, 4465-591 Le¸ ca do Balio, Portugal

Abstract In this contribution, Lemaitre’s ductile damage model is coupled with Hill’s orthotropic plasticity criterion. The coupling between damaging and material behaviour is accounted for within the framework of Continuum Damage Mechanics (CDM). The resulting constitutive equations are implemented in the Abaqus/Explicit code, for the prediction of fracture onset in sheet metal forming processes. The damage evolution law takes into account the important effect of micro-crack closure, which dramatically decreases the rate of damage growth under compressive paths. © 2006 Elsevier B.V. All rights reserved. Keywords: Finite elements; Damage; Sheet forming; Anisotropy

1. Introduction The great majority of sheet metal forming processes involves finite inelastic strains that are restricted by geometric instabilities due to necking and strain localization. Therefore, the ability to predict the formability limit is of paramount importance in order to avoid these unwanted phenomena, optimize the process, study the influence of each parameter on the necking occurrence and subsequently improve the press performance. Although forming limit diagrams, associated with finite element simulations, are currently performed by powerful commercial codes, when complex strain paths are involved these predictions often fail to give the right answer. Ductile failure is a result of internal degradation described throughout metallographic observations by the nucleation, growth and coalescence of voids and micro-cracks that may lead to macroscopic collapse. The fact that internal degradation (or damage) accompanies large plastic deformation suggests that these two dissipative processes, although different in nature, influence each other and should, therefore, be coupled at the constitutive level. Consequently the theory of Continuum Damage Mechanics (CDM) may provide a better insight to the physi-

cal phenomenon and play a significant role in the study of the formability and fracture onset in sheet metal forming processes. Here, Lemaitre’s ductile damage evolution law is coupled with Hill’s orthotropic plasticity criterion. The resulting constitutive equations are implemented and assessed for the prediction of fracture onset in sheet metal forming processes. 2. Constitutive modelling The accurate prediction of a macro-crack initiation is still a largely open question. Due to the complex interaction between various phenomena that precede failure, the modelling of the gradual internal deterioration, is not trivial. One of the possible tools to solve this problem is continuum damage mechanics, which represents a local approach to failure. In the framework of CDM the damage variable represents, in an average sense, the effective internal degradation, which is reflected on the local load carrying capacity of the material. The ultimate phase of the damage evolution is detected by a local criterion and corresponds to the failure of the representative volume element (RVE). 2.1. Hill’s orthotropic model

∗

Corresponding author. Tel.: +351 225081742; fax: +351 225081445. E-mail addresses: [email protected] (P. Teixeira), [email protected] (A.D. Santos), [email protected] (F.M. Andrade Pires), [email protected] (J.M.A. C´esar de S´a). URL: http://www.inegi.up.pt, http://www.fe.up.pt. 0924-0136/$ – see front matter © 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.jmatprotec.2006.04.059

When polycrystalline aggregates with random crystal orientation undergo manufacturing processes characterised by extreme straining along preferential directions (rolling is a typical example), the re-orientation of crystals according to the

P. Teixeira et al. / Journal of Materials Processing Technology 177 (2006) 278–281

direction of straining results in a final crystallographic arrangement whose macroscopic behaviour is clearly anisotropic. Hill’s criterion [1] has been introduced as an orthotropic extension of the standard von Mises criterion in order to model the anisotropy often found in formed steel. With σ ij denoting the stress tensor components on an orthonormal basis {e1 , e2 , e3 } whose vectors coincide with the principal axes of plastic orthotropy, the yield function associated with the Hill criterion can be cast in the following form: ¯ = (F (σ22 − σ33 )2 + G(σ33 − σ11 )2 + H(σ11 − σ22 )2 Φ(σ, σ) 1/2

2 2 2 +2Lσ23 + 2Mσ13 + 2Nσ12 )

− σ¯

(1)

where σ¯ is the relative yield stress which defines the size (state of hardening) of the yield surface in the space of stress components. The coefficients F, G, H, L, M, N are constants which are determined by testing the material in different orientations. They are defined as: F=

H ; r90

N=

1 (r0 + r90 )(2r45 + 1) 2 r90 (r0 + 1)

G=

1 ; r0 + 1

H = r0 G;

L = M = 1.5; (2)

where r0 , r45 and r90 are the anisotropic coefficients (Lankford’s r-values). Their values are the anisotropic material data given in terms of ratios of width strain to thickness strain along different directions. 2.2. Lemaitre’s ductile damage model

where σ is the stress tensor for the undamaged material. In addition, the damage variable, D, can assume values between 0 (undamaged state) and 1 (rupture). The evolution law for the internal variables can be derived from a potential of dissipation which is decomposed into plastic, ψp , and damage, ψd , components as
 −Y s+1 r ψ = ψp + ψd = Φ + (4) (1 − D)(s + 1) r for a process accounting for isotropic hardening and isotropic damage, in which r and s are material and temperature dependent properties and Φ and Y are, respectively, the yield function, Eq. (1), and the damage energy release rate, given by −1 [(1 + ν)σ+ : σ+ − νtr σ2 ] 2E(1 − D)2 −

hc [(1 + ν)σ− : σ− − ν−tr σ2 ] 2E(1 − hc D)2

a constant that describes the effect of partially micro-crack/void closure. The tensors σ + and σ − represent, respectively, the tensile and compressive components of σ, defined as: σ+ =

3 

σi ei ⊗ ei ;

3  σ− = − −σi ei ⊗ ei

i=1

(5)

where E and ν denote, respectively, the Young’s modulus and the Poisson’s ratio of the undamaged material and hc (0 ≤ hc ≤ 1) is

(6)

i=1

where σ i is the principal stresses and {e1 , e2 , e3 } is an orthornormal basis of vectors along the principal directions. The symbol   stands for the Macauley bracket,  x, if x ≥ 0 x = (7) 0, if x < 0 It is possible to show [2,4], that the damage energy release rate is a function of the square of the triaxiality ratio, which may be defined, at each point, as the ratio between the hydrostatic stress and the equivalent stress. This means that, if in the dissipation potential no distinction is made for the compressive or tensile hydrostatic states, the damage growth may be the same in both cases. One important feature of the above model is the fact that, as a result of the definition of the damage energy rate given in Eq. (5) there is a clear distinction between rates of damage growth observed for compressive and tensile stress states. The effect of crack closure is incorporated into the damage evolution equation by making use of the above tensile/compressive split of σ and the crack closure parameter hc . By the hypothesis of generalised normality, the plastic flow equation is obtained as ∂ψ ∂σ and the evolution law for the internal variables as ∂ψ α˙ = −γ˙ = γ˙ ∂σ¯
 1 −Y s ∂ψ ˙ = γ˙ D = −γ˙ ∂Y 1−D r ε˙ p = γ˙

By assuming homogenous distribution of micro-voids and the hypothesis of strain equivalence [2,3], the effective stress ˜ can be represented as: tensor, σ, σ σ˜ = (3) 1−D

Y =

279

(8)

(9) (10)

where α is the internal variable associated with isotropic hardening and γ˙ is the plastic consistency parameter, which subjected to the loading/unloading condition expressed by the complementarity law: γ˙ ≥ 0,

Φ ≤ 0,

˙ = 0. γΦ

(11)

2.3. Numerical integration scheme In explicit transient dynamic finite element schemes, the computational cost of an analysis depends crucially on the CPU time spent in Gauss point calculations. Therefore, the efficiency of numerical integration schemes for integration of the constitutive equations has a direct impact on the overall efficiency of the finite element framework. In this context, the use of more complex constitutive models can potentially incur a dramatic increase in analysis time. In the case of the present fully coupled damage model, this issue becomes particularly relevant. A possible alternative to avoid the problem and have a fast integration algorithm that estimates damage evolution with crack

280

P. Teixeira et al. / Journal of Materials Processing Technology 177 (2006) 278–281

closure effects is to decouple the above damage evolution, Eq. (10), from the plasticity equations. In this case, the decoupled state update algorithm follows the steps: (a) Given the incremental strain, ε, compute the updated stress, σ n+1 , and the incremental plastic multiplier γ, by using the conventional Hill’s model (possibly with kinematic hardening). (b) With the updated stress and incremental multiplier at hand, compute the new damage state, Dn+1 , by solving (for Dn+1 ) the Backward Euler discrete version of Eq. (10) which reads: Dn+1 = Dn +

γ 1 − Dn+1




−Yn+1 r

s (12)

Obviously, the damage variable does not affect the plasticity evolution law in the above-decoupled implementation. Damage growth does not cause plastic softening. Thus, if the decoupled approach is adopted, the hardening material properties of the plasticity model have to be calibrated accordingly. 3. Numerical example In this section, we present a numerical example aiming to illustrate the essential features of the model. This example considers the simulation of a sheet metal forming process. The problem consists of a circular sheet subjected to internal pressure. The sheet lies on a rigid cylindrical die and its edge is assumed clamped. The geometry and the finite element mesh are shown in Fig. 1. Due to the symmetry of the problem, only one quarter of the domain is considered in the finite element simulation.

Table 1 Process parameters Blank element type Integration points Number of elements

Solid (eight node) 1 588

Pressure evolution

Ramped (maximum value = 7.2 MPa)

Table 2 Material properties Young modulus, E (GPa) Poisson coefficient, ν Initial yield stress, σY0 (MPa)

210 0.3 175

Yield stress, σ Y (MPa) Damage data (exponent), s Damage data (denominator), r (MPa) Crack closure parameter, hc

488.35(0.001 + ε¯ )0.2427 1.0 1.25 0.2

A mesh with a layer of 588 eight-noded solid elements (C3D8R from the ABAQUS library) is used in the discretization of the metal sheet. The surface of the dies is discretized by three-noded rigid triangular elements (R3D3 from the ABAQUS library). The main process parameters are listed in Table 1. The rolling direction is assumed to be along one-axis (Fig. 1). The material properties employed in the numerical simulation are listed in Table 2. Three simulations were carried out with different anisotropy conditions: (a) the material is regarded as isotropic (r0 = 1.0, r45 = 1.0, r90 = 1.0); (b) with a higher r0 value (r0 = 3.0, r45 = 1.0, r90 = 1.0); (c) with a higher r90 value (r0 = 1.0, r45 = 1.0, r90 = 3.0). The final damage distribution for these three simulations is depicted in Fig. 2.

Fig. 1. Tool geometry.

Fig. 2. Damage contours at intermediate stages: Dmax = 0.11 in (a and b) and Dmax = 0.18 in (c).

P. Teixeira et al. / Journal of Materials Processing Technology 177 (2006) 278–281

281

for the centre of the specimen, and the two outside zones of irregular mesh. This is a typical problem associated with local softening models, as the one used here, that may suffer from mesh orientation dependence. This may be avoided or strongly reduced by adopting non-local or gradient models [4]. In this initial study, only the standard local model was implemented. The possibility of including, for example, gradient models in the ABAQUS environment is being studied at the moment. Finally, it is noted that in this particular example the crack closure effect does not play a significant role as the stress history is dominated by tensile stresses. 4. Conclusions

Fig. 3. Damage evolution (isotropic case).

It can be seen in Fig. 2 that, qualitatively, the damage contours obtained by the numerical analysis are in agreement with the expected. For the isotropic case, the damage contours assume a form of concentric circles around the centre of the specimen. This is also illustrated in Fig. 3 in which the damage evolution versus displacement, for the isotropic case and for some selected nodes, is represented. In the anisotropic cases, a variation of the anisotropy coefficients leads to different damage evolution for directions parallel and transverse to rolling direction; therefore, the damage contours assume a form of concentric ellipses around the centre of the specimen. This is caused by different levels of plastic strain in each direction. Nevertheless, the contour regularity, in both cases, is perturbed near the transition region between the regular mesh, used

A constitutive model for the prediction of damage on sheet metal forming processes was presented. The example used shows that the model was able to capture the damage evolution and predict the location of possible material failure. More complex deformation histories are needed to validate the indications given by the numerical test performed. Further improvements will be considered in future work, namely the possibility of inclusion of non-local damage models in the ABAQUS environment. References [1] R. Hill, The Mathematical Theory of Plasticity, Oxford University Press, London, 1950. [2] J. Lemaitre, A continuous damage mechanics model for ductile fracture, J. Engng. Mat. Tech., Trans. ASME 107 (1985) 83–89. [3] J. Lemaitre, J.L. Chaboche, Mechanics of Solid Materials, Cambridge University Press, 1990. [4] J.M.A. C´esar de S´a, P.M.A. Areias, C. Zheng, Damage modelling in metal forming problems using an implicit non-local gradient model, Comput. Methods Appl. Mech. Eng., in press.