Abstract
Probabilistic damage tolerance assessment is an essential method to evaluate the safety of aeroengine rotors. The stress intensity factors (SIFs) are the core parameters. The weight function method can calculate SIF efficiently. The available weight functions for corner cracks are suitable for cracks with universal stresses. However, when cracks are under bivariant stress distributions, the lack of the weight function database makes the damage tolerance assessment impractical. Therefore, this paper derives the point weight functions for corner cracks, which are suited for cracks in real-world components with two-directional stress distributions. The response surface method is used to build the surrogate model of the point weight functions. By integrating the weight functions with the stress distributions, SIFs can be calculated with high accuracy. In sum, 81% of the differences between the point weight function method and finite element results are less than 10%.
1. Introduction
Aeroengine rotors are typical life-limited parts [1] that operate in complex environments with remarkably high temperatures, pressures, and rotational speeds [2, 3]. The probabilistic damage tolerance assessment methodology is applied to address the threat of defects that may occur anywhere in the components [4–7]. The probability of failure of these rare defects is predicted based on a generalized stress-strength interference model, which requires the simulation of the crack propagation process under complex loads [8, 9]. The stress intensity factor (SIF) is the core parameter assessed during a crack propagation simulation, as it reflects the severity and stress concentration around a crack [10].
SIF solution research has been a central topic in fracture mechanics in recent decades. Various analytical and numerical methods for determining SIFs have been proposed. Representative SIF solutions of through cracks in two-dimensional geometries have been generated and collected in several well-known handbooks [11]. As to cracks in three-dimensional geometries, Newman and Raju [12, 13] have proposed SIF solutions for the corner, surface, and embedded cracks in three-dimensional geometries. Unfortunately, these solutions are expressed in simple load cases, including uniform remote tension and bending stress. Instead, finite element analysis or boundary element method is developed to obtain the SIFs of the cracks with arbitrary stresses. Nevertheless, the detailed three-dimensional numerical approach only calculates the singular SIF of the crack with specific loading and geometry with a complex preprocess to model and mesh the geometry. Thus, these numerical methods hardly meet the requirements of the highly efficient probabilistic damage tolerance assessment because millions of crack propagation processes need to be calculated [14–16].
To efficiently obtain the SIFs for any specific cracked body geometry that is applied with an arbitrary stress distribution, the weight function method was envisioned by Bueckner [17] and studied further by Rice [18]. The weight function, also called Green’s function of a crack surface, is determined by the geometry of the cracks in the plate. The SIFs are calculated efficiently by integrating Green’s function and the stress distribution on the crack surface. The standardized analytical weight function method was developed mainly concentrating on the through-thickness cracks in two-dimensional geometries [19, 20]. For elliptical embedded cracks and semielliptical surface cracks in three-dimensional geometries, the universal weight function was proposed [21–23]. These two weight functions are widely used to calculate SIFs with univariate stress distributions [24, 25].
Nevertheless, these weight functions are limited when the cracks are under complex bivariant stress distributions [26, 27]. The corner cracks are more likely to bear complex bivariant loadings because of the stress concentration, thermal loading, and machining or strengthening residual stress. Oore and Burns [28] proposed the point weight function method to calculate SIFs under bivariant stress distributions. Orynyak et al. then developed point weight functions for elliptical cracks [29, 30]. Finally, the point weight function method was applied on the corner cracks at holes for calculating SIFs under bivariant stress gradients and implemented in production software for design and life management [31–33]. However, few efforts have been devoted to explaining the details for establishing point weight functions for corner cracks in three-dimensional finite rectangular plates. The lack of the weight function database makes the probabilistic damage tolerance assessment or reliability analysis of corner cracks with arbitrary complex loadings impractical. Moreover, previous research ignored analyzing the accuracy of the point weight function method. To ensure the reliability of the derived point weight function database, the quantified error of the point weight function is supposed to be considered when calculating crack propagation process.
This study derives the point weight functions for corner cracks to calculate SIFs of corner defects under arbitrary stress distributions and help conduct the probabilistic damage tolerance assessment [25]. The challenge of the point weight function method is to derive the coefficients of the weight functions determined by the geometric parameters of a structure with cracks. The response surface method [34, 35] is an effective way to derive coefficients in the point weight functions. The errors induced by the fitting method are less than 1%. Then, the SIFs for corner cracks calculated by the point weight function method are validated. When bivariant stress distributions are applied onto the crack surface, the differences between the point weight function method and finite element results are less than 5 when and . Meanwhile, 81% of the relative errors are less than 10%. The errors of the SIFs calculated by the point weight function method may be induced by the response surface method and the accuracy of the reference SIFs. The errors induced by the response surface method are less than 1%. Moreover, the 1% error of the reference SIFs can cause a 1.25% deviation of the calculated SIFs.
This paper is divided into five sections. The theoretical formula of the point weight function is analyzed in Section 2. Section 3 introduces the derivation of the coefficients in point weight functions determined by multiple geometric parameters, including the length and thickness of a flat plate. In Section 4, the accuracy of the point weight function is quantified. Finally, the conclusions are summarized in Section 5.
2. Point Weight Function Theory for Bivariant Stress Field
Stress distributions of structural components are complicated such as aeroengines. These complex stress gradients are mainly induced by local stress concentrations or multiple loads, including centrifugal, thermal, and residual stresses. When the stresses vary in bivariant directions (vertical and along with the crack depth), the point weight function can help obtain accurate SIFs of cracks. The SIF considered in this paper is a mode I corner crack in rectangular geometry with finite length and width, as shown in Figure 1. The SIFs for two points along the corner cracks, and , are critical for the crack propagation calculation [31].
When the point weight function is known, the SIF, the core parameter used in probabilistic damage tolerance analysis, can be obtained by integrating the product of the load with the weight function for the crack surface. This method considers two variations compared with the universal weight function method, which only considers stress distribution changes in the direction.
The weight function is called Green’s function of a crack surface. The SIF is given by [27] where is the stress perpendicular to the crack plane in the corresponding uncracked body and is the point weight function. The integration is carried out for the crack dimensions, characterized by in the direction and in the direction. Both and are the radii of the corner crack. Only circular cracks are considered in this paper; thus, .
The point weight functions for quarter-circular corner cracks are, respectively, given by [26] where denotes the location where the point load is applied. The two points and do not exist on the crack surface. The locations of and are those symmetrical to the location of point load with respect to the - and -axes. is the point along the crack perimeter. and are the values of the -axes and -axes. represents the distance between the point along with the crack front to the ellipse center , represents the distance between the point in the crack surface to the ellipse center is the distance between and is the distance between and , and is the distance between and . The length parameters, , , and , correct the free surface effects for a surface crack by assuming a symmetrical stress distribution for an imaginary prolonged crack extending into the other semicircles. , , and are the critical coefficients of the point weight function, as determined by the crack size and geometry parameters.Π
Equation (2) can be expressed in cylindrical coordinates where and are the angles of and .
With the integration of stress distribution and point weight function, the SIF of along the surface crack can be obtained by
The SIFs for points and are concerned . For each point, by using three known SIFs and the associated three reference stress distributions (), six simultaneous equations are given as
The coefficients , , and can be determined by solving Equation (5). Finally, the point weight functions can be derived, and the SIFs of corner cracks under unknown stress distributions can be calculated by integrating the weight functions and stress distributions as Equation (4).
3. Derivation of the Point Weight Functions for Corner Cracks
3.1. Reference Stress Intensity Factors
, , , , , and are the critical coefficients of the point weight functions obtained by solving Equation (5). The construction of the equations leads to the acquisition of the reference SIFs under reference stress distributions. Three-dimensional finite elements are used to model the symmetric structure containing a quarter-circular corner crack to obtain the reference SIFs. Figure 1 shows the geometry and the coordinate system of the corner cracks. and are the thickness and length of the finite plate, respectively.
Finite element analysis is performed using ABAQUS version 6.14 with 20-node isoperimetric three-dimensional solid elements and reduced integration. As shown in Figure 2, three-dimensional prism elements with four midside nodes at the quarter points are used to model the square root singularity at the crack tip, and the separate crack tip nodal points are constrained and assumed to have the same displacement.
The SIFs for two points of corner cracks, and , shown in Figure 1, are critical for the crack propagation calculation [31]. It is worth mentioning that the reference solutions are determined as the SIFs taken at 2.5° inside the two points and . This approach aims to avoid the numerical ambiguities resulting from the corner point crack singularity associated with the calculation of the SIF directly at the surface.
The SIF is calculated from the -integral by using the virtual crack extension method. The analysis is conducted with a linear classic elastic material model with Young’s modulus (210 GPa) and Poisson’s ratio (0.361). The relationship between and SIF can be used to calculate the reference SIF where is the equation for the plane stress and is the formula for the plane strain.
The reference loads are applied directly onto the crack surface. Three types of loading are applied to each crack geometry, with the stress distribution given by where is a constant, in which the value of 400 MPa is assigned and is the crack radius whose value is 0.003 m in this paper.
3.2. Validation of the Finite Element Method
The accuracy of the finite element results needs to be validated because the reference SIFs are critical inputs for the point weight function method. Therefore, the SIFs calculated by ABAQUS are compared with FADD3D results [32] to validate the accuracy of the reference SIFs. The univariant stress distribution is given by
The bivariant stress distribution is given by
The bivariant stress distribution decreases to zero at both surface tips, while the univariant stress distribution is uniform in the -direction and applies a much larger resultant force on the crack plane than the bivariant stress field.
Figure 3 shows the SIF results along the crack tip calculated by FADD3D and ABAQUS. The circular crack is analyzed, and . The geometry of the finite plate is set as and . is the shape factor associated with the ratio, defined by
Notably, the predicted results exist sharp changes at the intersection of the crack perimeter with a free surface ( and tips) because of the complex nature of the singularity at that point. Therefore, the solution values for and tips are generally taken at the first node inside the free surface, 2.5° from the free surface.
MAPE is the average of the absolute relative errors of the actual and predicted values. The smaller value of MAPE indicates effective model prediction. where is the SIF results calculated by the weight function method, is the SIF results calculated by finite element method, and represents the number of the SIF values.
The maximum error between these two numerical methods is less than 3.12%. The average absolute relative error MAPE of universal stress distribution between ABAQUS and FADD3D is 1.37%, while the average absolute relative error MAPE of bivariant stress distribution is 1.80%. The validation ensures the accuracy of the reference SIFs in this paper.
3.3. Weight Function Derivation
The SIFs are determined by the geometry size parameters and . Supposing is infinite, and the other three variables are then regarded as critical parameters and normalized as . The coefficients of weight functions are also determined by the geometry size parameters , given by where , 0.85, 0.75, 0.6, 0.3, 0.12, and 0.06. is the same as .
The discrete data are calculated by finite element analysis. The response surface method [34, 35] is applied to obtain the relationship between the coefficients and the normalized size parameters. The theoretical formula is given by where is the highest order of polynomial; represent the undetermined coefficients; represent the normalized size parameters; and represent and , respectively; and is the error term. The undetermined coefficient is supposed to be obtained by the regression algorithm, which is denoted as . Thus, Equation (14) can be expressed by
3.4. Stress Fit Formula
The SIF denoted by can be expressed by integrating the stress and the weight function as Equation (4). In this case, the bivariate stress function is another critical input of the point weight function method. After a finite element analysis of rotor disks, the stress distribution is extracted. Then, the stress fit process can help obtain the stress distribution function, given by where and represent the coordinate value on the crack surface, and are the degrees of a polynomial function, and is the coefficient. The Cartesian coordinates are converted to cylindrical coordinates as
4. Validation of the Point Weight Functions
4.1. Comparison with the Universal Weight Function Method
The universal weight function derived by Glinka is suitable to calculate the stress intensity factor (SIF) for quarter-elliptical corner crack [23] in a finite thickness plate. As shown in Figure 4, are the plate’s thickness, width, and height. The universal weight function derived by Glinka is suitable when are infinite.
The universal weight function method takes univariant stress distribution into consideration. Thus, two types of stresses are applied to the crack surface, given by
Then, the SIFs are calculated using the point weight function method and universal weight function method [22, 23]. Meanwhile, the finite element method results calculated by ABAQUS are compared with these two weight function methods, as shown in Figure 5.
(a) Comparison of SIF results when cracks under stress distribution Equation (17)
(b) Comparison of SIF results when cracks under stress distribution Equation (18)
Figure 5 illustrates that the point weight function method and the universal weight function method fit the finite element method results well. The MAPE between the point weight function and the universal weight function method is less than 0.64% for point and 1.46% for point . Therefore, the accuracy of the point weight function method is verified under univariant stress distributions based on the comparison of the point weight function method results and the universal weight function method results.
4.2. Comparison with the Finite Element Results
4.2.1. Stress Distributions
To validate the SIF results calculated by the point weight function derived in Section 3, arbitrary stress distributions are applied onto the crack surface to replace the reference stresses (). Fourteen stress distributions along the and directions are applied onto the crack surface, as shown in Figure 6. The stress gradients are listed in Table 1. is assigned 400 MPa, and is assigned 0.003 m in this paper. Stresses 1-5 change in the direction as shown in Figure 6(a). Only these stress distributions are concerned by the universal weight function method. Stresses 6-9 change in the direction as shown in Figure 6(b), while stresses 10-14 are the bivariate gradients shown in Figures 6(c)–6(g). The point weight function method can calculate the SIFs when there are cracks under these three types of stresses. The varies associated with two dimensionless geometric ratios and . To obtain the weight function coefficients when and , 49 geometries are considered when , 0.85, 0.75, 0.6, 0.3, 0.12, and 0.06 and , 0.85, 0.75, 0.6, 0.3, 0.12, and 0.06.
(a) Stresses 1-5
(b) Stresses 6-9
(c) Stress 10
(d) Stress 11
(e) Stress 12
(f) Stress13
(g) Stress 14
4.2.2. Validation Results
The varies associated with two dimensionless geometric ratios and . To obtain the weight function coefficients when and , the surrogate model is established as Equation (10), and the surrogate model of the point weight function coefficients is shown in the appendix. The SIFs calculated by the obtained point weight functions with (surrogate model) are compared with finite element results to validate the accuracy of the weight function method. The SIF values are normalized by
These 49 geometries with different sizes are numbered as shown in Figure 7(a). The SIF results under stresses 1-14 are calculated, and Figures 7(a)–7(l) show the results of stresses 1, 5, 6, 7, 10, and 14. The normalized SIFs calculated by the finite element method for point B and point A of the 49 geometries are marked with hollow and solid points, respectively. With the surrogate model, the SIFs are calculated and shown as three-dimensional surfaces ( and ) in Figures 7(a), 7(c), 7(e), 7(g), 7(i), and 7(k). Figures 7(b), 7(d), 7(f), 7(h), 7(j), and 7(l) show the relationship between the finite element results and the weight function method results of the 49 geometries.
(a) 3D SIF results under stress 1
(b) 2D SIF results under stress 1
(c) 3D SIF results under stress 5
(d) 2D SIF results under stress 5
(e) 3D SIF results under stress 6
(f) 2D SIF results under stress 6
(g) 3D SIF results under stress 7
(h) 2D SIF results under stress 7
(i) 3D SIF results under stress 10
(j) 2D SIF results under stress 10
(k) 3D SIF results under stress 14
(l) 2D SIF results under stress 14
Results show that the SIFs are determined by the geometries with the same loading. When and are close to 1, the SIFs are 3-6 times larger than SIFs when and are close to 0. The validation result also shows that the SIFs calculated by the two methods are consistent when and are small.
The relative error between the weight function method and finite element method is defined as where is the SIF results calculated by the weight function method and is the SIF results calculated by a finite element method. Similarly, the absolute error between the weight function method and finite element method is defined as
The specific errors between the weight function method and finite element method results of these 14 stresses are shown in Table 2. The relative errors are larger than 100% when the absolute SIF value is small, less than 1 . Most of the differences between the weight function method and the finite element results are less than 5 when and except for stress 3, 10, 11, and 13 cases. The errors are related to the monotonicity of the stress distribution. Stress 1 can represent most stress concentration conditions that occur in rotor disks.
The cumulative distributions of and are shown in Figures 8(a) and 8(b). Results show that 81% of point and 54% of point are less than 10%. As shown in Figure 8(b), for point , 97% are less than 5 MPa while 93% for point . The absolute values of point are usually larger than point ; thus, the relative error of point is more minor.
(a) Relative errors
(b) Absolute error
In sum, the derived point weight functions can help obtain SIFs of corner cracks in disks under bivariant stress distributions. The relative error may increase when the absolute SIFs are less than 1, and 81% of the relative errors for point and 54% for point are less than 10%. Meanwhile, 90% of the absolute errors between the point weight function and finite element method results are less than 5 .
4.3. Error Analysis
4.3.1. Errors Induced by the Surrogate Model
Based on the finite element analysis mentioned in Section 3.1, Forty-nine geometries with different sizes (, 0.85, 0.75, 0.6, 0.3, 0.12, and 0.06) are established, and the associated discrete weight function coefficients are obtained by solving Equations (5). But to obtain the point weight function coefficients when and , the surrogate model is established by the response surface method. The surrogate model will induce errors in the SIFs.
The hexagonal polynomials () fit well, and the coefficient of determination (-square) value reaches 0.9995. The stress distributions 1-14 mentioned in Table 1 are applied on the crack surface, and the SIFs of the 49 geometries are calculated. The relative error of the SIFs calculated by discrete and surrogate model is given as where is the SIF results calculated by the surrogate model and is the SIF results calculated by discrete . The cumulative distribution of is shown in Figure 9. Results show that 93% is less than 1%. The response surface method helps obtain the surrogate model of weight function coefficients with relatively small errors.
4.3.2. Errors Induced by the Reference Stress Intensity Factors
The accuracy of the reference SIFs directly ensures the effectiveness of the as shown in Equations (5). It is assumed that the three reference solutions for point , , and deviate to the actual finite element results ranging from -20% to 20%. Thus, the three reference solutions are given as where is the offset scale and ranges from -20% to 20%.
The proportion of the new coefficients derived by the offset reference solutions are defined as where are the coefficients derived by the deviated reference solutions.
Similarly, the proportion of the SIFs calculated by the new coefficients are given as where are the SIFs calculated by the weight function method. Equations (23)–(25) also apply to point .
The show a multivariate linear correlation of the offset scale in Figure 10 because the are solved by Equation (5). When the reference SIFs have a 1% error, the will change to 5 times of the original value. The SIFs calculated by the point weight function also change, as shown in Figure 11. When the reference SIFs have a 1% error, the SIFs for points and have an error up to 1.25%. The coefficients of the point weight function are highly dependent on the reference SIFs. Thus, the precise reference SIFs are significant to the point weight function method.
(a) Point
(b) Point
The errors of the SIFs calculated by the point weight function method may be induced by the response surface method and the accuracy of the reference SIFs. First, the response surface fitting process of the discrete weight function coefficients causes the deviation of the SIFs. By fitting the coefficients with high order polynomial, the errors are less than 1%. Second, the accuracy of the reference SIFs calculated by the finite element method will influence the results of the weight function coefficients up to 5 times since these coefficients are calculated by solving simultaneous equations. A 1% error of the reference SIFs can cause a 1.25% deviation of the calculated SIFs.
In addition to these two factors, the selection of the types of the reference stress distributions influences the accuracy of the point weight function. Since the uniform and linear reference stresses are used in this paper, the SIF results of linear distributions are accurate. Further research is necessary to study the influence of the selection of the reference loadings. Moreover, the point weight function equations can approximate Green’s function instead of equaling to Green’s function completely. The SIF errors are inevitably induced because the weight functions deviate from Green’s function, which is the limitation of the weight function method itself.
In this case, during the crack propagation calculation and damage tolerance assessment, the SIF errors are supposed to be considered. Meanwhile, the development of the accurate finite element solutions is highly essential, which can provide the reference solutions to the point weight function method. Apart from the weight function method, a new surrogate model between the SIFs and the stress distributions of different geometries with cracks may be developed by other methods such as the machine learning method with enough SIF data.
5. Conclusions
This paper derives the point weight functions for corner cracks in finite plates to efficiently calculate SIFs when cracks are subjected to complex stress distributions induced by thermal, residual, and other types of stresses and conduct the probabilistic damage tolerance assessment. The conclusions are as follows: (1)Based on the point weight function for the cracked structure, the SIFs for corner cracks can be obtained by integrating the two-directional stress in the uncracked body. In most cases (more than 90%), when the stress distributions are close to that of the rotor disks, the differences between the point weight function method and finite element results are less than 5 when and . The relative error may increase when the absolute SIFs are less than 1 , and 81% of the relative errors for point and 54% for point are less than 10%(2)Obtaining the coefficients in the weight function, which are determined by the geometric parameters of the plates, is the critical process for the point weight function method. The response surface method is efficient for deriving the point weight function to reflect the relationship between geometric parameters and SIFs. The errors induced by the fitting method are less than 1%(3)The reference stress intensity factors under associated loadings are essential input for the point weight function method because the unknown point weight function coefficients are calculated by solving simultaneous equations. The finite element method can help obtain the reference stress intensity factors. The accuracy of the reference SIFs will influence the results of the weight function coefficients. A 1% error of the reference SIFs can cause five times of deviation of the coefficients and a 1.25% deviation of the calculated SIFs. Thus, the accurate SIFs are the basis of the point weight function
Appendix
The Polynomials of the Point Weight Function Coefficients of Corner Cracks
The coefficients for the normalized SIFs are shown in Table 3.
Nomenclature
| : | Stress intensity factor |
| : | Stress distribution |
| : | Weight function |
| : | Crack size |
| : | Crack size |
| : | Crack radius |
| : | Radius of point |
| : | Coefficients of point weight function |
| : | Length parameter of |
| : | Length parameter of |
| : | Length parameter of |
| : | Angle of |
| : | Angle of |
| : | Young’s modulus |
| : | Poisson’s ratio |
| : | Plate thickness |
| : | Plate length |
| : | Estimated response surface method coefficients |
| : | Response surface method coefficients |
| : | Normalized size parameter |
| : | Error term |
| : | -axis value |
| : | -axis value |
| : | Coefficients in stress function |
| : | integral value |
| : | Coefficient of the stress distribution |
| : | Relative error of fitting method |
| : | Relative error of weight function method |
| : | Absolute error of weight function method. |
Data Availability
All data are attached in Appendix in the manuscript.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
The work was funded by the National Natural Science Foundation of China and Civil Aviation Administration of China. The work was supported by the Innovation Team of Complex System Safety and Airworthiness of Aero-Engine from the Co-Innovation Centre for Advanced Aero-Engine of China.