Abstract

In the industrial field, the control of the structural integrity of a cracked metal part remains a capital maneuver in the prediction of failures in order to accurately define the functional lifetime of these parts. In this case, the Stress Intensity Factor (SIF) plays a central role to assess the damaging of a crack in a structure under a given loading. The purpose of this paper is to examine and check the analytical calculation of the SIF through an empirical approach of the form factor with those computing in mode I, using in the linear elastic domain, by the finite element method with ABAQUS software for two different cracked sections. The first defect studied in this work has a form of a circular section with a central crack under uniform tensile load. The second defect has a form of elliptical section with a central crack under uniform tensile load. The form factor function Y used, remains the main element for the analytical stress intensity factor determination. For each form of section, the obtained results were examined by comparing with the results given by the ABAQUS software. Indeed, after examination, the model was finally approved and validated for the circular section, on the other hand, it couldn’t be validated and adopted unfortunately in the elliptical case.

1. Introduction

Cracks are present in all structures, they can exist in the basic defect form in the material or may be induced during construction, these cracks are the main reason for the most failures that occur in structures and parts of machines in service, subjected to static or dynamic forces ¹. The purpose of fracture mechanics is to study and predict the cracks initiation and propagation in solids. The start of the study of brittle materials rupture began in 1920, with the work of Griffith, before reappearing in the 1950s and 1960s, when the discipline took off really with the new works of Irwin and Rice. As for the study of the rupture of ductile materials, it only begins at the late of 1960s and through the 1970s, with the fundamental works of Rice and Tracey ² and Gurson ³. Linear elastic mechanics is interested in the rupture of brittle materials. It is widely used by engineers because it allows the use of global energy criteria such as stress intensity factors FIC ⁴. Griffith ⁵ was interested in the problem of rupture in an elastic cracked medium from an energetic view point. He thus highlighted a variable called later the rate of energy restitution characterizing the fracture, and whose critical value is a characteristic of material ⁶. The first theoretical developments in the analysis of stress and strain fields nearly to a crack in elasticity. These studies, in particular by Irwin ⁷, allowed to define the FIC, characterizing the state of stress of the region in which the rupture occurs. The development of the finite element method made it possible to study numerical the mechanic of rupture, thus proposing more precisely solutions to more complex problems. Then appeared a multitude of methods allowing to calculate the stress intensity factors ⁸. Among these methods the method of the principle of superposition, extrapolation of displacements and the collocation method borders ⁹.

The main objective of this paper is to compute the stress intensity factor KI, in mode I, using in the linear elastic domain, by the finite element method with ABAQUS software. The first defect studied in this work has a form of a circular section with a central crack under uniform tensile load. The second defect has a form of elliptical section with a central crack under uniform tensile load.

2. Numerical Modeling

In brittle fracture mechanics, the Stress Intensity Factor K is the only significant parameter allowing to know the state of stresses and strain at any point of crack. The stress intensity factor is expressed in mode I in the following form ¹⁰:

(1)

Where Y is a dimensionless number called a form factor. It is a function of geometry and the crack length a and 𝜎 is the stress tensor. Then, the form factor is given by the following formula ¹⁰:

(2)

Furthermore, in this present study, were extracted the Stress Intensity Factor values of KI in mode I, simulating the behavior in the linear elastic domain. The ABAQUS software offers a way to evaluate the full contour by the classical finite element method (FEM) ¹¹, which usually requires to define explicitly the crack front.

The original formulation to calculate the integral J was defined by RICE ¹², in elastic linear area cracked plane (plane strain or plane stress state), in the absence of loading on crack (eq 3). This integral is independent of the chosen contour

(3)

Where W is the elastic strain energy density, T is the traction vector, u is the displacement vector and C is the counter clockwise contour beginning on the lower crack surface and ending on any point on the upper crack surface.

In the ABAQUS software, the special command (* CONTOUR INTEGRAL) ¹³ were accessed that is dedicated to calculate the integral J in the crack front given by the discretized formula of J with loading (eq 4):

(4)

Where np is the number of Gauss points, w_p is the weight integration, ξ_k is the coordinates of elements in local landmarks, A* is the area of the surface between the contour C₁ and C₂ and q is a smoothed function that can be chosen as q=1 for C₁ and q=0 for C₂.

The integral J is equivalent to the rate of the energy restitution G with Irwin formula ¹⁰. It was (the integral J) used to calculate the factor KI according to the following Rice equation ¹²:

(5)

With for plane stress and for plane strain.

To compute the KI Stress Intensity Factor numerically from the integral J, equation (5) was introduced in the ABAQUS code.

The numerical values of displacements are given by FEM method ¹¹ with the interpolation functions. Let the 8-node 2D isoparametric element shown in Figure (2-a). Using the interpolation functions N_i (ξ, η) corresponding to node i with (x_i, y_i) coordinates in the global reference and (ξ_i, η_i) in the local one ¹⁴.

(6)

The interpolation functions (shape functions) of the isoparametric element with 8 nodes are ¹⁴:

In this section, we will study the case of a defect with a circular section and for different crack lengths (a). The defect shape investigated is a circular section with a crack which located on the outer surface (Figure 3). The material is considered to be isotropic with a Young's modulus E and a Poisson's ratio ν. We apply a constant tensile load σ all around the section as shown in Figure 3.

The geometric characteristics of the defect used in this study are presented in the following Table 1.

The study of convergence is a very important phase for the validity of a result that determines which the best accuracy of the results with suitable computation time. The choice and generation of the mesh size have a great influence on the results, thus, it must take carefully consideration during the simulations.

In the other hand we will study in this section the case of a defect with an elliptic section defined by the two semi-axes with: a=10cm and b (varaiable). The defect shape investigated is an elliptic section with a central crack under a uniform tensile load σ (Figure 5). The material is considered to be isotropic with a Young's modulus E and a Poisson's ratio ν.

The geometric characteristics of the defect used in this study are present in the following Table 2.

In the same way as previously, we studied the convergence of the value of the integral J and of the stresses intensity factor according to the number of elements and for various values of the ratio b/a. as mentioned in Figure 6.

3. Results and Discussions

The numerical results of the integral J, of the stresses intensity factor KI and the form factor Y obtained using the ABAQUS software were presented in Table 3.

Then, the analytical expression of the form factor Y is obtained by the weight function method ¹⁵, and the coefficients of the polynomial are determined by identification with the results obtained by ABAQUS (Table 3). We are therefore looking for Y in the following form according to data of cracks handbook ¹⁶:

(7)

The form factor expression obtained in eq (7) contains six undetermined coefficients (C₁, C₂, C₃, C₄, C₅, C₆). To determine them, we use the ten values given in Table 3, which leads to have an overdetermined system of 10 equations with 6 unknowns as showen below:

(8)

The resolution of this system of equations (eq 8) was carried out using the Matlab software, and the obtained results for C₁, C₂, C₃, C₄, C₅ and C₆ are injected directly in eq 7, then we obtain the final expression of the form factor:

(9)

In order to check the accuracy of the obtained expression, we will recalculate the form factor with this expression for each value of a/R and we will compare them with the results already obtained by ABAQUS. All the obtained results are shown in the table below:

We note from the table that the analytical expression (eq 9) gives very close results and we also notice that the relative error is small with a maximum error of about 0.91%. Then to conclude, this analytical expression of the form factor was validated after the comparison of these results with those obtained with ABAQUS.

It should also be noted that the expression obtained for the form factor is only valid in the interval where the ratio a/R ∈ [0.1 - 1].

From these obtained results from the form factor, we will use an empirical formula given by equation 1 to calculate the analytical stresses intensity factor in an indirect way, and the results obtained thus will then be compared with those given by ABAQUS

After the comparison of the various values of (SIF) KI ABAQUS and KI analytical given in Table 5. It seen that in general, the agreement is quite good and the difference between the two quantities is quite smaller with a maximum reliance error of 0.90%. We can then say that this empirical approach to calculate the SIF via the form factor function is validated and gives fairly good results compared to other much longer and more complex methods.

As we did previously for the circular case, we will present in this section the variation of the integral J, the stresses intensity factor K and the form factor Y as a function of b/a. For each value of b/a, we took the result of the last two contours which remain stable and unchanged. The results are summarized in Table 6 below.

Then, by the same way in the precedent case (circular), the analytical expression of the form factor Y is obtained by the weight function method ¹⁵ and the coefficients of the polynomial are determined by identification with the results obtained by ABAQUS (Table. 6). We are therefore looking for Y in the following form according to data of cracks handbook ¹⁶:

(10)

The form factor expression obtained in eq (10) contains six undetermined coefficients (C₁, C₂, C₃, C₄, C₅, C₆). To determine them, we use the ten values given in Table 6, which leads to have an overdetermined system of 10 equations with 6 unknowns as showen below:

(11)

The resolution of this system of equations (eq 11) was carried out using the Matlab software, and the obtained results for C₁, C₂, C₃, C₄, C₅ and C₆ are injected directly in eq 10, then we obtain the final expression of the form factor:

(12)

In order to check the accuracy of the obtained expression, we will recalculate the form factor with this expression for each value of b/a and we will compare them with the results already obtained by ABAQUS. All the obtained results are shown in the table below:

We note from the table that the analytical expression (eq 12) gives very close results and we also notice that the relative error is small with a maximum error of about 0.16%. Then to conclude, this analytical expression of the form factor was validated after the comparison of these results with those obtained with ABAQUS.

As we did in the circular case, from these obtained results from the form factor, we will use an empirical formula given by equation 1 to calculate the analytical stresses intensity factor in an indirect way, and the results obtained thus will then be compared with those given by ABAQUS

In Contrary to what has just been seen in the circular case, in the elliptical one the empirical approach studied in this paper is not validated because it diverges and provided results very far from the numerical calculations made by ABAQUS with a relative error percentage included between 29 and 69%. In this case it is necessary to readjust the expression of adopted the form factor function in order to take into account the geometric specificities of this elliptical case.

4. Conclusion

Accurate and simple analytical solutions for the SIF evaluation of 2-D planar cracks propagation are necessary for safety control. In this paper, an empirical approach using form factor function was satisfactory examined for cracks in two simple geometries. A first determination made, by the finite element method using the ABAQUS software, of the Stress Intensity Factor in mode I and of the shape for two types of cracked sections under uniform normal loading: a circular section for different crack lengths and an elliptical section for several b/ ratios (major axis / minor axis). Then the numerical results obtained were compared to those calculated by the empirical approach. The accuracy of the empirical model is excellent solely in the circular case, on the other hand in the elliptical case, the results diverge with a fairly large error percentage, and consequently the model cannot be adopted.

References