A hybrid prognostics approach to estimate the residual useful life of a planetary gearbox with a local defect

Zhe Cheng1

1Science and Technology on Integrated Logistics Support Laboratory, National University of Defense Technology, Changsha, P. R. China

1Corresponding author

Journal of Vibroengineering, Vol. 17, Issue 2, 2015, p. 682-694.
Received 9 June 2014; received in revised form 16 July 2014; accepted 7 September 2014; published 31 March 2015

Copyright © 2015 JVE International Ltd.

Table of Contents Abstract Full-text  Download PDF References
Views 4
Reads 1
Downloads 11
since May, 2018

A hybrid prognostics approach for the monioring of a planetary gearbox with the local defect is presented. This hybrid method can predict the remaining useful life (RUL) of planetary gearbox with a fatigue crack. The method consists of a dynamical model for simulation data generation, a statistical algorithm for feature selection and weighting, and a modified grey model for RUL prediction. Experimental studies are conducted to validate and demonstrate the feasibility of the proposed method for RUL prediction of a cracked sun gear in planetary gearbox. And the validation has a promising result.

Keywords: planetary gearbox, defect, dynamical model, feature selection and weighing, residual useful life, prognosis.

1. Introduction

Engine and drive train failures are the leading causes of mishaps in rotorcrafts. As a critical component in main transmission gearbox of rotorcraft, the defect in planetary gear set is widely researched. One common defect is gear tooth crack caused by cyclic stress fatigue, which may lead to gear tooth breakage or complete failure of main transmission gearbox. However it is particularly challenging to detect and predict the remaining useful life (RUL) of these components, as the gear tooth fault of planetary gear set is located deep inside the main transmission gearbox.

As failure prognosis for planetary gearbox/rotorcraft transmission system has been studied widely, there are a lot of papers published in several key journals and at conferences. The common research papers in this area are data-driven approaches. Heimes [1] utilized advanced recurrent neural network architecture to estimate the RUL of the system. Peel [2] used a data driven approach of multiple neural network model to estimate the RUL left of an unspecified complex system. Orchard [3, 4] introduced a methodology for accurate and precise prediction of a failing component based on particle filtering and learning strategies. Swanson [5] used Kalman filters to track changes in relevant features like vibration levels, mode frequencies, or other waveform signature features, and then used the predicted states to estimate the future failure hazard, probability of survival, and RUL in an automated and objective methodology. Tse [6] used recurrent neural networks to forecast the rate of machine deterioration.

Physics-based methods have also been used to predict the failure time in mechanical systems. Ray [7] presents a stochastic model of fatigue crack damage in metallic materials to predict failure, and this method was validated by K-L analysis of test data. Li [8] presented a hybrid model-based method to predict the RUL of a gear with a fatigue crack. Experimental studies were conducted to validate and demonstrate the feasibility of the proposed method for prognosis of a cracked gear. Kacprzynski [9, 10] utilized stochastic, physics-of failure models controlled failures of single spur gear teeth and then on a helical gear contained within a drive train system. The method was validated with transitional run-to-failure data developed at Penn State ARL. Gu [11] presented a physics-of-failure based prognostics and health management approach for effective reliability prediction. The implementation procedure of this approach included failure modes, mechanisms, and effects analysis, data reduction and feature extraction from the life-cycle loads, damage accumulation, and assessment of uncertainty.

Generally, the methods aforementioned have their own shortages: 1) data-driven methods cannot provide sufficient guarantee about their prediction accuracy and definite physical meaning on the relationship between the index and damage; 2) physics-based approaches need high computational costs associated with the physical models and these methods are only suitable for off-line applications; 3) the data sets used for the prognosis of planetary gearbox include strong environmental noise components and many frequency components of other moving parts.

In view of the above analysis, a hybrid prognostics approach is presented for planetary gearbox. This method is composed of a dynamical model for vibration response analysis, a two-step statistical algorithm for feature selection and a modified grey model for RUL prediction. In comparison with the existing methods reported above, the proposed method harnesses the merits that the computational costs of analytical model are lower and the relationship between indicator and fault is easy to understand. This method has better prediction accuracy compared to the data driven method. Damage seeded tests and run-to-failure test are conducted to test the performance of the proposed method.

2. The dynamical model of planetary gearbox with defect

The model of the gear tooth failure can be used to analyze the dynamic features in order to provide new signal source and suitable tools for predicting these failures. In the models, gear tooth defect is commonly assessed by the variation of dynamic parameters due to defect, such as stiffness and damping. The finite element analysis (FEA) and the analytical model are the general approach to simulate the defect seeded in part [12-20]. Considering that FEA requires much computations time, this research will focus on the tooth stiffness reduction due to defect based on analytical model of planetary gear set.

2.1. The dynamical model of the healthy planetary gearbox

Due to the multiple components and its orbital motion of the planets, the epicyclic stage of transmission system is more complex than the ordinary gear train. As the internal ring gear is fixed to the top of the transmission casing and the planet gear support is very rigid, the central displacement of the ring gear and the planet gears in the planetary gearbox can be ignored to simplify the model. Then a lumped parameter dynamical formulation is employed to develop the dynamical model of a common planetary gearbox. Deduced from Lagrange equations, the dynamical differential equations of this planetary gearbox could be listed as:

I s θ ¨ s + i = 1 N D s p i + P s p i r b s = T D ,   I p i θ ¨ p i - ( D s p i - D r p i + P s p i - P r p i ) r b p i = 0 ,   I c + i = 1 N m p i r c 2 θ ¨ c - i = 1 N D s p i + D r p i + P s p i + P r p i r c c o s α = - T L ,

where θs, θpi and θc denote the rotation angle of sun-gear, planet gear and carrier respectively. TD and TL denote the driving torque and the loading torque respectively. r denotes radius and rb is the radius of basic circle. m is mass, I is the rotational inertia of subcomponents, for planet gear, sun gear and carrier, I=1/2mrb2. N is the number of planet gear and α is the mesh angle. The subscripts r, s, pi and c denote ring gear, sun gear, the ith planet gear and carrier. D and P is the adhesive engaging force and the elastic engaging force respectively.

Eq. (1) is positive semi-definite, nonlinear equations, which have N+ 2 degrees-of-freedom (dof’s). To translate the angle displacements of rigid bodies into the relative linear displacement, the relative displacement between the ith planet gear xspi and the sun gear, and the relative displacement between the carrier and the sun gear are defined as:

x s p i = θ s r b s - θ p i r b p i - θ c r c c o s α , x s c = θ s r b s - 2 θ c r c c o s α .  

To simplify the model, the manufacturing errors and the mesh errors in the model are ignored as can be referred in [15]. After that, the internal meshing side clearance and the external meshing side clearance are defined as 2bspi and 2brpi, the adhesive engaging force D and the elastic engaging force P can be presented as:

D s p i = C s p i x ˙ s p i ,                                                           D r p i = C r p i ( x ˙ s c - x ˙ s p i ) , P s p i = K s p i ( t ) f ( x s p i , b s p i ) ,                                           P r p i = K r p i t f x s c - x s p i , b r p i ,

where Kspi denotes the mesh stiffness between the planet gear and the sun gear, Krpi denotes the mesh stiffness between the internal ring gear and the planet gear; Cspi denotes the mesh damping between the planet gear and the sun gear; Crpi denotes the mesh damping between the ring gear and the planet gear. K(t) is time varying mesh stiffness, and f(x,b) is nonlinear clearance function, which is defined as:

f x , b = x - b , x > b , 0 ,   - b x b x + b , x < - b , ,

where b is the clearance constant.

After that, Eq. (2) and Eq. (3) are substituted into Eq. (1) and simplifying the equations, the dynamic model of 2K-H planetary gear set can be obtained. The equivalent mass is denoted as M, for the planet gear and the sun gear, M=I/rb2, for the carrier, Mc=Ic/rbc2+i=1N(mpi/cos2α).

2.2. The dynamical model of planetary gearbox with fatigue crack seeded

Fatigue crack is one of the common damages in planetary gearbox, and the cause of this damage is cyclic stressing of the gear tooth beyond its fatigue limit. Generally, the bending fatigue crack starts in the root section and it progresses until the tooth or part of it breaks off. The common fatigue crack occurs in the usual tooth root fillet section, so the tooth fatigue crack of sun gear is seeded in the model of planetary gear set in this section, as can be seen in Fig. 1.

Fig. 1. Gear-mesh stiffness variation of the gear pair with tooth crack



To simplify the model of damage, the shape of tooth root crack is approximated with straight lines, which are defined by the length l and the direction angle ς of the crack, as shown in Fig. 1. The level of crack is defined as s, which is determined by l and ς, s=s(lj,ς). For the mathematical tractability to solve the dynamical model, the crack size and its threshold value are defined as cj and cf,cj=2m'cotς-ljsinς,c~j=cj/cf,j= 1, 2, 3,…. And then sj=s(lj,ς)=s(c~j), sj [0 %, 100 %]. Referring [16, 17], the gear mesh stiffness of the gear pair is calculated by taking into account the geometric change due to the tooth crack. The details can be referred to [16].

The variation amplitude of gear mesh stiffness, which is caused by sun gear tooth crack, is defined as ΔK(fd,t), thus the time varying mesh stiffness of gear pair with tooth crack can be expressed as:

K s p i t = K 0 f 0 , t + Δ K f d , t ,

where K0f0,t is the mesh stiffness of the healthy gear pair, f0 is the mesh frequency, and fd is the meshing frequency of the gear pair with defect. fd is the function of the rotary frequency of sun gear fs, the rotary frequency of carrier fc and the number of planet gears N, fd=N(fs-fc). The time vary mesh stiffness is illustrated in Fig. 1(b). The dynamical model of planetary gearbox with sun gear tooth crack is acquired by substituting Eq. (5) into the model deduced above.

2.3. Simulation of the dynamical models

The parameters in the dynamical model above should be set up according to Table 1. Parameter s is varied as [0 %:5 %:100 %] for the defect seeded model. Matlab 7.1 is used to execute this simulation, and four-order Runge-Kutta method is selected as the solution for the models. The step size and the duration for the solution procedure are set up to 0.0001 second and 1 second respectively.

Table 1. Parameter value in the models

Parameter (Unit)
Parameter (Unit)
Modulus (mm)
r c (m)
Tooth number of sun gear
Driving torque (N⋅m)
Tooth number of planet gear
Loading torque (N⋅m)
Tooth number of ring gear
K s p i (N/m)
Number of planet gear
K r p i (N/m)
Tooth width (mm)
C s p i (N⋅s/m)
Pressure angle (Deg)
C r p i (N⋅s/m)
m s (kg)
b s p i , brpi (m)
m p i (kg)
I p i (m)
m r (kg)
I s (m)
r b s (m)
I c (m)
r b p i (m)

The dynamic responses of the dynamical model in both the time domain and the frequency domain are shown in Fig. 2, which is for a healthy planetary gearbox (s= 0 %) and 3 different levels of crack (s= 20 %, 40 %, 80 %). The healthy model is characterized by the dominance of the gear mesh frequency, denoted by 1X and its harmonics of 2X-6X, and etc. For the defect seeded models, the amplitude modulation in time domain of the gear mesh signal can be clearly observed, as can be seen in Fig. 2(a). As a consequence, many new frequency components appear at the left side of the dominant frequency (near 1/2X). The amplitude of the new frequency components increases with the growth of crack level, but the amplitude of dominant frequency and its harmonics decreases at the same time, as can be seen in Fig. 2(b). These features are very close to the experimental observations in the previous study literatures [16, 17, 21].

Fig. 2. Dynamical response of the models

a) Time domain

b) Frequency domain

2.4. Validation of the dynamical models

The test data sets in health condition and damage seeded condition are used to validate the dynamical model deduced in this paper. Above all, Root Mean Square (RMS) of each data set is calculated to reduce the effect of environment noise as the aid of the common algorithm TSA, and then the curves of error percentage of RMS between test data and simulation data of the dynamical models are depicted in Fig. 4. The real line is the error percentage between test data set and simulation data set in health condition, which is about 5 %, and the imaginal line is for damage seeded condition, which is about 15 %, as can be seen in Fig. 3. As the reason of that the damage seeded in the dynamical model above is simplified, which is different from the damage seeded in test. As a result, the imaginal line is above the real line obviously. Generally, the error percentage is low and the accuracy of dynamical model would be improved by considering more nonlinear factors of the damage seeded in the model.

Fig. 3. Validation results of dynamical models

3. Feature selection and weighting

From our literature review, 27 features have been used in different cases of gearbox condition monitoring. The features are organised into four groups and assigned with serial numbers.

1) Features in time domain are used most frequently in gearbox diagnosis [21]. They include root mean squared (RMS), crest factor (CF), energy ratio (ER), kurtosis, standard deviation, energy operator, absolute mean value, clearance factor and impulse factor, which are assigned with the serial numbers 1 to 9 respectively for the ease of identification in the process of feature selection.

2) There are many other traditional statistical feature parameters, such as FM0, FM4, NA4, M6A, M8A, NB4, NA4*, NB4*, M6A* and M8A*, generally used in planetary gearbox condition monitoring [22]. The serial numbers of these features are 10 to 19.

3) Other kinds of feature parameters based on frequency spectrum of vibration signal are widely used to detect and diagnose faults in rotorcraft power trains, such as mean frequency (MF), frequency centre (FC), root mean square frequency (RMSF) and standard deviation frequency (STDF) [23]. The serial numbers of these features are 20 to 23.

4) Other than the features presented above, there are some important features which have been validated in literatures, which are named Intra-Revolution Energy Variance (IREV) [24], spectrum kurtosis (SK) [25, 26], local spectrum kurtosis (LSK) [27], NSR [28]. The serial numbers of these features are 24 to 27.

The simulation data sets generated by the dynamical model above are utilized to calculate all the feature parameters listed above.

3.1. Feature selection

An optimal feature suitable for prognostics should have three merits: 1) Sensitivity, it means that the feature has a optimal classification distance; 2) Stability, which means the feature has the stable classification performance in different operation conditions (including loading and rotational speed). The target features will be selected from the 27 feature parameters above.

A statistic algrithm named two-sample Z-test is commonly used to measure the distance of two class samples [24, 27]. In this research, this algrithm is applied as sensitivity analysis algorithm, and then it is modified to analyze the stability and the relationality of feature parameters. The feature selection procedure can be described as follows:

1) Sensitivity analysis.

For the kth feature parameter, calculating the classification distance Zk of the healthy samples and the fault seeded samples under the identical operating condition by two samples Z-test procedure:

Z k = X - k , 1 - X - k , 2 S X k , 1 2 n 1 + S X k , 2 2 n 2 ,

where X1 and X2 are the healthy sample set and the fault seeded sample set, SX and X- are the standard deviation and the mean of X, respectively. And n is the number of samples for each sample set. k=1, 2,...,K, K is the number of feature parameters above. As the value of Zk increases, the sensitivity of the kth feature parameter will grow, that means the value of the kth feature parameter will have an obvious change while there is damage.

2) Stability analysis.

The classification distance of the kth feature parameter under the jth condition is defined as Zk,j, j=1, 2,..., J. Then, the similarity ratio of classification distance for different condition, which is named stability in this research, is calculated:

S k = S Z k / n z k m = 1 , n = 1 J Z k , m - Z k , n ,

where Zk is the classification distance vector of the kth feature parameter for different condition, Zk,m is the mth element of Zk, SZk is the standard deviation of Zk, J is the number of operation conditions, n is the element number of Zk. The numerator in Eq. (7) shows the general statistical characteristic of Zk, and the denominator in Eq. (7) denotes the statistical diversity of classification distance of the kth feature parameter for different operation conditions. So Sk indicates the stability of the kth feature parameter to different operation condition. As the value of Sk increases, the stability of the kth feature parameter will rise, which means the kth feature parameter will keep stable while the operation condition (e.g. rotary speed and loading) changes.

According to the feature selection algorithms in this section, the feature parameter suitable for prognosis can be selected from the potential feature set when their evaluation criterion is predefined.

3.2. Feature weighting

As the feature vector includes more than one feature parameters, the feature weighting is necessary for fusing the analysis results of all the parameters in feature vector to achieve more reliable results. In this research, the weighting of feature parameter is determined based on the sensitivity and stability performance of the feature parameter calculated above. And also note that the weighting sum of all the feature parameters is equal to 1. Thus the parameters in feature vector can be weighted as follows:

W 1 Z 1 + S 1 = W 2 Z 2 + S 2 = = R i Z i + S i = = R n Z n + S n , i = 1 N W i = 1 ,

where i is the serial number of the feature parameter, N is the number of feature parameter in the target vector.

4. RUL estimation based on modified grey model

4.1. Prognosis approach based on grey model

The simplest and most elementary form of grey model is GM(1,1) [31]. It is proposed as a forecasting model to extrapolate the original data sequence. The concept of the grey theory was first derived by Deng in 1986 [31]. It describes GM(1,1) model as approximations via least-square estimations at the first-order differential equation level. GM(1,1) has been widely used in the prediction of mechanical systems [31-38].

Vibration signal is an important information resource for monitoring rotary mechanicals. With the propagation of crack, the value of vibration indicators will increase generally. This process has two characteristics: 1) Xti 0 and 2) X(ti)X(ti-1), which show the data set has strong tendency. These two characteristics meet the need to build grey model (GM), so GM(1,1) can be used to predict the future condition of planetary gearbox. The prediction procedure is as bellows:

Define x(0)=x01,x02,,x0n as the feature condition series, x(1) is the accumulated generating operation (AGO) of x(0), x1k=i=1kx0i, k1,,n; z(1)is the mean generating operation of x(1),z1k=0.5x1k+0.5x1k-1,k2,,n. GM(1,1) is:

d x 1 ( t ) d t + a x   1 t = u .

The discrete model for prediction is:

x ( 1 ) k + 1 = x   1 0 - u / a e - a k + u a ,


a , u T = a ^ = B T B - 1 B T y N ,       B = - z 1 2 1 - z 1 n 1 ,       y N = x 0 2 , , x 0 n .

4.2. RUL estimation based on the modified grey model

To improve the prediction accuracy of the grey model in Eq. (10), a residual efficient ε is added to the grey input u of the grey model for the addition of novel information from the physics-based models:

d x ( 1 ) ( t ) d t + a x   1 t = u + ε .

Then the discrete form is:

x ( 1 ) k + 1 = x ( 1 ) 1 - u + ε a e - a k + u + ε a .

The residual efficient of feature Fi is defined as εi. After the feature parameters used for RUL prediction are determined, ε could be calculated based on the given information of simulation signal from dynamical model developed above:

ε = a x ^   1 k + 1 - a x ^ ( 1 ) 1 - u e - a k - u 1 - e - a k ,

where x^(1)1 and x^(1)k+1 are AGO of the time series of feature parameter F, and the elements of the time series are calculated based on simulation signal.

The one-step and multi-step prediction result of Fi is xi(1)k+1, which could be obtained based on test data or real data according to Eq. (13). After that the value of prediction results xi(1)k+1 are fused utilizing the feature weight vector Wi according to Eq. (15):

x ( 1 ) k + 1 = i = 1 N x i ( 1 ) k + 1 W i ,

where xi(1)k+1 is the prediction result of Fi, and Wi is the weight of Fi. Then x (1)k+1 in Eq. (15) is the expected value of future condition prediction.

In this research, RUL Nf is expressed as the stress cycle of gear crack tooth n. Nf is corresponding to x (1)k+1 according to Eq. (16):

N f = ξ x   ( 1 ) k + 1 ,

where ξ is determined based on the historical data of damage evolution test.

4.3. Description of the RUL estimation result

Logarithmic normal distribution (LND) is often used to describe the prognosis result of mechanical parts. In this research, the RUL probability distribution function of tooth fatigue crack in planetary gearbox is considered to belong to LND. The common form of LND is two-parameter-distribution, the density function is:

f t = 1 t σ 2 π e - 1 2 l n t - μ σ 2 ,

where t is time, σ and μ are the standard deviation and the mean, respectively.

As there are a lot of uncertainty issues in equipment operation, the accuracy of prediction result will decrease quickly while the prediction time extends, and then the variance of RUL distribution will grow. By iteration of modified grey model, multi-step prediction will be carried. Variance σ could be expressed as the function of prediction life step n, σ=σ(n)=σen-1. The RUL distribution is obtained by substituting the expressions of variance and expectation into Eq. (17):

f n = 1 n σ n 2 π e - 1 2 l n n - N f σ n 2 .

Eq. (15) and Eq. (18) are combined, and then the prognosis result of tooth root fatigue crack in planetary gear sets is generated.

The relationship of RUL N and failure time T is correlated with operation profile in test/field, such as the rotary speed variation.

To evaluate the prognosis results, the prediction error η is defined as the difference between the actual failure time Taf and the predicted failure time Tpf:

η = 100 × T a f - T p f T a f - T c % ,

where Taf is the actual failure time, Tpf is the predicted failure time and Tc is the current time.

5. Validation with test data

5.1. Test rig and the tooth crack propagation test

As can be seen in Fig. 4(a), the test rig provides torque and speed reduction through three stages, which simulates the main transmission system of general rotorcraft. The second stage is a planetary gearbox, of which the internal gear rig is splined to the top of gearbox casing.

To validate the prognostics approach in this research, a number of experiments have been conducted. The vibration signal generated by the planetary gearbox is picked up by an accelerometer bolted on the top of planetary gearbox casing and the electrical signal is transferred to the data acquisition system. The sampling frequency fs is set up to 10 kHz. The signal is low-pass filtered at 2 kHz through a 4th order filter to limit aliasing distortion and retain waveform integrity of the effective signal. Data is stored in a PC for post processing. A number of 204,800 data points has been acquired in tests corresponding to a time-history length of 20 s.

In this research, an incipient crack is seeded in the sun gear, as can be seen in Fig. 4(b), and a run-to-failure test is carried out. In the test, the interval of sampling is 10 min. For every 50 minutes, the test is paused and the actual size of fatigue crack in propagation is checked. When the fatigue crack level reaches the threshold value, the test is ended. The threshold value of crack level is set up according to the experience of the prior test. The total run time of this crack seeded test is 736 min. From the test, the run-to-failure transitional data of fatigue crack propagation in sun gear of planetary gearbox is acquired.

Fig. 4. Configuration of the test rig and the sun-gear with tooth root crack

a) The test rig

b) The sun-gear with tooth root crack

5.2. Validation results with test data

Based on the simulation data of the dynamical models, feature selection and weighting has been carried out before the test data validation. Only a few valueble features are kept and used, as can be seen in Fig. 5. As there is no prior knowledge of setting the threshold, a median value is used as selection criteria. In this research, [20, 20] is selected to be the threshold vector to illustrate the proposed methodology.

Fig. 5. Feature selection results



In the feature selection procedure, the feature parameters are calculated based on the simulation data sets of the dynamical models with diferent level damge seeded. The procedure of the feature selection could be clearly dipicted in Fig. 5. The first step is the sensetivity analysis of feature parameters, the threshold value for this step is set up to 20, and the feature parameters #10, #11, #12, #13, #14, #16, #17, #18, #19, #22, #24, #25, #26, #27 are in the second step for stability analysis, as shown in Fig. 5(a). In the second step, the threshold value is set up to 20, and the selection result of this step is the feature parameters #25, #26 and #27, as shown in Fig. 5(a).

The feature parameters #25, #26 and #27 have been selected to be the remaining features, as shown in Fig. 5(b). These features are SK, LSK and NSR, and they are labeled F1, F2 and F3 in this paper for convenience. Following the procedure in Eq. (9), the weights of these selected features are 0.2752, 0.2013 and 0.5235 respectively.

To validate the proposed RUL prediction approach in this research, the test data sets X1~X4 at time 360 min and 460 min are used for analysis. 50 feature samples of each test data set are extracted. The RUL for the test data at t= 360 min and t= 460 min is predicted. The prediction result of test data is recorded in Table 2, where F0 denotes the healthy case, Ft denotes the failure case, FA is the prediction result at t= 360 min, FB is the prediction result at t= 460 min. The RUL at t= 360 min is 681 min with the confidence of 95 %. And the RUL at t= 460 min is 705 min with the confidence of 95 %, while the real RUL of this test is 736 min.

According to Eq. (17), the prediction accuracy at t= 360 min is 85.37 % and at t= 460 min is 88.77 %. It could be seen from these prediction results that, as the prediction time reaches failure time, the accuracy of prediction result will be improved, in other words, the expected RUL of prediction result will approximate the real RUL, while the variance of prediction result decreases quickly. Compared to the traditional GM(1,1), the modified GM(1,1) has better accuracy and lower prediction error, which can be seen in Table 2.

Table 2. Comparation of the prediction results for the modified GM and the traditional model

Model type
Prediction start time (min)
Failure time for 95 % confidence (min)
Real failure time (min)
Prediction error (min)
Accuracy (%)
Traditional GM(1,1)
Modified GM(1,1)

6. Conclusions

A hybrid RUL esitimation methodology for planetary gearbox is devised by integrating dynamical models for damage simulation signal generation, a two-step statistical algorithm for feature selection and weighting, and a modified grey model for prediction, to forecast the RUL of a cracked sun-gear in planetary gearbox. The proposed method is validated with the sun gear run-to-failure test. The experimental study results indicate that this prognostics approach has high accuracy to predict RUL of sun gear tooth fatigue crack in planetary gearbox. However, the RUL prediction of this approach is a little conservative. It is suggested that the method should improve the model parameter of dynamical model and GM(1,1) for certain applications.


  1. Heimes F. Recurrent neural networks for remaining useful life estimation. International Conference on Prognostics and Health Management, CA, USA, 2008.
  2. Peel L. Data driven prognostics using a kalman filter ensemble of neural network models. International Conference on Prognostics and Health Management, CA, USA, 2008.
  3. Orchard M. A Particle filtering framework for failure prognosis. Proceedings of World Tribology Congress, Washington, D.C., USA, 2005.
  4. Orchard M. Advances in uncertainty representation and management for particle filtering applied to prognostics. International Conference on Prognostics and Health Management, CA, USA, 2008.
  5. Swanson D. A general prognostic tracking algorithm for predictive maintenance. IEEE Aerospace Conference, Big Sky, 2001, p. 2971-2977.
  6. Tse P., Atherton D. Prediction of machine deterioration using vibration based fault trends and recurrent neural networks. Journal of Vibration and Acoustics, Vol. 121, Issue 3, 1999, p. 355-362.
  7. Ray A. Stochastic modeling of fatigue crack damage for risk analysis and remaining life prediction. Journal of Dynamic System, Measurement and Control, Vol. 121, Issue 3, 1999, p. 386-393.
  8. Li C., Lee H. Gear fatigue crack prognosis using embedded model, gear dynamic model and fracture mechanics. Mechanical Systems and Signal Processing, Vol. 19, Issue 4, 2005, p. 836-846.
  9. Kacprzynski G., Sarlashkar A., Roemer M., etc. Predicting remaining life by fusing the physics of failure modeling with diagnostics. Journal of Materials, Vol. 56, Issue 3, 2004, p. 29-35.
  10. Kacprzynski G., Maynard K., Roemer M. Enhancement of physics-of-failure prognostic models with system level features. Aerospace Conference Proceedings, 2002.
  11. Gu J., Pecht M. Prognostics and health management using physics-of-failure. Proceedings of the Annual Reliability and Maintainability Symposium, 2008.
  12. Pimsarn M., Kazerounian K. Efficient evaluation of spur gear tooth mesh load using pseudo-interference stiffness estimation method. Mechanism and Machine Theory, Vol. 37, Issue 8, 2002, p. 769-786.
  13. Wang J. Numerical and experimental analysis of spur gears in mesh. Ph.D., Curtin University of Technology, 2003.
  14. Ambarisha V., Parker R. Nonlinear dynamics of planetary gears using analytical and finite element models. Journal of Sound and Vibration, Vol. 302, Issue 3, 2007, p. 577-595.
  15. Chaari F., Fakhfakh T., Haddar M. Dynamic analysis of a planetary gear failure caused by tooth pitting and cracking. Journal of Failure Analysis and Prevention, Vol. 6, Issue 2, 2006, p. 39-44.
  16. Chaari F., Baccar W. Effect of spalling or tooth breakage on gearmesh stiffness and dynamic response of a one-stage spur gear transmission. European Journal of Mechanics A/Solids, Vol. 27, Issue 4, 2008, p. 691-705.
  17. Hbaieb R., Chaari F., Fakhfakh T., Haddar M. Influence of eccentricity, profile error and tooth pitting on helical planetary gear vibration. Machine Dynamics Problems, Vol. 29, Issue 1, 2005, p. 5-32.
  18. Walha L., Louati J., Fakhfakh T., Haddar M. Dynamic response of two stages gear system damaged by teeth defects. Machine Dynamics Problems, Vol. 29, Issue 1, 2005, p. 107-124.
  19. Choy F., Polychuk V., Zakarajsek J., Handschuh R., Townsend D. Analysis of the effects of surface pitting and wear on the vibrations of a gear transmission system. Tribology, Vol. 29, Issue 1, 1996, p. 77-83.
  20. Yuksel C., Kahraman A. Dynamic tooth loads of planetary gear sets having tooth profile wear. Mechanism and Machine Theory, Vol. 39, Issue 7, 2004, p. 695-715.
  21. Dempsey P., Lewicki D., Le D. Investigation of current methods to identify helicopter gear health. Aerospace Conference, Big Sky, Montana, 2007.
  22. Samuel P., Pines D. A review of vibration-based techniques for helicopter transmission diagnostics. Journal of Sound and Vibration, Vol. 282, Issue 1, 2005, p. 475-508.
  23. Yang B., Kim K. Application of dempster-shafer theory in fault diagnosis of induction motors using vibration and current signals. Mechanical Systems and Signal Processing, Vol. 20, Issue 2, 2006, p. 403-420.
  24. Wu B., Saxena A. An approach to fault diagnosis of helicopter planetary gears. IEEE Autotestcon, 2004, p. 475-481.
  25. Barszcz T., Randall R. Application of spectral kurtosis for detection of a tooth crack in the planetary gear of a wind turbine. Mechanical Systems and Signal Processing, Vol. 23, Issue 4, 2009, p. 1352-1365.
  26. Antoni J., Randall R. The spectral kurtosis: application to the vibratory surveillance and diagnostics of rotating machines. Mechanical Systems and Signal Processing, Vol. 20, Issue 2, 2006, p. 308-331.
  27. Cheng Z., Hu N., Gao J. An approach to detect damage quantitatively for planetary gear sets based on physical models. Advanced Science Letters, Vol. 4, Issue 4, 2011, p. 1695-1701.
  28. Cheng Z., Hu N. Application of GRA in damage detection and mode identification of planetary gear set. Journal of Grey System, Vol. 23, Issue 4, 2011, p. 315-326.
  29. Wen K. Grey systems: Modeling and prediction. Yang’s Scientific Research Institute, Tucson, 2004.
  30. Wen K., Changchien S., Yeh C., Wang C., Lin H. Apply MATLAB in Grey System Theory. CHWA Publisher, Taipei, 2006.
  31. Deng J. Introduction to grey system theory. Journal of Grey System, Vol. 1, Issue 1, 1989, p. 1-24.
  32. Deng J. Control problems of grey system. System Control Letters, Vol. 1, Issue 5, 1982, p. 285-294.
  33. Deng J. A Novel grey model GM(1,1|τ,r): Generalizing GM(1,1). Journal of Grey System, Vol. 13, Issue 1, 2001, p. 1-8.
  34. Lin S., Wu S. An Intelligent Web-based GRA/Cointegration Analysis for Systematic Risk. International Journal of Computers, Vol. 4, Issue 4, 2010, p. 223-234.
  35. Wang Y., Liao R., Chen W., Du L. A on-line forecast method for predicting the gas-in-oil concentrations in oill-filled transformer. International Conference on Condition Monitoring and Diagnosis, Changwon, Korea, 2006.
  36. Anjali C., Nirmal K. Applying grey theory prediction model on the dga data of the transformer oil and using it for fault diagnosis. Transactions on Power Systems, Vol. 4, Issue 1, 2009, p. 43-52.
  37. Cempel C. Forecasting the global and partial system condition by means of multidimensional condition monitoring methods. Journal of Theoretical and Applied Mechanics, Vol. 46, Issue 3, 2008, p. 777-797.
  38. Cempel C. Decomposition of the symptom observation matrix and grey forecasting in vibration condition monitoring of machines. International Journal of Applied Mathematics and Computation Science, Vol. 18, Issue 4, 2008, p. 569-579.
  39. Chen Z., Yang Y., Hu Z. Detecting and predicting early faults of complex rotating machinery based on cyclostationary time series model. Journal of Vibration and Acoustics, Vol. 128, Issue 5, 2006, p. 666-671.
  40. Cheng Z., Hu N., Zhang X. Crack level estimation approach for planetary gearbox based on simulation signal and GRA. Journal of Sound and Vibration, Vol. 331, Issue 26, 2012, p. 5853-5863.