Evaluation of fracture toughness and its scatter in the DBTT region of different types of pressure vessel steels

Engineering Failure Analysis 18 (2011) 172–185

Contents lists available at ScienceDirect

Engineering Failure Analysis journal homepage: www.elsevier.com/locate/engfailanal

Evaluation of fracture toughness and its scatter in the DBTT region of different types of pressure vessel steels M.K. Samal a,⇑, J.K. Chakravartty b, M. Seidenfuss c, E. Roos c a

Reactor Safety Division, Bhabha Atomic Research Centre, Trombay, Mumbai 400 085, India Mechanical Metallurgy Section, Bhabha Atomic Research Centre, Trombay, Mumbai 400 085, India c Institut für Materialprüfung Werkstoffkunde und Festigkeitslehre, Universität Stuttgart, Germany b

a r t i c l e

i n f o

Article history: Received 14 July 2010 Received in revised form 24 August 2010 Accepted 24 August 2010 Available online 15 September 2010 Keywords: Fracture toughness master curve Ductile-to-brittle transition Nonlocal formulation Damage mechanics Rousselier’s model

a b s t r a c t Safety analysis of critical components and the overall plant is a major and important task in the field of mechanical engineering. Therefore, it is essential to know the allowable loads and the corresponding failure behaviour (e.g., crack initiation, growth and instability) of the component. Most of the safety critical components, such as pressure vessels of nuclear reactors, are made of ferritic grade low-alloy steels. When these components operate in the ductile-to-brittle transition (DBTT) regime of the material (due to irradiation embrittlement and other material aging and degradation mechanisms taking place), the two mechanisms that compete with each other in the failure process are: ductile and cleavage fracture. From fracture mechanics point of view, fracture toughness of the material must be adequate to prevent the failure. However, there is considerable scatter observed in the fracture toughness data and the analyst must account for this in the safety analysis. Master curve approach according to ASTM E-1921 standard is popularly employed for this purpose. However, it requires the data for transition temperature T0, which is dependent upon the specimen geometry and loading configurations employed in the laboratory tests. Its transferability to safety analysis of components is questionable. In this work, a combined model for ductile and cleavage fracture is used to predict the fracture toughness scatter and its variation with temperature in the DBTT range. It is demonstrated that the above data for fracture toughness can be predicted once we know the material stress–strain data at different temperatures and a single set of Weibull statistics parameters for cleavage fracture. Extensive experimental investigations have been carried out on two types of pressure vessel steels in the DBTT region using different kinds of specimens to validate the predictions of the model. Ó 2010 Elsevier Ltd. All rights reserved.

1. Introduction and motivation for current research Prevention of failure of pressurized and high-energy components and systems has been an important issue in the design of all types of power and process plants. Instead of the traditional fracture mechanics based approaches, it is now possible to go for the detailed modelling of different fracture or material degradation processes using local approaches [1–3]. There are two classes of local approaches, i.e., uncoupled and coupled types. In the uncoupled type of models for ductile fracture (e.g., Rice and Tracey’s model [1]), the material damage is evaluated from the local stress and strain fields in a post-processing

⇑ Corresponding author. Present address: Department of Mechanical Engineering, Ohio State University, Columbus, OH 43201, United States. Tel.: +91 22 25591523/+1 614 6883288; fax: +91 22 25505151. E-mail addresses: [email protected], [email protected] (M.K. Samal). 1350-6307/$ - see front matter Ó 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.engfailanal.2010.08.018

fracture toughness KJc / MPa√m

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

300

173

Pf= 5% Pf=50% Pf=95%

200

100

-100

-80

-60

-40

-20

temperature / °C Fig. 1. Scatter in the experimental fracture toughness data for 22NiMoCr3-7 steel.

exercise. In the coupled type of models for ductile fracture (popularly known as continuum damage models), the effect of material damage is usually considered as an internal variable in the expressions representing the constitutive behaviour of engineering materials (e.g., Rousselier’s [2] and Gurson–Tvergaard–Needleman’s [3–5] model). These local damage models have been used by many researchers to predict the load-deformation and fracture resistance behaviour of different types of specimens and components [6–9]. On the other hand, cleavage fracture data (e.g., fracture toughness) tend to be highly scattered. So, in addition to the existence of a critical stress for initiation of cleavage fracture, there is also an inherent statistics attached to this behaviour. The typical fracture toughness variation with temperature in the DBTT region for a pressure vessel steel material (DIN22NiMoCr3-7) is shown in Fig. 1. This data has been obtained from fracture tests conducted with standard 1T compact tension (CT) specimen at different temperatures. The geometry of such a specimen is shown in Fig. 2. The dimensions for a 1T standard CT specimen are: width W = 50 mm, height H = 60 mm, thickness B = 25 mm. The specimens are usually side-grooved in order to induce plane strain condition during the test and with a 20% side-groove; the net thickness BN is 20 mm. The initial crack length to width ratio (i.e., a0/W) is usually designed to be maintained in the range of 0.5–0.55. One can see from Fig. 1 that the scatter and the mean value of the fracture toughness increase as the temperature increases. At the lower end of temperature scale, the fracture process is purely of cleavage nature, whereas at the higher end of temperature scale, there is some amount of stable ductile crack growth before the initiation of unstable cleavage fracture process. The mean value of fracture toughness increases because of improvement in ductility with temperature. The statistical scatter also increases because of interplay of two different mechanisms (i.e., stable ductile crack growth followed by the unstable cleavage fracture event). In the lower shelf temperature, only one mechanism (i.e., cleavage fracture) is dominant and hence the scatter is less. In literature, prediction of the statistical scatter in the cleavage fracture process is popularly captured by applying Beremin’s model [10] which uses Weibull’s statistical model [11] along with weakest link theory.

B F W

Bn

VLL

a0

F Fig. 2. Standard compact tension specimen used in the fracture toughness test.

174

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

Fig. 3. Maximum principal stress distribution (in MPa) in a typical CT specimen which is used for calculation of Weibull stress as a loading parameter.

The Weibull stress model originally proposed by the Beremin group [10] provides a framework to quantify the relationship between macro and micro-scale driving forces for cleavage fracture in a fully three dimensional continuum. Beremin’s group introduced the scalar Weibull stress rw as a loading parameter (which governs the probability of cleavage fracture for the component and loading under consideration). This parameter is computed by integrating a weighted value of the maximum principal (tensile) stress over the fracture process zone as

rw

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Xn Vi m ¼ ðriI Þm i¼1 V0

ð1Þ

where riI is the maximum principal stress in the material point ‘i’, Vi is the volume associated with the material point (due to finite discretization) and V0 is a characteristic or reference volume in order to render the weight attached to the maximum principal stress in Eq. (1), a dimensionless quantity. The typical principal stress distribution obtained from a numerical (i.e., finite element) analysis of the CT specimen in 2D plane strain domain is shown in Fig. 3. One can see the dominance of the principal stresses in the crack tip region and their rapid decay with distance from the crack tip. This parameter, i.e., the Weibull stress rw serves the role of a crack tip loading parameter and the corresponding probability of cleavage fracture can be calculated as

  m  rw Pf ¼ 1  exp 

ru

ð2Þ

The model described above adopts a two-parameter Weibull description for calculation of the cumulative failure probability. The two parameters being denoted as: m and ru, where m denotes the Weibull modulus (shape parameter) which quantifies the statistical scatter and ru is a scale parameter which sets the value of rw at 63.2% failure probability. The model assumes the Weibull parameters to be inherent properties of the material that describe the formation of a certain distribution of metallurgical scale cracks once plastic deformation occurs as the precursor to cleavage. These parameters are usually determined by combined experiment and numerical simulation of the specimens at purely cleavage region using maximum likelihood method [11]. Calibration of Weibull parameters for the Beremin’s model is an important aspect of cleavage fracture study. Gao et al. [12] proposed a procedure to calibrate the parameters of the Weibull stress model using fracture toughness data obtained from two sets of test specimens that exhibit different constraint levels at fracture. They also introduced a threshold Weibull stress value corresponding to the minimum value of fracture toughness of the material. Using the three-parameter Weibull stress model with parameters calibrated according to the proposed procedure, Gao et al. [13] predicted the distributions of measured fracture toughness values in various specimens of a A515-70 pressure vessel steel, including surface crack specimens subject to different combinations of bending and tension loading conditions. Recent efforts on the Weibull stress model have focused on estimation of effects of temperature, loading rate, etc. on the values of these model parameters. Gao et al. [14,15] studied effects of loading rate on the Weibull stress model for the A515-70 steel and found that the Weibull shape parameter ‘m’ is rate independent over the range of low-to-moderate loading rates. The dependency of the Weibull scale parameter ru on the loading rate characterizes the effect of loading rate on fracture toughness. On the other hand, in order to predict the fracture probability in the DBT region using a numerical model, Beremin’s model for cleavage fracture needs to be coupled with the models for ductile fracture as both the processes are competing with each other in this regime. When simulating the crack-tip stress field at a lower temperature (when stress gradients are large) using finite element (FE) method, one needs to use a very fine mesh, which cannot be used in case of classical (or local) damage mechanics model as the mesh size is pre-defined for a material (which is typically of the order of 0.2 mm). This order of FE mesh is too coarse to correctly capture the low temperature crack-tip stress field and hence the Weibull parameters are

175

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

experiment simulation (Pf =0.05) simulation (Pf =0.50) simulation (Pf =0.95) stable crack growth

350 300 250

0.8

0.6

200 0.4

150 100

0.2

stable crack growth / mm

fracture toughness KJc / MPa√m

1.0 400

50 0

0.0 -100

-80

-60

-40

-20

temperature / °C Fig. 4. Comparison of predicted fracture toughness variation obtained with elastic–plastic analysis with experimental data.

unable to predict the fracture toughness transition curve when the classical damage mechanics or elasto-plastic models are combined with the Beremin’s model. This point is demonstrated here with an example. Fig. 4 shows the results of analysis of combined elastic–plastic and Beremin’s model in the DBTT regime of the pressure vessel steel DIN22NiMoCr3-7. The same 1T CT specimens as shown in Fig. 2 is used for the analysis. 3D iso-parametric 20-noded brick elements were used in the analysis. The FE mesh of the CT specimen is shown in Fig. 5. Due to symmetry, one-quarter of the specimen is modelled along with the imposition of symmetric boundary conditions. The material stress–strain data for different temperatures were obtained by testing standard round tensile specimen and these are shown in Fig. 6 for various temperatures of tests. The Weibull parameters employed in the analysis are: m = 34.6 and ru = 2005 MPa. These data of Weibull parameters adequately describe the pure cleavage fracture process in the lower end of the temperature scale as seen in Fig. 4. The correspondence between the numerical prediction and experimental data are excellent at this temperature. However, the numerical model is quite unable to predict the mean values and scatter of fracture toughness in the higher end of temperature scale. The problem lies in the inability to model small amounts of ductile crack growth (before the onset of unstable cleavage fracture) in the FE analysis when local damage or elastic–plastic models are used. Again, the minimum amount of stable crack growth that can be simulated is of the order of one element size (e.g., of the order of 0.2–0.4 mm). In our experiments, we observed that the average stable crack growth (before onset of instability) is of the order of 0.2 mm at 20 °C and hence to simulate very small amounts of crack growth (of the order of one-tenth of 0.2 mm), we require mesh sizes of the order of 0.02 mm. Clearly, the local damage models cannot be used to handle this order of mesh sizes and hence a mesh-independent damage mechanics model must be employed for the analysis. Although, many researchers have used an empirical variation of Weibull parameters with temperature to predict the fracture toughness transition curve along with elasto-plastic analysis for calculation of Weibull stress [16,17], such an assumption lacks a sound physical basis. From the experimentalists’ point of view, the master curve approach is now-a-days being increasingly employed for prediction of the temperature dependence of fracture toughness of ferritic steels in ductile-to-brittle transition regime [18]. The approach is based on the ASTM E-1921 standard. This approach is based on empirical assumptions and considers the variation and scatter in the fracture toughness in the transition region as a property of ferritic steels and evaluates fracture

Z Y X Fig. 5. 3D FE mesh of the CT specimen with consideration of side-groove geometry (one-quarter model considering symmetry).

176

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185 1600

true stress / MPa

1200

800 0

- 20 C 0 - 60 C 0 -100 C

400

0 0.0

0.5

1.0

1.5

2.0

true plastic strain Fig. 6. Material stress–strain curve at various temperatures obtained from tensile tests of 22NiMoCr3-7 steel.

toughness as a function of temperature for any given fracture probability using the weakest link theory. The nature of variation has been obtained from wide-scale experimental data and suitable functions for their variation with temperature have been proposed. Wallin and his group [19–21] have pioneered the development and application of this approach in structural safety analysis. However, the application of this method requires the data for transition temperature T0, which is dependent upon the specimen geometry and loading configurations employed in the laboratory tests. Its transferability to safety analysis of components is questionable. In an earlier work [22], the authors used a combined model (i.e., the nonlocal Rousselier’s model for ductile fracture and Beremin’s model for cleavage fracture) to predict the fracture toughness scatter and its variation with temperature in the DBTT range for a typical pressure vessel steel. The fracture toughness results of compact tension specimen tests at the DBTT regime was compared with the results of the new model. In this work, we explore the applicability of this model to other specimen types (e.g., three point bend specimens with crack at the centre and loaded like a simply supported beam with load at the centre). From the FE analysis point of view, we also consider the side-groove in the current finite element discretization (which was ignored in the previous analysis), hence modelling the geometry more accurately. The stress distribution with such an analysis will be more accurate compared to the case when side-grooves were neglected for simplicity. As the Weibull stress calculation depends upon the accurate computation of crack-tip stress field, this is very important as accuracy of the predicted cleavage fracture probability depends upon the accuracy of Weibull stress magnitude. It is also important to demonstrate that the new model is able to predict the cleavage failure probabilities for other types of ferritic steels. For this purpose, extensive experimental investigations have been carried out on two different types of pressure vessel steels in the DBTT region using different kinds of specimens (compact tension and three point bend) to validate the predictions of the model for different types of specimens and materials. The fracture experiments have also been conducted at different temperatures (20, 60 and 100 °C for material-1, 10, 40, 75 and 110 °C for material-2) for the two materials in order to demonstrate that the model can predict the fracture toughness variation over the whole temperature range in the DBTT region. It is demonstrated that both the scatter and mean values for fracture toughness can be predicted once we know the material stress–strain data at different temperatures and a single set of Weibull statistics parameters for cleavage fracture for the material. The paper is divided into five sections. The mesh-independent nonlocal model and its finite element implementation which is necessary for simulation of sub-millimetre scale crack growth is described briefly in Sections 2 and 3. The experimental investigations on the fracture toughness variation using different types of specimens and materials are detailed in Section 4. The results of the combined ductile and cleavage fracture model used in the analysis of these specimens are discussed in Section 4. The numerical predictions have been always compared with the experimental data in order to demonstrate their validity. Section 5 ends the discussion with a conclusion and scope for future research. 2. Nonlocal regularization of material damage As discussed in the previous section, the results of analysis using a local damage model are dependent on size of discretization used in the numerical treatment and hence, a very fine mesh cannot be employed in the crack path due to restriction in the choice of element size. Due to the onset of damage, material softening takes place and these changes the nature of the governing differential equation (e.g., loss of ellipticity for elasto-static problems and loss of hyperbolicity for elasto-dynamic problems [23,24]). When the microscopic aspect of material damage is considered, it can be realized that damage development in a micros-structure is not strictly a point-function. It depends upon the state of stress, strain and damage of surrounding regions also due to the underlying micros-structure (i.e., grain boundary, dislocations, dislocation sources and obstacles, etc.). Moreover, there are phenomena involving motion of dislocations which drive the plastic deformation and damage accumulation due to dislocation pile-up at obstacles, etc. which are not explicitly represented in continuum

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

177

mechanics formulations. Hence, we need a regularization technique where the damage development in a material point will be coupled to the state of damage of surrounding points in a region with a characteristic dimension, which depends upon the material. Such a scheme is outlined in Fig. 7. Fig. 7 shows the importance of consideration of influence of surrounding points in a characteristic volume X on the damage development at a material point. The weight given to the points nearer to the point under consideration should be higher and it diminishes exponential with distance from the above point. Such a weightage scheme is shown in Fig. 7 which is the Gaussian weight function W. The increment of the nonlocal variable in a material _ is mathematically defined as a weighted average of the increpoint ~ x, i.e. the increment of nonlocal void volume fraction d, _ ment of the local void volume fraction f in a domain X [Fig. 7], i.e.,

_ ~ dð xÞ ¼

1 Wð~ xÞ

Z

Wð~ y; ~ xÞf_ ð~ yÞdXð~ yÞ

ð3Þ

X

where ~ y is the position vector of the infinitesimally small volume dX. In this work, the regularization treatment is performed on the damage rate and not on the damage. This has been done in order to make it convenient for use in the subsequent incremental nonlinear FE formulation of the problem. The function Wð~ y; ~ xÞ is the Gaussian weight function given by

Wð~ y; ~ xÞ ¼

1 8p3=2 l

3

exp 

j~ x ~ y j2 4l

! ð4Þ

2

and the normalisation factor is the integral of the Gaussian weight function, i.e.,

Wð~ xÞ ¼

Z

Wð~ y; ~ xÞdX

ð5Þ

X

The length parameter l determines the size of the volume, which effectively contributes to the nonlocal quantity and is related to the scale of the microstructure. The above integral nonlocal kernel holds the property that the local continuum is retrieved if l ? 0. By expanding f_ ð~ yÞ in Taylor’s series around the point x, we get

@ f_ 1 @ 2 f_ 1 @ 3 f_ f_ ð~ yÞ ¼ f_ ð~ xÞ þ ðyi  xi Þ þ ðyi  xi Þðyj  xj Þ þ ðy  xi Þðyj  xj Þðyk  xk Þ @xi 2! @xi @xj 3! @xi @xj @xk i 1 @ 4 f_ ðy  xi Þðyj  xj Þðyk  xk Þðyl  xl Þ þ    þ 4! @xi @xj @xk @xl i

ð6Þ

where Einstein’s summation convention applies to the indices i, j, k and l. Substituting Eq. (6) into Eq. (3), we get

_ ~ dð xÞ ¼ f_ ð~ xÞ þ C length r2 f_ ð~ xÞ þ dlength r4 f_ ð~ xÞ þ    2

where the Laplacian operator is defined by r ¼

P

ð7Þ 2

@ i @x2 i

4

2 2

and r ¼ ðr Þ , etc. The gradient parameters Clength and dlength have

the dimensions of length to an even power and the odd derivatives of Eq. (6) vanish due to the nature of the Gaussian weight function of Eq. (4) when substituted in Eq. (3). Taking the Laplacian of Eq. (7), we obtain

Fig. 7. Regularization of the internal variable ‘f’ through a Gaussian weighted integral in a characteristic volume X.

178

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

_ ~ xÞ ¼ r2 f_ ð~ xÞ þ C length r4 f_ ð~ xÞ þ dlength r6 f_ ð~ xÞ þ    r2 dð

ð8Þ 4_

xÞ and other higher order terms Substracting Clength x Eq. (7) from Eq. (8), we obtain Eq. (9) when the terms containing r f ð~ are neglected.

d_  f_  C length r2 d_ ¼ 0

ð9Þ

The above Eq. (9) is the diffusion equation for damage and the increment of the nonlocal variable ‘d’ is linked to the increment of local void volume fraction ‘f’ though a characteristic length parameter ‘Clength’ and the Laplacian of increment of nonlocal damage ‘d’. This is an implicit description of damage diffusion and it needs to be solved along with the mechanical equilibrium equation as discussed in the following section. 3. Finite element implementation of nonlocal damage model For numerical simulations, a nonlocal form of the material yield surface has been constructed from the classical Rousselier’s model [2] as shown in Eq. (10) below. In this yield function, the void volume ‘f’ is replaced by the nonlocal material damage ‘d’ as

/¼

req 1d

þ Drk d exp



p ð1  dÞrk



 Rðeeq Þ ¼ 0

ð10Þ

The equilibrium equation in the continuum to be solved along with the damage diffusion Eq. (9) is written as

r  r þ fb ¼ 0

ð11Þ

The associated boundary conditions are

r  njCf ¼ fm

ð12Þ

ujCu ¼ u0

ð13Þ

rd_  njCd ¼ 0

ð14Þ

where rij is the Cauchy stress tensor and fb is the body force per unit volume and fm is the applied mechanical traction at the surface. Eq. (12) is the traction or force boundary condition, njCf is the normal to the boundary Cf, Eq. (13) is the geometric or essential displacement boundary condition and Eq. (14) is the Neumann or force boundary condition for the damage degree of freedom and njCd is the normal to the boundary Cd of the domain Xt+Dt. In our analysis, we employ an incremental procedure and use the updated Lagrangian formulation to express the equilibrium configuration of the body. Assuming additive decomposition of total strain increment into elastic e_ e and plastic e_ p parts, we can write

e_ ¼ e_ e þ e_ p

ð15Þ

The yield function can be written in terms of mean hydrostatic and deviatoric parts of stress tensor and other field variables as

/ðp; q; Ha ; dÞ ¼ 0

ð16Þ

where p and q are the hydrostatic pressure and von Mises equivalent stress respectively. Ha is internal state variable such as hardening, I is the Kronecker–Delta or second order identity tensor, s is the deviatoric part of stress tensor r. The increase in void volume fraction (due to combined void nucleation and growth process) during plastic deformation in the ductile fracture process can be written as a function F, i.e.,

Fðp; q; Ha ; dÞ ¼ f_ ¼ f_ growth þ f_ nucleation ¼ ð1  dÞe_ p : I þ Aðeeq Þe_ eq

ð17Þ

where e_ p is the increment in plastic strain tensor ep, I is the tensor equivalent to Kronecker–Delta function, e_ eq is the increment in equivalent plastic strain of the matrix material and A is the void nucleation constant which obeys Gaussian distribution (as a function of eeq) and can be expressed as

Aðeeq Þ ¼

fN 1 pffiffiffiffiffiffiffi exp  2 sN 2p



eeq  eN sN

2 !

ð18Þ

The constitutive behaviour of the new material model can be finally obtained in the following forms (i.e., increment of stress tensor and increment of void fraction potential) as functions of increment of total strain tensor and damage variable

@ r ¼ C ep : @ e þ C ed : @d

ð19Þ

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

179

and

@F ¼ C de : @ e þ C dd : @d

ð20Þ

where the material tangent stiffness matrices can be defined as

C ep ¼ 2G

  q 0 4 q 3 nn  2Gmql nI  kmpn In G 1  m J þ Kð1  m ÞII þ  qn pl qtr 3 qtr 2

C ed ¼ 2Gmqd n  kmpd I

ð21Þ

C de ¼ C d11 ðmpl I þ mpn nÞ þ C d12 ðmql I þ mqn nÞ þ C d13 n þ C d14 I C dd ¼ C d11 mpd þ C d12 mqd þ C d15 In the above expression, J is the fourth order unit tensor and J 0 ¼ J  II3. The details of the derivation and the co-efficients of Eq. (21) can be found in Ref. [23]. For implementation of the above model in a finite element framework, the weak forms of the governing equations are expressed in the updated Lagrangian setting as tþDt

R¼

Z

tþDt

tþDt

tþDt X

fbi dui d

Xþ

Z

tþDt

tþDt C f

t si dui dtþDt Cf

ð22Þ

and

Z

2_ _ d_  f_  C ddð length r dÞ  dX ¼ 0

ð23Þ

X

respectively. Expressing the generalized displacement, strain and damage vectors at any material point inside the finite ele_ as [Fig. 8] ^ and b d) ment in terms of the generalized nodal variables (u

^; u ¼ Nu u b _ d_ ¼ N d d;

^ ¼ Bu u b rd_ ¼ Bd d_

tþDt e

ð24Þ

where Nu and Nd are the shape function matrices for the displacement and nonlocal damage variable, Bu and Bd are matrices containing the derivatives of the shape function Nu and Nd respectively. Expanding Eqs. (22) and (23), using consistent material matrices from Eqs. (21) and substituting the expressions for u;tþDt e; d_ and rd_ from Eq. (24), we get the following equations

Z  Z  Z  Z  b ^ ÞT ^þ ^þ dðDu BTu C ep Bu  dX Du BTNL tþDt rBNL  dX Du BTu tþDt rBNL  dX þ BTu C ed Nd  dX d_ X X X X  Z Z T T  Nu fb dX  Nu ts dC ¼ 0 X

 Z  Z  Z  Z b b b _ X ^ NTd Nd  dX d_  N Td C de Bu  dX u NTd C dd Nd  dX d_ þ BTd C length Bd  dX d_ þ NTd tþDt dd X X X X X  Z Z _ X ¼0  NTd tþDt f_ dX þ BTd C length tþDt rdd

b d d_ T

ð25Þ

C

Z

X

ð26Þ

X

Fig. 8. An arbitrary region mapped into four-noded master element with damage ‘d’ as an extra nodal degree of freedom in addition to displacements u and v.

180

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

where ts is the surface traction and is given as ts ¼ rnjCf (Cauchy’s traction law), tþDt r is the matrix and tþDt r is the vector containing Cauchy’s stress components at time t + Dt and BNL is the nonlinear strain–displacement transformation matrix. ^ ÞT and d b For arbitrary value of dðDu d_ T , the terms inside the brackets of Eqs. (25) and (26) should be zero and hence these equations reduce to the set of finite element algebraic equations, which can be written in convenient (matrix) form as

A



K uu þ K NL

K ud

K du

K dd



^ Du ^ Dd



( ¼A

fmext  fmint

)! ð27Þ

fdint

where A is the assembly operator which is used to assemble to element stiffness matrices and the matrices and force vectors of all the elements in the domain X. The domain and boundary of each element are represented by Xe and Ce respectively. Hence the element level stiffness and force vectors can be written as

Z

Xe

BTu C ep Bu  dX;

Z

BTNL tþDt rBNL  dX Z Z K ud ¼ BTu C ed N d  dX; K du ¼  NTd C de Bu  dX Xe Xe Z Z Z K dd ¼ NTd Nd  dX  NTd C dd Nd  dX þ BTd C length Bd  dX Xe Xe Xe Z Z Z fmext ¼ NTu fb dX þ NTu t s dC; fmint ¼ BTu tþDt r  dX Xe Ce Xe Z Z Z fdint ¼ NTd tþDt d_  dX  NTd tþDt f_  dX þ BTd C length rtþDt d_  dX K uu ¼

Xe

K NL ¼

Xe

Xe

ð28Þ

Xe

where the left superscript ‘t + Dt’ refers to the quantities at the end of current time step [and the values are at the previous iteration process (i  1)] of the incremental nonlinear finite element analysis. When the iteration process converges, the values of the variables at iteration steps i and (i  1)also converge. The assembled global FE equations are solved for the global degrees of freedom when we specify the required boundary and loading conditions. It may be noted that the assembled (i.e., global) internal damaged force vectors of the elements becomes a null vector. The stiffness terms Kud, Kdu and Kdd in the element stiffness matrix are contributions of the new nonlocal formulation. The above model is used along with the Beremin’s model [Eq. (2)] for cleavage fracture to simulate the failure probability in the DBTT region. Further details of this model and its implementation in a finite element program can be found in Ref. [24].

4. Results and discussion In order to demonstrate the efficacy of the nonlocal damage model as discussed in the previous section, we provide several examples in this section and discuss the results. The first example is the demonstration of mesh-independent nature of the solutions. The standard 1T compact tension specimen is simulated again with the new model with different sizes of discretization in the crack propagation region. The material of the CT specimen is DIN 22NiMoCr3-7. The results of predicted fracture resistance (J–R) curve (at room temperature) for different mesh sizes employed in the analysis are shown in Fig. 9 along with the experimental data. It can be observed that the results of the new model are mesh-independent and they are also very close to the experimental data. There are other advantages of the nonlocal model. Hence, the model can be used along with a very fine mesh for prediction of failure probability in the DBTT region. The results of such analyses are discussed in the following paragraphs.

1600

experiment 0.4mm mesh 0.2mm mesh 0.1mm mesh 0.02mm mesh

J-resistance / N.mm

-1

1400 1200 1000 800 600 400 200 0 0

1

2

3

crack growth / mm Fig. 9. Predicted fracture resistance curve at room temperature with different mesh sizes in the crack path.

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

fracture toughness KJc / MPa√m

400 350 300 250

181

experiment simulation (without crack growth), Pf =0.05 simulation (without crack growth), Pf =0.50 simulation (without crack growth), Pf =0.95 simulation (with crack growth), Pf =0.05 simulation (with crack growth), Pf =0.50 simulation (with crack growth), Pf =0.95

200 150

master curve (Pf =0.95)

100 50

master curve (Pf =0.05)

0 -100

-80

-60

-40

-20

temperature / °C Fig. 10. Predicted fracture toughness scatter and its variation with temperature by the combined nonlocal ductile damage and cleavage fracture model.

We consider here again the problem of prediction of fracture toughness variation with temperature in the DBTT region for the same material (i.e., DIN22NiMoCr3-7) as discussed earlier in the context of employing the local damage model for analysis. The standard 1T CT specimen has been analysed with the new nonlocal formulation at different temperatures and the results are shown in Fig. 10. It is now very clear that the nonlocal model along with Beremin’s formulation is able to predict the fracture toughness scatter and mean values across the whole DBTT temperature range very accurately. This is due to the ability to use very fine mesh and predict the micron-scale stable crack growth before unstable fracture by cleavage. The results are also compared with the predictions of the master curve which is an empirical equation. However, one needs to know the value of transition temperature T0 in order to plot such a curve. This value of transition temperature in the master curve is also dependent upon the specimen type, geometry and loading conditions and it is also different in case of actual components with varying geometry, loading and crack configurations. It is not possible to know the value of T0 in advance for the component geometry of interest, which is a major limitation of this empirical master curve approach. However, a conservative estimate of T0 is usually used for safety evaluation which may not be optimum. These nonlocal models are based on physically and microscopically motivated damage processes and as these are based on material constitutive formulations; they do not have such limitations as discussed above. Hence, the model should be able to predict the fracture toughness variation across the DBTT regime for different specimen types, material types and loading configurations which will be demonstrated in the following paragraphs. Another specimen type (i.e., a standard single-edged cracked specimen loaded in three point bending as shown in Fig. 11) is considered here for the demonstration. The results of the nonlocal model are shown in Fig. 12. It is again observed that the nonlocal model is able to predict the scatter in fracture toughness very accurately irrespective of specimen type and geometry. In order to demonstrate the applicability of the model to different materials, a different pressure vessel steel (referred here as material-2) is considered. This material is similar to the standard pressure vessel steel of grade SA508Cl.2. Standard axi-symmetric tensile specimen of 6 mm diameter have been prepared from this material and tested at different temperatures in order to evaluate the material stress–strain curve which will be an input for our analysis along with the Webull parameters. The material stress–strain data for this material is presented in Fig. 13. In order to compare the results of analysis with the master curve approach as discussed earlier, the fracture toughness transition temperature T0 is evaluated following the standard ASTM procedure [18]. The load–displacement data directly obtained from tests of standard 1T CT specimens at different sub-zero temperatures are shown in Fig. 14. As can be observed from Fig. 14, the magnitude of load

a

F

L = 4W Fig. 11. Dimensions of the three point bend specimen.

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

fracture toughness KJc / MPa√m

182

400

experiment simulation (P =0.05) simulation (P =0.50)

300

simulation (P =0.95)

200

master curve (P =0.95)

100 master curve (P =0.05)

0 -100

-80

-60

-40

-20

temperature / °C Fig. 12. Comparison of predicted fracture toughness variation for the standard 1T three point bend specimen with experiment and master curve.

Fig. 13. Material stress–strain data at different temperatures for the material-2.

80

o -40 C

70

o -10 C

Load (kN)

60

o -75 C

50 40

o -110 C

30 20 10 0 0

2000

4000

6000

8000

10000

Displacement (micron) Fig. 14. Experimental load–displacement data of the CT specimens for the material-2.

at cleavage fracture is lowest for the test temperature of 110 °C whereas it is highest for the test at 10 °C. There is considerable scatter in both fracture load and fracture ductility data at all temperatures. However, the mean values and scatter of both fracture load and fracture ductility increase with increasing test temperature. This will be clearer when the data is presented in terms of fracture toughness KJC, which is a combination of load and ductility. This data is presented in Fig. 15. Note the variation in scatter and mean values of fracture toughness, which is in line with our earlier argument. The fracture transition temperature T0 is an important parameter in construction of master curve as discussed earlier. This parameter can be evaluated following standard ASTM approach and it should be independent

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

183

1000

800

KJC, MPa m0.5

o

To = -80 C 600

KJC (0.95)

K

400

KJC (0.5)

, M = 30 -40 C

KJC (0.05)

200 KJC (lower bound)

0 -150

-100

-50

0

50

100

o

Temperature ( C) Fig. 15. Determination of fracture toughness transition temperature T0 following ASTM standard.

-50 o

To = 74 C

T0 (oC)

-75

-100

-125

-150 -120

-110

-100

-90

-80

-70

-60

-50

-40

-30

o

Temperature ( C) Fig. 16. Variation of fracture toughness transition temperature T0 for data at different temperatures.

of test temperature. This means, one should get same value of T0 irrespective of which data of fracture toughness is used (i.e., data at 110, 75 and 40 as shown in Fig. 15). However, we obtain a slight variation in the magnitude of T0 when different data sets were used and its variation with temperature is shown in Fig. 16. An average value of 74 °C is obtained when all the data set are considered as shown in Fig. 16 with a lower bound of 84.9 °C and upper bound of 70.3 °C. Analysis is carried out with the mesh-independent damage model and the probability of failure is obtained at different temperatures. The predicted data from numerical simulation is directly compared with the experimental failure probability at all the temperatures in Figs. 17–19. An excellent correspondence between experiment

Probability of failure (Pf)

1.0

Test temperature = -110 deg. C FE analysis exp.

0.8

0.6

0.4

0.2

0.0 0

20

40

60

80

100

0.5

KJC (MPa.m ) Fig. 17. Probability of cleavage failure of CT specimens when tested at 110 °C (experiment vs. FE analysis).

184

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

and analysis has also been obtained for this material, thus validating the new approach. The data of fracture toughness variation and its scatter with temperature has been plotted along with the experimental data and data from master curve approach in Fig. 20. It clearly demonstrates the ability of the model for a reliable safety analysis. The local damage models are clearly inadequate for this purpose of safety analysis due to their inability to handle varying mesh sizes.

Probability of failure (Pf)

1.0

Test temperature = -75 deg. C FE analysis exp.

0.8

0.6

0.4

0.2

0.0 0

50

100

150

200

0.5

KJC (MPa.m ) Fig. 18. Probability of cleavage failure of CT specimens when tested at 75 °C (experiment vs. FE analysis).

Probability of failure (Pf)

1.0

Test temperature = -40 deg. C FE analysis exp.

0.8

0.6

0.4

0.2

0.0 0

50

100

150

200

250

300

0.5

KJC (MPa.m ) Fig. 19. Probability of cleavage failure of CT specimens when tested at 40 °C (experiment vs. FE analysis).

500 450

exp. data FE simulation (P =0.05) f

FE simulation (P = 0.5)

400

master curve (P =0.5) f

f

FE simulation (P =0.95) f

KJC MPa.m

0.5

350 master curve (P =0.95)

300

f

250 200 150 master curve (Pf=0.05,T =-70.3 C)

100 50 0 -140

master curve (Pf=0.05, T =-84.9 C)

-120

-100

-80

-60

-40

-20

0

20

o

Temperature ( C) Fig. 20. Variation of fracture toughness in the DBT region (experiment, FE simulation and prediction of master curve).

M.K. Samal et al. / Engineering Failure Analysis 18 (2011) 172–185

185

5. Conclusions In this work, a scheme for nonlocal regularization of the damage mechanics constitutive equations was presented. It was demonstrated that such a scheme is necessary in order to carry the safety analysis task in situations involving simulation of sub-millimetre crack growths, small specimens, steep stress gradients, etc. Though local damage models have been used effectively to predict the ductile crack growth behaviour of different types of specimens and components, the issue of using a constant mesh size limits its applications in many other situations. With the new model, the analysts will have an effective tool for efficient and reliable safety analysis. The future task will involve application of the model to predict the statistics of fracture toughness variation due to other potential degradation mechanisms such as irradiations embrittlement and other material aging and degradation mechanisms. References [1] Rice JR, Tracey DM. On the ductile enlargement of voids in triaxial stress fields. J Mech Phys Solids 1969;17():201–17. [2] Rousselier G. Ductile fracture models and their potential in local approach of fracture. Nucl Eng Des 1987;105:97–111. [3] Gurson AL. Continuum theory of ductile rupture by void nucleation and growth: part-I: yield criteria and flow rules for porous ductile media. Trans ASME J Eng Mater Technol 1977;99:2–15. [4] Tvergaard V, Needleman A. Analysis of cup-cone fracture in a round tensile bar. Acta Metall 1984;32(1):157–69. [5] Needleman A, Tvergaard V. An analysis of ductile rupture in notched bars. J Mech Phys Solids 1984;32:461–90. [6] Kussmaul K, Eisele U, Seidenfuss M. On the applicability of local approach models for the determination of the failure behaviour of steels with different toughness. In: Mehta, et al., editors. Fatigue and fracture mechanics in pressure vessel and piping. vol. 304. New York: ASME; 1995. p. 17–25. [7] Pitard-Bouet J-M, Seidenfuss M, Bethmont M, Kussmaul K. Experimental investigations on the shallow crack effect on the 10 MnMoNi 55 steel and computational analysis in the upper shelf by means of the global and local approaches. Nucl Eng Des 1999;190:171–90. [8] Pavankumar TV, Samal MK, Chattopadhyay J, Dutta BK, Kushwaha HS, Roos E, et al. Transferability of fracture parameters from specimens to component level. Int J Pres Ves Pip 2005;82(5):386–99. [9] Eberle A, Klingbeil D, Schicker J. The calculation of dynamic JR-curves from the finite element analysis of a Charpy test using a rate-dependent damage model. Nucl Eng Des 2000;198(1–2):75–87. [10] Beremin FM. A local criterion for cleavage fracture of a nuclear pressure vessel steel. Metall Trans A 1983;14A:2277–87. [11] Khalili A, Kromp K. Statistical properties of Weibull distribution. J Mater Sci 1991;26:6741–52. [12] Gao X, Ruriggieri C, Dodds RH. Calibration of Weibull stress parameters using fracture toughness data. Int J Fract 1998;92:175–200. [13] Gao X, Dodds RH, Tregoning RL, Joyce JA, Link RE. A Weibull stress model to predict cleavage fracture in plates containing surface cracks. Fatigue Fract Eng Mater Struct 1999;22:481–93. [14] Gao X, Dodds RH, Tregoning RL, Joyce JA. Weibull stress model for cleavage fracture under high-rate loading. Fatigue Fract Eng Mater Struct 2001;24:551–64. [15] Gao X, Dodds RH. Loading rate effects on parameters of the Weibull stress model for ferritic steels. Eng Fract Mech 2005;72:2416–25. [16] Seebich HP. Doctoral dissertation. Institut für Materialprüfung Werkstoffkunde und Festigkeitslehre, Universität Stuttgart, Germany; 2007. [17] Petti J, Dodds RH. Calibration of the Weibull stress scale parameter ru using the master curve. Eng Fract Mech 2005;72:91–120. [18] ASTM-E 1921-05. Standard test method for determination of reference temperature for ferritic steels in the transition range. American Society for Testing and Materials; 2005. [19] Wallin K. The scatter in KJC-results. Eng Fract Mech 1984;19:1085–93. [20] Wallin K. Fracture toughness transition curve shape for ferritic structural steels. In: Keoh ST, Lee KH, editors. Fracture of engineering materials and structures. Elsevier Applied Science; 1991. p. 83–8. [21] Wallin K. Master curve method: a new concept for brittle fracture. Int J Mater Product Technol 1999;14:342–54. [22] Samal MK, Seidenfuss M, Roos E, Dutta BK, Kushwaha HS. Experimental and numerical investigation of ductile-to-brittle transition in a pressure vessel steel. Mater Sci Eng 2008;A496:25–35. [23] Peerlings RHJ, De Borst R, Brekelmans WAM, Geers MGD. Localisation issues in local and nonlocal continuum approaches to fracture. Eur J Mech ASolids 2002;21:175–89. [24] Samal MK, Seidenfuss M, Roos E, Dutta BK, Kushwaha HS. Finite element formulation of a new nonlocal damage model. Finite Elements Anal Des 2008;44(6–7):358–71.