In fracture mechanics, the energy release rate, G {\displaystyle G} , is the rate at which energy is transformed as a material undergoes fracture. Mathematically, the energy release rate is expressed as the decrease in total potential energy per increase in fracture surface area, and is thus expressed in terms of energy per unit area. Various energy balances can be constructed relating the energy released during fracture to the energy of the resulting new surface, as well as other dissipative processes such as plasticity and heat generation. The energy release rate is central to the field of fracture mechanics when solving problems and estimating material properties related to fracture and fatigue.
Definition
The energy release rate G {\displaystyle G} is defined as the instantaneous loss of total potential energy Π {\displaystyle \Pi } per unit crack growth area s {\displaystyle s} ,
G ≡ − ∂ Π ∂ s , {\displaystyle G\equiv -{\frac {\partial \Pi }{\partial s}},}
where the total potential energy is written in terms of the total strain energy Ω {\displaystyle \Omega } , surface traction t {\displaystyle \mathbf {t} } , displacement u {\displaystyle \mathbf {u} } , and body force b {\displaystyle \mathbf {b} } by
Π = Ω − { ∫ S t t ⋅ u d S + ∫ V b ⋅ u d V } . {\displaystyle \Pi =\Omega -\left\{\int _{{\mathcal {S}}_{t}}\mathbf {t} \cdot \mathbf {u} \,dS+\int _{\mathcal {V}}\mathbf {b} \cdot \mathbf {u} \,dV\right\}.}
The first integral is over the surface S t {\displaystyle S_{t}} of the material, and the second is over its volume V {\displaystyle V} . The figure on the right shows the plot of an external force P {\displaystyle P} vs. the load-point displacement q {\displaystyle q} , in which the area under the curve is the strain energy. The white area between the curve and the P {\displaystyle P} -axis is referred to as the complementary energy. In the case of a linearly-elastic material, P ( q ) {\displaystyle P(q)} is a straight line and the strain energy is equal to the complementary energy.
Prescribed displacement In the case of prescribed displacement, the strain energy can be expressed in terms of the specified displacement and the crack surface Ω ( q , s ) {\displaystyle \Omega (q,s)} , and the change in this strain energy is only affected by the change in fracture surface area: δ Ω = ( ∂ Ω / ∂ s ) δ s {\displaystyle \delta \Omega =(\partial \Omega /\partial s)\delta s} . Correspondingly, the energy release rate in this case is expressed as
G = − ∂ Ω ∂ s | q . {\displaystyle G=-\left.{\frac {\partial \Omega }{\partial s}}\right|_{q}.}
Here is where one can accurately refer to G {\displaystyle G} as the strain energy release rate.
Prescribed loads When the load is prescribed instead of the displacement, the strain energy needs to be modified as Ω ( q ( P , s ) , s ) {\displaystyle \Omega (q(P,s),s)} . The energy release rate is then computed as
G = − ∂ ∂ s | P ( Ω − P q ) . {\displaystyle G=-\left.{\frac {\partial }{\partial s}}\right|_{P}\left(\Omega -Pq\right).}
If the material is linearly-elastic, then Ω = P q / 2 {\displaystyle \Omega =Pq/2} and one may instead write
G = ∂ Ω ∂ s | P . {\displaystyle G=\left.{\frac {\partial \Omega }{\partial s}}\right|_{P}.}
G in two-dimensional cases In the cases of two-dimensional problems, the change in crack growth area is simply the change in crack length times the thickness of the specimen. Namely, ∂ s = B ∂ a {\displaystyle \partial s=B\partial a} . Therefore, the equation for computing G {\displaystyle G} can be modified for the 2D case:
Prescribed Displacement: G = − 1 B ∂ Ω ∂ a | q . {\displaystyle G=-\left.{\frac {1}{B}}{\frac {\partial \Omega }{\partial a}}\right|_{q}.}
Prescribed Load: G = − 1 B ∂ ∂ a | P ( Ω − P q ) . {\displaystyle G=-\left.{\frac {1}{B}}{\frac {\partial }{\partial a}}\right|_{P}\left(\Omega -Pq\right).}
Prescribed Load, Linear Elastic: G = 1 B ∂ Ω ∂ a | P . {\displaystyle G=\left.{\frac {1}{B}}{\frac {\partial \Omega }{\partial a}}\right|_{P}.}
One can refer to the example calculations embedded in the next section for further information. Sometimes, the strain energy is written using U = Ω / B {\displaystyle U=\Omega /B} , an energy-per-unit thickness. This gives
Prescribed Displacement: G = − ∂ U ∂ a | q . {\displaystyle G=-\left.{\frac {\partial U}{\partial a}}\right|_{q}.}
Prescribed Load: G = − ∂ ∂ a | P ( U − P q B ) . {\displaystyle G=-\left.{\frac {\partial }{\partial a}}\right|_{P}\left(U-{\frac {Pq}{B}}\right).}
Prescribed Load, Linear Elastic: G = ∂ U ∂ a | P . {\displaystyle G=\left.{\frac {\partial U}{\partial a}}\right|_{P}.}
Relation to stress intensity factors The energy release rate is directly related to the stress intensity factor associated with a given two-dimensional loading mode (Mode-I, Mode-II, or Mode-III) when the crack grows straight ahead. This is applicable to cracks under plane stress, plane strain, and antiplane shear. For Mode-I, the energy release rate G {\displaystyle G} rate is related to the Mode-I stress intensity factor K I {\displaystyle K_{I}} for a linearly-elastic material by
G = K I 2 E ′ , {\displaystyle G={\frac {K_{I}^{2}}{E'}},}
where E ′ {\displaystyle E'} is related to Young's modulus E {\displaystyle E} and Poisson's ratio ν {\displaystyle \nu } depending on whether the material is under plane stress or plane strain:
E ′ = { E , p l a n e s t r e s s , E 1 − ν 2 , p l a n e s t r a i n . {\displaystyle E'={\begin{cases}E,&\mathrm {plane~stress} ,\\\\{\dfrac {E}{1-\nu ^{2}}},&\mathrm {plane~strain} .\end{cases}}}
For Mode-II, the energy release rate is similarly written as
G = K I I 2 E ′ . {\displaystyle G={\frac {K_{II}^{2}}{E'}}.}
For Mode-III (antiplane shear), the energy release rate now is a function of the shear modulus μ {\displaystyle \mu } ,
G = K I I I 2 2 μ . {\displaystyle G={\frac {K_{III}^{2}}{2\mu }}.}
For an arbitrary combination of all loading modes, these linear elastic solutions may be superposed as
G = K I 2 E ′ + K I I 2 E ′ + K I I I 2 2 μ . {\displaystyle G={\frac {K_{I}^{2}}{E'}}+{\frac {K_{II}^{2}}{E'}}+{\frac {K_{III}^{2}}{2\mu }}.}
Relation to fracture toughness Crack growth is initiated when the energy release rate overcomes a critical value G c {\displaystyle G_{c}} , which is a material property,
G ≥ G c , {\displaystyle G\geq G_{c},}
Under Mode-I loading, the critical energy release rate G c {\displaystyle G_{c}} is then related to the Mode-I fracture toughness K I C {\displaystyle K_{IC}} , another material property, by
G c = K I C 2 E ′ . {\displaystyle G_{c}={\frac {K_{IC}^{2}}{E'}}.}
Calculating G There are a variety of methods available for calculating the energy release rate given material properties, specimen geometry, and loading conditions. Some are dependent on certain criteria being satisfied, such as the material being entirely elastic or even linearly-elastic, and/or that the crack must grow straight ahead. The only method presented that works arbitrarily is that using the total potential energy. If two methods are both applicable, they should yield identical energy release rates.
Total potential energy The only method to calculate G {\displaystyle G} for arbitrary conditions is to calculate the total potential energy and differentiate it with respect to the crack surface area. This is typically done by:
calculating the stress field resulting from the loading, calculating the strain energy in the material resulting from the stress field, calculating the work done by the external loads, all in terms of the crack surface area.
Compliance method If the material is linearly elastic, the computation of its energy release rate can be much simplified. In this case, the Load vs. Load-point Displacement curve is linear with a positive slope, and the displacement per unit force applied is defined as the compliance, C {\displaystyle C}
C = q P . {\displaystyle C={\frac {q}{P}}.}
The corresponding strain energy Ω {\displaystyle \Omega } (area under the curve) is equal to
Ω = 1 2 P q = 1 2 q 2 C = 1 2 P 2 C . {\displaystyle \Omega ={\frac {1}{2}}Pq={\frac {1}{2}}{\frac {q^{2}}{C}}={\frac {1}{2}}P^{2}C.}
Using the compliance method, one can show that the energy release rate for both cases of prescribed load and displacement come out to be
G = 1 2 P 2 ∂ C ∂ s . {\displaystyle G={\frac {1}{2}}P^{2}{\frac {\partial C}{\partial s}}.}
Multiple specimen methods for nonlinear materials
In the case of prescribed displacement, holding the crack length fixed, the energy release rate can be computed by
G = − ∫ 0 q ∂ P ∂ s d q , {\displaystyle G=-\int _{0}^{q}{\frac {\partial P}{\partial s}}\,dq,}
while in the case of prescribed load,
G = ∫ 0 P ∂ q ∂ s d P . {\displaystyle G=\int _{0}^{P}{\frac {\partial q}{\partial s}}\,dP.}
As one can see, in both cases, the energy release rate G {\displaystyle G} times the change in surface d s {\displaystyle ds} returns the area between curves, which indicates the energy dissipated for the new surface area as illustrated in the right figure
G d s = − d s ∫ 0 q ∂ P ∂ s d q = d s ∫ 0 P ∂ q ∂ s d P . {\displaystyle Gds=-ds\int _{0}^{q}{\frac {\partial P}{\partial s}}\,dq=ds\int _{0}^{P}{\frac {\partial q}{\partial s}}\,dP.}
Crack closure integral Since the energy release rate is defined as the negative derivative of the total potential energy with respect to crack surface growth, the energy release rate may be written as the difference between the potential energy before and after the crack grows. After some careful derivation, this leads one to the crack closure integral
G = lim Δ s → 0 − 1 Δ s ∫ Δ s 1 2 t i 0 ( Δ u i + − Δ u i − ) d S , {\displaystyle G=\lim _{\Delta s\to 0}-{\frac {1}{\Delta s}}\int _{\Delta s}{\frac {1}{2}}\,t_{i}^{0}\left(\Delta u_{i}^{+}-\Delta u_{i}^{-}\right)\,dS,}
where Δ s {\displaystyle \Delta s} is the new fracture surface area, t i 0 {\displaystyle t_{i}^{0}} are the components of the traction released on the top fracture surface as the crack grows, Δ u i + − Δ u i − {\displaystyle \Delta u_{i}^{+}-\Delta u_{i}^{-}} are the components of the crack opening displacement (the difference in displacement increments between the top and bottom crack surfaces), and the integral is over the surface of the material S {\displaystyle S} . The crack closure integral is valid only for elastic materials but is still valid for cracks that grow in any direction. Nevertheless, for a two-dimensional crack that does indeed grow straight ahead, the crack closure integral simplifies to
G = lim Δ a → 0 1 Δ a ∫ 0 Δ a σ i 2 ( x 1 , 0 ) u i ( Δ a − x 1 , π ) d x 1 , {\displaystyle G=\lim _{\Delta a\to 0}{\frac {1}{\Delta a}}\int _{0}^{\Delta a}\sigma _{i2}(x_{1},0)u_{i}(\Delta a-x_{1},\pi )\,dx_{1},}
where Δ a {\displaystyle \Delta a} is the new crack length, and the displacement components are written as a function of the polar coordinates r = Δ a − x 1 {\displaystyle r=\Delta a-x_{1}} and θ = π {\displaystyle \theta =\pi } .
J-integral In certain situations, the energy release rate G {\displaystyle G} can be calculated using the J-integral, i.e. G = J {\displaystyle G=J} , using
J = ∫ Γ ( W n 1 − t i ∂ u i ∂ x 1 ) d Γ , {\displaystyle J=\int _{\Gamma }\left(Wn_{1}-t_{i}\,{\frac {\partial u_{i}}{\partial x_{1}}}\right)\,d\Gamma ,}
where W {\displaystyle W} is the elastic strain energy density, n 1 {\displaystyle n_{1}} is the x 1 {\displaystyle x_{1}} component of the unit vector normal to Γ {\displaystyle \Gamma } , the curve used for the line integral, t i {\displaystyle t_{i}} are the components of the traction vector t = σ ⋅ n {\displaystyle \mathbf {t} ={\boldsymbol {\sigma }}\cdot \mathbf {n} } , where σ {\displaystyle {\boldsymbol {\sigma }}} is the stress tensor, and u i {\displaystyle u_{i}} are the components of the displacement vector. This integral is zero over a simple closed path and is path independent, allowing any simple path starting and ending on the crack faces to be used to calculate J {\displaystyle J} . In order to equate the energy release rate to the J-integral, G = J {\displaystyle G=J} , the following conditions must be met:
the crack must be growing straight ahead, and the deformation near the crack (enclosed by Γ {\displaystyle \Gamma } ) must be elastic (not plastic). The J-integral may be calculated with these conditions violated, but then G ≠ J {\displaystyle G\neq J} . When they are not violated, one can then relate the energy release rate and the J-integral to the elastic moduli and the stress intensity factors using
G = J = K I 2 E ′ + K I I 2 E ′ + K I I I 2 2 μ . {\displaystyle G=J={\frac {K_{I}^{2}}{E'}}+{\frac {K_{II}^{2}}{E'}}+{\frac {K_{III}^{2}}{2\mu }}.}
Computational methods in fracture mechanics A handful of methods exist for calculating G {\displaystyle G} with finite elements. Although a direct calculation of the J-integral is possible (using the strains and stresses outputted by FEA), approximate approaches for some type of crack growth exist and provide reasonable accuracy with straightforward calculations. This section will elaborate on some relatively simple methods for fracture analysis utilizing numerical simulations.
Nodal release method If the crack is growing straight, the energy release rate can be decomposed as a sum of 3 terms G i {\displaystyle G_{i}} associated with the energy in each 3 modes. As a result, the Nodal Release method (NR) can be used to determine G i {\displaystyle G_{i}} from FEA results. The energy release rate is calculated at the nodes of the finite element mesh for the crack at an initial length and extended by a small distance Δ a {\displaystyle \Delta a} . First, we calculate the displacement variation at the node of interest Δ u → = u → ( t + 1 ) − u → ( t ) {\displaystyle \Delta {\vec {u}}={\vec {u}}^{(t+1)}-{\vec {u}}^{(t)}} (before and after the crack tip node is released). Secondly, we keep track of the nodal force F → {\displaystyle {\vec {F}}} outputted by FEA. Finally, we can find each components of G {\displaystyle G} using the following formulas:
G 1 NR = 1 Δ a F 2 Δ u 2 2 {\displaystyle G_{1}^{\text{NR}}={\frac {1}{\Delta a}}F_{2}{\frac {\Delta u_{2}}{2}}}
G 2 NR = 1 Δ a F 1 Δ u 1 2 {\displaystyle G_{2}^{\text{NR}}={\frac {1}{\Delta a}}F_{1}{\frac {\Delta u_{1}}{2}}}
G 3 NR = 1 Δ a F 3 Δ u 3 2 {\displaystyle G_{3}^{\text{NR}}={\frac {1}{\Delta a}}F_{3}{\frac {\Delta u_{3}}{2}}} Where Δ a {\displaystyle \Delta a} is the width of the element bounding the crack tip. The accuracy of the method highly depends on the mesh refinement, both because the displacement and forces depend on it, and because G = lim Δ a → 0 G NR {\displaystyle G=\lim _{\Delta a\to 0}G^{\text{NR}}} . Note that the equations above are derived using the crack closure integral. If the energy release rate exceeds a critical value, the crack will grow. In this case, a new FEA simulation is performed (for the next time step) where the node at the crack tip is released. For a bounded substrate, we may simply stop enforcing fixed Dirichlet boundary conditions at the crack tip node of the previous time step (i.e. displacements are no longer restrained). For a symmetric crack, we would need to update the geometry of the domain with a longer crack opening (and therefore generate a new mesh).
Modified crack closure integral Similar to the Nodal Release Method, the Modified Crack Closure Integral (MCCI) is a method for calculating the energy release rate utilizing FEA nodal displacements ( u i j ) {\displaystyle (u_{i}^{j})} and forces ( F i j ) {\displaystyle (F_{i}^{j})} . Where i {\displaystyle i} represents the direction corresponding to the Cartesian basis vectors with origin at the crack tip, and j {\displaystyle j} represents the nodal index. MCCI is more computationally efficient than the nodal release method because it only requires one analysis for each increment of crack growth. A necessary condition for the MCCI method is uniform element length ( Δ a ) {\displaystyle (\Delta a)} along the crack face in the x 1 {\displaystyle x_{1}} -direction. Additionally, this method requires sufficient discretization such that over the length of one element stress fields are self-similar. This implies that K ( a + Δ a ) ≈ K ( a ) {\displaystyle K(a+\Delta a)\approx K(a)} as the crack propagates. Below are examples of the MCCI method with two types of common finite elements.
4-node elements
The 4-node square linear elements seen in Figure 2 have a distance between nodes j {\displaystyle j} and j + 1 {\displaystyle j+1} equal to Δ a . {\displaystyle \Delta a.} Consider a crack with its tip located at node j . {\displaystyle j.} Similar to the nodal release method, if the crack were to propagate one element length along the line of symmetry (parallel to the x 1 {\displaystyle x_{1}} -axis) the crack opening displacement would be the displacement at the previous crack tip, i.e. u j {\displaystyle {\boldsymbol {u^{j}}}} and the force at the new crack tip ( j + 1 ) {\displaystyle (j+1)} would be F j + 1 . {\displaystyle {\boldsymbol {F}}^{j+1}.} Since the crack growth is assumed to be self-similar the displacement at node j {\displaystyle j} after the crack propagates is equal to the displacement at node j − 1 {\displaystyle j-1} before the crack propagates. This same concept can be applied to the forces at node j + 1 {\displaystyle j+1} and j .