Finite element computation of fatigue growth rates for mode I cracks subjected to welding residual stresses

Engineering Fracture Mechanics 78 (2011) 2505–2520

Contents lists available at ScienceDirect

Engineering Fracture Mechanics journal homepage: www.elsevier.com/locate/engfracmech

Finite element computation of fatigue growth rates for mode I cracks subjected to welding residual stresses Chin-Hyung Lee, Kyong-Ho Chang ⇑ Department of Civil & Environmental Engineering, Chung-Ang University, 221, Heukseok-Dong, Dongjak-Ku, Seoul 156-756, Republic of Korea

a r t i c l e

i n f o

Article history: Received 17 August 2010 Received in revised form 1 May 2011 Accepted 12 June 2011

Keywords: Welds Welding residual stress Fatigue crack growth rate Finite element analysis Modified J-integral definition

a b s t r a c t This paper presents a finite element modeling procedure for predicting fatigue crack growth rate in butt welds subject to mode I loading condition. Sequentially coupled three-dimensional thermal–mechanical finite element model to simulate welding residual stress was first developed. The weld-induced residual stress effect on the fatigue crack growth rate was then modeled by calculating the stress intensity factor due to the residual stress field based on the superposition rule of the linear elastic fracture mechanics. The results demonstrated the significance of the residual stresses in assessment of the fatigue crack growth rate in the welds. Ó 2011 Elsevier Ltd. All rights reserved.

1. Introduction Most of the steel structures in engineering practice today are fabricated by welding. The advantage of welding as joining process includes high joint efficiency, simple set up and low fabrication cost. However, due to the intense concentration of heating in localized zone and subsequent cooling during welding, highly non-uniform temperature distributions occur across the welds and base metal. The temperature gradients cause the local plastic strains to be formed, thus leading to unavoidable residual stresses there [1–3]. The presence of residual stresses in welded structures can significantly affect the fatigue behavior during external cyclic loading [4,5]. Tensile residual stresses are generally detrimental to the fatigue life by increasing the rate of fatigue crack growth. The microstructure and hardness changes produced by welding process also have adverse effects on fatigue crack growth rate. Nevertheless, welding residual stress has been identified as the most influential factor [6,7]. Accurate prediction and efficient estimation of the residual stresses are therefore essential for structural integrity and fatigue life assessment of the welded part. However, accurately predicting weld-induced residual stresses is a very difficult task because of the complexity of welding process which includes localized heating, temperature dependence of material properties and moving heat source, etc. Although there are a number of techniques for predicting welding residual stresses, with the availability of material properties, the finite element method is convenient and useful [8–13]. It can be used to simulate welding temperature field, welding residual stress field and welding deformation. The superposition rule of the linear elastic fracture mechanics is frequently employed to calculate fatigue crack growth rate in welding residual stress field [14]. The key task is to determine the stress intensity factor resulting from the weldinduced residual stresses using either the weight function method [15] or the finite element method. The weight function ⇑ Corresponding author. Tel.: +82 2 820 5337; fax: +82 2 823 5339. E-mail address: ifi[email protected] (K.-H. Chang). 0013-7944/$ - see front matter Ó 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.engfracmech.2011.06.006

2506

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

Nomenclature c ct

specific heat parameter to reflect stress increment due to the dependence of physical and mechanical material properties on temperature E Young’s modulus h convection heat transfer coefficient I arc current J magnitude of the J-integral Kx, Ky, Kz thermal conductivity effective stress intensity factor Keff Kapplied stress intensity factor due to external applied loading Kresidual stress intensity factor resulting from residual stresses M, P, R number of Gaussian points per element in x, y and z direction, respectively m, C material constants of the Paris’ equation ni unit vector drawn normal to the outside of C N number of elements in V Nk bi-linear shape function Nx, Ny, Nz direction cosines of the outward drawn normal to the boundary q boundary heat flux qj sufficiently smooth function Q rate of moving heat generation per unit volume Q(t) distributed heat flux QA heat input from the welding arc QM energy induced by high temperature melt droplets nodal values of qj for kth node Qjk r(t) radial coordinate with the origin at the arc center on the surface of work piece r0 arc beam radius R stress ratio Reff effective stress ratio S1, S2 path which surrounds the crack tip (Fig. 9) T temperature T0 room temperature uj displacement vector U arc voltage Uik nodal displacements V volume enclosed by S1 and S2 (Fig. 9) wm, wp, wr weights for Gaussian points W mechanical strain energy density Wtotal total strain energy density W p jinitial state plastic work done before initial state x1 coordinate along local crack-driving direction Xik nodal coordinates b parameter equal to unity for plane stress conditions e emissivity e0 initial strain e0ijk nodal values of initial strain for kth node ee, ep em de dee dep deth ds dT dr d1i DKeff

g gm @ gm/@xi

elastic strain and plastic strain mechanical strain which is the sum of elastic and plastic strains total strain increment elastic strain increment plastic strain increment thermal strain increment path length along C (Fig. 8) temperature increment stress increment Kronecker’s delta range of effective stress intensity factor arc efficiency factor natural coordinates inverse Jacobian matrix

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

r rij q [Dd]

C P

2507

Stefan–Boltzmann constant stress tensor density elastic–plastic material matrix path which surrounds the crack tip (Fig. 8) potential energy

method has been successfully used by several research works [16,17]. However, most weight functions were limited to simple geometries or require finite width corrections. Moreover, some weight functions are in very complicated forms that the calculation process involves solving complex integral equations. On the other hand, the finite element method is a more general and robust tool for complex geometries and loading conditions. It has been successfully used in the fatigue crack growth problem [18–21]. As mentioned above, recently, finite element simulations to predict welding temperature field, welding residual stress field and welding deformation have been employed more frequently, and finite element methods to simulate fatigue crack growth in a residual stress field have been developed. Nevertheless, very few works have been done on the development of an efficient and reliable simulation method in order to include the two aspects in the same package due to the complexity of the processes involved. Servetti and Zhang [20] compared the fatigue crack growth rates in welding residual stress field calculated from the different empirical fatigue crack growth laws such as the Walker equation, the Harter T-method and the NASGRO equation, which are based on the Paris’ law, using the finite element method. They have shown that the Harter T-method and the NASGRO equation give better predictions compared to the Walker equation. However, in their work, measured welding residual stresses were input into the finite element model without simulating the welding process. Barsom and his coworker [18,19] developed finite element models for the analysis of fatigue crack growth in welds incorporating the residual stresses from the welding simulation based on the linear elastic fracture mechanics, but these models were confined to two-dimensional problems. This study aimed to present a finite element modeling procedure for predicting fatigue crack growth rate in butt welds subject to mode I loading condition, and focused on the quantitative assessment of the influence of weld-induced residual stress on the fatigue crack growth rate. Three-dimensional thermal–mechanical finite element model was first developed in order to accurately predict the weld-induced residual stresses. Then, the rate of fatigue crack growth in the welds subject to the applied mechanical stress in conjunction with the residual stress was predicted by calculating the stress intensity factors resulting from the residual stress field using the modified J-integral definition [22], giving path independent J-integral values for mode I crack in a three-dimensional residual stress bearing body, and the applied mechanical stress field only employing the conventional J-integral definition, respectively, based on the superposition rule of the linear elastic fracture mechanics. The fatigue crack growth rate due only to the applied mechanical stress was also computed for comparison. The scope of this paper has been confined to prediction of the fatigue crack growth rate for the given crack configuration; hence the modeling of the crack extension is not practical concern of this work. 2. Finite element modeling of welding residual stress It is generally understood that welding involves highly complex phenomena, arising from the interactions between heat transfer, metallurgical transformation and mechanical fields. Complex numerical approaches are then needed to accurately model the welding process. However, as far as welding residual stress modeling is concerned, numerical procedures can be significantly simplified, as discussed in the forthcoming [8–13]. In this study, three-dimensional finite element modeling of welding residual stress was performed using the in-house finite element code [23], which has been extensively verified against experiments and numerical results found in the literature [24–27]. The process of welding is a coupled thermal–mechanical process. Since the thermal field has a strong influence on the stress field with little inverse influence, sequentially coupled analysis was carried out to starting with the transient heat transfer analysis. The heat transfer analysis is based on the heat conduction formulation with the moving heat source. The resulting temperature history solutions were then employed as the thermal loading for thermal stress evolution in the subsequent transient thermal–mechanical analysis. Welding residual stresses are the final state of the thermal stresses after all welding passes are over and the work piece is cooled down to the ambient temperature. Theoretical considerations of the heat transfer and thermal–mechanical analyses are briefly described below. 2.1. Heat transfer analysis The spatial and temporal temperature distribution T(x, y, z, t) during welding satisfies the following governing partial differential equation for the three-dimensional transient heat conduction with internal heat generation and considering q, K and c as functions of temperature only.

2508

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

      @ @T @ @T @ @T @T þ þ þ Q ¼ qc Kx Ky Kz @x @x @y @y @z @z @t

ð1Þ

The general solution is obtained by applying the following initial and boundary conditions:

Tðx; y; z; 0Þ ¼ T 0 ðx; y; zÞ

ð2Þ

  @T @T @T Nx þ K y Ny þ K z Nz þ q þ hðT  T 0 Þ þ reðT 4  T 40 Þ ¼ 0 Kx @x @y @z

ð3Þ

where r = 5.67  108 J/(m2 K4 s), e = 0.2 [3] and T0 = 20 °C. According to the nature of arc welding, the heat input to the work piece can be divided into two portions. One is the heat of the welding arc, and the other is that of the melt droplets. The heat of the welding arc is modeled by a surface heat source with a Gaussian distribution, and that of the melt droplets is modeled by a volumetric heat source with uniform density. At any time t, a point lying on the surface of the work piece within the r0 receive the Q(t) according to the following equation:

QðtÞ ¼

"  2 # 3Q A rðtÞ exp  r0 pr20

ð4Þ

where

Q A ¼ gIU  Q M

ð5Þ

The heat of the welding arc was assumed to be 40% of the total heat input, and the heat of the melt droplets 60% of the total heat input [28]. The arc efficiency factor is assumed as 0.85 for the flux cored arc welding process used in the present analysis. The heat flux was applied during the time variation that corresponded to the approach and passing of the welding torch. To take account of the heat transfer due to fluid flow in the weld pool, the thermal conductivity for molten metal was assumed to increase linearly between the solidus temperature and 3000 K by a factor of three as suggested in [3]. The liquid-to-solid phase transformation effects of the weld pool were modeled by taking into account the latent heat of fusion. The latent heat, the heat energy that the system stores and releases during the phase change, was accounted for assuming the value 270 J/g between the solidus temperature 1450 °C and the liquidus temperature 1500 °C [29]. However, the latent heat during solid-state phase transformation was not considered, since the influence of the metallurgical phase transformation [3] on welding residual stress is insignificant in the mild carbon steel used in this work [29]. In the thermal analysis, the process of sequential weld filler deposition was simulated using a consistent filler activation/ deactivation scheme. This scheme keeps track of the movement of the weld torch and updates the status of weld filler (deposited or not). For a weld filler that is not yet deposited at a given time, a value for thermal conductivity equivalent to that of air is assigned. This process is called filler deactivation. The deactivated fillers can be regarded as virtual fillers that are not actually present. After the weld filler is deposited, it is reactivated and the thermal conductivity is made to change from air value to that of the material used. Through such a consistent filler activation/deactivation process, the effects of sequential weld metal deposition can be properly modeled. 2.2. Thermal–mechanical analysis The second step of the residual stress analysis involves the use of the thermal histories predicted by the previous heat transfer analysis as an input (thermal loading) for the thermal–mechanical analysis. During the welding process, as mentioned earlier, solid-state phase transformation has a negligible effect on the evolution of residual stresses. Therefore, addictive strain decomposition can be used to decompose the differential form of the total strain into three components as follows:

fdeg ¼ fdee g þ fdep g þ fdeth g

ð6Þ

The elastic strain increment is calculated using the isotropic Hook’s law with temperature-dependent Young’s modulus and Poisson’s ratio. The thermal strain increment is computed using the coefficient of thermal expansion. For the plastic strain increment, a rate-independent elastic–plastic constitutive equation is considered with the Von Mises yield criterion, temperature-dependent mechanical properties and linear isotropic hardening rule. The incremental form of stress–strain relation can be written as

fdrg ¼ ½Dd fdeg  fct gdT

ð7Þ

where [Dd] is divided into ½Ded  for the elastic range and ½Dpd  for the plastic range. The filler activation/deactivation scheme was also used in the thermal–mechanical analysis. The mechanical aspect of the scheme is to change the stiffness of the weld filler. For the weld filler which the welding torch has not yet approached, a severely reduced material stiffness is assigned [30]. At the time of application of weld filler deposition, it is reactivated and the temperature-dependent mechanical properties of the material are assigned with no record of strain history to bring it into existence without incurring strain incompatibilities.

2509

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

Since the thermal elastic–plastic analysis is a nonlinear problem, the incremental calculation technique was employed in solving the problem. The full Newton–Raphson iterative solution technique [31] was used for obtaining the solution. 2.3. Verification In order to confirm the accuracy of the finite element analysis method described in the previous sections, a specimen was constructed using double ‘V’ butt joint welding, with a length, width and thickness of L = 600 mm, W = 600 mm, t = 20 mm, respectively. Details on the materials used and the preparation of the specimen are given in [32]. Residual stress measurements were carried out on the two axis strain gauges with the layering technique. Measuring stresses by strain gauges using the layering technique can obtain the residual stresses on the surface of the structure to be evaluated. The detailed procedure for measuring stresses can be found in [24]. Three-dimensional thermo-mechanical finite element analysis for the experimental model using the finite element method was also performed employing the same welding conditions and process parameters as those used in the fabrication of the specimen. The modeling procedure is similar to that given in the following section. Fig. 1a and b portrays the longitudinal (rx) and transverse (ry) residual stresses at a cross-section of the weld piece perpendicular to the welding direction, respectively. The symbols represent the results of the experimental value measured by the strain gauges, the solid curve representing the results of the finite element analysis calculated along the experimental measurement location. The results are normalized by the respective yield stress (rY) of the weld metal and base metal. It can be seen that the residual stress distributions computed by the finite element analysis show very good agreement with those determined by the experiment. Therefore, the finite element method used here is considered appropriate for analyzing the residual stresses induced by the welding process. 2.4. Finite element model Finite element simulation of a butt joint welding process was conducted using the sequentially coupled heat transfer and thermal–mechanical analyses. Two 400 mm  200 mm  6 mm plates with a single V-groove joint between them, as shown in Fig. 2a, were assumed to be welded by single-pass. The welding parameters chosen for this analysis were as follows:

Longitudinal stress (σx/σY)

(a)

Analysis Experiment

1.2 1 0.8 0.6 0.4 0.2 0 -0.2 -0.4 -0.6

0

50

100

150

200

250

300

Distance (mm)

Transverse stress (σy/σY)

(b)

Analysis Experiment

0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1 -0.2

0

50

100

150

200

250

300

Distance (mm) Fig. 1. Comparison of the predicted with the measured residual stresses: (a) longitudinal residual stress distribution and (b) transverse residual stress distribution.

2510

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

(a) Not to scale

(b)

Z

X

Y

Fig. 2. Analysis model: (a) dimensions of welded plate and shape of groove and (b) finite element model.

welding method, flux cored arc welding process; welding current, 240 A; welding voltage, 30 V; and welding speed, 5 mm/s, respectively. The three-dimensional finite element mesh model used in the simulation is shown in Fig. 2b with eight-noded isoparametric solid elements. Four layers are used to discretize the computational domain. Note that the finite element model only represents one half of the weld piece due to the symmetry with respect to the weld centerline. In order to facilitate data mapping between the heat transfer and thermal–mechanical models, the same finite element mesh refinement was used except for the element type and applied boundary conditions. For the heat transfer model, the element type is one which has single degree of freedom, temperature, on its each node. For the thermal–mechanical model, the element type is the other with three translational degrees of freedom at each node. A fine grid was adopted in the welds and its vicinity in order to apply heat flux more accurately when the moving heat source passed the welds at specific time steps. Element size increases progressively with distance from the weld centerline. Mesh sensitivity study has been conducted to examine the dependence of finite element mesh size on the accuracy of the analysis results. As a result, the present finite element mesh with the smallest element size of 0.5 mm (transverse)  1.5 mm (thickness)  10.0 mm (longitudinal) was considered to give sufficiently accurate results using a reasonable amount of computer time and memory. As the plate was not clamped during welding, no boundary conditions except those to prevent rigid body motion of the weld piece were applied. SM400, one of the most widely used structural steels in Korea, has been used as the material of choice. The chemical composition of the material is presented in Table 1. In the finite element simulation, temperature-dependent thermo-physical and mechanical properties of the material were incorporated. Fig. 3 shows the temperature-dependent thermo-physical properties such as thermal conductivity, specific heat and density [26]. It should be noted that units in Fig. 3 are organized so that they can be shown on one graph for clarity. Temperature-dependent thermo-mechanical properties such as Young’s modulus, thermal expansion coefficient, Poisson’s ratio and yield stress are given in Fig. 4 [26]. Both the yield stress and elastic modulus are reduced to 5.0 MPa and 5.0 GPa, respectively, at the melting temperature to simulate low strength at high temperatures [18]. For the weld metal, autogenous weldment was assumed, i.e. the temperature-dependent thermophysical and mechanical properties of the weld metal are the same as those of the base meal [33].

Table 1 Chemical composition of the material used (mass%). Material

C

Si

Mn

P

S

Ni

Cr

Mo

Cu

V

Al

N

SM400

0.14

0.045

1.14

0.014

0.022

0.01

0.01

0.002

0.04

0.002

0.007

0.0049

2511

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

10 Density -3 3 ρ (g/mm ) x 10

ρ, χ, λ

8

Specific heat c ( / g · ) x 10 -1

6

4 Heat conductivity -2 λ ( J / mm · s · ) x 10

2

0

0

200

400

600

800

1000

1200

1400

Temperature ( ) Fig. 3. Temperature-dependent thermo-physical constants of the material [26].

Linear expansion coeff. α=1.2×10-5(1/ Poisson's ratio ν=0.3

400

)

300 150 200

σy

100

E

100

0

50

0

200

400

600

800

1000

1200

1400

Young's modulus, E (GPa)

Yield stress, σy (MPa)

200

0

Temperature ( ) Fig. 4. Temperature-dependent thermo-mechanical properties of the material [26].

The work hardening behavior induced by the thermal cyclic loading during welding in and around the welds should be carefully taken into account when a numerical method is utilized to accurately predict welding residual stresses. In this study, a linear strain hardening of the base metal and the weld metal was assumed with the rate of 500 MPa for the temperature range 20–700 °C and 20 MPa for the temperatures above 1000 °C [3]. A linear transition between the hardening rate at 700 °C and that at 1000 °C was assumed. 2.5. Results and discussion Results are first prepared for the transient thermal cycles associated with the welding. Temperature histories are predicted for the finite element model at the welds. Fig. 5 shows the thermal cycles at different locations remote from the welding start position. The curve with hollow marks represents the temperature history at the location which is 100 mm away from the start position (1/4 position), the broken curve represents the temperature history at the location 200 mm away from the start position (2/4 position), and curve without any mark represents the temperature history at the location 300 mm away from the start position (3/4 position). From the figure, it can be observed that the maximum temperature at the weld pool is 2000 °C or so. This result corresponds well with the welding process in practice. Furthermore, it is clear that the temperature histories at the three locations are almost identical. Therefore, it can be concluded that the temperature field is steady when the welding torch moves along the joint. Results are next presented for the weld-induced residual stresses. From the stress analysis all the stress and strain components can be obtainable throughout the plate. Here, we will discuss only the relevant data. In this discussion, the profiles of longitudinal (acting parallel to the weld line – rx), transverse (acting perpendicular to the weld line, i.e. acting along the plate width – ry) and through-thickness (acting normal to the weld line, i.e. acting along the plate thickness – rz) residual stresses are predicted. Fig. 6a and b shows the longitudinal, transverse and through-thickness residual stress profiles at a cross-section of the weld piece perpendicular to the weld line and half way through the weld length. The stress profiles predicted are within the top layer of 1.5 mm in thickness and within the bottom layer of 1.5 mm in thickness, respectively, and plotted as a function of distance from the weld centerline. From the simulated results, it can be shown that the maximum longitudinal

2512

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520 1/4 position 2/4 position 3/4 position

2400

Temperature (ºC)

2000 1600 1200 800 400 0

0

50

100

150

200

250

300

Time (s) Fig. 5. Thermal cycles at several locations remote from the welding start position during the welding.

(a) 600

Transverse Longitudinal Through-thickness

Residual stress (MPa)

500 400 300 200 100 0 -100 -200

0

25

50

75

100

125

150

175

200

Distance from the weld centerline (mm)

(b) 600

Transverse Longitudinal Through-thickness

Residual Stress (MPa)

500 400 300 200 100 0 -100 -200

0

25

50

75

100

125

150

175

200

Distance from the weld centreline (mm) Fig. 6. Residual stress profiles at a cross-section of the weld piece perpendicular to the weld line: (a) top surface and (b) bottom surface.

tensile stress is imposed at the welds due to the resistance of the material against contraction attributed to cooling after welding. The longitudinal stresses level out in compression away from the welds for self-equilibrating purpose. On the other hand, the transverse stresses are wholly tensile and level out to zero. Like the longitudinal stress components, the highest tensile stress arises at the welds. The through-thickness stress components have fluctuating profiles that vary between tensile and compressive. It is worth noting that the longitudinal tensile residual stresses in and around the welds exceed the room temperature yield stress, which is also noted in experimental work such as [34]. The strain hardening could be considered to be the main factor to explain the higher residual stresses [26]. The results also reveal that the residual stress distributions are similar through the thickness, i.e. notable through-thickness variation of the residual stresses is not produced for the given finite element model. This is mostly attributed to the fact that for the case of thin-walled plate, the heat deposited during butt welding is high enough to result in a uniform temperature increase through the thickness at the welds;

2513

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520 Transverse Longitudinal Through-thickness

600

Residual stress (MPa)

500 400 300 200 100 0 -100 -200 -300 -400

0

50

100

150

200

250

300

350

400

Distance along the weld length (mm) Fig. 7. Residual stress profiles at a cross-section of the weld piece along the weld length.

hence the longitudinal shrinkage during cooling after welding, which creates the residual stresses, is almost the same through the thickness. Fig. 7 depicts the residual stress profiles at a cross section of the weld piece along the weld length beside the weld centerline for depth of 1.5 mm from the top surface. It is shown that the stress distributions are symmetrical with respect to the mid-length plane except for the region near the welding start and end points of which stresses have no influence on the further calculation of the J-integral values for the through-thickness center crack in welds [35]. It should be noted that the longitudinal residual stresses are primarily dominated by work piece’s restraint in the longitudinal direction, and their typical distribution along the weld length is given in the results. The transverse residual stress distribution, however, can be attributed to both longitudinal and transverse restraints [36]. The results demonstrate the contribution of the longitudinal shrinkage to the variation of the transverse residual stresses along the weld length. 3. Finite element analysis of fatigue crack growth rate in butt welds As described earlier, in the analysis of fatigue crack growth rate in welding residual stress field, the key task is to determine the stress intensity factor due to the weld-induced residual stresses. There are a few approaches for evaluating crack tip stress intensity factor by finite element method, such as the crack tip displacement extrapolation [37], the modified virtual crack closure technique [20,38] and the J-integral [39]. In the displacement extrapolation method, the displacement and stress data of the finite element analysis at the vicinity of the crack tip are assumed to obey their asymptotic behavior and the crack tip parameters are computed. This method is simple but it is not easy to guarantee the accuracy. The modified virtual crack closure technique (also known as virtual crack closure integral method) is another way to compute the energy release rate. In this method, energy that is required to close the crack for one finite element length is calculated by multiplying the nodal reaction force perpendicular to the crack growth path at the crack tip node and the opening displacement at the node immediately behind the crack tip. The advantages of the technique are in its simplicities in computing the energy release rate such that only the nodal reaction force and opening displacement are used in the computation and that using their appropriate vector components, the energy release rate is decomposed into the modes I–III components [40]. The Jintegral approach can also obtain the energy release rate through the domain integral method. The domain integral method utilizes the subroutines in the finite element program to perform integrations based on finite element geometries. Therefore, the domain integral approach is well suited for fracture mechanics analysis using the finite element method. These methods are widely used practices and easily implemented in commercial finite element packages. In this study, the J-integral technique was adopted to calculate the stress intensity factor resulting from the residual stresses. However, when utilizing finite element technique using the J-integral to evaluate the energy release rate for a crack containing residual stresses, significant problems arise, i.e. for purely mechanical loading, the values of the J-integral determined by the domain integral technique are path independent (domain independent) outside the crack tip, but the J-integral for a crack in a residual stress field no longer possesses the path independency [35]. Some researchers have reported path dependent values of the J-integral when analyzing a crack containing residual stresses [41,42]. These works demonstrated the need to determine path independent Jintegral values when dealing with a crack which is subject to residual stresses. This study employed the modified J-integral definition [22] which yields path independent values for mode I crack in a three-dimensional residual stress field in order to obtain the stress intensity factor due to the weld-induced residual stresses. 3.1. Path independent J-integral in a three-dimensional residual stress field The J-integral has emerged over the last years as one of the leading parameters to characterize crack-propagation in solids. The J-integral derived by Rice [39] was proposed as the strain energy release rate at the crack tip in linear or nonlinear solid and written as

2514

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

x2

ni ds x1

0

Fig. 8. Typical contour for evaluation of the J-integral.

X2

ni X1

S1

S2 V Vt

X3 Fig. 9. Volume V formed by S1 and S2.

J¼

Q Z   d @uj nids Wd1i  rij ¼ @x1 da C

ð8Þ

where

W¼

Z

em ij

0

rijdemij

ð9Þ

Although the J-integral is a potential parameter to characterize the crack tip driving force of cracked bodies made of elastic and elastic–plastic materials, it is difficult to compute the unit normal vector to the line contour in a finite element analysis, and path independence is not valid in the vicinity of the crack tip undoubtedly due to the numerical solution being unable to model the singular behavior [43]. Therefore, alternative forms of computing this integral have been proposed [44]. Through the use of divergence theorem, Rice’s line integral can be transformed to an equivalent area or domain integral in two dimensions. In three dimensions, an analogous volume integral expression which is suited for evaluation in three-dimensional finite element analysis has been developed by Li et al. [45] using the virtual crack extension method [46] and Eshelby’s energy momentum tensor:

J¼

Q Z   @qj d @u rik i  Wdjk dV ¼ @x @x da j k V

ð10Þ

where qj is non-zero within V and zero on the boundary of V [47]. The J-integral in Eq. (10) can be reduced for the case in which the crack is extended in the direction of its plane as follow [48]:

J¼

Z  V

rik

 @ui @q1  Wd1k dV @x1 @xk

ð11Þ

For a purely mechanical load, in the absence of a body force and crack face traction, the J-integral is path independent and can be estimated from fields remote from the crack tip. But the J-integral leads to path dependent values in the presence of residual stress [35]. Therefore, it is necessary to introduce a modified path independent J-integral definition considering residual stress. Appling the divergence theorem to Eq. (8) and generalizing the results for three dimensions subject to in-plane displacement yields

ðJÞvolume ¼  where

  Z  Z  @q @ui @q1 @W @ eij dV  q dV W 1  rij  rij @x1 @x1 @x1 @xj @x1 1 V V

ð12Þ

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520 m @ em @W @W @ eij ij ¼ ¼ rij @x1 @ em @x @x1 1 ij

2515

ð13Þ

A residual stress field in welds is a self balancing stress distribution caused by internal plastic strains due to non-uniform heat expansion. Such residual stress problems can then be treated as initial strain problem [35]. When initial strain is absent, the second term in the integral of Eq. (12) is zero. Otherwise, em ij is then written as

emij ¼ eeij þ epij  e0ij

ð14Þ

Inserting Eqs. (13) and (14) into Eq. (12) leads to a path independent J-integral definition considering residual stress and initial strain as follow:

ðJÞvolume ¼ 

 Z Z  @q @ui @q1 dV þ W 1  rij @x1 @x1 @xj V V

rij

@ e0ij @x1

q1 dV

ð15Þ

Here W must be adjusted to account for the initial strain energy density, i.e.

W ¼ W total  W p jinitial state

ð16Þ

In order to incorporate initial plastic deformation, Lei et al. [35] suggested the definition of the ‘initial state’ of the material. At this initial state, the process simulations that give the desired residual stress distribution are completed, but boundary conditions have not yet been changed to introduce the crack into the material. Hereafter, the value of the J-integral determined using the Eq. (15) will be referred as Jmod, whereas that calculated by the standard J-integral definition will be labeled Junmod. It should be noted that the modified J-integral is equivalent to the stress intensity factor, K, under small scale yielding conditions and provides the intensity of the crack tip stress field for elastic–plastic conditions [35]. 3.2. Finite element procedure For the chosen eight-noded isoparametric solid element (1 6 gi 6 1, i = 1–3), the coordinates and the displacements of an arbitrary point referred to a global reference system are given by

xi ¼

8 X

Nk X ik

k¼1

ui ¼

8 X

Nk U ik

i ¼ 1—3

ð17Þ

k¼1

The expression of q1 is written in a similar way,

q1 ¼

8 X

Nk Q 1k

ð18Þ

k¼1

where Q1k = 1 if the kth node is on S1 and Q1k = 0 if the kth node is on S2 [49] and the gradient of q1 is given by 8 X 3 @q1 X @N k @ gm ¼ Q @xi @ gm @xi 1k k¼1 m¼1

ð19Þ

In a similar way, the initial stain in an element can be expressed as

e0ij ¼

8 X

Nk e0ijk

ð20Þ

k¼1

The gradient of

@ e0ij @xi

¼

e0ij is given by

8 X 3 X @Nk @ gm 0 e @ gm @xi ijk k¼1 m¼1

ð21Þ

Employing Gaussian integration, the discretized form for the volume integral becomes

J ðvolumeÞ ¼

N X M X P X R X n¼1 m¼1 p¼1 r¼1

"(

) #  @ e0ij @uj @q1 rij  Wd1i þ rij q jJj @x1 @xi @x1 1

wm wp wr

ð22Þ

nm gp fr

3.3. Influence of weld-induced residual stresses on the fatigue crack growth rate A quantitative assessment of the weld-induced residual stress effect on the fatigue crack growth rate can be made by the principle of superposition of the linear elastic fracture mechanics based on the Keff as

K eff ¼ K applied þ K residual

ð23Þ

2516

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

The rate of fatigue crack growth can then be described in the Paris’ equation, which is employed due to the availability of the material constants, in terms of the DKeff given by

da=dN ¼ CðDK eff Þm

ð24Þ

Under cyclic loading, the DKeff and the Reff are calculated as

    min max min DK eff ¼ K max applied þ K residual  K applied þ K residual ¼ K applied  K applied ¼ DK applied

Reff

K min eff ¼ max ¼ K eff

K min applied þ K residual

ð25Þ

!

K max applied þ K residual

ð26Þ

As shown in the equations, the range of effective stress intensity factor does not change when the residual stresses are superimposed, whereas the effective stress ratio will continually change as the crack propagates through the superimposed residual stress field for a constant stress ratio and hence will cause a mean stress alteration at the crack tip. Therefore, the analysis of fatigue crack growth in a residual stress field should incorporate the effective stress ratio in addition to the stress range [19]:

da CðDK eff Þm ¼ dN ð1  Reff Þ

ð27Þ

It is obvious from Eqs. (25)–(27) that the fatigue crack growth rate in welds can be evaluated if the solutions for the stress intensity factor Kapplied and Kresidual are known, i.e. the total stress intensity factor for the combined residual and applied mechanical stresses does not need to be calculated. 3.4. Prediction of fatigue crack growth rate in butt welds under mode I loading Three-dimensional elastic–plastic finite element analysis has been performed to calculate the rate of fatigue crack growth in butt welds under mode I loading condition in accordance with Eq. (27). The same geometry and material as in the residual stress analysis were used, except for the applied boundary conditions. As the previous analysis, an elastic–plastic constitutive material model with the Von Mises yield criterion and linear isotropic hardening behavior was employed. The material properties at room temperature given in Fig. 4 were used for the finite element analysis. The material constants of the Paris’ equation was assumed as m = 3, C = 3  1013 mm/cycle as recommended by Maddox [4] for wide range of structural steels. Due to the symmetry, only a quarter of the welded plate was modeled [22] employing three-dimensional eight-noded isoparametric solid elements with in-plane mesh design shown in Fig. 10. The finite element mesh used here has more refinement near the symmetry plane compared to that used in the residual stress analysis, but otherwise has the same features, i.e. four layers of elements were used to compute the domain integral through the crack front. The element width, i.e. the path size used in calculation of the J-integral values, near the crack tip region is determined to be 0.1 mm by a mesh convergence study. Two different element size (0.2 and 0.1 mm) of mesh refinement has been considered, all containing a crack length (a) of 2.5 mm. Remote mode I tension loading has been applied, and Fig. 11 shows the J-integral values due to the mechanical tensile stress of 80 MPa with no residual stresses under the linear elastic condition versus the domain size measured directly ahead of the crack tip. Detailed description of the finite element analysis conditions will be given later in this paper. The J-integral values under mechanical loading only have been determined by the standard J-integral definition. Superimposed are the results of the weight function approach [50] from which stress intensity factor can be calculated using the linear elastic fracture mechanics solution. Confidence in the present finite element analysis is gained from the path independence of the J-integral values, and the result for element size of 0.1 mm is founded to be in close agreement (within 3%) with that of the linear elastic fracture mechanics solution. The residual stresses and plastic strains existing after the finite element simulation of the welding process detailed in Section 2.5 were introduced as initial conditions, or mapped, into the finite element model. After an initial equilibrium step (no additional loading applied) required due to the discrepancy of the mesh density between the two finite element models was invoked, a center through-thickness fatigue crack in welds perpendicular to the weld line was introduced into the model by releasing the coincident crack surface nodes, which allows for mode I loading condition. The fatigue crack length is 2.5 mm which is considered as representative of a conservative yet reasonable weld flaw acceptance criteria [51]. The release of nodes was performed such that the crack had a parallel geometry through the thickness. Indeed, the presence of similar residual stresses through the plate thickness as explained in Section 2.5 suggests that an actual crack will be likely to grow in a parallel manner through the thickness. The fatigue crack was assumed to grow in a direction normal to the applied loading giving pure mode I loading condition [42]. As stated before, for the evaluation of the fatigue crack growth rate in welds, the stress intensity factors Kapplied and Kresidual should be known. As small scale yielding conditions hold for the case considered here, which was satisfied with the plastic zone size (the size of yielded elements) near the crack tip that did not exceed 0.2 times the width of the analysis model [52] under the

2517

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

200 mm

a = 2.5 mm Fig. 10. Finite element mesh with crack-tip detail.

LEFM J_0.1mm J_0.2mm

1 0.9 0.8

J (N/mm)

0.7 0.6 0.5 0.4 0.3 0.2 0.1 0

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

Distance from the crack tip (mm) Fig. 11. Comparison of the J-integral values for the crack subject to mechanical loading only under the linear elastic fracture mechanics assumption as a function of the element size.

applied mechanical stresses superimposed on the residual stress field, the stress intensity factors can be calculated via the relation.

sffiffiffiffiffi EJ K¼ b

ð28Þ

2518

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

40

J_unmod J_mod

35

J (N/mm)

30 25 20 15 10 5 0

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

Distance from the crack tip (l/a) Fig. 12. Comparison of Jmod and Junmod for the residual stress field.

Therefore, the respective J-integral values for the residual stress field and the applied mechanical loading alone were first computed using Eqs. (15) and (11), respectively. The initial state of Eqs. (14)–(16) was chosen as the post welding condition, i.e. the state after the residual stresses and plastic strains are input and the equilibrium step is performed, prior to introduction of the crack. The J-integral values are computed for various paths by the procedure described in this work. Comparison of Jmod and Junmod for the residual stress field are given in Fig. 12, which are plotted as a function of distance from the crack tip (l) normalized by the crack length. Average values of the J-integral along the entire crack front are reported since the present finite element model produces similar results through the crack front. It can be observed that the J-integral values obtained from the modified J-integral definition show good path independence except for those calculated on the paths near the crack tip at which high stress gradients exist. On the other hand, due to the presence of residual stress, the path dependence of the J-integral values computed by the standard J-integral definition can be clearly seen in the figure. Fig. 13 shows the J-integral values for the different levels of mechanical loads only which correspond to the stress ranges considered here under R = 0. Prescribed displacements were applied to the finite element model which contains no residual stresses along one edge to generate applied mechanical mode I loads. All nodes at the edge are constrained to move together in loading direction. The loading of the model was carried out in increments up to the mechanical tensile stress of 160 MPa. As before, the average J-integral values through the crack front are reported. It is clear that the J-integral values are path independent for all loading levels and the path dependence near the crack tip does not occur due to the relatively low stress gradients compared to the residual stresses. Moreover, it can be found that the J-integral value for the mechanical loading of 80 MPa is larger than that calculated from the linear elastic fracture mechanics (see Fig. 11) owing to the yielding of the material at the crack tip region. After conversion of the J-integral values into the respective stress intensity factors, the fatigue crack growth rate was calculated using Eq. (27). Fig. 14 portrays the predicted fatigue crack growth rates in the welds for the constant stress ratio (R = 0) versus the effective stress intensity factor ranges corresponding to the applied stress ranges. For comparison, the fatigue crack growth rates due only to the applied mechanical stress ranges are also given in the figure. It is immediately apparent that the rate of fatigue crack growth for the case with the residual stresses is higher than that for the case without the residual stresses. This is because the longitudinal tensile residual stresses in the welds which, in this case, act as the crack

80MPa 100MPa 120MPa 140MPa 160MPa

10 9 8

J (N/mm)

7 6 5 4 3 2 1 0

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

Distance from the crack tip (l/a) Fig. 13. The J-integral values for the different levels of mechanical loads only.

2519

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520

50

Applied

da/dN (m/cycle)×E-12

Residual+Applied

40 30 20 10 0 10

15

20

25

30

35

40

K (MPa m) Fig. 14. Predicted fatigue crack growth rates with respect to the effective stress intensity factor ranges under the constant stress ratio.

driving force accelerate the rate of fatigue crack growth in conjunction with the applied mechanical tensile stresses. The figure also shows that the difference of fatigue crack growth rate between the two cases is increased drastically with increasing the stress range applied, i.e. the fatigue crack growth rate due to the combined residual stress and applied stress range of 160 MPa is increased over three times that due only to the applied stress range. This implies that for the given geometry and applied stress ranges, the residual stresses could result in a reduction in the fatigue life of the structure up to 70% compared to that of the cracked body without the residual stresses. 4. Conclusions This study has described a finite element modeling procedure capable of predicting fatigue crack growth rate in butt welds subject to mode I loading condition. For the analysis of a crack in welds, residual stress analysis and fracture mechanics analysis must be implemented consecutively. Therefore, first, sequentially coupled three-dimensional thermal–mechanical finite element model has been developed to simulate the residual stress states in butt welds of thin-walled steel plate. In the finite element model, temperature-dependent material properties, work hardening behavior of the material welded, and weld filler variation with time were taken into account. After the simulation, a center through-thickness fatigue crack in welds perpendicular to the weld line was introduced into the model, which satisfied the mode I condition. The influence of weld-induced residual stresses on the fatigue crack growth rate was then explored by calculating the stress intensity factors due to the residual stress field using the modified J-integral definition and the applied mechanical stress field only employing the conventional J-integral definition, respectively, based on the superposition rule of the linear elastic fracture mechanics. Numerical results have shown that notable through-thickness variation of the residual stresses is not produced for the thin-walled plate model. Therefore, the J-integral values for a crack in the residual stress field along the crack front are similar; hence the average J-integral value through the thickness can be employed to calculate the fatigue crack growth rate. However, the through-thickness variation of the residual stresses in welds of thick-walled plate is considerable; thus influencing the J-integral values for a crack in the welds along the crack front significantly, to which the finite element procedure for predicting fatigue crack growth rate in welds presented in this work can also be applicable, at which the future works are aimed. Moreover, the results have demonstrated that welding residual stresses severely affect the fatigue crack growth rate in the welds; thus substantially decreasing the fatigue life of the structure, i.e. for the given geometry and applied stress ranges, the residual stresses could result in a reduction in the fatigue life of the structure up to 70% compared to that of the cracked body without the residual stresses. References [1] Hibbitt HD, Marcal PV. A numerical thermo-mechanical model for the welding and subsequent loading of fabrication structure. Comput Struct 1973;3:1145–74. [2] Karlsson RI, Josefson BL. Three-dimensional finite element analysis of temperatures and stresses in a single-pass butt-welded pipe. J Press Vessel Technol 1990;112:76–84. [3] Taljat B, Radhakrishnan B, Zacharia T. Numerical analysis of GTA welding process with emphasis on post-solidification phase transformation effects on the residual stresses. Mater Sci Engng A 1998;246:45–54. [4] Maddox SJ. Fatigue strength of welded structures. Abington Publishing; 1991. [5] Withers PJ. Residual stress and its role in failure. Rep Prog Phys 2007;70:2211–64. [6] Bussu G, Irving PE. The role of residual stress and heat affected zone properties on fatigue crack propagation in friction stir welded 2024-T351 aluminium joints. Int J Fatigue 2003;25:77–88. [7] Ghidini T, Dalle DC. Fatigue crack propagation assessment based on residual stresses obtained through cut-compliance technique. Fatigue Fract Engng Mater Struct 2006;30:214–22. [8] Lindgren L-E. Finite element modelling and simulation of welding – Part 1. Increased complexity. J Therm Stresses 2001;24:141–92. [9] Lindgren L-E. Finite element modelling and simulation of welding – Part 2. Improved material modeling. J Therm Stresses 2001;24:195–231.

2520 [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52]

C.-H. Lee, K.-H. Chang / Engineering Fracture Mechanics 78 (2011) 2505–2520 Lindgren L-E. Finite element modelling and simulation of welding – Part 3. Efficiency and integration. J Therm Stresses 2001;24:305–34. Goldak J, Akhlagi M. Computational welding mechanics. Springer; 2005. Lindgren L-E. Numerical modelling of welding. Comput Methods Appl Mech Eng 2006;195:6710–36. Lindgren L-E. Computational welding mechanics (thermomechanical and microstructural simulations). Woodhead Publishing; 2007. Itoh YZ, Suruga S, Kashiwaya H. Prediction of fatigue crack growth rate in welding residual stress field. Engng Fract Mech 1989;33:397–407. Anderson TL. Fracture mechanics fundamentals and applications. CRC Press; 1995. Glinka G, Shen G. Universal features of weight functions for cracks in mode I. Engng Fract Mech 1991;40:1135–46. Pouget G, Reynolds AP. Residual stress and microstructure effects on fatigue crack growth in AA2025 friction stir welds. Int J Fatigue 2008;30:463–72. Barsoum Z. Residual stress analysis and fatigue of multi-pass welded tubular structures. Engng Fail Anal 2008;15:863–74. Barsoum Z, Barsoum I. Residual stress effects on fatigue life of welded structures using LEFM. Engng Fail Anal 2009;16:449–67. Servetti G, Zhang X. Predicting fatigue crack growth rate in a welded butt joint: the role of effective R ratio in accounting for residual stress effect. Engng Fract Mech 2009;76:1589–602. Labeas G, Diamantakos I. Numerical investigation of through crack behaviour under welding residual stresses. Engng Fract Mech 2009;76:1691–702. Chang KH, Lee CH. Residual stresses and fracture mechanics analysis of a crack in welds of high strength steels. Engng Fract Mech 2007;74:980–94. Lee CH. A study on the mechanical characteristics of high strength steel for the application to the steel bridge. Ph.D. Thesis, Chung-Ang University, Korea; 2005. Lee CH, Chang KH. Numerical investigation of the residual stresses in strength-mismatched dissimilar steel butt welds. J Strain Anal Eng 2008;43:55–66. Lee CH, Chang KH. Three-dimensional finite element simulation of residual stresses in circumferential welds of steel pipe including pipe diameter effects. Mater Sci Engng A 2008;487:210–8. Lee CH, Chang KH, Lee CY. Comparative study of welding residual stresses in carbon and stainless steel butt welds. Proc Inst Mech Engng – Part B: J Engng Manuf 2008;222:1685–94. Chang KH, Lee CH. Finite element analysis of the residual stresses in T-joint fillet welds made of similar and dissimilar steels. Int J Adv Manuf Technol 2009;41:250–8. Pardo E, Weckman DC. Prediction of weld pool and reinforcement dimensions of GMA welds using a finite element model. Metall Mater Trans B 1989;20:937–47. Deng D, Liang W, Murakawa H. Determination of welding deformation in fillet-welded joint by means of numerical simulation and comparison with experimental measurements. J Mater Process Technol 2007;183:219–25. Lindgren L-E, Hedblom E. Modelling of addition of filler material in large deformation analysis of multipass welding. Commun Numer Methods Engng 2001;17:647–57. Bathe KJ. Finite element procedures. Prentice Hall; 1996. Lee CH, Chang KH, Park HC, Lee JH. Characteristics of mechanical properties at elevated temperatures and residual stresses in welded joint of SM570TMC steel. J KSSC 2006;18:395–403. Teng TL, Chang PH, Tseng WC. Effect of welding sequence on residual stresses. Comput Struct 2003;81:273–86. Paradowska A, Price JWH, Ibrahim R, Finlayson T. A neutron diffraction study of residual stress due to welding. J Mater Process Technol 2005;164– 165:1099–105. Lei Y, O’dowd NP, Webster GA. Fracture mechanics analysis of a crack in a residual stress field. Int J Fract 2000;106:195–216. Dong P, Zhang J. Residual stresses in strength-mismatched welds and implications on fracture behavior. Engng Fract Mech 1999;64:485–505. Bao R, Zhang X, Yahaya NA. Evaluating stress intensity factors due to weld residual stresses by the weight function and finite element methods. Engng Fract Mech 2010;77:2550–66. Rybicki EF, Kanninen MF. A finite element calculation of stress intensity factors by a modified crack closure integral. Engng Fract Mech 1977;9:931–8. Rice JR. A path independent integral and the approximations analysis of strain concentration by notches and crack. J Appl Mech 1968;35:379–86. Okada H, Kawai H, Araki K. A virtual crack closure-integral method (VCCM) to compute the energy release rates and stress intensity factors based on quadratic tetrahedral finite elements. Engng Fract Mech 2008;75:4466–85. Hou YC, Pan J. A fracture parameter for welded structures with residual stresses. Comput Mech 1998;22:281–8. Pavier MJ, Poussard CGC, Smith DJ. Effect of residual stress around cold worked holes on fracture under superimposed mechanical load. Engng Fract Mech 1999;63:751–73. Woo KS, Hong CH, Shin YS. An extended equivalent domain integral method for mixed mode fracture problems by the p-version of FEM. Int J Numer Methods Engng 1998;42:857–84. Nikishkov GP, Atluri SN. An equivalent domain integral method for computing crack-tip integral parameters in non-elastic, thermal mechanical fracture. Engng Fract Mech 1987;26:851–67. Li FZ, Shih CF, Needleman A. A comparison of methods for calculating energy release rates. Engng Fract Mech 1985;21:405–21. Parks DM. The virtual crack extension method for nonlinear material behavior. Comput Methods Appl Mech Eng 1977;12:353–64. Shih CF, Moran B, Nakamura T. Energy release rate along a three-dimensional crack front in a thermally stressed body. Int J Fract 1986;30:79–102. Haddi A, Weichert D. On the computation of the J-integral for three-dimensional geometries in inhomogeneous materials. Comput Mater Sci 1996;5:143–50. Delorenzi HG. Energy release rate calculations by the finite element method. Engng Fract Mech 1985;21:129–43. Wang QZ. Some simple Mode-I SIF expressions of finite width strip with a center crack derived by using an approximate weight function. Engng Fract Mech 1998;60:37–45. Chi WM. Prediction of steel connection failure using computational fracture mechanics. Ph.D. Thesis, Stanford University, US; 1999. Ren XB, Zhang ZL, Nyhus B. Effect of residual stresses on the crack-tip constraint in a modified boundary layer model. Int J Solids Struct 2009;46:2629–41.