Advances in Materials Science and Engineering
Volume 2021 (2021), Article ID 9987079, 10 pages
https://doi.org/10.1155/2021/9987079

Research Article

Three-Dimensional Stress Intensity Factor Correction Method Based on Thickness Size Effect

Xiaochen Zhang,¹ Hanwen Li,¹ Weiying Meng,¹ Tian Zhang,¹ Zhenpeng Liu,² and Xiao Lv³

¹School of Mechanical Engineering, Shenyang Jianzhu University, Mailing Address: No. 25 Hunnan Middle Road, Hunnan District, Shenyang 110168, China
²Shenyang Jianzhu University, Mailing Address: No. 25 Hunnan Middle Road, Hunnan District, Shenyang 110168, China
³Shanghai MicroPort Medical (Group) Co, Ltd., ZhangDong Rd, ZJ Hi-Tech Park, Shanghai 201203, China

Correspondence should be addressed to Weiying Meng; mengweiying025@163.com

Received 30 March 2021; Revised 5 July 2021; Accepted 12 August 2021; Published 18 August 2021

Academic Editor: Georgios I. Giannopoulos

Copyright © 2021 Xiaochen Zhang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

This paper presents the fatigue crack growth test of 7N01 aluminum alloy specimens with different thicknesses under three-point bending. Moreover, the effect of specimen thickness on the fatigue crack growth life was analyzed. The influence mechanism of the thickness size effect on the stress intensity factor was explored. A three-dimensional stress intensity factor calculation model based on the thickness size effect was proposed by combining the finite element method and the interaction integral method. The results show that the fatigue crack growth life of 7N01 aluminum alloy decreases with the increase of specimen thickness. The numerical solutions of the stress intensity factor (the proposed model) are higher than the analytical solutions (the traditional formula) in the steady-state crack growth stage. The difference between numerical solutions and analytical solutions increases with the increase of crack growth length and specimen thickness. Compared with the analytical solutions, the relationship curve between crack length “a” and cycle number “N” from the numerical solutions is closer to the experimental data. This result proves that the proposed calculation model can correct the stress intensity factor with the thickness size effect.

1. Introduction

From the perspective of fracture mechanics, the solution of stress intensity factor is an important subject for evaluating fracture mechanical properties of materials [1, 2]. In engineering practice, the fracture mechanics analysis of the crack propagation process is essentially a three-dimensional problem. The three-dimensional crack is one of the most common and dangerous defects in the mechanical structure. Therefore, accurate stress analysis of the cracked structure is the key to ensure the safety and reliability of the structure [3–6].

The stress intensity factor is an important parameter to characterize the strength of the stress field near the crack tip of an elastic body under external force. The solution of the stress intensity factor is often related to the geometry of the specimen, the crack growth length, and the load type [7, 8]. The main methods to calculate the stress intensity factor are the analytical method and numerical method. The analytical method, as a traditional method, is generally used in practical engineering, and its accuracy is limited by the geometric structure and static load conditions. The numerical method is widely concerned at home and abroad because it is not affected by the geometric structure or load conditions, which mainly includes the finite element method [9], meshless method [10], and boundary element method [11].

The study on the stress intensity factor has always been a hot issue at home and abroad. Nehar et al. [12] simulated the crack growth at the interface of brittle materials and biomaterials by the extended finite element method and solved the stress intensity factor by the displacement jump method and verified its accuracy by comparing with the theoretical value. Zhou and Xue [13] used the stress extrapolation method based on the extended finite element to calculate the stress intensity factor of the crack in a plate. It was found that the model has high accuracy and is easy to implement for the evaluation of pure type I cracks. Ai et al. [14] used FRANC3D software to solve the stress intensity factor and verified it through fatigue life prediction. Its basic principle is still based on the finite element method, and the mesh needs to be divided repeatedly and encrypted. Raviraj et al. [15] conducted three-dimensional finite element simulation for stress intensity factors of CT specimens with different compositions or different length-width ratios and verified the validity of the model by comparing with the test data. At present, the solution of stress intensity factors is mostly confined to two-dimensional problems in engineering applications, and the effect of thickness size on its calculation results is always ignored. However, the existing three-dimensional finite element model cannot reasonably describe the change of the total energy caused by microcracks. Because the stress intensity factor at the crack tip has a strong three-dimensional effect, if we want to accurately evaluate the fracture mechanical properties of materials and reasonably predict their fatigue crack growth life, we need to start from the three-dimensional point of view, reasonably describe the total energy change of the system, and correct the calculation results of the stress intensity factor with the thickness size effect.

In this paper, a three-point bending SE (B) specimen with different thicknesses is designed for the fatigue crack growth test of 7N01 aluminum alloy used in a high-speed train (the initial material that has not undergone service), and the influence of specimen thickness on the fatigue crack growth life is discussed. The finite element method is used to establish a three-dimensional stress intensity factor calculation model, and the interaction integral method is used to solve the stress intensity factor. The influence of thickness size effect on the stress intensity factor is analyzed.

2. Test Cases

The fatigue crack growth test of 7N01 aluminum alloy used in the high-speed train was carried out in the air environment of the laboratory. The three-point bending SE (B) specimen is designed according to the standard recommended in GB/T 4161‐2007, the standard test method for determination of plane-strain fracture toughness [16], as shown in Figure 1. The length of the specimen is 140 mm, the notch at the front edge of the crack is prefabricated by the wire cutting method, and the depth of the notch is 7 mm. In order to study the fatigue crack growth behavior of materials and the influence of the thickness size effect on the stress intensity factor at the crack tip, according to the relevant provisions in GB/T 6398-2017 [17] on the thickness design size range for the three-point bending SE (B) specimen, the specimen is designed in two thickness sizes. The specific dimensions are 30 mm in width W, 15 mm in thickness B and 30 mm in width W, 10 mm in thickness B. The material studied in this paper comes from the proximate matter of the sleeper beam for the high-speed railway, and the thickness values of 10 and 15 mm were selected based on base metal thickness from practical engineering.

Figure 1: Standard three-point bending SE (B) specimen.

The Instron 8800 electrohydraulic servo fatigue tester is used as the test equipment, as shown in Figure 2. The main parameters of the equipment are as follows. The maximum load is ±100 KN; the load measurement accuracy is ±0.005% of the full scale for the load sensor or ±0.5% of the displayed value (whichever is greater); the strain measurement accuracy is ±0.005% of the full scale for the extensometer or ±0.25% of the displayed value (whichever is greater); the conversion accuracy of the control channel is 0.003%, the system resolution of the whole machine is 19 bits, and the dynamic measurement error within 100 Hz is less than 0.5%. The test is carried out at the same stress level, the span S is 120 mm, the stress ratio R is 0.1, and the test frequency is 10 Hz. The test results are shown in Figure 3. According to the test results, the a-N curve and curve are drawn, as shown in Figure 4. The comparison of the test data curve for the SE (B) specimen with different thicknesses is obtained.

Figure 2: Instron 8800 electrohydraulic servo fatigue tester.

Figure 3: Loading device and crack growth length test.

(a)

(b)

Figure 4: Comparison diagram on the test data curve for specimens with different thicknesses: (a) a-N curve and (b) curve.

As shown in Figure 4(a), as the design thickness of the three-point bending SE (B) specimen increases, the fatigue crack propagation life of the specimen decreases. The fatigue crack growth life of the specimen with the thickness of 10 mm is 1.33 times that of the specimen with the thickness of 15 mm. As shown in Figure 4(b), in the steady-state crack growth stage, the fatigue crack growth rate of the specimen with the thickness of 10 mm is slightly lower than that of the specimen with the thickness of 15 mm. The test results show that the thickness size effect has a certain effect on the fatigue crack growth life and fatigue crack growth behavior. When the fracture mechanical properties of materials are evaluated, different thicknesses of the specimen design often affect the accuracy of the test results.

In addition, the expression of stress intensity factor amplitude for the SE (B) specimen is as follows [17]:

where

is the load range; W is the width of the specimen; B is the thickness of the specimen; is the ratio of the crack length to the specimen width, . According to the expression of stress intensity factor amplitude , it can be known that the thickness size is not a parameter that affects the calculation result of the stress intensity factor; that is, when the specimen has the same width, crack growth length, and stress level, the calculation results of the stress intensity factor are exactly the same. Therefore, the influence of the thickness size effect on the stress intensity factor is neglected, which will lead to inaccurate calculation results.

3. Finite Element Analysis of the Three-Dimensional Stress Intensity Factor

Due to the complexity of the structure and stress state in the practical engineering, it is difficult to obtain the theoretical solution of the stress-strain field at the crack tip. Therefore, not all stress intensity factor calculation methods are suitable for the three-dimensional case. The finite element method is relatively accurate and not restricted by the crack geometry and load form, so the finite element method is more suitable for three-dimensional problems. Meanwhile, in the theory of linear elastic fracture mechanics, for the description of the change of the system total energy caused by the preexisting microcracks, the total energy release from the interaction integral method is more reasonable than the energy release rate.

Therefore, in order to explore the influence of the thickness size effect on the stress intensity factor, ANSYS finite element software was used to establish the finite element calculation model of the three-dimensional SE (B) specimen, and the total energy change of the system was reasonably described by combining the interaction integral method. Finally, a three-dimensional stress intensity factor calculation model is constructed based on the thickness size effect to realize the solution of the stress intensity factor at the crack tip.

3.1. Solving the Stress Intensity Factor

The interaction integral method is used to calculate the stress intensity factor of the three-dimensional finite element model at the crack tip. The interaction integral method is to separate and obtain the stress intensity factors of type I in the real field by establishing the auxiliary field at the crack tip. Among them, the auxiliary field at the crack tip is an arbitrary displacement field and stress field satisfying the equilibrium conditions and geometric relations, while the real field at the crack tip represents the actual displacement field and stress field [18]. Rice used the J-integral to calculate the stress intensity factor and proved that the value of the stress intensity factor is independent of the selected integration loop [19]. According to the superposition principle of stress, strain, and displacement in linear elastic materials, the real field and the auxiliary field are substituted by superposition into the loop integral defined by the J-integral, and the expression is obtained as follows [20]:

By separating the J-integral corresponding to the auxiliary field and the real field, the stress intensity factor from the interaction integral method can be obtained. The solution expression is as follows:

where represents the stress intensity factors of type I in the real field; represents the stress intensity factors of type I in the auxiliary field; and E is the elastic modulus of the material.

3.2. Establishing the Three-Dimensional Finite Element Model

A three-dimensional finite element model is established for the three-point bending SE (B) specimen with different thicknesses. The elastic modulus E is 69 GPa, and Poisson’s ratio V is 0.3 in the model. The three-dimensional model of the three-point bending SE (B) specimen is obtained, as shown in Figure 5.

Figure 5: Three-dimensional model of the three-point bending SE (B) specimen.

Because the stress and strain at the crack tip are infinite, that is to say, a singularity appears at the crack tip. Therefore, the whole structure is divided into two parts for finite element meshing. Among them, one part is the area including the crack, that is, the crack body; the other part is the part other than the crack body, that is, the noncrack body. The two-part mesh is connected by a multipoint constraint equation. For the crack body, the degenerated 20-node isoparametric singular element is used for the simulation. It is performed by degrading a certain surface of the common 20-node isoparametric element into a line (the front line of the crack) and moving 4 intermediate nodes adjacent to the line onto 1/4 side length of the line, which constitutes a degenerated 20-node isoparametric singular element [21].

The meshing of the crack tip and the processing of intermediate nodes are shown in Figure 6. The element constructed in this way can meet the singularity requirements of fracture mechanics and can be applied to the calculation of linear elastic and elastic-plastic fracture mechanics parameters. In order to reduce the calculation scale, an 8-node solid element is used to simulate the noncracked body. After setting the boundary conditions and applying the load, the solver of ANSYS software can be used to obtain the stress-strain field of the crack tip under the real loading and the stress-strain field of designable auxiliary. By substituting the solution results into equations (3) and (4), the stress intensity factor at the crack tip solved by the interaction integral method can be obtained as follows.

Figure 6: (a) Crack tip meshing. (b) Processing of intermediate nodes at the crack tip.

In the process of solving the stress intensity factor at the crack tip by the interaction integral method, the nodes at the crack tip are integrated by winding according to the path of the surrounding nodes, and each winding integral corresponds to the solution of one stress intensity factor. In order to ensure the accuracy of the calculation, at least 10 values of stress intensity factors are output for each node along the thickness direction of the crack front edge, and a more stable solution is selected as the output result of the node. Finally, the mean value of the stress intensity factors at all nodes of the crack front is calculated, which is taken as the final calculation result of the stress intensity factors at the crack length.

4. Thickness Size Effect of the Three-Dimensional Stress Intensity Factor

4.1. Comparison and Analysis of Calculation Results of the Stress Intensity Factor

The analytical solution of the stress intensity factor is obtained by equations (1) and (2), and the relationship curve between the stress intensity factor and fatigue crack growth length is established based on the analytical solution of the stress intensity factor and the numerical solution of the three-dimensional stress intensity factor obtained by the finite element method. The comparison between the analytical and numerical solutions of the curve of stress intensity factors K and fatigue crack growth length a for SE (B) specimens with different thicknesses is shown in Figure 7.

(a)

(b)

Figure 7: K-a relation curve of the SE (B) specimen with different thicknesses: (a) K-a relation curve of the specimen with thickness 10 mm; (b) K-a relation curve of the specimen with thickness 15 mm.

As shown in Figure 7, the numerical solution of the stress intensity factor is larger than the analytical solution in each stage of steady-state crack growth. As the crack growth length increases, the gap between the numerical solution and the analytical solution also increases. This is mainly because with the increase of crack growth length, the net cross-sectional area carried by the sample decreases, and large plastic deformation occurs at the crack tip, which results in the crack closure effect. As shown in Figure 7(a), the numerical solution and analytical solution of the stress intensity factor for the specimen with the thickness of 10 mm are relatively close, and the difference is between 2.531% and 6.414%. With the increase of specimen thickness, the gap between the numerical solution and analytical solution of the stress intensity factor increases. As shown in Figure 7(b), the difference between the numerical solution and analytical solution of the stress intensity factor for the specimen with the thickness of 15 mm is 3.961%∼8.908%. The comparison between K-a relation curves of different thickness specimens shows that the thickness size effect has a certain influence on the solution of the stress intensity factor, and the value of the stress intensity factor obtained by the analytical method is relatively conservative.

4.2. Correction of the Three-Dimensional Stress Intensity Factor

In order to verify the calculation accuracy of the three-dimensional stress intensity factor obtained from the finite element model proposed in this paper, the Paris formula is theoretically derived as follows:

According to the crack growth increment measured by the test and the stress intensity factor range corresponding to the crack length, the corresponding number of crack growth cycles is calculated by equation (5), which is as follows:

By substituting the numerical solution of the three-dimensional stress intensity factor calculated by the finite element method into equation (5), a new set of material constants C and M can be obtained. Then, according to the above derivation process based on the Paris formula, the numerical solution of the a-N relation curve of the stress intensity factor can be fitted. Similarly, the analytical solution of the a-N relation curve of the stress intensity factor is fitted. The two new a-N curves were compared with the test data.

The comparison of the a-N curve fitted by different methods is shown in Figure 8. In Figure 8, the fatigue crack growth life in the a-N relation curve obtained by the analytical solution is higher than that obtained by the numerical solution and that from the test data, so it is dangerous to evaluate the material based on the curve fitting results obtained by the analytical solution. The a-N curve obtained by the numerical solution is closer to the test data, and the difference between the fatigue crack growth life for specimens with different thicknesses and test data is less than 2.578%. The proposed solution method of the three-dimensional stress intensity factor by considering the thickness size effect has a certain correction effect for evaluating the fracture mechanical properties of materials. As the thickness of the specimen increases, the accuracy of the analytical solution of the stress intensity factor decreases, and the correction effect of the proposed solution method of the three-dimensional stress intensity factor becomes more significant.

(a)

(b)

Figure 8: Comparison diagram on the a-N curve fitted by different methods. (a) Comparison diagram on the a-N relation curve of the specimen with thickness 10 mm. (b) Comparison diagram on the a-N relation curve of the specimen with thickness 15 mm.

4.3. Effect of Thickness Size on the Stress Intensity Factor

In order to further verify the effect of the thickness size on the solution of the stress intensity factor, the distribution state of the stress intensity factor at each node of the crack tip is truly presented. The numerical and analytical solutions of stress intensity factors under different fatigue crack growth lengths are obtained, and the data distribution along the thickness direction is plotted. When the fatigue crack growth length is 12 mm, 15 mm, and 18 mm, the distribution of the numerical solution and analytical solution of the stress intensity factor along the thickness direction is shown in Figures 9–11, where the abscissa in the figure is the ratio of the length in the thickness direction of the specimen to the total thickness.

(a)

(b)

Figure 9: Distribution of numerical and analytical solutions for the stress intensity factor along the thickness direction, a = 12 mm. (a) Numerical and analytical solutions for the specimen with thickness 10 mm. (b) Numerical and analytical solutions for the specimen with thickness 15 mm.

(a)

(b)

Figure 10: Distribution of numerical and analytical solutions for the stress intensity factor along the thickness direction, a = 15 mm. (a) Numerical and analytical solutions for the specimen with thickness 10 mm. (b) Numerical and analytical solutions for the specimen with thickness 15 mm.

(a)

(b)

Figure 11: Distribution of numerical and analytical solutions for the stress intensity factor along the thickness direction, a = 18 mm. (a) Numerical and analytical solutions for the specimen with thickness 10 mm. (b) Numerical and analytical solutions for the specimen with thickness 15 mm.

As shown in Figures 9–11, the numerical solutions of the stress intensity factor are the smallest at the edge of the thickness direction, which are almost the same as the analytical solution, where the influence of the thickness size effect on the stress intensity factor can be ignored. The closer to the center of the thickness direction, the more significant the influence of the thickness size effect on the stress intensity factor and the greater the gap between the numerical solution and analytical solution of the stress intensity factor which will reach the maximum at the center. The influence of the thickness size effect is not considered in the solution process of the analytical method, the three-dimensional components with a certain thickness are equivalent to the plane stress state or plane strain state, so the values read at each node along the thickness direction have little difference. An analytical solution of the stress intensity factor is used to evaluate the fatigue crack growth behavior and fatigue crack growth life of the material, which will produce certain errors, and the evaluation results are not accurate enough. By comparing the stress intensity factors of different thicknesses of the SE (B) specimens in different crack growth stages, it can be found that, with the increase of crack growth length, the gap between the numerical solutions at the center and those at the edge for the specimen thickness direction increases. The difference between the numerical solution and analytical solution of the specimen with the thickness of 15 mm is higher than that of the specimen with the thickness of 10 mm at each stage of crack growth.

The essence of the thickness size effect is the effect of the plastic zone at the crack tip on the stress intensity factor. When a plastic zone exists at the crack tip, the premise of crack growth is that the driving force at the crack tip can promote the cracking of the plastic zone. The plastic zone at the crack tip has the function of absorbing plastic deformation, which can improve the closure degree of the crack tip and reduce the stress concentration at the crack tip, which hinders the crack growth. The larger the size of the plastic zone at the crack tip, the more significant the ability to block the crack growth and the higher the fatigue crack growth life. Therefore, the stress intensity factor decreases with the increase of the size of the plastic zone at the crack tip. When the thickness of the specimen is thin, the crack tip is in the plane stress state; with the increase of specimen thickness, it begins to transition from the plane stress state to the plane strain state; when the specimen is over a certain thickness, the crack tip is in the plane strain state. Similarly, the distribution rule of the stress intensity factor along the thickness direction of the specimen can be understood as the mutual transformation process of the plane strain state and plane stress state. At the center of the specimen thickness direction, the plastic zone of the crack tip is greatly constrained, and the size of the plastic zone is small. At this time, the crack tip is in the plane strain state. From the center of the thickness direction to the edge of the specimen, the constraint from the specimen thickness direction center is gradually reduced at the plastic zone of the crack tip, and the size of the plastic zone is gradually increased. The crack tip is in the plane stress state near the surface. In summary, under the same objective conditions, the thinner the thickness of the specimen, the larger the size of the plastic zone at the crack tip, the more obvious the hindrance of the crack growth, and the smaller the stress intensity factor at the crack tip. Conversely, the thicker the specimen, the greater the stress intensity factors at the crack tip.

5. Conclusion

(1)The fatigue crack growth test results of SE (B) specimens with different thicknesses show that, as the thickness of the specimen increases, the fatigue crack growth life decreases.(2)The comparison between the numerical solution and analytical solution of the stress intensity factor shows that the numerical solution is higher than the analytical solution in all stages of steady-state crack growth, and the gap between the numerical solution and the analytical solution increases with the increase of the specimen thickness. Among them, the difference range of the specimen with thickness 10 mm is 2.531%–6.414%; the difference range of the specimen with thickness 15 mm is 3.961%–8.908%. Because the influence of the thickness size effect on the stress intensity factor is not considered, the value of the stress intensity factor obtained by the analytical method is relatively conservative.(3)Based on the Paris formula, the theoretical derivation is carried out to obtain the numerical solution, and the a-N relation curve obtained by the numerical solution and analytical solution is compared with the test data. The results show that the a-N curve obtained by the numerical solution is closer to the test data, and the analytical method is dangerous for the evaluation of the fracture mechanical properties of materials. The three-dimensional calculation model proposed in this paper can correct the stress intensity factor with thickness size effect.(4)The distribution rule of the numerical and analytical solutions of the stress intensity factor along the thickness direction of the specimen shows that the difference between the numerical solution and analytical solution increases with the increase of crack growth length. Moreover, the thicker the specimen is, the more obvious the difference between the numerical solution and analytical solution is.(5)The essence of the thickness size effect is the effect of the plastic zone at the crack tip on the stress intensity factor. The thinner the specimen thickness, the larger the size of the plastic zone at the crack tip, the stronger the ability of hindering crack growth, and the smaller the stress intensity factor at the crack tip. Conversely, the thicker the specimen, the greater the stress intensity factor at the crack tip.

Data Availability

The figure and table data used to support the findings of this study are included within the article, and the article permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This research was subsidized by the National Natural Science Foundation of China (no. 52005352), the Key Laboratory of Vibration and Control of Aero-Propulsion System, Ministry of Education, Northeastern University (no. VCAME202007), the “Seedling Cultivation” Project for Young Scientific and Technological Talents of Liaoning Education Department (no. lnqn201908), and the Doctoral Research Startup Foundation of Liaoning Province (no. 2019-BS-198).

References

U. Zerbst, K. Mädler, and H. Hintze, “Fracture mechanics in railway applications--an overview,” Engineering Fracture Mechanics, vol. 72, no. 2, pp. 163–194, 2005.

U. Zerbst, M. Vormwald, C. Andersch, K. Mädlerc, and M. Pfuffa, “The development of a damage tolerance concept for railway components and its demonstration for a railway axle,” Engineering Fracture Mechanics, vol. 72, no. 2, pp. 209–239, 2005.

J. C. Newman, E. L. Anagnostou, and D. Rusk, “Fatigue and crack-growth analyses on 7075-T651 aluminum alloy coupons under constant and variable- amplitude loading,” International Journal of Fatigue, vol. 20, no. 4, pp. 142–153, 2013.

S. K. Jha and K. S. Ravichandran, “Effect of mean stress and aging on fatigue crack growth in a metastable titanium alloy, Ti-10V-2Fe-3Al,” Metallurgical and Materials Transactions A, vol. 31, no. 3, pp. 703–714, 2000.

Y. M. Shao, J. Liu, and C. K. Mechefske, “Driving axle housing failure analysis of a mining dump truck based on the load spectrum,” Engineering Failure Analysis, vol. 18, no. 3, pp. 1049–1057, 2011.

T. Ghidini and C. D. Donne, “Fatigue life predictions using fracture mechanics methods,” Engineering Fracture Mechanics, vol. 76, no. 1, pp. 134–148, 2009.

S. Hasunuma, T. Ogawa, and M. Kikuchi, “Fracture mechanics approach to small fatigue crack growth and coalescence under low cycle fatigue,” Journal of the Society of Materials Science Japan, vol. 60, no. 12, pp. 1060–1065, 2011.

S. S. You and D. Fricke, “Advances of virtual testing and hybrid simulation in automotive performance and durability evaluation,” SAE International Journal of Materials and Manufacturing, vol. 4, no. 1, pp. 98–110, 2011.

R. Branco and F. V. Antunes, “Finite element modeling and analysis of crack shape evolution in mode-I fatigue middle cracked tension specimens,” Engineering Fracture Mechanics, vol. 75, no. 10, pp. 3020–3037, 2008.

N. Sageresan and R. Drathir, “Crack propagation in concrete using meshless method,” CMES Computer Modeling in Engineering & Sciences, vol. 32, no. 2, pp. 103–112, 2008.

P. H. Wen and M. H. Aliabadi, “Dual boundary element method for modeling curved crack paths,” International Journal of Fracture, vol. 176, no. 1, pp. 127–133, 2012.

K. C. Nehar, B. E. Hachi, F. Cazes, and M. Haboussi, “Evalution of stress intensity factors for bi-material interface cracks using displacement jump methods,” Acta Mechanica Sinica, vol. 33, no. 6, pp. 1–14, 2017.

B. Zhou and S. Xue, “Displacement extrapolation method for stress intensity factor based on extended finite method,” Mechanics in Engineering, vol. 39, no. 4, pp. 371–378, 2017.

S. Ai, M. Yu, X. Cheng, and J. Wang, “Analysis and application of three-dimensional crack growth based on Franc3D,” Journal of Mechanical Strength, vol. 40, no. 1, pp. 251–254, 2018.

M. K. Raviraj, C. M. Sharanaprabhu, and G. C. Mohankumar, “Effect of specimen crack lengths on stress intensity factor for Al6061-TiC composites using experimental and 3D numerical methods,” International Journal of Structural Integrity, vol. 8, no. 5, pp. 506–515, 2017.

Gb/T. 4161-2007, Metallic Materials—Determination of Plane-Strain Fracture Toughness, General Administration of Quality Supervision, Inspection and Quarantine of P.R.C., Beijing, China, 2007.

Gb/T. 6398-2017, Metallic Materials-Fatigue Testing-Fatigue Crack Growth method, General Administration of Quality Supervision, Inspection and Quarantine of P.R.C., Beijing, China, 2017.

D. Wang, X. Gao, and J. Shi, “Study on relation M-integral and change of total potential energy in three-dimensional solids,” Journal of Water Resources and Architectural Engineering, vol. 7, no. 1, pp. 36–38, 2009.

J. R. Rice, “A path independent integral and the approximate analysis of strain concentration by notches and cracks,” Journal of Applied Mechanics, vol. 35, no. 2, pp. 379–386, 1967.

J. Gong, S. Zhang, L. He et al., “Computation of the stress intensity factor based on the interaction integral method,” Journal of Nanchang Hangkong University: Natural Science, vol. 29, no. 1, pp. 42–48, 2015.

X. B. Lin and R. A. Smith, “Fatigue growth analysis of interacting and coalescing surface defects,” International Journal of Fracture, vol. 85, no. 3, pp. 283–299, 1997.