A fracture criterion for fracture simulation of ductile metals based on micro-mechanisms

A fracture criterion for fracture simulation of ductile metals based on micro-mechanisms

Accepted Manuscript A fracture criterion for fracture simulation of ductile metals based on micromechanisms Shen Yan, Xianzhong Zhao PII: DOI: Referen...

2MB Sizes 0 Downloads 108 Views

Accepted Manuscript A fracture criterion for fracture simulation of ductile metals based on micromechanisms Shen Yan, Xianzhong Zhao PII: DOI: Reference:

S0167-8442(17)30399-3 https://doi.org/10.1016/j.tafmec.2018.02.005 TAFMEC 1998

To appear in:

Theoretical and Applied Fracture Mechanics

Received Date: Revised Date: Accepted Date:

27 August 2017 6 January 2018 7 February 2018

Please cite this article as: S. Yan, X. Zhao, A fracture criterion for fracture simulation of ductile metals based on micro-mechanisms, Theoretical and Applied Fracture Mechanics (2018), doi: https://doi.org/10.1016/j.tafmec. 2018.02.005

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

A fracture criterion for fracture simulation of ductile metals based on micro-mechanisms a

Shen Yan , Xianzhong Zhao

a,b

a. Department of Structural Engineering, College of Civil Engineering, Tongji University, Shanghai, China

b. State Key Laboratory of Disaster Reduction in Civil Engineering, Tongji University, Shanghai, China Abstract: A correct representation of the fracture limit for material subject to various stress states is essential for the accurate prediction of fracture behaviour. In this paper, a comprehensive fracture criterion for ductile metals is developed. The proposed criterion is composed of two competing and correlated sub-criterions corresponding to the stress states where the shear fracture and tensile fracture modes can be captured. Especially, a new shear fracture criterion is developed based on microscopic studies on the void nucleation and evolution behavior. According to this criterion, fracture in low positive and negative triaxiality regime is triaxiality sensitive and Lode parameter dependent. The tension fracture criterion is derived from the Rice and Tracey void growth equation, thus it predicts the ductility loss due to stress triaxiality in an exponential function. The overall fracture criterion has five free parameters including the strain hardening exponent, and is proved to be capable of predicting fracture correctly in almost all types of specimen tests in the literatures. A calibration strategy by using four specimen tests is proposed, and specimen tests are performed to calibrate and validate the fracture criterion for Q345 steel. Numerical simulation of the tests are performed, and the experimental and numerical results agree well, which validates the proposed fracture criterion and shows its application in fracture s modelling of large-scale engineering problems.

Keywords: Fracture criterion; Shear fracture; Micro-mechanisms; Experiments; Numerical simulations.

1

1. Introduction

The three-dimensional stress state in isotropic materials can be uniquely identified by three principal stresses (σ1, σ2, σ3) or three coordinate-free invariants of the stress tensor. In plasticity theories, a combination of {I1, J2, θ} is frequently used; I1 is the first invariant of the stress tensor, J2 is the second invariant of the deviatoric stress tensor, and θ is the Lode angle related to the normalized third deviatoric stress invariant. Classic metal plasticity theory assumes that hydrostatic stress σh, which is one third of I1, has no or negligible effect on the plastic yielding. The well-known J-2 plasticity proposed by von Mises is among those simplest criteria, and suggests that the plastic yielding is driven by J2 only. When the research on metal materials is extended to study the fracture behavior of materials, however, the significant influence from σh came into play, which has been ascertained by numerous experimental evidences [1-4]. It was observed that material ductility decreases with increasing stress triaxiality T, which is the normalized hydrostatic stress and is defined as σh/σeq where σeq is the equivalent von Mises stress. However, ductility is hardly affected by θ or its equivalent parameters such as the Lode angle parameter or the Lode parameter. These experimental findings were annotated by many theoretical or empirical fracture criteria and models [5-7]. Many of them follow the void growth equations developed by Rice and Tracey [8], and predict the ductility reduction due to the stress triaxiality by using an exponential function [1, 2, 9]. Taking the influence of σh into account, these criteria and models enable much better prediction and simulation of metal fracture behavior than the simplistic uniaxial critical strain criterion or the maximum shear stress criterion. According to the above T-based criteria and models, one may draw a conclusion that the fracture strain in simple shear is greater than that in simple tension due to the absence of hydrostatic tension. Nevertheless, recent experiments [3, 10-13] showed that sometimes the simple shear (or torsion) fracture strain can be lower than the simple tension fracture strain. This suggests that these criteria and models may no longer be applicable in low positive or negative 2

triaxiality regime where a fracture mode of shear other than hydrostatic tension dominates. Wierzbicki and Xue [14] and Bai and Wierzbicki [15] assumed the loss of ductility under shear condition comes from the change of θ, which varies from 0 in simple tension to π/6 in simple shear, and modified the high-triaxiality exponential fracture models by introducing a correction term (or correction terms) of θ or its equivalent parameters. The effect of θ can be taken away from their modified models by giving a constant value. If θ is taken as π/6, the models shall be reduced to T-based models derived from the Rice–Tracey void growth equation [8], implying that fracture under simple shear and plane strain tension (both stress states have a Lode angle of θ=π/6) can be characterized by the same void growth rule. This somehow deviates from the fracture mechanism at microscopic scale because limited void growth is observed in low triaxiality fracture [12]. Hence, these modified models consider the fracture behavior in an overall viewpoint without distinguishing the shear fracture from tension fracture. Similar views were also shared by Wilkins et al. [16] and Xue [17, 18]. By assuming separable dependences of the fracture locus on T and θ (or their equivalent parameters), they investigated the dependences in a phenomenal fashion and developed their models in continuum approaches. Research at the BMW R&D Center and MATFEM Co. in Munich [19] proposed a fracture model referred to as CrachFEM from a very different perspective. CrachFEM distinguishes shear fracture mode from ductile fracture mode, and the one with the lowest calculated fracture strain determines the fracture limit. Although acknowledging the existence of two different fracture modes and mechanisms, CrachFEM is phenomenal in nature because the formulas used to calculate both the shear and the ductile fracture strains are phenomenological. Here, the authors do not imply micro-mechanism-based models should have evident advantages over phenomenological models. Indeed, some phenomenological models have excellent performance in predicting and simulating material fracture in engineering structures. However, models with solid microscopic mechanism bases may be more attractive and promising because the macroscopic fracture is resulted from the damage at microscopic scale. At microscopic scale, ductile fractures develop due to the nucleation, growth and coalescence of microscopic voids 3

that initiate at inclusions and second phase particles [20]. In hydrostatic tension, it is already well acknowledged that the fracture properties are controlled by the growth of voids: the growing voids reach a critical size, relative to their spacing, and a local plastic instability develops between voids to result in failure. With respect to shear fracture, experimental observations showed shallow and elongated dimples at fracture surfaces, suggesting void elongation with limited growth and shear localization between voids [12]. Hence, most of the studies relate shear fracture to the localized deformation in a thin void-containing band [16, 19, 21]. Nonetheless, it is not well addressed what is the microscopic mechanism that causes the failure of the shear band and finally the fracture. This probably explains why the current dominant literatures on shear fracture models are limited to macroscopic phenomenological models. This present study is intended to develop a comprehensive fracture criterion for ductile metals. In order to capture the metal fracture in both shear and tension modes, the criterion consists two competing while multiplying branches. The key difference to the CrachFEM model is this proposed model computes the strain for fracture based on microscopic mechanisms. It should be noted that high-triaxiality tension fracture is already well understood with micro damaged based fracture criteria, e.g., Rice–Tracey void growth equation. The focus of the study is the shear fracture criterion. Both stress triaxiality and Lode angle play significant roles in shear damage and fracture and the effects of both parameters are taken into account. Here, the Lode parameter L instead of the Lode angle θ or the normalized third deviatoric stress invariantξis used to distinguish different shear stress states in three dimensions. Thus, by combining the tension and shear branches, this proposed criterion covers the whole stress state space, constituting a significant difference from the above mentioned micromechanics-based fracture models which can only account for the tension fracture mode. The proposed criterion is validated by a large number of test results reported in the literatures. Specimen tests of Q345 steel, are performed to calibrate the free parameters for this material. Numerical simulation of the tests are also carried out to test the performance of this proposed fracture criterion in fracture simulation in engineering problems. 4

2. Characterization of the stress state

Let σij be the Cauchy stress tensor and σ1, σ2 and σ3 be its principal stress values, thus the three-dimensional stress state in isotropic materials can be uniquely identified by three principal stresses, as shown in Fig. 1. In the space of principal stresses, a circular cylindrical coordinate system can be defined by setting the z axis passing through the origin and having equal angles with respect to the principal stress axes. An arbitrary stress state at a point P can be identified by z, r and θ (θ is normally referred as the Lode angle). Each point on the z axis corresponds to a hydrostatic stress state, and the plane with z=0 is called the π plane which corresponds to a deviatoric stress state. All three components (z, r, θ) can be related to the invariant of the stress tensor or the stress deviator tensor, as expressed in Eq. 1 ~ Eq. 3.

I1

z

(1)

3

r  2J2

cos  3  

(2)

27 2

J3

 3J 2 

(3)

3 2

where I1 is the first invariant of the stress tensor σij, and J2 and J3 are the second and the third invariant of the stress deviator tensor sij. They can be calculated by

I1  tr  

(4)

J2 

1 S  : S  2

(5)

J3 

1 S  S  : S  3

(6)

and the hydrostatic stress σh and the equivalent von Mises stress σeq can be calculated by

h 

I1 3

(7)

 eq  3J 2

(8)

5

Stress triaxiality is defined as the ratio of hydrostatic stress to the equivalent stress T

h  eq

(9)

and thus it is clear that, as shown in Fig. 1, there exist an infinite number of stress states for a given stress triaxiality, each corresponding to a point on the surface of a cone with coordinate axis of z as the axis. Hence, the Lode angle is important in distinguishing those different stress states. In this study, Lode parameter L defined in Eq. 10 is adopted to describe the Lode angle dependence; L ranges from –1 to 1 and relates to θ uniquely through L= L

2 2   1   3 1   3

tan(θ–π/6). (10)

3. Shear fracture criterion with pressure and Lode dependence

3.1 Assumptions Metal fractures due to nucleation, growth and coalescence of microscopic voids, whereas their importance is different. In materials where the second phase particles and inclusion are well bonded to the matrix, void nucleation is often the critical step, and fracture occurs soon after the void forms. When void nucleation occurs with little difficulty, the fracture properties are controlled by the growth and coalescence of voids; the growing voids reach a critical size, and a local plastic instability developed between voids results in failure [20]. In fact, the nucleation difficulty decreases as the hydrostatic stress increases and thus can also be significantly different even in the same metal material [22, 23]. Void nucleates easily under high stress triaxiality, which constitutes part of the reason why void growth and coalescence controls tension fracture. In low positive and negative triaxiality regime where shear fracture mode dominates, void nucleation plays clearly a more important role. The importance of void nucleation can be demonstrated by the fracture process of the shear band formed at the outer ring in a uniaxial tensile test along 45° from the tensile axis; as previously introduced, the formation and fracture of the shear band is believed to be the major mechanism of shear failure. In a uniaxial tensile test, the penny-shaped flaw

6

resulted from the growth and further coalescence of void due to triaxial stress state at neck produces shear deformation bands at 45° from the tensile axis. This concentration of strain provides sufficient plasticity to nucleate voids in the smaller more numerous particles. Since the small particles are closely spaced, an instability occurs soon after these smaller voids form, resulting in total fracture of the specimen and the cup and cone appearance of the matching surfaces [20]. Ponte Castañeda and Zaidman [24] and Danas and Ponte Castañeda [25] showed that shear band formation is resulted from the shape change and abrupt collapse of voids along a compressive direction (with small, but finite porosity) in shear. The collapse of voids can dramatically soften the response of the porous material, and in turn lead to the localized deformation inside a thin band leading to unloading outside the band. The shallow small elongated shear dimples observed at shear fracture surfaces [12] are proof of the void-collapse-induced shear band. Connecting the dots, we believed that the shear fracture mechanism at microscopic scale can be described as follows. After the nucleation of voids initiating at inclusions and second phase particles, the voids experience shape change other than dramatic volume expansion, and the abrupt collapse of voids leads to the formation of a thin localized void-containing deformation band. Inside the shear band, voids nucleate at smaller more numerous particles driven by remarkable plasticity, causing the failure of the shear band and the final fracture of the metal material. Therefore, a shear fracture criterion is proposed by hypothesizing that the strain to fracture εfs at low positive and negative triaxiality regime can be calculated as

 fs  c1 nc  c2 cl

(11)

where εnc is the equivalent plastic strain that causes the void nucleation, and εcl is the critical strain that leads to the collapse of voids and the formation of the localized shear band; c1 and c2 are corresponding coefficients. 3.2 Void nucleation strain Most of the expressions proposed for void nucleation prediction use the same fundamental assumption. It is assumed 7

that void formation takes place at a critical local stress; the local stress is resulted from the internal stresses built up as a result of the inhomogeneity in deformation between the inclusions and the matrix [26]. A number of models for estimating void nucleation stress have been proposed. For the case of small particles (<1μm), models incorporating dislocation–particle interaction are provided by Brown and Stobbs [27] and Goods and Brown [26]; while for larger (≈10μm) and widely spaced (>50μm) inclusions – the metal materials in structural applications fall into this category, models based on continuum theory is certainly more appropriate [23]. The most widely used continuum model is due to Argon, et al. [22]. They analyzed the solution of problems involving non deformable spherical inclusions in a perfectly plastic or elastic matrix, and found that the critical local stress σc that causes initiation of voids is approximately equal to the sum of the hydrostatic stress and the effective (von Mises) stress

 c   h   eq

(12)

The Beremin research group in France [23] applied the Argon et al. criterion to experimental data for a carbon manganese steel, but argued that the a semi-empirical relationship (Eq. 13) gave better predictions of void nucleation at MnS inclusion that were elongated in the rolling direction

 c   h  k   eq   ys 

(13)

where σys is the yield strength at given temperature while k is a function of the inclusion shape and the loading direction. The Beremin model provides more general applicability than the Argon model by introducing a free parameter k, however, the basic concept of these two models is identical, namely, the nucleation strain decreases as the hydrostatic stress increases. In this study, the Argon model which is the most widely used is taken as an approximate estimation for the void nucleation at second phase particles and inclusions, thus the equivalent stress to void nucleation is found to be dependent on the stress triaxiality, as expressed by

 eq =

c

(14)

T 1

Hollomon's equation, a power law relationship between the stress σ and the amount of plastic strain εp (i.e., σ=Kεpn 8

where K is the strength index and n the strain hardening exponent), can be applied to describe the plasticity of most metals. The equivalent plastic strain to void nucleation εnc can be calculated as

 nc   eq =

c 1

T  1 n

(15)

3.3 Void collapse strain The initial void volume fraction in a ductile metal is usually very small, such that the evolution of each void after its nucleation is almost independent [28]. If the nucleated void is initially spherical and periodically arrayed, a macroscopic material point can be represented by a unit cell that contains a discrete spherical cavity placed at its center and is exposed to stresses acting parallel to the principal geometrical axes, as schematically illustrated in Fig. 2. Hence, the microstructural damage of a material point can be investigated by studying this unit cell; by employing finite element (FE) numerical approaches, this methodology is generally referred to as the “FE unit cell based micromechanical study”. Several such studies have been conducted to investigate the void evolution rules under different stress states; most of these studies focused on hydrostatic tension [28, 29], and showed an increasing growth speed of void as the stress triaxiality increases. The unit cell study performed in this paper, however, looks at the low positive and slightly negative triaxiality regime to develop a criterion by which the void collapse is predicted. Another typical application of the cell model studies is to capture mechanical features of the interaction between a series of cell elements forming the fracture process zone of the crack tip and the interaction between the fracture process and the plastic dissipation in the background material [30-32]. Due to the three-fold symmetry of a cubic cell, only 1/8 of a cell is modeled in this study. The model is developed in the general purpose code ABAQUS/Standard, and 8-node linear brick finite elements with reduced integration are used in all computations. As shown in Fig. 3, the length of the 1/8 cell (Li0) is 59 μm and the radius of the void (R0) is 20 μm. Hence, the initial void volume fraction (f0) is 2%, which is a relatively high value such that the effect of void evolution will be more pronounced. The matrix material of the cell model is modelled by isotropic von Mises plasticity with a 9

strain hardening stress–strain relationship. The strength index K equals to 950 and the strain hardening exponent n equals to 0.22. The Young’s modulus and Poisson’s ratio are E = 205000 MPa and ν = 0.3, respectively. The boundary conditions of the model are as follows. The inner planes of symmetry of the cell model are fixed against normal displacements, i.e., boundary condition of Ui =0 is imposed for plane xi=0. All outer planes of the cell model have to behave only as rigid movable planes in coordinate directions during any type of loading, simulating the lattice periodicity and the point wise symmetry with respect to every void as a center. This feature is realized by imposing a uniform normal displacement Ui for the face at xi= Li0 through an equation-type constraint. The cell model is loaded with macroscopic stress σi through pressure loads at the outer planes. Since both the stress triaxiality and the Lode parameter have notable influences on the void evolution, a broad spectrum of T and L is investigated with the values of T range from –0.3 to 0.3 with an interval of 0.1, and L range from –1 to 1 with an interval of 0.215. Different stress states are realized by imposing a consistent relation of the applied stress σi. According to Eq. 4, 5, 7, 8 and 9, the stress triaxiality is a function of the three principle stresses:

T=

2  1   2   3  3  1   2    2   3    3   1  2

2

(16)

2

Solve σ2 from Eq. 10, and get:

2 

1 1 1  L 1  1  L   3 2 2

(17)

Substitute σ2 back into Eq. 16, one get:  3  L 2  9T 2  3  L2   32  2  9  L2   18T 2  3  L2  1 3   3  L 2  9T 2 3  L2  12  0      

(18)

Solve Eq. 17, the third principle stress σ3 is a function of the first principle stress σ1 under any given T and L:

3 

18T 3  L2  L2  27T 2  9T 2 L2  9

3  L 

2

 9T 2  3  L2 

1

(19)

and σ2 can be calculated from Eq. 17 and therefore is also a function of σ1. Fig. 4 shows the deformation of one eighth of the initially spherical void in low positive and slightly negative 10

triaxiality regime which is what this study focuses. Voids do not grow as they were observed in prior-mentioned high-triaxiality unit cell studies; instead, the voids have a notable decrease in size and an abrupt collapse in most cases. This agrees well with those experimental observations of shallow and elongated dimples at shear fracture surfaces [12]. Moreover, both T and L play significant roles. The stress triaxiality controls the pace of void size decrease: as the stress triaxiality decreases from 0.3 to –0.3, the rate of void shrinkage remarkably accelerates. The Lode parameter determines the void shape: for the negative value of L = –1 the void becomes prolate elliptical with one elongated principal direction, whereas for the positive value of L = 1 oblate elliptical void with two elongated principal directions is observed. In order to quantitatively determine the void shrinkage behavior against strain development, the void volume fraction f and the macroscopic equivalent strain εeq are calculated from Eq. 20 to Eq. 25 [29]. V   L10  U1  L20  U 2  L30  U3  Ve  V0 1  f 0 

 L  U k  i  ln  i 0  L  U k 1 i  i0

k 

 eq

E

k 

h

(21)

V0 V 1  f0   e V V

f  1

 i

3 1  2v 

(20)



(22)

   

1   1 k    2 k  3

(23)



2



1   2 k    3 k  3



2





1

2 2 1   3 k   1 k   3 

(24)

k

 eq    eq i 

(25)

i 0

where, V0 and V is the initial and current volume of this 1/8 cell cubic, f0 and f the initial and current void volume fraction, Ui and Ui (k) (i =1, 2, 3) the xi -directional displacement of the outer plane with the location xi = Li0, at the current state and at the end of the k th incremental step, respectively. ΔVe is an approximate correction term for the elastic dilatational change in cell volume due to the imposed hydrostatic stress σh; Δεi (k) (i =1, 2, 3) and Δεeq(k) is the strain 11

increment in the xi direction and the equivalent strain increment in the k th incremental step. Fig. 5 presents the evolution of normalized void volume fraction V/V0 under different stress states, from which the effects of both stress triaxiality and Lode parameter can be clearly observed. The speed of void shrinkage accelerates with the decrease of stress triaxiality and the increase of Lode parameter. The shrinkage of void softens the response of porous material and can in turn lead to the formation of localized shear band outside which unloading occurs. Hence, if the critical strain εcl is taken as the strain when the void vanished completely, it can be predicted by Eq. 26 and is clearly triaxiality sensitive and Lode parameter dependent.

 cl  p1  p2 e p T  p L 3

(26)

4

where p1 to p4 are free parameters. In this unit cell study, a combination of p1=0.53, p2=0.48, p3=2.37, p4=1.50 gives a good prediction of the void collapse strain, as shown in Fig. 6. 3.4 Fracture criterion Substitute Eq. 15 and Eq. 26 back to Eq. 11, the shear fracture strain εfs can be calculated as

 fs  c1 nc  c2 cl  c1

c

T  1

1 n





 c2 p1  p2 e p3T  p4 L  c1 c c2 p2 

p1  e p3T  p4 L p2 1

T  1 n

(27)

Setting p1/ p2=q1, c1εcc2 p2=1/q2, p3=q3, p4=q4, Eq. 27 becomes

 fs 

q1  e q3T  q4 L

(28)

1

q2 T  1 n

Free parameters q1 to q4 can be calibrated with sufficient specimen test data; specifically, at least four types of specimen tests representing different stress states are required to obtain four sets of {εfs, T, L} relations. However, we prefer not to conduct so many tests for ease of operation. Shear specimen tests are much more difficult to perform than conventional tensile specimen tests. In order to convert the tension loading applied at specimen ends into desired shear loading within the gauge section, shear specimens (such as the butterfly specimen designed by Bao and Wierzbicki [3], the hat specimen by El-Magd et al. [10] and the modified Lindholm specimen by Graham et al. [33]) all have relatively 12

more complex shapes. Some other studies employed torsional testing system to generate shear fracture [13, 33], whereas the commonly adopted notched tube specimen is also difficult to prepare. Moreover, the fracture strains εfs and the corresponding stress states {T, L} can hardly be obtained from the test data directly. Therefore, numerical simulations of the specimen tests are often required to identify the fracture initiation location where all components of the stress and strain tensors are extracted to calculate the local strain and to determine the stress states [3]. The numerical modelling process is complicated and time-consuming due to the significant nonlinearity of the localized plastic deformation in the FE analysis. Therefore, Eq. 28 is desired to be simplified based on existing experimental evidences and mathematical treatment to reduce the number of free parameters. Many experimental studies showed that for ductile metals, the strain leading to shear fracture is less than that to simple tension fracture [10-13]. In the case of plane stress, the ductility under pure shear (T=0, L=0) is the local minima [3, 34]. Hence, take the partial derivative and the second derivative of εfs, with respect to T, we have 1 1  dL  1 1  q3T  q4 L eq3T  q4 L  q3  q4  q1 T  1 n    q1 T  1 n  q2  e dT  n   2 2 q1 T  1 n



 T

s f T 0



0

(29)

T 0

 2 fs T 2

0

(30)

T 0

In plane stress space, the condition σ3 = 0 uniquely relates the triaxiality T and the Lode parameter L through

1  27  1    L  3 tan  arccos   T  T 2      3  6   2  3

(31)

thus the derivative and the second derivative of L over dT at T=0 in the plane stress space are dL dT

 T 0

3 3 2

(32)

d 2L 0 dT 2 T  0

(33)

Substitue Eq. 32 and Eq. 33 back into Eq. 29 and Eq. 30, we have 13

q2  nq3 

q3 

3 3 nq4  1 2

3 3 1 q4   1 2 n

(34)

(35)

and Eq. 28 becomes

  s f

nq3 

3 3 nq4  1  eq3T  q4 L 2 1

q1  T  1 n

(36)

Like the exponent of 1.5T in the Rice–Tracey void growth equation [8], q3T–q4L reflects the microstructural characteristic of the void evolution, such that q3 and q4 can be regarded as constants regardless of material types. In FE unit cell study, we have q3=p3=2.37 and q4=p4=1.50 for ideal material. In real metal material, however, there are inclusions and second phase particles inside the voids after their nucleation, thus the shapes of the voids cannot be arbitrarily altered, nor can the voids completely vanish. But we find that, with appropriate free parameters of p1 to p4, Eq. 26 can always be used to represent the influences of the stress triaxiality and the Lode parameter on the critical strain when the void size reduces to a certain value. For instance, if the critical state of shear band formation due to void collapse and unit cell softening is reached when the void has shrunken by 15% of its initial size, Eq. 26 also gives a good prediction of the critical strain under different stress states whereas a bigger p3 of 4.96 and smaller p4 of 0.94 should be used instead. The micro-mechanism reason behind the mathematical fitting reason to adopt a smaller q4 (or p4) is that, the effect of the Lode parameter, i.e., controlling the void shape, in real metals is indeed smaller than what was shown in the FE unit cell study due to the existence of inclusions inside the voids. Therefore, the value of q4 must be much smaller than the 1.50 obtained in previous FE unit cell study of an idealized hollow spherical void to reflect the actual microstructural responses. Since q4 is very small in value and thus has marginal influence on the fracture strain prediction, it is taken as a constant for different metal materials to reduce the number of free parameters, while q3 remains a free parameter to keep flexibility and versatility of application in various metallographically different metals. In this study, q4 is set as 0.4 (a small positive value) to make a concise form of the fracture strain calculation formula. 14

By replacing q1 and q2 with A and B, Eq. 36 and Eq. 35 become: An  n  1  e AT 0.4 L

 fs 

B  T  1

1 n

, with A 

1 2 n

(37)

4. Fracture criterion based on different fracture mechanisms

In high triaxiality regime, many widely accepted fracture models are derived from the Rice and Tracey void growth equation, such as those by Hancock and Mackenzie [1], and Kanvinde and Deierlein [9], as previously introduced. They all describe the ductility loss due to stress triaxiality in exponential functions. In this study, Eq. 38 is adopted to calculate the strain to high-triaxiality tension fracture εft.

 ft  Ce DT

(38)

Comparing with the those based on the original Rice–Tracey void growth equation, Eq. 38 introduces another free parameter (i.e., D) other than adopting a constant coefficient of 1.5 to ensure better precision for a variety of metallographically different metal materials. Therefore, an overall fracture criterion can be defined by distinguishing the shear fracture mode and the tension fracture mode, and taking the lowest calculated fracture strain as the fracture limit εf. According to Eq. 12 or Eq. 14, no void nucleates when T≤ –1, such that no microscopic damage is caused in the material and fracture shall not be induced. Therefore,  

f  



 min   ,  s f

T  1 t f



(39)

T  1

where εfs and εft can be calculated through Eq. 37 and Eq. 38, respectively. Bao [35] reported fifteen different specimen tests on 2024–T351 aluminum alloy, covering a wide range of stress triaxiality from –0.3 to 0.9. The test results can be used to validate the fracture criterion given by Eq. 37 to Eq. 39. A summary of the specimens is presented in Fig. 7. Test # 1 to 3 are typical tensile specimen tests in high triaxiality regime where the fracture strain is hardly affected by the Lode parameter, therefore, their equivalent strains to fracture εf 15

are plotted against the stress triaxiality T alone; test #1 is a standard axisymmetric tensile specimen test, giving the true stress–strain curve and a hardening exponent of n=0153. Test # 4 to 15 represent almost exactly plane stress states, thus their results are also displayed on the εf–T plane while the Lode parameter L can be calculated through Eq. 31. The validation strategy is to use a necessary amount of tests to calibrate the free parameters A to D, and use additional tests to check the accuracy of the calibrated criterion. The two different branches of the fracture criterion can be calibrated separately. The two free parameters {C, D} characterizing tension fracture in Eq. 38 will be determined from Test # 3 and 15. These two types of specimens, circumferentially notched round specimen and rectangular solid bar specimen, are commonly adopted while representing different stress triaxiality, thus are chosen in this calibration. Substituting the εf–T data from Test # 3 and 15 ({T=0.93, εf=0.16} and {T=0.37, εf=0.35}) into Eq. 38 and solving a set of nonlinear algebraic equation, one gets C=0.59 and D=1.40. The predicted fracture strain using these two values in tension fracture mode is shown in Fig. 7 by the thick dotted line. Apparently, the curve passes through the test points # 3 and 15; other test points closely surround the curve. A different set of tests, #6 an upsetting cylinder test and #10 a pure shear specimen test, are selected for calibration of the shear fracture branch of the fracture criterion. Substituting the test results ({T=–0.23, L=0.63, εf=0.38} and {T=0.01, L=–0.03, εf=0.21}) into Eq. 37 and solving a set of nonlinear algebraic equation, the following values of the free parameters are obtained: A=7.33 and B=6.08, and the shear fracture locus is shown in Fig. 7 by a thin solid line. Alternatively, one can reduce the number of necessary tests for shear branch calibration down to one by introducing a simple presumption-in plane stress space, the shear branch and the tension branch is divided by the simple axisymmetrical tension (T=1/3, L=–1). This presumption is based on an evident experimental phenomenon that tension fracture always occurs with stress triaxiality higher than 1/3. For the calibrated tension branch, strain to fracture equals to 0.37 when T=1/3, thus substituting {T=1/3, L=–1, εf=0.37} and the test results of Test #10 ({T=0.01, L=–0.03, εf=0.21}) into Eq. 37, one gets: A=6.52 and B=5.49. The corresponding shear fracture locus is also shown in Fig. 7 by 16

the thick solid line, which is observed to be very close to the thin-solid-line fracture locus calibrated by two specimen tests. Other test points lie closely to both curves, indicating an excellent performance of the proposed shear fracture criterion. The fracture locus predicted by the second calibration strategy provides satisfactory precision while reducing half of the workload, thus are recommended in future application. It is to noted that, following the presumption of a fixed dividing point between the shear branch and the tension branch in the plane stress space, one can of course calibrate the shear branch first by two specimen tests (e.g. Test #6 and 10) and then calibrate the tension branch by only one test (Test #3). However, we prefer not to do so for two reasons. On one hand, the tensile specimen tests are easy to perform, and the corresponding numerical simulation are also of less difficulty. On the other hand, numerical simulation of the upsetting cylinder test (Test #6) is complex despite the simple specimen geometry, because trail calculations are needed to find an appropriate friction coefficient between the loading platen and the end surfaces of the specimen [3]. The proposed fracture criterion predicts fracture of metallic materials correctly in all types of specimen tests. Using the same test results reported by Bao [35], Wierzbicki et al. [34] calibrated and validated eight fracture criteria, as shown in Fig. 8, including the constant equivalent strain criterion, the maximum shear stress criterion, the fracture forming limit diagram (FFLD) [36], the Johnson–Cook model [2], the Xue–Wierzbicki fracture criterion [14], the Wilkins model [16], the CrachFEM fracture models [19], and the Cockcroft–Latham model [37]. It is observed that the criterion proposed in this study shows much better performance than most of the above criteria which only work well in certain ranges of stress triaxiality. Fig. 9 shows the 3D locus of the calibrated fracture criterion in (εf, T, L) space, in which the shear branch is calibrated by Test #10 and the dividing point at simple axisymmetrical tension, i.e., the second calibration strategy. In negative triaxiality regime, the locus predicts an increase of ductility with decreasing stress triaxiality; when the stress triaxiality is smaller than –1/3, the ductility increase is rather drastic. In an axisymmetrical compression case of {T=–0.4, L=1}, the fracture strain amounts up to 1.03, which is indeed very high and almost unattainable in standard tests for this 17

material (2024–T351 aluminum alloy). The predicted drastic increase of fracture strain coincides well with the trend of the test results, as well as the empirical fracture locus reported by Bao and Wierzbicki [3] who assume a “cut-off” value of the stress triaxiality at a negative stress triaxiality of T=–1/3 below which no fracture would occur. In positive triaxiality regime, the locus shows an increasing ductility as the triaxiality increases in the low positive triaxiality regime, followed by an abrupt transition to the standard high triaxiality regime with the opposite trend. At the intermediate levels of triaxiality, ranging from 1/3 (L=–1) to 0.57 (L=1) for this calibrated material, there is a competition between shear fracture and tension fracture. We also note that for fixed, small values of stress triaxiality, the ductility decreases as the Lode parameter increases from –1 to +1. This influence of the Lode parameter agrees well with the conclusion given by Danas and Ponte Castañeda [25].

5. Parameter calibration for Q345 steel and fracture simulation

5.1 Specimen tests Fracture is an important mode of failure in steel buildings. During the Northbridge Earthquake in 1994, fractures of beam–column connections were observed in over 150 steel buildings [38]. In the following year, the Kobe Earthquake caused unexpectedly considerable destruction of steel buildings; among 337 seriously damaged steel buildings that were assessed as unrepairable, cracks and fractures were regarded as common failure modes [39]. Therefore, accurately calibrated fracture criterion of the constructional steel materials brings great benefits in evaluating the structural performances and thus reducing the loss from natural or man-made disasters. In mainland China, Q345 steel with a nominal yield strength of 345 MPa is the most commonly adopted mild constructional steel. Hence, four types of specimens made from this metal material were tested to calibrate the free parameters following the calibration strategy proposed in this paper. a. smooth round bar (SRB) specimens, as shown in Fig. 10a. The specimens were used to obtain the true stress–strain curve and the strain hardening exponent n. 18

b. notched round bar (NRB) specimens, as shown in Fig. 10b. Specimens NRB10, NRB5 and NRB0.5 had circumferential notch radius of 10 mm, 5 mm and 0.5 mm, respectively. These specimens represent various high triaxiality stress states. c. rectangular solid bar (RSB) specimens, as shown in Fig. 10c. The specimens represent a medium triaxiality stress state. d. shear (S) specimens, as shown in Fig. 10d. These specimens were developed with a greatly weakened central o

o

o

region, and specimens S0, S30 and S60 featured an angle of 0 , 30 and 60 , respectively, between the loading direction and the line joining the centers of the circles at the edges of the central region. In order to ensure highly localized plastic strain and increase the possibility of crack formation in the central region, two arc channels were machined on both the top and the bottom surfaces of the central region. Preliminary FE simulations revealed the necessity of releasing the transverse restraint in the central region in generating a desired shear stress state, therefore, a hole was drilled at each end of the specimen. The shear specimen was loaded with the aid of the loading accessories specifically designed for this test program, as shown in Fig. 10e. Specimens NRB5, RSB and S0 were used to calibrate the free parameters, while the others helped to validate the fracture criterion that had been calibrated. All specimens were cut from a 1200 mm long, 440 mm wide and 20 mm thick block of Q345 steel. Fig. 10f shows the testing setup. Each specimen was pulled at the ends using an MTS–810 Material Testing System. The tests were carried out at room temperature under displacement-control mode at quasi-static loading rates, following the ASTM E8/E8M standard [40]. The instantaneous gauge length was measured (original gauge length is 50 mm), and the pulling force was recorded. For each type of specimen, a total of three specimens were tested and an average force–deformation curve was obtained by means of simple averaging. Test results of the SRB specimens were reported in a previous paper [41], in which the true stress–strain relation was determined, as shown in Fig. 11. The strain hardening exponent is equal to the strain to necking, i.e., n=0.145. 19

Fig. 12 presents the force–deformation curves obtained in the NRB, RSB and shear tests, among which the curves of the NRB tests were previously reported in [41] and are redrawn here (Fig. 12a-12c). For NRB5, RSB and S0 test, there is an approximately linear drop of the pulling force in each curve, and the beginning of the drop is taken as the fracture initiation in this study. The corresponding gauge deformation is denoted as the gauge deformation to fracture uf, which amounted to 1.94mm, 14.9mm and 3.31mm, respectively, for NRB5, RSB and S0 test. Due to material hardening, load could still increase even if the structure starts to crack. However, cracks were observed to grow very rapidly in all tests, such that it is reasonable to take the beginning of the drop as an approximate indication of the fracture initiation since those two stages are very close to each other. Special attention was paid to the shear specimens representing shear fracture mode. Taking specimen S0 for instance, plastic deformation was highly localized at its central region with fracture initiated at the center, as shown in Fig. 13. The crack propagated rapidly along a straight line to both the upper and lower edges, causing a slightly tilted cut-through crack. 5.2 Parameter validations For all tests, plastic deformations and stress states around the fracture initiation locations are too complex to be accurately calculated from the testing data, therefore, numerical simulations of the NRB5, RSB and S0 tests were carried out to obtain all individual components of the stress and strain tensors at fracture locations. The FE models are developed using the general purpose code ABAQUS/Standard. As shown in Fig. 14, 4-node axisymmetrical elements are used in computation of NRB5 specimen; 8-node linear brick elements with reduced integration and a small number of 4-node linear tetrahedron elements are used for RSB specimen, and 8-node linear brick elements with reduced integration are used for S0 specimen. It is to be noted that non-local effect exists when this proposed model is used under some complicated non-uniform stress-strain state. Therefore, modern theories of fracture generally accept the need for a material length scale, i.e., a parameter with the dimensions of length which is incorporated in some way or other into the theoretical model [42]. 20

Some models use physical lengths related to the microstructure of the material, such as the grain size or inclusion spacing, while others put forward continuum-mechanics based models to determine the length scales. The former approach was adopted in a previous research [41], and the characteristic length of the Q345 steels is determined as 0.25 mm. Therefore, a mesh size of approximately 0.25 mm is adopted for the necking regions in NRB5 and RSB specimens and the central region of the S0 specimen. The material is modeled by isotropic von Mises plasticity with a strain hardening stress–strain relationship shown in Fig. 11; no damage and fracture properties of the material are considered. The boundary conditions of the FE models are shown in Fig. 14. For the S0 specimen, the two pins were modeled as rigid bodies with a fixed boundary applied at one pin and a displacement loading applied on the other; surface-to-surface frictionless contacts are introduced to model the interaction between the pins and the holes. For all tests, the correlation between the experimental and numerical results is almost perfect in terms of the pulling force versus gauge deformation responses, as shown in Fig. 15, thus it is reasonable to study the fracture initiation behavior from the numerical simulations. The locations of fracture initiation are determined from experimental observation: in the tensile tests (i.e., NRB5 and RSB test) fractures initiate at the inner core of the necking cross sections, while fracture in the S0 test initiates from the outer surface of the localized central region’s center. Numerical results have confirmed that all the above locations have the highest equivalent strain developed. Therefore, the fracture strain εf is determined as the equivalent strain corresponding to the gauge deformation to fracture uf, and found to be 0.34, 1.16 and 0.86 for NRB5, RSB and S0 test, respectively. The stress state during the cause of loading process is not constant, as shown in Fig. 16 and Fig. 17, thus the average stress triaxiality Tav and the average Lode parameter Lav are defined to characterize the equivalent stress state:

Tav  Lav  where

1

f 1

f

f



0



f

0

T d

(40)

Ld

(41)

is the equivalent strain. In NRB5 and RSB tests, Tav amount to 0.86 and 0.44, respectively. In S0 test, Tav is 21

equal to 0.02, and Lav is –0.07. Substituting the two sets of Tav–ε f relations from RSB and NRB5 tests, i.e., {T=0.44, εf=1.16} and {T=0.86, εf=0.34}, into Eq. 38 and solving a set of nonlinear algebraic equation, we get C=4.25 and D=2.95. Hence, the tension branch of the fracture criterion for Q345 steel is:

 ft  4.25e2.95T

(42)

From Eq. 42, the fracture strain at T=1/3 is 1.59. Therefore, substituting {T=1/3, L=–1, εf=1.59} and {T=0.02, L=–0.07, εf=0.89} obtained from the S0 test into Eq. 37 and solving a set of nonlinear algebraic equation, the free parameters in the shear branch of the fracture criterion are obtained as A=7.08 and B=1.38, and we have:

 fs 

0.17  e7.08T  0.4 L 1.38 T  1

(43)

6.9

Eq. 39, 42 and 43 constitute a fracture locus for Q345 steel. The fracture locus can be projected into a εf–T plane in the plane stress space, as shown in Fig. 18. It is to be noted that the specimen tests were carried out at a quasi-static loading rate, thus the above fracture locus is only applicable for static or quasi-static loading scenarios. 5.3 Numerical simulation of the specimen fracture The calibrated fracture criterion can be embedded into FE programs to predict or simulate the fracture phenomena in engineering problems. In most engineering cases, however, the local stress state around the fractured region changes over the loading history due to the presence of large-scale yielding. Therefore, the fracture criterion is adapted in an incremental form to take into account the current stress state in every increments: 1

d  s

d t 

1

 fs 1



t f

  

  

B  T  1 n An  n  1  e AT 0.4 L

 

(44)

1   C  e DT

(45)

where Δd s and Δd t are the incremental damage caused by shear and tension loading, respectively. In each increment, the shear fracture strain εfs and the tension fracture strain εft are calculated according to the current stress state at the 22

material point, and the one with a smaller value determines the loading mode, under which the corresponding incremental damage is caused. To the authors’ knowledge, there is still limited knowledge about the interaction between the shear damage and the tension damage so far, therefore, a linear accumulative total damage is assumed in this study and thus fracture is predicted to occur when

 d   d s

s

t

1

(46)

t

When it comes to simulating fractures and cracks by finite element or meshfree methods, the last two decades have seen a rapid development of various methods including the interelement separation models [43], the extended finite element method [44] and the cohesive segments method [45]. A cracking-particle method [46-48] has recently been proposed to present a meshfree method that does not require an explicit crack representation. This method has been proved to be a robust and efficient approach for modeling discrete cracks and has been applied to several two- and three-dimensional problems in statics and dynamics, where the method does not show any “mesh” orientation bias. Peridynamics [49] has attracted broad interests for researchers in computational solid mechanics since no techniques such as smoothing the normals of the crack surfaces in order to avoid erratic crack paths is necessary in this method; more recently, a dual-horizon peridynamics formulation has been proposed [50, 51]. In this study, the fracture in specimen tests are simulated using the pseudo element deleting method introduced by Yan et al. [41] in an implicit time integration solver ABAQUS/Standard. This method deletes element indirectly by means of reducing Young’s Modulus of this element down to zero, such that fractures and cracks are simulated by elements with no stress but entitled with infinite deformability. More advanced models based on discrete fracture approaches were for instance proposed by Areias et al. [52, 53]. Therefore, Eq. 44 to 46 are embedded into the general purpose code ABAQUS/Standard as a user subroutine USDFLD, the flow chart of which is shown in the Appendix A, to simulate the NRB, RSB and shear tests with special focus on the fracture phenomenon. FE models similar to that in Section 5.2 are developed, and the calibrated free 23

parameters are adopted. FE solvers including ABAQUS/Standard do not allow non-positive values for material properties, therefore, Young’s Modulus of the elements to be deleted is assigned with a sufficiently small value (1 MPa in this study). Fig. 19, Fig. 20 and Fig. 21 show the comparisons between the tests and FE results, in terms of force versus gauge deformation responses, of the NRB, RSB and shear tests respectively. It is observed that numerical simulations give excellent predictions on the fracture initiation as well as the crack propagation. As shear fracture is the major focus of this study, Fig. 22 presents a more detailed look into the fracture initiation and crack propagation in the shear specimen, herein taking specimen S0 as an example. Fracture initiates at the outer surface of the localized central region, and crack propagates quickly along the direction of the plastic deformation band and towards the inner core of the central region as well, forming a slightly tilted cut-through crack in the end. The fracture process predicted by the numerical simulation, as well as the tilting angle of the fracture surface, coincides well with the experimental observation (see Fig. 13). Therefore, the fracture criterion proposed in this study may be used for future application in large-scale engineering problems.

6. Conclusion

This present study develops a comprehensive fracture criterion for ductile metals. To capture the most significant feature of fracture behavior including both shear fracture and tension fracture, the proposed criterion is composed of two competing while cooperating branches corresponding to different fracture modes of shear and tension. A shear fracture criterion is newly developed based on the significant roles of void nucleation and void-collapse-induced shear band formation in shear fracture. Argon void nucleation formula is adopted to predict the critical strain to nucleation, and a series of finite-element unit cell based micromechanical study is carried out to determine the void collapse strain. According to this criterion, fracture in low positive and negative triaxiality regime is triaxiality sensitive and Lode parameter dependent. The tension fracture criterion is derived from the Rice and Tracey 24

void growth equation, thus it predicts the ductility loss due to stress triaxiality in an exponential function; the Lode parameter has no influence on the fracture in high triaxiality regime. The overall fracture criterion comprises five free parameters including the strain hardening exponent n, and is validated by a set of fifteen tests on 2024-T351 aluminum alloy reported by Bao [35]. A calibration strategy of using four specimens including smooth and notched round bar, rectangular solid bar and pure shear specimen is suggested. The proposed criterion is proved to be capable of predicting fracture correctly in all types of specimen tests. Following the calibration strategy, specimen tests are performed to calibrate and validate the fracture criterion for Q345 steel. Numerical simulation of the tests are performed, and the correlation of experimental and numerical results in terms of the force–deformation responses is very good. Numerical results also give excellent predictions on the fracture initiation as well as the crack propagation, indicating that the proposed fracture criterion can be applied to fracture simulation in large-scale structural engineering problems.

Acknowledgement

The work presented in this paper was funded by the National Natural Science Foundation of China (No. 51678432). The authors would like to thank Xiao-jing Cai (Engineering Mechanics Experimental Center in Shanghai Jiao Tong University) and Jin-tai Liu (Department of Structural Engineering in Tongji University) for their assistance in specimen tests.

Appendix A. Flow chart of the USDFLD subroutine in Abaqus/Standard

Fig. A-1 References [1] J.W. Hancock and A.C. Mackenzie, On the mechanisms of ductile failure in high-strength steels subjected to multi-axial

stress-states, J MECH PHYS SOLIDS 24 (1976) 147-169.

25

[2] G.R. Johnson and W.H. Cook, Fracture characteristics of three metals subjected to various strains, strain rates, temperatures and

pressures, ENG FRACT MECH 21 (1985) 31-48.

[3] Y.B. Bao and T. Wierzbicki, On fracture locus in the equivalent strain and stress triaxiality space, INT J MECH SCI 46 (2004)

81-98.

[4] L. Driemeier, M. Brünig, G. Micheli and M. Alves, Experiments on stress-triaxiality dependence of material behavior of

aluminum alloys, MECH MATER 42 (2010) 207-217.

[5] P. Brozzo, B. DeLuca and R. Rendina, A new method for the prediction of formability limits in metal sheets, Sheet metal forming

and formability, in Proceedings of the 7th Biennial Conference of the International Deep Drawing Research Group, Amsterdam,

Netherlands, 1972, pp. 9-13.

[6] D.M. Norris, J.E. Reaugh, B. Moran and D.F.

ui ones, A plastic-strain, mean-stress criterion for ductile fracture, J ENG

MATER-T ASME 100 (1978) 279-286.

[7] M. Oyane, T. Sato, K. Okimoto and S. Shima, Criteria for ductile fracture and their applications, J MECH WORK TECHNOL 4

(1980) 65-81.

[8] J.R. Rice and D.M. Tracey, On the ductile enlargement of voids in triaxial stress fields, J MECH PHYS SOLIDS 17 (1969)

201-217.

[9] A.M. Kanvinde and G.G. Deierlein, Void growth model and stress modified critical strain model to predict ductile fracture in

structural steels, J STRUCT ENG 132 (2006) 1907-1918.

[10] E. El-Magd, H. Gese, R. Tham, H. Hooputra and H. Werner, Fracture criteria for automobile crashworthiness simulation of

wrought aluminium alloy components, MATERIALWISS WERKST 32 (2001) 712-724.

[11] T. Wierzbicki, Y.B. Bao and Y.L. Bai, A new experimental technique for constructing a fracture envelope of metals under

multi-axial loading, in Proceedings of the 2005 SEM Annual Conference and Exposition on Experimental and Applied

Mechanics, Portland, USA, 2005, pp. 1295-1303.

26

[12] I. Barsoum and J. Faleskog, Rupture mechanisms in combined tension and shear - Experiments, INT J SOLIDS STRUCT 44 (2007) 1768-1786.

[13] T. Coppola, L. Cortese and P. Folgarait, The effect of stress invariants on ductile fracture limit in steels, ENG FRACT MECH

76 (2009) 1288-1302.

[14] T. Wierzbicki and L. Xue, On the effect of the third invariant of the stress deviator on ductile fracture, Technical report, Impact

and Crashworthiness Laboratory, Massachusetts Institute of Technology, Cambridge, USA, 2005.

[15] Y.L. Bai and T. Wierzbicki, A new model of metal plasticity and fracture with pressure and Lode dependence, INT J

PLASTICITY 24 (2008) 1071-1096.

[16] M.L. Wilkins, R.D. Streit and J.E. Reaugh, Cumulative-strain-damage model of ductile fracture: simulation and prediction of

engineering fracture tests, ucrl-53058, Technical report, Lawrence Livermore National Lab., Livermore, USA, 1980.

[17] L. Xue, Damage accumulation and fracture initiation in uncracked ductile solids subject to triaxial loading, INT J SOLIDS

STRUCT 44 (2007) 5163-5181.

[18] L. Xue, Stress based fracture envelope for damage plastic solids, ENG FRACT MECH 76 (2009) 419-438.

[19] H. Hooputra, H. Gese, H. Dell and H. Werner, A comprehensive failure model for crashworthiness simulation of aluminium

extrusions, INT J CRASHWORTHINES 9 (2004) 449-464. [20] T.L. Anderson, Fracture Mechanics - Fundamentals and Applications, second ed., CRC press, New York, 1997. [21] I. Barsoum and J. Faleskog, Micromechanical analysis on the influence of the Lode parameter on void growth and coalescence,

INT J SOLIDS STRUCT 48 (2011) 925-938.

[22] A.S. Argon, J. Im and R. Safoglu, Cavity formation from inclusions in ductile fracture, Metallurgical Transactions A 6 (1975)

825-837.

[23] F.M. Beremin, Cavity formation from inclusions in ductile fracture of A508 steel, Metallurgical Transactions A 12 (1981)

723-731.

27

[24] P. Ponte Castañeda and M. Zaidman, Constitutive models for porous materials with evolving microstructure, J MECH PHYS

SOLIDS 42 (1994) 1459-1497.

[25] K. Danas and P. Ponte Castañeda, Influence of the Lode parameter and the stress triaxiality on the failure of elasto-plastic

porous materials, INT J SOLIDS STRUCT 49 (2012) 1325-1342.

[26] S.H. Goods and L.M. Brown, The nucleation of cavities by plastic deformation, Acta Metallurgica 27 (1979) 1-15.

[27] L.M. Brown and W.M. Stobbs, The work-hardening of copper-silica v. equilibrium plastic relaxation by secondary dislocations,

PHILOS MAG 34 (1976) 351-372.

[28] N. Bonora, D. Gentile, A. Pirondi and G. Newaz, Ductile damage evolution under triaxial state of stress: theory and experiments,

INT J PLASTICITY 21 (2005) 981-1007.

[29] M. Kuna and D.Z. Sun, Three-dimensional cell model analyses of void growth in ductile materials, INT J FRACTURE 81 (1996)

235-258. [30] L. Xia and C.F. Shih, Ductile crack growth - I. A numerical study using computational cells with microstructurally-based length scales, J MECH PHYS SOLIDS 43 (1995) 233-259.

[31] C. Ruggieri, T.L. Panontin and R.H. Dodds, Numerical modeling of ductile crack growth in 3-D using computational cell

elements, INT J FRACTURE 82 (1996) 67-95. [32] J. Faleskog, X.S. Gao and C. Shih, Cell model for nonlinear fracture analysis - I. Micromechanics calibration, INT J FRACTURE 89 (1998) 355-373.

[33] S.M. Graham, T.T. Zhang, X.S. Gao and M. Hayden, Development of a combined tension–torsion experiment for calibration of ductile fracture models under conditions of low triaxiality, INT J MECH SCI 54 (2012) 172-181.

[34] T. Wierzbicki, Y.B. Bao, Y.W. Lee and Y.L. Bai, Calibration and evaluation of seven fracture models, INT J MECH SCI 47

(2005) 719-743.

[35] Y.B. Bao, Prediction of ductile crack formation in uncracked bodies, PhD Thesis, Massachusetts Institute of Technology,

28

Cambridge, USA, 2003.

[36] S.P. Keeler and W.A. Backofen, Plastic instability and fracture in sheets stretched over rigid punches, TRANS ASM 56 (1963)

25-48.

[37] M.G. Cockcroft and D.J. Latham, Ductility and the workability of metals, J I MET 96 (1968) 33-39.

[38] S.A. Mahin, Lessons from damage to steel buildings during the Northridge earthquake, ENG STRUCT 20 (1998) 261-270.

[39] Architectural Institute of Japan, Damage and Lessons of Steel Structures in the Hyogoken–Nanbu Earthquake, Technical report, Tokyo, Japan, 1995.

[40] ASTM international, ASTM E8/E8M - 13a Standard Test Methods for Tension Testing of Metallic Materials, West

Conshohocken, PA, United States.

[41] S. Yan, X. Zhao and A. Wu, Ductile fracture simulation of constructional steels based on yield-to-fracture stress-strain

relationship and micro-mechanism based fracture criterion, J STRUCT ENG (2017) (Accepted).

[42] D. Taylor. The theory of critical distances, ENG FRACT MECH 75 (2008) 1696-1705.

[43] X.P. Xu and A. Needleman, Numerical simulations of dynamic crack growth along an interface, INT J FRACTURE 74 (1996)

289-324.

[44] N. Moës, J. Dolbow and T. Belytschko, A finite element method for crack growth without remeshing, INT J NUMER METH

ENG 46 (1999) 131-150.

[45] J.J.C. Remmers, R.D. Borst and A. Needleman, A cohesive segments method for the simulation of crack growth, COMPUT

MECH, 31 (2003) 69-77.

[46] T. Rabczuk and T. Belytschko, Cracking particles: a simplified meshfree method for arbitrary evolving cracks, INT J NUMER

METH ENG 64 (2004) 2316-2343.

[47] T. Rabczuk and T. Belytschko, A three-dimensional large deformation meshfree method for arbitrary evolving cracks,

COMPUT METHOD APPL M 29 (2007) 2777-2799.

29

[48] T. Rabczuk, G. Zi, S. Bordas and H Nguyen-Xuan, A simple and robust three-dimensional cracking-particle method without

enrichment. COMPUT METHOD APPL M 199 (2010) 2437-2455.

[49] S.A. Silling, Reformulation of elasticity theory for discontinuities and long-range forces, J MECH PHYS SOLIDS 48 (2000)

175-209.

[50] H. Ren, X. Zhuang and T. Rabczuk, Dual-horizon peridynamics, INT J NUMER METH ENG 108 (2016) 1451-1476.

[51] H. Ren, X. Zhuang and T. Rabczuk, Dual-horizon peridynamics: A stable solution to varying horizons. COMPUT METHOD

APPL M 318(Suppl. C) (2017) 762-782.

[52] P. Areias, T. Rabczuk and J.C. de Sá, A novel two-stage discrete crack method based on the screened Poisson equation and local

mesh refinement, COMPUT MECH 58 (2016) 1003-1018.

[53] P. Areias, M.A. Msekh and T. Rabczuk, Damage and fracture algorithm using the screened Poisson equation and local

remeshing, ENG FRACT MECH 158 (2016) 116-143.

30

z σ2

π plane

N

θ

(σ1, σ2, σ3) P (z, r, θ) r

σ1

O

σ3

Fig. 1 Stress state in the space of principal stress. Void

x2

Matrix

x1 x3 Macrocosmic

Microcosmic

Fig. 2 A macroscopic material point and a microstructure-based unit cell. x3

L0 = 59 μm

x2

x1

R0 = 20 μm

Fig. 3 Finite element model of 1/8 of a unit cell. L = -1

L=0

L=1

deformed

T = 0.3

T=0

T = -0.3

initial

Fig. 4 Evolution of the 1/8 initially spherical void in low positive and negative triaxiality regime. 31

Void volume fraction - f / f0

2.0

L = -1 T = -0.3 0 0.3 1.5

1

0

1.0

0.5

0.0 0.0

0.2

0.4

0.6

0.8

1.0

1.2

Strain - εeq

Fig. 5 Evolution of the void volume fraction against the equivalent strain under different stress state. Unit cell study Proposed formula 4 3

εcl 2 1 0 0.3

0.2

0.1

0.0

-0.1

T

-0.2

-0.3 1.0

0.5

0.0

-0.5

-1.0

L

Fig. 6 Void collapse strain in the space of the stress triaxiality and the Lode parameter.

9

Fraction strain - εf

0.8

Tension branch Shear branch calibrated by Test # 10 and 6 Shear branch calibrated by Test #10 and the dividing point Experiment (plane stress) Experiment (axial symmetry) Round, smooth 1 2-3 Round, notch Flat-grooved 4 5-8 Upsetting Cylinder Round, notch (compression) 9 10-11 Flat dog-bone shear specimen Plate with a circular hole 12 Dog-bone 13 Pipe 14 Solid bar 15 Dividing point (T=1/3, L=-1)

13 5

1

6 7

0.4 8

11

12

15 14

10

2

3

4 full points used for calibration

-0.4

-0.2

0.0

0.2

0.4

0.6

0.8

1.0

Stress triaxiality - T

Fig. 7 Prediction of the proposed fracture criterion compared with a set of test points reported by Bao [35].

32

0.4

0.8

Ductile

0.4

Fraction strain - εf

0.8

Fraction strain - εf

0.4

Fraction strain - εf

0.8

Shear 0.2 0.4

0.6

0.8

1.0

-0.4 -0.2 0.0

0.8

1.0

-0.4 -0.2 0.0

0.8

0.4

0.2 0.4

0.6

0.8

1.0

0.8

1.0

Stress triaxiality - T

(c)

Fraction strain - εf

0.4

0.6

(b)

Fraction strain - εf

(a)

0.8

0.2 0.4

Stress triaxiality - T

Stress triaxiality - T

0.8

0.4

Fraction strain - εf

-0.4 -0.2 0.0

C1=-0.07, C2=1.02, C3=-1.62

C1=0.13, C2=0.13, C3=-1.5 -0.4 -0.2 0.0

0.2 0.4

0.6

0.8

1.0

Stress triaxiality - T

-0.4

-0.2

0.6

(e)

0.8

1.0

-0.4 -0.2 0.0

0.2 0.4

0.6

Stress triaxiality - T

(f)

Fraction strain - εf

0.4

0.2 0.4

Stress triaxiality - T

(d)

0.8

-0.4 -0.2 0.0

0.0 0.2 0.4 0.6 Stress triaxiality - T

0.8

1.0

(g) Fig. 8 Calibration and comparision of seven fracture criteria by Wierzbicki et al. (a) Xue–Wierzbicki model; (b) CrachFEM; (c) maximum shear stress; (d) fracture forming limit diagram; (e) Wilkins model; (f) Johnson–Cook fracture model; (g) Cockcroft-Latham model. (redrawn from [34])

33

0.8

εf 0.4 1.0

0.0 -0.4 -0.2

0.5 0.0 0.0 L

0.2

T

0.4 -0.5

0.6 0.8 1.0

Fig. 9 Fracture locus in the space of stress triaxiality and Lode parameter.

34

-1.0

t=5

Φ16

Φ16

Φ16

25

t=5

200

200 70

200 70

R=10 for NRB10 R=5 for NRB5

Φ10

A Φ10

Φ10

15

NRB0.5

R=0.5

A

(a)

(b)

2.5

(c)

12

30o for S30 60o for S60

4

45° 45°

B

B 200

A

t = 3.25

S0 A

S30 / S60

20

0.75 30

B-B

(d) (e) (f) Fig. 10 Specimens of Q345 steel and test setup. (a) SRB specimens; (b) NRB specimens; (c) RSB specimens; (d) shear specimens; (e) loading accessories assembled on the shear specimens; (f) test setup.

35

1000

True stress (MPa)

800

Linea

r pos

t- n e

s tr cking

ain h

arde

ning

600 Necking point 400

200 εnk=0.145 0 0.0

0.2

0.4 0.6 True strain

0.8

1.0

Fig. 11 Ture stress–strain curves of Q345 steel [41]. 60

Force (kN)

50 40 30 20

NRB10-1 NRB10-2 NRB10-3 average curve

10 0

0

1

2 Displacement (mm)

3

4

(a) 60 NRB5-1 NRB5-2 NRB5-3 average curve

Force (kN)

50 40 30 20 10 0

0

1

2 Displacement (mm)

(b)

36

3

4

80

Force (kN)

60

40 NRB0.5-1 NRB0.5-2 NRB0.5-3 average curve

20

0 0.0

0.5

1.0 1.5 Displacement (mm)

2.0

(c) 50

Force (kN)

40 30 20 RSB-1 RSB-2 RSB-3 average curve

10 0

0

5

10 15 Displacement (mm)

20

(d) 4

Force (kN)

3

2

S0-1 S0-2 S0-3 average curve

1

0 0

1

2 3 Deformation (mm)

(e)

37

4

5

4

Force (kN)

3

2

S30-1 S30-2 S30-3 average curve

1

0 0.0

0.5

1.0 1.5 2.0 Deformation (mm)

2.5

3.0

(f) 4

Force (kN)

3

2

1

0 0.0

S60-1 S60-2 S60-3 average curve 0.5

1.0 1.5 Displacement (mm)

2.0

(g) Fig. 12. Force–deformation curves obtained in the specimen tests. (a) NRB10; (b) NRB5; (c) NRB0.5; (d) RSB; (e) S0; (f) S30; (g)S60.

Undeformed

Plastic flow, pre-fracture

Totally fractured

Fig. 13 Plastic deformation and fracture within the localised shear region of the S0 specimen.

38

Displacement

Displacement

Translational constraint

Translational constraint

(a)

(b)

Displacement

A A

Translational constraint Rotational constraint

(c) Fig. 14 FE models of the specimen tests. (a) NRB5; (b) RSB; (c) S0. 60

NRB5 - Test NRB5 - FE RSB - Test RSB - FE S0 - Test S0 -FE

Force (kN)

40

20 full points denote initation of fracture

0 0

5

10 Deformation (mm)

15

20

Fig. 15 Comparison between experimetal and numerical results for load-deformation curve.

39

1.2 NRB5 RSB Stress triaxiality - T

0.9

0.86

0.6 0.44 0.3

0.34

0.0 0.0

1.16 0.5 1.0 Equivalent strain - 

1.5

Stress triaxiality - T & Lode parameter - L

Fig. 16 Stress triaxiality at the fracture initiation locations in the NRB5 and RSB tests. 0.2 T L 0.1 0.02 0.0 -0.07 -0.1

0.86

-0.2 0.0

0.3

0.6

0.9

1.2

Equivalent strain -ε

2.0

1.5

Fracture strain - εf

Fig. 17 Stress triaxiality and Lode parameter at the fracture ininitation location in the S0 test.

1.0

0.5

-0.6

-0.4

-0.2

0.0 0.0

0.2

Stress triaxiality -T

Fig. 18 Fracture locus of Q345 steel in the plane stress space.

40

0.4

0.6

80

NRB10 - Test NRB10 - FE NRB5 - Test NRB5 - FE NRB0.5 - Test NRB0.5 - FE

Force (kN)

60

40

20

0 0

1

2 Deformation (mm)

3

4

Fig. 19 Comparison between experimetal and numerical force–deforamtion responses in the NRB tests. 50 RSB - Test RSB - FE

Force (kN)

40

30

20

10

0 0

5

10 Deformation (mm)

15

20

Fig. 20 Comparison between experimetal and numerical force–deforamtion responses in the RSB tests. 5 S0 - Test S0 - FE S30 - Test S30 - FE S60 - Test S60 - FE

Force (kN)

4

3

2

1

0 0

1

2 3 Deformation (mm)

4

5

Fig. 21 Comparisons between the experimetal and numerical force–deforamtion responses in the shear tests.

41

3.18 mm

3.27mm

3.31mm

3.42mm

3.64mm

Outer layer

Gauge deformation

Middle layer

Undeformed shape

Legend 1000 MPa

Inner layer

500 MPa

0 MPa

Fig. 22. Fracture initiation and crack propagation in the central region of S0 specimen.

Start of iterative step

Acquire state variables ( 、 、 、 ) End of iterative step

Yes No Acquire stress and strain )

Update state variables ( , etc.)

Yes No Calculate

Yes No

No Calculate



E = 1.0 MPa

No Yes Calculate

Calculate

Fig. A-1. Flow chart of the USDFLD subroutine in Abaqus/Standard

42

Yes

Highlights  

A new shear fracture criterion is proposed based on microscopic studies on the void nucleation and evolution behavior. A two-branch fracture criterion is established to capture the metal fracture in both shear and tension modes, and shows better performance than many other criterions.  A calibration strategy is proposed and specimen tests are performed to calibrate the fracture criterion for Q345 steel.

43