Numerical and experimental analysis of the notch effect on fatigue behavior of polymethylmethacrylate metal based on strain energy density method and the extended finite element method

Using the Williams’ series expansions [23, 24], the stress relationships can be expressed nearness of the notch tip as follow: 4σθθσrrτrθ=∑i=12λirλi−1aifi,θ(θ)fi,r(θ)fi,rθ(θ)$$\left\{ {\matrix{ {{\sigma _{\theta \theta }}} \cr {{\sigma _{rr}}} \cr {{\tau _{r\theta }}} \cr } } \right\} = \sum\limits_{i = 1}^2 {{\lambda _i}{r^{{\lambda _i} - 1}}{a_i}} \left\{ {\matrix{ {{f_{i,\theta }}(\theta )} \cr {{f_{i,r}}(\theta )} \cr {{f_{i,r\theta }}(\theta )} \cr } } \right\}$$ where λ_i (i = 1,2) are characteristic values that are dependent on the notch opening angle (Fig. 1) [8]; r is the distance from notch tip; and a_i (i = 1, 2) are defined as [24–26]: a1=K1Nλ12π[(1+λ1)+χ1(1−λ1)],a2=K2Nλ22π[(1+λ2)+χ2(1−λ2)]$$\matrix{ {{a_1}} \hfill & = \hfill & {{{K_1^N} \over {{\lambda _1}\sqrt {2\pi } [(1 + {\lambda _1}) + {\chi _1}(1 - {\lambda _1})]}},} \hfill \cr {{a_2}} \hfill & = \hfill & {{{K_2^N} \over {{\lambda _2}\sqrt {2\pi } [(1 + {\lambda _2}) + {\chi _2}(1 - {\lambda _2})]}}} \hfill \cr } $$ And f₁(θ), f₂(θ) represents stress functions that correspond to symmetric (Mode I) and antsymmetric (Mode II) components, respectively, and their expressions are given in [8, 24, 25], K1N,K2N$$K_1^N,\,K_2^N$$ are the notch stress intensity factors (NSIFs) that correspond to the opening (Mode I, i = 1) and sliding (Mode II, i = 2) modes, respectively [25]: K1N=2πlimr→0r1−λ1σθθ(θ→0),K2N=2πlimr→0r1−λ2τrθ(θ→0)$$\matrix{ {K_1^N} \hfill & = \hfill & {\sqrt {2\pi } \mathop {\lim }\limits_{r \to 0} {r^{1 - {\lambda _1}}}{\sigma _{\theta \theta }}(\theta \to 0),} \hfill \cr {K_2^N} \hfill & = \hfill & {\sqrt {2\pi } \mathop {\lim }\limits_{r \to 0} {r^{1 - {\lambda _2}}}{\tau _{r\theta }}(\theta \to 0)} \hfill \cr } $$ The expression of the strain energy density in bidimensional problems (W) can be written at any point as [8]: 5W=12E[σθθ2+σrr2+σzz2−2vσθθσrr+σθθσzz+σrrσzz+2(1+v)τrθ2]$$\matrix{ W \hfill & = \hfill & {{1 \over {2E}}[\sigma _{\theta \theta }^2 + \sigma _{rr}^2 + \sigma _{zz}^2} \hfill \cr {} \hfill & {} \hfill & { - 2v\left( {{\sigma _{\theta \theta }}{\sigma _{rr}} + {\sigma _{\theta \theta }}{\sigma _{zz}} + {\sigma _{rr}}{\sigma _{zz}}} \right)} \hfill \cr {} \hfill & {} \hfill & { + 2(1 + v)\tau _{r\theta }^2]} \hfill \cr } $$ where E is Young’s modulus, and the values of stress σ_zz equals zero for plane stress condition and σ_zz = ν(σ_θθ + σ_rr) for plane strain condition.

Substituting the stress components Eq. (4) into Eq. (5), and integrating the strain energy density over a control area on a boundary defined by a distance of the radius R₀, the averaged SED (W−)$$(\mathop W\limits^ - )$$ expression can be calculated: 6W−=1A0∫0R0∫−γ+γW.rdrdθ$$\mathop W\limits^ - = {1 \over {{A_0}}}\int_0^{{R_0}} {\int_{ - \gamma }^{ + \gamma } {W.rdrd\theta } } $$ with A₀ as a zone of the control area;

The elastic strain energy density per cycle, W_e, can be formulated in terms of the stress amplitude applied and elastic modulus, E [27–29]: 7ΔWe=(Δσ)22E$${\rm{\Delta }}{W_e} = {{{{({\rm{\Delta }}\sigma )}^2}} \over {2E}}$$ In Mi et al. [30] the total strain energy density was determined by its relationship with fatigue life using the following equation without taking a plasticization [31]: 8ΔWT=A2NfB$$\Delta {W_T} = A{\left( {2{N_f}} \right)^B}$$ where (2N_f) is the fatigue life; A is the strain energy density coefficient; and B is the strain energy density index.

For fatigue problems [9–11], the combined mode I and II can easily be evaluated by means of the SED averaged over a control volume, considered as a material property and assumed to be a circular sector of radius R_c, as shown in Figure 5, according to Lazzarin and Zambardi [8]. It becomes a circle with radius R_c in the case of cracks or V-notches in mode I or mode mixed I + II (Fig. 6). For static loading, the radius R_c can be evaluated according to the following expression: 9Rc=(1+v)(5−8v)4πKIcσt2$${R_c} = {{(1 + v)(5 - 8v)} \over {4\pi }}{\left( {{{{K_{Ic}}} \over {{\sigma _t}}}} \right)^2}$$ Under fatigue limit conditions, R_c depends on smooth specimen fatigue limit Δσ₀ and on the threshold behavior ΔK_th, [6, 32] the relationship becomes: 10Rc=(1+v)(5−8v)4πΔKthΔσ02$${R_c} = {{(1 + v)(5 - 8v)} \over {4\pi }}{\left( {{{\Delta {K_{th}}} \over {\Delta {\sigma _0}}}} \right)^2}$$ When the ν = 0.3 equation (10) gives a critical radius: 11Rc≈0.27ΔKthΔσ02$${R_c} \approx 0.27{\left( {{{\Delta {K_{th}}} \over {\Delta {\sigma _0}}}} \right)^2}$$

Fig. 4.

Fig. 5.

Fig. 6.

El Haddad, Topper, and Smith [33] provide a parameter to link the stress intensity threshold term ΔK_th and fatigue limit Δσ₀ with critical crack length l0=ΔKthΔσ02⋅1π<=>l0⋅π=ΔKthΔσ02$${l_0} = {\left( {{{\Delta {K_{th}}} \over {\Delta {\sigma _0}}}} \right)^2} \cdot {1 \over \pi } < = > {l_0} \cdot \pi = {\left( {{{\Delta {K_{th}}} \over {\Delta {\sigma _0}}}} \right)^2}$$ And the critical crack length can be obtained by the following expression: KIc=σYπl0=>l0=KIcσY2⋅1π, from which l0=0.222mmΔKthΔσ02=0.697mm≅0.70mm$$\matrix{ {{K_{Ic}} = {\sigma _Y}\sqrt {\pi {l_0}} = > {l_0} = {{\left( {{{{K_{Ic}}} \over {{\sigma _Y}}}} \right)}^2} \cdot {1 \over \pi },} \hfill \cr {\quad {\rm{\;from which\;}}{l_0} = 0.222{\rm{mm}}} \hfill \cr {{{\left( {{{\Delta {K_{th}}} \over {\Delta {\sigma _0}}}} \right)}^2} = 0.697{\rm{mm}} \cong 0.70{\rm{mm}}} \hfill \cr } $$ Hence, the radius of the control volume Rc=0.27ΔKthΔσ02=0.189mm$${R_c} = 0.27{\left( {{{\Delta {K_{th}}} \over {\Delta {\sigma _0}}}} \right)^2} = 0.189{\rm{mm}}$$ In the case of a combined bending and shearing action on a component weakened by a V-notch and under conditions of linear elasticity, the SED averaged over the control volume, which embraces the notch tip, and can be calculated by means of the following expression [4, 27, 34]: 12ΔW−1=cwe1EΔK1NRC1−λ12;ΔW−2=cwe2EΔK2NRC1−λ22$$\matrix{ {\Delta {{\mathop W\limits^ - }_1} = {c_w}{{{e_1}} \over E}{{\left[ {{{\Delta K_1^N} \over {R_C^{1 - {\lambda _1}}}}} \right]}^2};} \hfill \cr {\Delta {{\mathop W\limits^ - }_2} = {c_w}{{{e_2}} \over E}{{\left[ {{{\Delta K_2^N} \over {R_C^{1 - {\lambda _2}}}}} \right]}^2}} \hfill \cr } $$ where ΔK₁ and ΔK₂ represent the Mode I and Mode II NSIF ranges, and R_c is the critical radius of the control volume related to Mode I and Mode II loadings. The terms e₁ and e₂ indicate two parameters that summarize the dependence on the V-notch geometry. Lazzarin and Zambardi [8] propose the following approximate formulas for ν = 0.3 under plane strain condition: 13e1=−5.373×10−6(2α)2+6.151×10−4(2α)+0.133,e2=4.809×10−6(2α)2−2.346×10−3(2α)+0.34.$$\matrix{ {{e_1}} \hfill & = \hfill & { - 5.373 \times {{10}^{ - 6}}{{(2\alpha )}^2}} \hfill \cr {} \hfill & {} \hfill & { + 6.151 \times {{10}^{ - 4}}(2\alpha ) + 0.133,} \hfill \cr {{e_2}} \hfill & = \hfill & {4.809 \times {{10}^{ - 6}}{{(2\alpha )}^2}} \hfill \cr {} \hfill & {} \hfill & { - 2.346 \times {{10}^{ - 3}}(2\alpha ) + 0.34.} \hfill \cr } $$ The parameters λ₁ and λ₂ are the eigenvalues of the Williams’ stress field solution for the N-SIF ΔK₁ and ΔK₂ for modes I and II, respectively, and (2α) is the opening angle. The evaluation of ΔK₁ and ΔK₂ is obtained numerically by the XFEM method from the following relationships: 14ΔK1=K1max−K1minΔK2=K2max−K2min$$\left\{ {\matrix{ {\Delta {K_1} = {K_{1max}} - {{\rm{K}}_{1min}}} \hfill \cr {\Delta {K_2} = {K_{2max}} - {{\rm{K}}_{2min}}} \hfill \cr } } \right.$$ The parameter c_w is the weighing parameter that takes into account the different loading ratio, and it is equal to 1.0 for R = 0 and 0.5 for R = -1.

In a highly stressed region, all stress and strain components are correlated to mode I and mode II (NSIF). The variation of total strain energy density averaged ΔW_T over a control volume surrounding the notch tip that can be evaluated from the following closed-form expression in a defined zone A₀ [2, 5, 30, 35]: 15ΔWT=ΔW1+ΔW2$${\rm{\Delta }}{W_T} = {\rm{\Delta }}{W_1} + {\rm{\Delta }}{W_2}$$ In numerical computation, the averaged SED can be evaluated directly from the FE results, once the size of the control volume is properly defined, ΔW_T, by summation of the strainenergies W_FEM,i calculated for each i-th finite element belonging to the control volume V [36]: 16ΔWT=cw∑VWFEM,iV$${\rm{\Delta }}{W_T} = {c_w}{{\sum\nolimits_V {{W_{FEM,i}}} } \over V}$$

In the case of blunt notches, the control volume is assumed to have crescent shape, with R_c being its maximum width as measured along the notch bisector line (Fig. 7b) as stated by Lazzarin and Berto [35]. Under mixed-mode conditions, the maximum principal stress σ_Max occurs on a notch edge point determined by the polar angle φ, out of the bisector line (Fig. 7b). The control volume is no longer centered with respect to the notch bisector, but rigidly rotated with respect to the notch bisector and centered on the point where the maximum principal stress reaches its maximum value (Fig. 7) [37–39]. This rotation is shown in Figure 7a where the control area is drawn for a U-shaped notch both under mode I loading (Fig. 7b) and mixed-mode loading (Fig. 7a). Thus, the procedure developed for mode I was extended to mixed I/II mode conditions.

Fig. 7.

Dealing with blunt notches under fatigue loading, the following expression can be used to evaluate the strain energy for notches in mode I, 17ΔW1=F(2α)H2α,RcRΔσtip2E$${\rm{\Delta }}{W_1} = F(2\alpha )H\left( {2\alpha ,{{{R_c}} \over R}} \right){{{\rm{\Delta }}\sigma _{tip}^2} \over E}$$ where F and H values are given in Berto and Lazzerin [21] as the functions of the openings angles, Poisson’s ratios, and the ratio of R_c to the notch tip radius; F(2α) can be expressed by the following relationship [21]: 18F(2α)=q−1q21−λ12π1+w~12$$F(2\alpha ) = {\left( {{{q - 1} \over q}} \right)^{2\left( {1 - {\lambda _1}} \right)}}{\left[ {{{\sqrt {2\pi } } \over {1 + {{\mathop w\limits^\~ }_1}}}} \right]^2}$$

To account for the shear contribution, Lazzarin et al. [40] dealt with a torsion problem in which they developed an expression for the strain energy density over volume control for the U-notch and the blunt V-notch. from this relation we used it in the shear action that contributes to the mode –II of sliding mode, this formula is expressed as follows [22]: 19ΔW2=ω22α,RcRΔτmax22G$${\rm{\Delta }}{W_2} = {\omega _2}\left( {2\alpha ,{{Rc} \over R}} \right){{{\rm{\Delta }}\tau _{max}^2} \over {2G}}$$ where, G is the shear modulus, and the value of ω₂ is taken as the same as the value reported in the torsion, and is a function of the notch opening angles and the ratio of R_c/R.

The presence of a geometric discontinuity defined by a notch in the specimen undergoing fatigue loading leads to crack initiation after a certain number of loading cycles, followed by its propagation. Numerically, the XFEM is used to represent the discontinuities independent of the mesh. The discontinuities can be modeled by enriching all discontinuous elements using enrichment functions that satisfy the discontinuous behavior and adding additional nodal degrees of freedom. It enables us to control the initiation and propagation of fatigue cracks without the need to remesh the structure.

Crack initiation is localized at the point where the principal stress criterion σ_I reaches its maximum value (see Fig. 12), triggering a zone surrounding the crack tip defined by the control volume, which follows the cracked tip during crack propagation (Fig. 9).

Fig. 9.

Fig. 10.

Fig. 11.

Fig. 12.

Based on the local approach of strain energy density which consists in the first step of defining the volume control zone where it rules the elaboration of the zone damage and then leads to the initiation and propagation of cracks. This zone is located around the notch tip and includes a quantity of elastic deformation energy. If this energy reaches the value of the critical energy defined by W_c, the structure undergoes a crack initiation and rupture. A calculation is made of the total local energy, which takes into account two quantities of energy linking to the bending and shearing energies ΔW_T = ΔW₁ + ΔW₂ normalized by the critical energy W_c described by the equations above.

Figure 10 shows the distribution of the ASED normalized by the critical energy in dimensionless units as a function of the fatigue life for the two types of notches, V and U. The results obtained are compared with those from the experiment. At low fatigue cycles, the energy values reach maximum values for both types of notches. The trend follows a decay towards a stabilization of the energy at a high life cycle. In terms of amplitude, the local energy of the sharp notch is higher than the energy values of the blunt notch, and this is due to the attenuation of the stress concentrations at the notch level. The presence of notches in a structure always leads to the amplification of energy concentrations and the reduction of the structure’s lifetime. The introduction of the fatigue behavior curve of the smooth specimen (lack of notches) in Figure 10 (blunt U notch) shows that the lifetime of these specimens is longer (more cycles) than the specimens with notches and with lower energy quantities. This implies that a high concentration of energy is present in specimens with notches.

Figure 11 shows approximations of the mean local energy SED by power equations of the experimental values for the different notches. Table 2 gives the approximation equations of all the notches with correlations coefficients. For a sharp V-notch, the critical volume becomes a circular sector of radius Rc centered at the notch tip (Fig. 6). For a blunt notch under mixed-mode loading, the critical volume is no longer centered at the notch tip, but rather at the point where the principal stress reaches its maximum value along the border of the notch (Fig. 7).

It was fundamentally assumed that the crescentshaped volume rotates rigidly under mixed mode, with no change in shape and size (Fig. 7) [3]. The results obtained are in good agreement with the experimental results (Fig. 10).

Figure 12 shows the area where the first maximum principal stress (σ_Imax) is located for both U and V notches. For the V notch, the stress is located at the notch tip (Fig. 12a), whereas for the U notch the stress is located far from the notch tip, at the point where the notch lip begins (Fig. 12b).

Another interest in this study is to analyze the variation of cyclic loading on the evolution of the average strain density energy under the effect of the presence of notches. Figures 13 and 14 show the distribution of the average strain energy normalized by critical energy via fatigue life: the high loads develop high-amplitude average SEDs. For low angles in V notches, the intensity of the average SED is greater than for more open angles. The lifetime is reduced compared to the more open angles. The values obtained by numerical calculation are in good approximation compared with those of the experimental ones.

Fig. 13.

Fig. 14.

For U-shaped notches with a significant radius, the high loading leads to high mean SED values compared to low-radius notches and a longer lifetime compared to sharp notches. Sharp V-notches develop higher mean SED values than blunt notches. This suggests that the loading amplitude plays an influential role on the ASED and the fatigue life of structures subjected to bending fatigue.

Several simulations were performed for a sharp notch V140° at different mesh element numbers. Figure 15 shows convergence to the reference result from a mesh element number of 3224 elements. As the mesh element number increases, the results tend to deviate from the reference value.

Fig. 15.

A statistical calculation was performed to compare the experimental and numerical results. The results of the calculation showed a good correlation between the two sets of results. The standard deviation of the experimental test results was slightly lower than that of the simulation test results. The minimum and maximum values of the two data sets were also very close (see Table 3).