Abstract

Parameters identification of cracked rotors has been gaining importance in recent years, but it is still a great challenge to determine the crack parameters including crack location, depth, and angle for operating rotors. This work proposes a new method to identify crack parameters in a rotor-bearing system based on a Kriging surrogate model and an improved nondominated sorting genetic algorithm-III (NSGA-III). A rotor-bearing system with a breathing crack is established by the finite element method and the superharmonic components are used as index to detect the cracks, the Kriging surrogate model between crack parameters and the superharmonic component amplitudes of the vibration response for rotors are constructed, and an improved NSGA-III is proposed to obtain the optimal crack parameters. Numerical experiments show that the proposed method can identify the crack location, depth, and angle accurately and efficiently for operating rotors.

1. Introduction

Rotors are one of the most important parts of rotating machinery, which are widely used in power stations, aircraft engines, motors, and so on. Transverse fatigue crack is easily produced on rotors under bending excitation, and slant fatigue crack is generated under torsional excitation due to the harsh operation conditions. The crack propagation speed might be extremely slow at the initial stage, as the crack has propagated to a certain depth, the crack propagation speed increases sharply and the shaft may malfunction within a few hours of operation, causing a catastrophic accident [1]. In order to reduce the maintenance costs of equipment and avoid the occurrence of serious accidents, it is of enormous significance to continuously monitor the health of rotors during operation.

Many researchers have studied crack monitoring for rotors based on vibration signals [26], the methods can be split into two groups depending on whether the rotors are rotating or not. The crack parameter identification methods under nonoperating state of rotors are usually based on modal parameters of the cracked rotors. Dong et al. [7] established a rotor system based on wavelet finite element method and determined the position and depth of the crack by the intersection of contours of the first three natural frequencies. Haji and Olutunde Oyadiji [8] put forward the concept of orthogonal natural frequencies (ONFs) and achieved the location of cracks by calculating the ONFs at different positions of a disc. Zapico-Valle et al. [9] and Yu et al. [10] used the modal parameters of rotors corresponding to different crack positions and depths as the neural network input to identify rotor crack parameters. Xiang et al. [11] estimated the crack parameters by minimizing the errors of natural frequencies between simulation and experimental results. But the modal parameters used in crack identification cannot be easily obtained in operating conditions, and it is more meaningful to identify the crack parameters for rotors under operating conditions.

Generally, crack parameters identification methods for operating rotors need to test the dynamic response of the rotor. Prabhakar et al. [12] detected rotor crack during start-up and shutdown of a rotor-bearing system qualitatively with continuous wavelet transform. Singh and Tiwari [13] detected the slope discontinuity in the elastic line of the shaft due to cracks and achieved localization of multiple cracks in a stepped shaft. Saravanan and Sekhar [14] identified the crack in a rotor-bearing system utilizing the concept of operational deflection shape (ODS) and kurtosis of vibration response. Zhang et al. [15] detected the singularity of the ODS in rotors to identify the crack location and used the approximate waveform capacity dimension (AWDC) to make the singularity more prominent. Ramesh Babu and Sekhar [16] proposed a concept called amplitude deviation curve method based on ODS and realized the two cracks localization in a rotating rotor. Lu et al. [17] proposed a multicrack location method based on proper orthogonal decomposition (POD) using fractal dimension (FD) and gapped smoothing method (GSM). Lees et al. [18] reviewed the rotor crack parameter identification method based on equivalent crack forces, which treated the effects of cracks as equivalent forces exerted on a healthy rotor system. Some researchers used equivalent crack forces method to obtain the crack position and depth of rotors, such as Sekhar [1921], Markert et al. [22], Pennacchi and et al. [23]. The second group of studies concentrates on crack detection and location for operating rotors, but it is difficult to quantitatively identify the shape and depth of the crack.

As the crack angle affects the stress, the speed of crack propagation, and fatigue life of rotors, some researchers studied the dynamic behavior of the rotor due to the slant crack [2428]. However, there is no relevant research on crack parameters identification if crack angle is considered for operating rotors at present, transverse crack and slant crack have similar stiffness, and crack detection results could be misleading. Focusing on this challenge, the problem of crack parameters identification is transformed into multiobjective optimization problem. Then, NSGA-III [29, 30] is improved to obtain the optimal crack location, depth, and angle at the same time, and the Kriging surrogate model [31, 32] is applied to increase the optimization speed.

In this paper, a new crack parameters identification method is proposed for operating rotors using Kriging surrogate model and an improved NSGA-III. In Section 2, the model of a rotor-bearing system with breathing cracks is constructed by the finite element method, and the influence of different crack parameters on the amplitude of superharmonic components is discussed. In Section 3, the Kriging surrogate model between crack parameters and the superharmonic component amplitudes is constructed. In Section 4, using an improved NSGA-III to predict the optimal crack parameters based on the Kriging surrogate model by minimizing the multiobjective function values related to the superharmonic components amplitudes. Numerical experiments show that the proposed method can identify the crack location, depth, and angle for operating rotors accurately and efficiently.

2. Superharmonic Components Analysis of a Rotor-Bearing System with a Breathing Crack

2.1. Equations of Motion of a Cracked Two-Disc Rotor-Bearing System

As shown in Figure 1, a two-disc rotor-bearing system is established by the finite element method, and the physical parameters are given in Table 1. The shaft is discretized into 60 Timoshenko beam elements with six degrees of freedom at each node, the two discs are assumed as rigid bodies with three translational and three rotational inertias, and the gyroscopic effects of which are considered. The ball bearings are considered as linear springs and dampers.

By assembling the cracked shaft element, the uncracked shaft element, the rigid discs, and the bearings, the equation of motion of the cracked rotor system in a fixed coordinate system can be written aswhere and are the mass matrices of the shaft and the two discs, respectively, is the shaft Rayleigh damping matrix, is the damping matrix of the bearings, is the rotational speed of the shaft, is the shaft gyroscopic matrix, is the gyroscopic matrix of the two discs, is the displacement matrix of the node, is the stiffness matrix of the shaft, is the stiffness matrix of the bearing, is the stiffness matrix of the crack element, is the excitation due to static unbalance of discs, is the gravitational force, and is the external excitation during operation.

Taking into account the bending and torsional coupling excitation caused by the eccentricity of the discs, the excitation acted upon the discs is

When rotors are rotating at a constant speed , the elements in can be expressed as [33]where is the torsional angle, is the unbalance orientation angle of the disc, and is the disc eccentricity.

2.2. Model of the Shaft Element with Breathing Crack

Figure 2(a) shows a cracked shaft element of length L and radius R, a is the crack depth, is the crack width, are resultant forces acted upon the 12 nodes, is the global coordinate system, and is the local coordinate system. is the distance between the center of the crack and the left end of the element, and is the crack angle between the cross section and the shaft centerline. The most probable angle is 45° for a slant crack and 90° for a transverse crack in practical engineering [25, 26, 34], and so the crack angles studied in this paper are 45° and 90°.

The flexibility matrix of the uncracked element and the additional flexibility matrix of the cracked element can be derived based on SERR theory [4], and the detailed expressions can be obtained from [24].

is the total flexibility matrix of the crack element, given by

The stiffness matrix of the cracked element and uncracked element can be represented aswhere is the transformation matrix written as

So as to simulate the breathing behavior of the crack, the stiffness of the crack element is calculated based on the classical theory of crack closure line (CCL) [24], which is a hypothetical line separating the closed and open parts of the crack (Figure 2). The crack front is subdivided into m parts, the opening mode stress intensity factor is calculated by Equation (7), and a positive corresponds to the open crack state and a negative to the closed state. Thus, the location of the CCL is determined, and the stiffness matrix of the crack element can be obtained as

The Newmark method [35] is utilized to solve the equations numerically, setting the Newmark parameters  = 0.5 and  = 0.25 to ensure numerical stability.

2.3. Superharmonic Components Analysis of the Cracked Rotor System

The superharmonic components in subcritical speed region can be used as index to detect the cracks in the rotor-bearing systems [36, 37]. The vertical vibration response of the rotor-bearing system is collected at four measuring points (illustrated as A, B, C, and D in Figure 1), and the rotating speed is 1/3 of the first critical speed (840 r/min) of the system. The 1X, 2X, and 3X superharmonic components amplitudes of the vibration response at the measuring point A are obtained through fast Fourier transform (FFT).

Figure 3 depicts the amplitude curves of superharmonic components of the vertical vibration response with different crack depths as the crack is located at the 21st element of the shaft. It shows that the 2X amplitude increases with increasing crack depths, the 1X and 2X amplitudes resulted from transverse crack are greater than those from slant crack, and the 1X, 2X, and 3X amplitudes by transverse crack are close to those by slant crack as the crack depth is small.

Figure 4 depicts the amplitudes curves of superharmonic components of the vertical vibration response with different crack location as the crack depth is 0.25 cm. It shows that the 1X, 2X, and 3X amplitudes generated by transverse crack are greater than those by slant crack and the 1X, 2X, and 3X amplitudes generated by transverse crack are close to those by slant crack as the crack is near the bearing.

As the 1X and 2X amplitudes of the vertical vibration response for the cracked rotor-bearing system are more sensitive to crack parameters, then the 1X and 2X amplitudes are selected as the index vectors for crack identification, and the superharmonic features of the four measuring points are utilized to avoid the multisolution problem in the process of crack parameters identification.

3. Construction of the Kriging Surrogate Model

The problem of crack parameters identification can be transformed to an optimization problem. Using optimization methods to estimate the optimal parameters via minimizing the objective function relating to output deviations such as dynamical response and modal parameters [38, 39]. However, the optimization process requires a large number of iterative steps, and each iterative step needs to repeat calculation of complex finite element models, making the problem very complicated. Thus, it is necessary to establish an uncomplicated relationship between the crack parameters and the superharmonic features for crack identification.

The Kriging surrogate model provides explicit functions to represent the relationships between inputs and outputs of a linear or nonlinear system [39], which is a statistics-based interpolation method and is not affected by random error [31, 32]. Construction of the Kriging model requires small amount of samples to obtain high estimation accuracy, significantly reducing calculation time.

Given initial crack parameters samples and the corresponding 1X and 2X amplitudes arewhere is the component of sampling crack parameters ( is the crack position, is the crack depth, and is the crack angle), represented as a three-dimensional variable vector, , , and .where is the set of the output vector, is the superharmonic components amplitudes at the measurement point, and and .

And their relationship can be written aswhere is a vector of a linear combination of selected functions, is the matrix of regression coefficients, is the stochastic process with a variance of and a mean of zero, and is the dimension of the vector. The covariance matrix is expressed by

is the Gauss correlation function matrix and is the correlation parameter, given by

When existing an untested crack parameter , the corresponding superharmonic component amplitudes can be estimated as the linear predictor, which is shown aswhere is a coefficient vector and is the column of the matrix .

The correlation vector can be written aswhere is the initial crack samples and is the new sample.

The predicted error isin which

To keep the predictor unbiased, it is required that

Under this condition, the mean squared error of Equation (13) is

The Lagrangian function for the problem of minimizing the mean squared error with the constraint of Equation (17) is

The gradient of Equation (19) with respect to is

According to the first order necessary conditions for optimality

Substituting the above equation into Equation (10) gives

Hence, the relation between crack parameters and the corresponding 1X and 2X amplitudes values has been deduced.

After the Kriging surrogate model is established, the squared multiple correlations (SMC) and the empirical integrated squared error (EISE) criterion [40] can be used to evaluate the preciseness of the Kriging surrogate model.where is the set of the vector of the Kriging surrogate model, is the component of the superharmonic component amplitudes calculated by finite element analysis, is the mean of the true values, and is the length of the vector . This paper adopts the multipoint adding strategy to update the Kriging surrogate model if the predicted superharmonic components amplitudes cannot meet and [41].

4. Crack Parameters Identification Based on an Improved NSGA-III

After the Kriging surrogate model between crack parameters and the superharmonic components amplitudes is established, crack parameters identification problem is transformed into a multiobjective optimization problem, which can be stated as follows:withwhere and are the superharmonic component amplitudes at the measurement point. is calculated via the Kriging surrogate model and is measured on actual cracked rotor, and and .

It is a great challenge for the classical multiobjective optimization algorithms when the optimization problem involves four or more objectives. To get the global optimal solution, the NSGA-III recently proposed by Deb and Jain [29, 30] may be a good choice to optimize 8 objective functions at the same time. NSGA-III is an evolutionary many-objective optimization algorithm based on reference point, which can ensure the diversity of population members via adaptively updating a number of well-spread reference points.

The crack angles considered are 45° or 90°, and the crack is located in the element in this paper. If the NSGA-III is directly used for crack parameters identification in the iterative process, the crack parameters of the population members are mostly decimal numbers, causing crack parameters identification results to be inaccurate. In addition, the decimal crack parameters reduce the convergence rate of the objective function and increase the optimization time.

Therefore, the NSGA-III needs to be improved for crack identification of the rotor-bearing system. Preprocessing of population members is added before members’ selection in every iterative process to obtain crack parameters accurately and efficiently, and the improved part of the NSGA-III is as follows.(1)Preprocessing of crack position:where is the crack position obtained by preprocessing of at the generation and is the rounding function.(2)Preprocessing of crack angle:where is the crack angle obtained after preprocessing of at the generation. The preprocessing of the crack angle can make the crack angle close to 45° or 90° change to 45° or 90°.

Replace and with and to form a new population.

There are many Pareto optimal solutions by the improved NSGA-III. As it is very difficult to have a set of crack parameters that minimize all the objective function values, an index is required to select the optimal crack parameters in the Pareto optimal solutions.

By defining as the selection index of the crack parameters identification results, can be obtained by weighing each objective function value:where is the weight coefficient, and and in this paper.

The corresponding identification result of the minimum value of is taken as the final crack parameters.

5. Numerical Investigation

5.1. Identification Results Based on a Kriging Surrogate Model

Numerical experiments are carried out to verify the effectiveness of the method. The detailed parameters of the rotor-bearing system are given in Table 1, 300 sets of samples of crack parameters are generated by Latin hypercube sampling [42], and then corresponding vertical vibration response of the rotor-bearing system assuming a 5% random white noise is collected at four measuring points. Then, the Kriging surrogate model between crack parameters and corresponding 1X and 2X amplitudes is established.

The improved NSGA-III is employed to estimate the crack parameters found on the Kriging surrogate model by minimizing the multiobjective function values given in Equation (25), and the initial parameters of the improved NSGA-III are as follows: the population size, maximum number of iterations, crossover percentage, and mutation percentage are 150, 100, 0.5, and 0.5, respectively.

The flowchart of the method for identifying crack parameters based on the improved NSGA-III and Kriging surrogate model is shown in Figure 5.

This paper considers eight cases (Table 2) to evaluate the performance of the proposed method and selects case 1 to analyze the identification process of crack parameters for rotors in detail.

Figure 6 illustrates the mean value of objective function varies with evolutionary generations, and it shows that the mean value of objective function optimized by the improved NSGA-III has already converged in the 30th generation, while the mean value of the objective function with the original NSGA-III not only has higher values but also has a slow convergence speed.

The crack parameters results identified with the improved NSGA-III are compared with results with the original NSGA-III (Figure 7), and it can be seen that the crack parameters obtained by the improved NSGA-III are closer to the actual crack parameters in Pareto optimal sets.

Figures 6 and 7 indicate that the improved NSGA-III can obtain Pareto optimal solutions accurately and efficiently.

The values of the objective function corresponding to each group of crack parameters in Pareto optimal solution are shown in Figure 8, and it is necessary to select a best set of solutions based on the index in Equation (29).

As shown in Figure 9, the selection index achieves the minimum value in the 90th population member, and the corresponding crack parameter is . Compared with other crack parameters in Pareto optimal solutions, the 90th set of parameters is closer to the actual values (marked as blue circle in Figure 7), which indicates that it is reasonable to select the optimal crack parameters by the index .

The identified crack parameters are listed in Table 2, and it is worth noting that results are very accurate except for case 3, as the crack depth and crack position in case 3 are near the boundary (crack location is close to the 55th element and crack depth is close to 0.1 cm), so the 1X and 2X amplitudes predicted by the Kriging surrogate model are less accurate. The accuracy of crack parameter identification results is improved by increasing the number of samples in the boundary area.

5.2. Comparison with the Results by Artificial Neural Networks

Artificial neural network (ANN) models are widely used in engineering for regression and classification, which have been proved to be effective in identifying the crack location and depth using modal parameters of static rotors [9, 10]. In order to explain the accuracy and efficiency of the proposed method in identifying the crack parameters, ANN is used to construct the relationship between superharmonic components amplitudes and crack parameters of a rotating rotor-bearing system.

The configuration of ANN used in this paper is backpropagation algorithm (BPA), and the input and target output for training are expressed, respectively, as and , where , is the superharmonic components matrix in Equation (9) and , is the crack parameters matrix in Equation (10). A single hidden layer containing 20 neurons is designed for training process, and the training progress of ANN uses Levenberg–Marquardt backpropagation algorithm.

The 300 sets of samples same as that in Section 5.1 are used to construct the neural network, assuming the 5% random white noise. Table 3 shows crack parameters identified by ANN, and it is imprecise compared with identified results obtained by the proposed method (Table 2) in this paper.

In addition, the crack parameters identification results obtained by the training neural network with 600 sets of samples are shown in Table 4, and it can be found that increasing the samples number trained by neural network can improve the identification accuracy of crack parameters.

However, increasing the number of samples requires more time for finite element calculation, and identified results are still imprecise compared with that obtained by the proposed method (Table 2) in this paper.

Based on the comparison result, it is clear to draw that the proposed method based on Kriging surrogate model and an improved NSGA-III can identify the crack parameters of operating rotors more accurately with fewer samples than ANN.

6. Conclusions

A new crack parameters identification method based on a Kriging surrogate model and an improved NSGA-III is presented in the paper. The numerical results clearly show that the proposed approach can identify the crack location, depth, and angle for operating rotors accurately and efficiently.(1)The 1X and 2X superharmonic components amplitudes of a rotor-bearing system in 1/3 subcritical speed region can be selected as the index vectors for crack parameters identification.(2)Crack parameters and their corresponding 1X and 2X amplitudes of vibration response for rotors are used to establish the Kriging surrogate model, which avoids complicated finite element calculation, demanding a large number of iterative steps.(3)NSGA-III is improved by adding the preprocessing of population members to avoid the disadvantages of the original NSGA-III for crack parameters identification, which can improve the accuracy of crack parameters identification and speed up the convergence of the multiobjective function values.(4)Compared with ANN, the proposed method can identify crack parameters for operating rotors more accurately with fewer samples.

It should be noted that only one crack is considered in this paper. There may be multiple cracks in rotors, and identification of multiple crack parameters is more complicated and requires further study.

Data Availability

All data generated or analyzed during this study are included in this published article.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This study was supported by the National Natural Science Foundation of China (no. 51875482).