Bi-periodically correlated random processes as a model for gear pair vibration

Abstract

The model of gear pair vibration in the form of bi-periodically correlated random processes (BPCRP) that describes its stochastic recurrence with two periods is proposed. Particular cases of this model are considered. It is shown that BPCRP model allows one to analyses unequally the mean and the covariance function of the additive and multiplicative components. There are considered technologies for the estimation of the Fourier coefficients of the mean and the covariance functions.

1. Introduction

2. Modeling of gear pair vibration

( 1 ) where W is a some load and   t  is an angular position of the gear. The terms x m   and x e   describe the contact properties of the gears, while terms x 1   and x 2   are caused by manufacturing error. It is supposed that each term x i   , i  1,2 is periodic with a rotation period Pi  1/f i

of the corresponding gear. There are three periodic terms in ( 1 ), namely x e   W x m   , x e  x 1   andx e  x 2   , which are periodic functions with period Pm  1/f m , P1 and P2 . The model of the cyclostationarity offered in [12, 18] was obtained by introducing the fluctuations of the angular position of the gears as a some random variable. The mean function of this random process includes the harmonic with frequencies f m , f 1 and f 2 . The covariance function consists of three different kinds of harmonics, in that, the harmonics with frequencies that are a linear combination of the rotation frequencies kf 1 lf 2 , the harmonics of the mesh frequency nf m and the harmonics with frequencies that are a linear combination of the mesh frequency and the rotation frequencies, i.e. nf m kf i . The first and the second order nonstationarities have been substantiated by the processing of vibration signals measured on the gear systems [12, 18], and the quantities that describe the structure of the cyclostationarity estimated by means of synchronous averaging were proposed to be used for fault detection.

In [20–22], after applying the synchronous averaging with the period P1 or P2 , the vibration signal is expressed as

M g t   Al 1 al t  cos2f lt bl t  l  , l 1 ( 2 ) ( 3 ) here M is the number of gear mesh harmonics, Al and l are the amplitude and the phase of the l th harmonic respectively. The modulation effects are described by the functions 1 al t  and bl t  , which are periodic with the considered rotation period. These functions are closely approximate to the signal’s deterministic component corresponding to one revolution of the selected gear.

In [19, 23] the gear vibration signal is modeled as

x    x 1   x 2   x 1,2   x c   n   , where x 1   and x 2   describe the deterministic periodic oscillations generated by the rotation of the output and input wheels respectively, x 1,2   is a component with period P12 r1P1 r2P2 , x c   is the second order cyclostationary process with period P12 and n   is a fluctuation component. The deterministic part of the signal ( 3 ) can be extracted by means of synchronous averaging with common period P12 of the shafts as far as it is possible [19]. Paper is dedicated to the development of the cyclostationary models of gear vibrations considered in the literature, their concretization, and the elaboration on this basis of other estimation techniques for the analysis of the modulation effects occurring in the vibration signals as the faults originate.

3. BPCRP representation

The efficiency of vibration signal processing for machinery condition monitoring can be explained mostly by their possibility to develop modulations caused by the appearance of faults. The modulation effects in the vibration model as a periodically correlated random processes (PCRP), which describe the stochastic recurrence with one period can be explained by the jointly stationary random processes k t  in their harmonic representation [8, 9, 24]: ik 2 t  t   k t e P1 ,

k Z  t    P2 e ik 2P t

k 1 , k Z where Z is a set of integer numbers and P1 is a period of the rotations for one of the wheels. Following this equation, we concludes that the modulation of the signals of two stochastic rhythms provided by the rotation of two wheels can be explained as where the harmonic of frequency 2 /P1 and its multiples are modulated for this once by PCRP with il 2 t period P2 : k(P2 ) t   kl t e P2 .

l Z

Then, for the random process (4), we have:

 t    kl t e i klt ,

k ,l Z where kl t  are jointly stationary random processes and kl k 2 /P1  l 2 /P2  . The random processes presented by series ( 5 ) are bi-periodically correlated random processes (BPCRP) [9, 25, 26]. As we can see from ( 5 ), a BPCRP is a sum of the amplitude and phase modulated harmonics. Here frequencies kl are the linear combination of the two main frequencies 10 k 2 /P1  and 01 l 2 /P2  . The modulating processes have the mathematical expectations mkl E kl t  which are the Fourier coefficients of the mean: ( 4 ) ( 5 ) ( 6 ) ( 7 ) ( 8 ) ( 9 )

It follows from (8) that

where are the cross-spectral densities of the modulating processes pq t  . The functions ( 8 ) and ( 9 ) are respectively the covariance and spectral components [9, 25, 26].

For the covariance function R t ,  E  t  t   ,  t   t  m t  , we have m t  E  t    mkle i klt .

k ,l Z R t ,  E  t    Rkl ( )e i klt ,

k ,l Z Rkl     rp k ,q l ,p ,qe i pq ,

p ,qZ f kl    1 

 Rkl  e id .

2  f kl     f p k ,q l ,p ,q   pq  ,

p ,qZ f pqkl    1 

 rpqkl  e id , 2  where and rpqkl   E pq t kl t   , pq t  pq t  mpq are the cross-covariance functions of the PCRP processes, and the “¯” signifies complex conjugation. Thus, the cross-covariance functions of the modulating processes defines the Fourier coefficients of the covariance function ( 7 ) in which the numbers are shifted by k and l . It follows from ( 8 ) those cross-correlations of modulating processes kl t  with different numbers k and l lead to bi-periodical non-stationarity of the second order. As the result of these correlations, it is appear the correlation of the spectral components, which can be characterized by the appropriate Fourier transformation of ( 8 ): rpq   E pq t pq t   : R 00     rpq  e i pq .

p ,qZ

This covariance function of the stationary approximation for the BPCRP is averaged BPCRP covariance function.

The zeroth spectral component f 00     f pq   pq  , ( 10 )

p ,qZ is a power spectral density of the stationary approximation for the BPCRP. It defines the spectral decomposition of the averaged in time instantaneous power R 0,t  for the oscillations.

We should note that the covariance and the spectral components are the total characteristics of the amplitude and the phase modulation of the BPCRP carrier harmonics. The zeroth spectral component, as can be seen from ( 10 ), is a sum of the power spectral densities of the modulating processes pq t  shifted by pq . The components f kl   explained in ( 9 ) are a sum of the shifted cross-spectral densities for modulating processes. Their numbers differs by k and l . Proceeding from the abovementioned, it is possible to conclude that the zeroth spectral function f 00   describes the spectral composition of the oscillations and the non-zeroth functions f kl   . It explains the correlations of the harmonics of this composition for the components with frequencies shifted by kl k 2 /P1  l 2 /P2  . When modulating processes of the corresponding numbers are mutually correlated, than these correlations are not equal to zero.

4. Method for statistical analysis

The time synchronous averaging (TSA) method was one of the early techniques used for the analysis of hidden periodicities [27, 28]. If the hidden periodicity is presented and modeled as a PCRP, then such technology was used for evaluation of its mean and covariance function [9, 25, 26]. It is so-called the coherent method [29, 30]. Synchronous averaging was also used for analysis of the vibration signals, which are characterized by the recurrence of two or more periods [3, 7, 9, 13, 18, 22]. We consider below its application for the estimation of BPCRP characteristics.

The coherent statistics of the BPCRP mean function have the form

mˆ l  1 Nl 1

  t nPl  , N l n 0 ( 11 ) where Pl is one of the non-stationarity periods and N l is the number of averaged periods. The mathematical expression of ( 11 ) for l  1 is equal to where

Emˆ t   1

If P1 nP2 and n is a natural number, then s N1 l P1 /P2   1 and Emˆ 1 t  m t  , i.e. formula ( 11 ) is the unbiased estimator of the BPCRP mean function. In other cases, formula ( 11 ) is a biased estimator of the mean additive component with periodP1 . The bias value depends on the ratio P1 /P2 and tends to zero as N 1   .

Using (11), we can form the formulae

mˆ k 0  1 P1mˆ 1 t e ik 2P1 tdt , mˆ 0l  1 P2mˆ 2 t e ik 2P2 tdt ,

P1 0 P2 0 which, in the general case, are the asymptotically unbiased estimators of the Fourier coefficients of the mean additive components.

It is easily see that unbiased estimators of the BPCRP mean function and its Fourier coefficients can be obtained using synchronous averaging with common period P : mˆ t   1 N1 t nP  , mˆ kl  1 Pmˆ t e i kltdt . ( 12 )

N n 0 P 0

Here N is the number of realization periods P which are averaged.

Taking into account ( 11 ), we can form the coherent estimators of the covariance function and its

Fourier coefficients:

( 13 ) ( 14 ) (17) Rˆ t ,   1 N

  t nP  mˆ t nP   t  nP  mˆ t  nP  , N n 0

Rˆkl    1 PRˆ t , e i kltdt .

P 0 Using the synchronous averaging of the BPCRP samples over one of the periods P2 in the form where L , r  1,2 are the numbers of the highest harmonics. The coefficients of the polynomials are r determined by the formulae mˆ kl  1 T

  t e kltdt ,

T T 1 T Rˆkl      t  mˆ t   t   mˆ t  e kltdt , (18)

T T where T is the length of signal realization. The number of harmonics to be taken into account in ( 15 ) and (16) can be obtained on the basis of the results of experimental data processing by means of the coherent method or stationary spectral estimation.

In the general case, employing formulae (17) and (18) leads to an increase of the additional errors caused by leakage effects. These effects are absent as T NP . Formulae (17) and (18) can then be rewritten in the form of ( 11 ) and ( 13 ). Indeed, mˆ kl 

Similarly,

Rˆkl    1 Pe i klt  1 N1  t nP  mˆ t nP   t  nP  mˆ t  nP dt .

P 0 N n 0 

The discrete estimators for the Fourier coefficients of the mean and covariance functions can be formed by substituting the integral transformations (17) and (18) by integral sums: 1 K 1 mˆ kl   t e i klnh ,

K n 0 Rˆkl rh   1 K 1

  nh  mˆ nh   n r h  mˆ n r h e i klnh .

K n 0 HereP1  M 1  1h , P2  M 2  1h and T Kh , where K rN M 1  1 nN M 2  1 .

To avoid the aliasing effects of the first and the second kinds [32], it is recommended to choose the sampling interval h in accordance with the inequalities h  Pi , h  Pi , i  1,2

2L1  1 4L2  1

If these inequalities are satisfied, the expressions ( 15 ) and (16) can be considered as the interpolation formulae for the estimators. We should note that in the case of T NP the component estimators coincide with the estimators determined by the least squares (LS) method [9, 31, 32]. However, using the LS method allows one to avoid the leakage errors in general case. These errors can be significant in cases when the values of the rotation frequency and/or their combinations are close. To construct the LS estimators for the mean and the covariance function we rewrite the series ( 6 ) and ( 7 ) in the form

M 1 m t  m  mlc coslt mls sinlt  , 0 l 1

M 2 R t ,  R 0     Rlc  coslt Rls  sinlt  l 1

, 1 1 where ml ml 1l2  2 mlc imls  , Rl   Rl 1l2    2 Rlc   iRls   ,

2 2 m0 m00 , R 0   R 00   , l  1l j j Pj , l 1  1,L1 , l 2  1,L2 . and N 1  2L1 L1  1 , N 2  2L2 L2  1 . The LS estimators for the Fourier coefficients of mean and covariance function are defined as the quantities which provide the minimum values of the quadratic functions

T   M1  F1 mˆ 0 ,mˆ 1c ,...,mˆ Mc1 ,mˆ 1s ,...,mˆ Ms1   0  t   mˆ 0  1 l mˆ lc coslt mˆ ls sinlt  dt , , 2 F2 Rˆ0   ,Rˆ1c   ,...,RˆMc 2   ,Rˆ1s   ,...,RˆMs 2   

T   M2     t ,   Rˆ0     Rˆlc  coslt Rˆls  sinlt  dt , 0   l 1  2 where  t ,    t  mˆ t   t   mˆ t   . They are the solutions of the system equations which represent the necessary conditions for the existence of the minimum of functionals (19) and (20): F1  0 , F1  0 , F1  0 ,r  1,M 1 , mˆ 0 mˆ rc mˆ rs (19) (20) (21) BFˆ20  0 BFˆr2c  0 BFˆr2s  0 , r  1,M 2 . (22)

The lag-dependent vanishing of the covariance function is the enough sufficient condition of the mean square consistency of the Fourier coefficients for the mean function. It also can be indicator of the asymptotic unbiasness of the estimators of covariance component. This condition is also sufficient for the consistency of mean square of the covariance component estimators for Gaussian BPCRP. For similar purposes, the series procedures were introduced, the latest of them are self-adaptive noise cancellation [33] and spectral method [34]. The best result are obtained using time synchronous averaging, however it requires a separate operations including individual resampling in each considered case.

5. Conclusions

The advantage of the LS estimators is the absence of the leakage effect. The possible bias of the LS estimators can be caused only by the previous inexact estimation of the mean function. When the realization length increases, values of the component estimators and the variances for the LS are quickly drawing together. So, the LS estimation can be rated as the preferable technique for statistical processing of the PCRP experimental time series.

6. References

References

[1] W. A. Gardner , Introduction to Random Processes with Applications to Signals and Systems , Macmillan, New York, 1985 .

[2] W. A. Gardner , Cyclostationarity in Communications and Signal Processing , IEEE Press, New York, 1994 .

[3] A. Napolitano , Generalizations of Cyclostationary Signal Processing: Spectral Analysis and Applications . IEEE Press, 2012 . doi: 10 .1002/9781118437926.

[4] A. Napolitano , Cyclostationary Processes and Time Series: Theory, Applications, and Generalizations. Elsevier, Academic Press, 2019 . 10 .1016/C2017-0-04240-4.

[5] E. G. Gladyshev , “ Periodically and Almost-Periodically Correlated Random Processes with a Continuous Time Parameter,” Theory Prob. Appl. , vol. 8 , pp. 173 - 177 , 1963 .

[6] Y. Dragan , I. Javorskyj , Rhythmics of Sea Waving and Underwater Acoustic Signals . Naukova Dumka, Kyjiv, 1980 (in Russian).

[7] Y. Dragan , V. Rozhkov , I. Javorskyj , The Methods of Probabilistic Analysis of Oceanological Rhythmics . Gidrometeostat, Leningrad, 1987 (in Russian).

[8] H. L. Hurd , A. Miamee , Periodically Correlated Random Sequences: Spectral Theory and Practice . Wiley, New York, 2007 . doi: 10 .1002/9780470182833.

[9] I. Javorskyi , Mathematical models and analysis of stochastic oscillations. Lviv, Karpenko Physico-Mechanical Institute of NAS of Ukraine, 2013 (in Ukraine) .

[10] V. Y. Mykhailyshyn , I. M. Javorskyi , Y. T. Vasylyna , O. P. Drabych and I. Yu . Isaev, “ Probabilistic models and statistical methods for the analysis of vibrational signals in the problems of diagnostics of machines and structures ,” Mater. Sci., vol. 33 , pp. 655 - 672 , 1997 .

[11] A. C. McCormick , A. K. Nandi , “ Cyclostationarity in rotating machine vibrations,” Mech. Syst. Signal Process ., vol. 12 , no. 2 , pp. 225 - 242 , 1998 . doi: 10 .1006/mssp. 1997 . 0148 .

[12] C. Capdessus , M. Sidahmed , J. L. Lacoume , “ Cyclostationary processes: Application in gear fault early diagnostics,” Mech. Syst. Signal Process ., vol. 14 , no. 3 , pp. 371 - 385 , 2000 . doi: 10 .1006/mssp. 1999 . 1260 .

[13] J. Antoni , F. Bonnardot , A. Raad , M. El Badaoui, “ Cyclostationary modeling of rotating machine vibration signals,” Mech. Syst. Signal Process ., vol. 18 , pp. 1285 - 1314 , 2004 . doi: 10 .1016/S0888- 3270 ( 03 ) 00088 - 8 .

[14] J. Antoni , “ Cyclostationarity by examples,” Mech. Syst. Signal Process ., vol. 23 , pp. 987 - 1036 , 2009 . doi: 10 .1016/j.ymssp. 2008 . 10 .010.

[15] R. B. Randall , J. Antoni, “ Rolling element bearing diagnostics - A tutorial,” Mech. Syst. Signal Process ., vol. 25 , no. 2 , pp. 485 - 520 , 2011 . doi: 10 .1016/j.ymssp. 2010 . 07 .017.