The prediction of curvilinear fatigue and fracture is a topic for ongoing research. A mathematical study has provided few curved crack path propagation solutions [1, 2], but experimental results can be accurately reproduced in computer ways. The main objective of fracture mechanics is to determine whether or not a structure will fail based on a crack. Crack analysis should start from field stresses in the crack tip, evaluated by the equivalent stress intensity factor (Keq). By comparing the Keq value to a peculiar characteristic of each material called fracture toughness or threshold stress intensity factor, it is possible to determine whether or not a crack member would fail when subjected to static or fatigue loading, respectively. In general, cracks’ initiation and propagation must be associated with stress intensity factors (SIFs) in a complicated state [3,4,5,6,7]. The prediction of curvilinear fatigue and fractures remains a topic of ongoing research. Few curved crack path propagation solutions were obtained by mathematical analysis, but computational approaches can accurately reproduce experimental results [8,9,10]. Numerical simulations typically calculate the energy derivatives in the first order concerning crack length or the equivalent SIFs. Meanwhile, new approaches and methods in several areas of study have been suggested and developed rapidly, including the finite element method (FEM), discrete element method (DEM) [11,12,13], element free Galerkin (EFG) method [14], extended finite element method (XFEM) [15,16,17,18], cohesive element method [19, 20], and phase-field method [21]. Nevertheless, as a versatile method to simulate crack propagation, the FEM is still prevalent. In general, the propagation of an internal crack in the engineering practice is essential because of its significant effect on the quality and stability of engineering structures. Thus, it is important to predict the crack propagation path for safety and reliability of engineering structures. Therefore, in many industries, the accurate estimation of the crack path and fatigue life estimation is of primary importance in terms of reliability. In various applications, such as aerospace manufacturing and the aviation industry, experimental studies are necessary for fatigue analysis. However, due to high costs, an accurate computational approach is required for crack propagation analysis to prognosticate the direction of crack growth and fatigue life in both static and dynamic loading conditions [22,23,24]. The failure is related to (a) the presence of flaws such as interfaces and cracks and (b) the nature of loads that fluctuate. Cracks tend to initiate and propagate when subjected to fluctuating loads until the structure does not bear the load contributing to complete failure. These cracks are considered to be fatigue cracks, and the expected fatigue life is one of the most significant parameters to determine the structure's safety. It can be computed by adding the number of loading cycles needed to propagate the fatigue crack and lead to failure. The calculation of the growth rate of the crack is generally based on the relationship within the range of SIFs and the crack geometry. One way to predict SIFs is using computational methods such as the XFEM. The XFEM introduced by Belytschko and Black in 1999 [25] was widely used in novel investigations. It is based on the conventional finite element (FE) framework. It uses a novel displacement feature allowing discontinuities to happen, to defeat the need to re-mesh continuously during the crack tip development process [26]. Extended work was offered to promote practical models to estimate the fatigue failure and fatigue growth to subdue fatigue failures. Several experimental tests are proposed [27, 28], but the performance of the experimental procedures is usually time-consuming and costly. An effective way to reduce the laboratory work, time, and costs is to incorporate a simulation methodology to including the numerical analysis. Many fatigue crack problems identified to date by the literature use different computational approaches in simulation for simplistic and complicated geometries in 2D and 3D [10, 22, 29,30,31,32,33,34,35]. In this paper, the ANSYS APDL 19.2 is used to precisely predict the mixed-mode stress intensity and correlated fatigue response for seven specimens of the modified compact tension specimen (MCTS). In particular, three methods were widely used to illustrate the fatigue estimation of different specimens: (a) the fracture mechanics based on Paris and Erdogan [36]; (b) the method of strain-life independently proposed by Coffin [37]; and (c) the method of stress-life proposed by Wöhler [38]. The first approach was employed in this study for predicting fatigue life by which the crack tip can be described separately by the SIFs. The primary motivation for this study was to perform a proper comparison of the predicted results for crack propagation parameters using the author's 2D software and the 3D ANSYS software as an alternate tool for simulating fatigue crack growth subjected to mixed-mode loading. The computational time to simulate the crack propagation in the developed 2D program is approximately one-fifth compared with the ANSYS modeling for the same geometry and mesh density.
There are three types of cracks that ANSYS can model: arbitrary, semi-elliptical, and pre-meshed. The Separating Morphing and Adaptive Remeshing Technology (SMART) uses the crack front in the pre-meshed crack approach, and the stress intensity factor is the failure criteria. The rendered node sets will distribute to the front of the crack, as well as to the top and bottom of the crack. By employing SIFs as the failure criteria, the front of the crack is represented in the pre-meshed crack process used by the “Smart Crack Growth” simulation tool. SMART automatically updates the mesh from crack geometry modifications based on crack propagation on each solution level, eliminating the need for long pre-processing sessions. The sphere of influence method around the crack tip that goes through the thickness can be used to improve the mesh around the crack tip. The geometric regions to be described are the crack tip, the crack's top surface, and the crack's bottom surface. Instead, each of these regions is associated with a node-set to be used for analysis. The ANSYS MECHANICAL APDL 19.2 software considers mixed-mode loading where the maximum circumferential stress criterion is implemented. The following is the formula for the direction of crack propagation [32, 33]:
Crack propagation angle according to the second mode of SIFs
In ANSYS Mechanical, the simulation of crack propagation is restricted to region II of the typical fatigue crack growth graph based on Paris’ law that can be represented as:
From Eq. (2), for a crack increment, the number of fatigue life cycles may be predicted as:
The equivalent range for the stress intensity factor formula can be found as follows [39]:
The developed program is a Visual Fortran-based simulation for analyzing two dimensional crack propagation processes under LEFM conditions. This program uses the FE approach to estimate the propagation of quasi-static fatigue cracks in 2D components, taking into account the fracture's mechanical boundary conditions. The adaptive mesh FE method includes four main components for crack trajectory simulation: the mesh optimization technique, the crack criterion, the direction criterion, and the crack propagation strategy. The mesh refinement can be controlled by the characteristic scale of each element, which the error estimator predicts. The step by step analysis is disturbed, and a new FE model is constructed when the error reaches a stated, cumulative error at some point throughout the simulation. Within the initial boundary conditions, the system refines the mesh as needed. After the new mesh is generated, the solution variables (displacement, stresses, strains) are mapped from the previous mesh into the newly created mesh. The analysis is then repeated and continued until the errors exceed the pre-decided number. SIFs are employed in LEFM as fracture criteria. The maximum circumferential stress theory [40], maximum energy release theory [41], and minimum strain-energy density theory [42] are just a few of the methodologies used to predict a crack's propagation direction. The maximum circumferential stress theory has been used in the developed program as a direction criterion.
The most successful engineering implementation of fracture mechanics is crack growth laws based on stress intensity factor ranges [43]. Irwin [44] defined the stress intensity factor, abbreviated as K, as a mechanical parameter that represents the state of stress at a crack point. For various common crack configurations, many analytical solutions for the SIFs were developed. On the other hand, these analytical solutions are restricted to simple crack geometries and loading conditions [2]. Small-sized yield requirements were implemented [45], which established a J-integral methodology for investigating the nonlinear behavior of the material [46]. The equivalent domain integral methodology replaces the finite-size domain with a divergence theorem by integration along the contour for FEM. The integral independent contour of the trajectory is as follows [46]:
For 2D structures, the contour integral i is replaced by area integral (Figure 2B). Then, Eq. (5) is modified as:
The step by step procedure of the developed program software is illustrated in Figure 3. The details of these steps are explained in the literature [29, 45, 47,48,49,50]. The adaptive mesh refinement was explained in detail by Bashiri and Alshoaibi [45].
Computational procedure of the fatigue crack growth developed program.
In LEFM simulation, the J-integral by definition takes into account translational mechanical energy balance in the tip region of the crack along x axis. Eq. (6) simplifies the computation of mode I and II of SIFs but they do not compute
The second method uses the relation between SIFs and
Where
Adaptive mesh refinement is an optimization method used in the field of FE mesh. The a posteriori error estimation used in this scheme is derived from the previous mesh generation. The stress error norm, also known as relative stress norm error, is used to estimate the mesh refinement error. The ratio of element norm stress error to the average standard stress error of the entire region was calculated using h-type adaptive mesh optimization. The mesh size for each element in this scheme is provided as:
However, the average norm stress error over the whole domain is represented as:
This means each element with
Depending on the user numbers of refinement, the existing mesh will be considered as the new background mesh and the advancing front method will be replicated. The element type was triangular with six nodes.
As illustrated in Figure 4, a four holes’ rectangular plate of 10 mm diameter and a 6 mm long edge crack is analyzed. The dimensions of the plate are 100 mm × 100 mm × 10 mm. The plate is supported on the bottom and is subjected to a tension force of 10 MPa. The plate is made of aluminium 7075-T6, whose properties are shown in Table 1.
Cracked plate with four holes and one edge crack.
Mechanical properties of aluminium 7075-T6 [51]
Modulus of elasticity, E | 72 GPa |
Poisson's ratio, |
0.33 |
Yield strength, |
469 MPa |
Ultimate strength, |
538 MPa |
Fracture toughness, |
3,288.76 MPa
|
The geometry is fixed in the x and y directions at the bottom face in the developed program, while in Ansys the geometry is fixed in the x, y, and z directions.
The meshed geometry created in Ansys workbench has 14,737 nodes and 7,177 elements, while the meshed geometry generated in the developed program has 3,933 nodes and 1,911 elements, as represented in Figure 5.
Initial generated mesh
The predicted crack growth paths in both software is shown in Figure 6, and these were also in good agreement with that path predicted by Knowles and Sternberg [51] using the fast multipole method (FMM) boundary element method.
Crack growth path:
The maximum and minimum principal stresses were compared side by side in both programs, as represented in Figure 7. In Figure 7, the crack starts to grow in a straight line as there is no influence of any holes. As it approaches the hole, the path propagates towards the hole. The crack deviated its direction because the first upper hole was not close enough to attract the crack to sink in the hole. It continued to propagate in a straight line up until it reached near the second upper hole, close enough to the crack trajectory so that the crack sinks into this hole. The predicted estimation of the mixed-mode stress intensity parameters
Crack growth trajectory
Figure 9 shows the comparison for the stress contour distribution in y direction between the developed program results, Ansys results, and the results obtained by Knowles and Sternberg [51] using FMM boundary element method.
Comparison for the stress contour distribution in y direction
A more complex mixed-mode crack growth problem was investigated to illustrate the performance of the developed program to represent the propagation of cracks from beginning notches in a three holes single edge cracked plate and precisely prognosticate the values of SIFs. The geometry of three holes single edge cracked plate, and its final adaptive mesh, is presented in Figure 10. This case was investigated for the first time by Ingraffea and Grigoriu [53]. The results revealed that small changes in the crack path's initial position significantly impacted the results. As a benchmark for numerical models, this example is commonly used to assess the system's predictive abilities. As indicated in Table 2, the geometry analyzed has 2L = 508 mm, width = 203.2 mm, and thickness = 12.7 mm, with two alternative configurations depending on the initial crack length and position. The initial meshes for this geometry, generated by the developed program and Ansys before starting the simulation, are shown in Figure 11.
Geometry and boundary condition of the three holes’ single edge crack plate.
Configurations of the three holes’ single edge cracked plate
1 | 25.4 | 152.4 |
2 | 38.1 | 127 |
Initial mesh of the three holes’ single edge cracked plate in
The length of the crack in the first configurations is a
Mechanical properties of the three holes single edge cracked plate [22]
Modulus of elasticity, |
205 GPa |
Poisson's ratio, |
0.3 |
Yield strength, |
516 MPa |
Fracture toughness, |
730 MPa
|
Threshold SIF, Δ |
80 MPa
|
Paris’ law coefficient, |
1.2 × 10−11 |
Paris law exponent, |
3 |
Specimens are constrained to all degrees of freedom from the left pin position, similar to the experimental specimen. In contrast, the correct pin only has a first degree of freedom of translation UX. For the first configuration of specimen number one, the simulated crack propagation using the developed program and Ansys software is shown in Figure 12 with an identical crack growth trajectory. The crack growth trajectory goes between the lower and middle holes in this specimen and approaches the middle hole from the right side. It implies that the shear stress intensity factor (KII) surrounding the holes has significantly increased, forcing the crack step-sizes to be decreased.
Predicted crack growth path for specimen 1,
The crack trajectory using the developed program and Ansys ends up in excellent agreement with experimental data presented in the literature [47, 48], numerical results using the coupled extended meshfree-smoothed meshfree method presented by Knowles and Sternberg [51], and XFEM results obtained by Zhang and Tabiei [21]. Numerical results obtained by Kanth et al. [14] using a polygonal XFEM with numerical integration as shown in Figure 13 A, B, C, D, E, and F, respectively.
Specimen 1, crack growth trajectory
The results of this simulation for dimensionless stress factors by using the developed program and Ansys are shown in Figure 14, which are almost identical to each other. As can be seen in this figure, the SIFs increase in magnitude as the crack approaches the hole, demonstrating the hole's influence on the crack trajectory. For the second configuration of specimen number two, the simulated crack propagation using the developed program and Ansys software is shown in Figure 15 with an almost identical crack growth trajectory. The crack in this specimen propagates above the lower hole to the middle hole's left side.
Dimensionless SIFs for specimen 1. SIFs, stress intensity factors.
Predicted crack growth path for specimen 2,
Figure 16 depicts the direction of crack growth for specimen 2, which was compared to the experimental crack trajectory performed in the literature [47, 48] and the numerical findings utilizing the coupled expanded meshfree-smoothed mesh-free technique provided by Knowles and Sternberg [51], as illustrated in Figure 16 A, B, C, and D respectively.
Specimen 1, crack growth trajectory,
The predicted values of the two modes of SIFs in the developed program and Ansys are depicted in Figure 17 with excellent agreement. As observed in the figure, the lower hole has no influence on the crack growth direction, since the initial crack was far from the lower hole, while the middle hole has a stronger influence, attracting the crack to sink into the hole.
Dimensionless SIFs for specimen 2. SIFs, stress intensity factors.
Regarding the simulation results of two- and three-dimensional analysis using the developed program and Ansys software, a satisfactory agreement was obtained between numerical and experimental crack propagation growth trajectories for the simulated specimens. However, for the first geometry there was a differential between the 2D and 3D
The three-dimensional simulation is preferable to the two-dimensional simulation to visualize better stresses, strains, deformation, and many other characteristics, despite the two-dimensional representation being sufficient.