Cyclostationary modeling for local fault diagnosis of planetary gear vibration signals

Cyclostationary modeling for local fault diagnosis of planetary gear vibration signals

Journal Pre-proof Cyclostationary modeling for local fault diagnosis of planetary gear vibration signals Ruo-Bin Sun, Zhi-Bo Yang, Konstantinos Grylli...

10MB Sizes 0 Downloads 35 Views

Journal Pre-proof Cyclostationary modeling for local fault diagnosis of planetary gear vibration signals Ruo-Bin Sun, Zhi-Bo Yang, Konstantinos Gryllias, Xue-Feng Chen

PII:

S0022-460X(20)30006-7

DOI:

https://doi.org/10.1016/j.jsv.2020.115175

Reference:

YJSVI 115175

To appear in:

Journal of Sound and Vibration

Received Date: 20 June 2019 Revised Date:

29 December 2019

Accepted Date: 4 January 2020

Please cite this article as: R.-B. Sun, Z.-B. Yang, K. Gryllias, X.-F. Chen, Cyclostationary modeling for local fault diagnosis of planetary gear vibration signals, Journal of Sound and Vibration (2020), doi: https://doi.org/10.1016/j.jsv.2020.115175. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier Ltd.

Cyclostationary modeling for local fault diagnosis of planetary gear vibration signals Ruo-Bin Suna,c , Zhi-Bo Yanga,b,∗, Konstantinos Grylliasc,d , Xue-Feng Chena,b a

School of Mechanical Engineering, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, China b State Key Laboratory for Manufacturing Systems Engineering, Xi’an Jiaotong University, Xi’an, Shaanxi 710049, China c Department of Mechanical Engineering, KU Leuven, Celestijnenlaan 300, Leuven 3001, Belgium d Dynamics of Mechanical and Mechatronic Systems, Flanders Make, Belgium

Abstract Owing to the rotation and reciprocation of mechanical equipment, their vibration signals inherently exhibit cyclic stationary characteristics. This paper provides a phenomenological signal model of the planetary gear vibrations under the cyclostationarity paradigm. The model is helpful to give comprehensive demodulation information which is the limit of classic spectral analysis of deterministic signals, so that the diagnosis of local faults in planetary gear sets can be more explicit. First, the effects of rotating speed fluctuation on periodic gear meshing vibrations are analyzed. It is shown that the stationary speed fluctuation will turn the periodicity of signals into second-order cyclostationarity. Next, an impulsive oscillation sequence with interval timing jitter is used to describe the dynamic response of the meshing stiffness change caused by local defects. Meanwhile, the influence of a ∗

Corresponding author Email address: [email protected] (Zhi-Bo Yang)

Preprint submitted to Journal of Sound and Vibration

January 6, 2020

vibration transfer path is also considered in the model. Based on the above analysis, the following procedure is utilized to analyze planetary gearbox vibrations: remove the periodic components by the self-adaptive noise cancellation method and then, apply the second-order frequency descriptor to reveal fault features of the signals. Both the simulation and experiment are used to validate the correctness of model derivation and the superiority of the cyclostationary analysis. Keywords: cyclostationary modeling, self-adaptive noise cancellation, cyclic spectral coherence, planetary gear fault diagnosis 1. Introduction Almost all machines have a periodic form of motion because of their design. Consequently, their vibration signals are usually conjectured to have periodic or cyclostationary characteristics [1]. Different from the periodicity of deterministic signals, cyclostationarity characterizes the non-stationary nature of stochastic signal process, which is more suitable and precise in describing the randomness of machinery vibrations. A signal is said to exhibit cyclostationarity if a cascade of linear and non-linear transformations produces a periodic component [2]. Under this broad definition, cyclostationary analysis is intended to reveal the harmonic components in random interference. More importantly, the non-linear transformations, such as the square of signals, enable the cyclostationary analysis to further extract periodic signals from the multiplicative noise [3]. Historically, the concept of cyclostationarity, which broke through the usual stationary assumption of stochastic signals, was first proposed in communication research [4]. In recent years,

2

the study of cyclostationary analysis methods has aroused widespread interest among the mechanical fault diagnosis research group. Some remarkable works have been carried out on the vibration fault diagnosis of rolling element bearings [5–8], gears [9, 10], and reciprocating machinery [11]. A series of outstanding works by Antoni et al. [12, 13] promotes the wide application of the cyclostationary model in mechanical vibration analysis. A thorough understanding of the structure of cyclostationary statistics is a prerequisite for fault diagnosis of a specific engineering object. That is to say, a cyclostationary modeling of vibration signals of the diagnostic target should be built in advance. The authors of [14] exhibited the modulation mechanism in rolling element bearing vibrations based on the cyclostationary analysis and demonstrated the model with two different types of bearing faults. The authors of [15] analyzed the dynamic behavior of hydroelectric turbines using cyclostationary modeling, and the model helped extrapolate small stationary load strain samples on turbine blades. The authors of [16] compared differences in cyclostationary models of internal-combustion engines, gearboxes, and rolling element bearings, which differ in the ratio between periodic and purely high-order cyclostationary parts. The planetary gear sets are important transmission components and are widely used in heavy-duty and large transmission ratios. However, the cyclostationary modeling of planetary gear vibration signals is still indigent. Our paper tries to fill this gap and make a more accurate diagnosis of local faults in planetary gear sets through cyclostationary analysis. It is worth mentioning that there is a certain similarity between the phenomenological model of planetary gear sets and of rolling bearings because of similar forms of motion. However, the

3

gear sets contain more harmonic vibration components, which also present cyclostationary features. Besides, the specific fault feature frequencies of the planetary sets and bearings will not affect the diagnosis. Cyclostationary modeling, also known as a phenomenological modeling, is a type of signal behavior description method. Despite the lack of exact expressions of dynamic equations, that is, differential equations of mass, damping, stiffness, and displacement, it uses specific statistical descriptors to directly specify the vibration response of healthy and defective systems. The planetary gear set dynamic models [17–19] mainly provide the modal information for optimizing the mechanical design. As for a component with relatively complex structure and motion relations, considering various physical factors in the dynamic model is a challenge. Generally, the defects are modeled as changes in parameters of dynamic equations. However, the system responses calculated from modified equations often have considerable differences with actual vibrations of a defective system [20]. This weakens the guidance significance of models for fault signal processing. On the contrary, phenomenological models derive and summarize the fault feature expression in their diagnostic plot according to the analysis of kinematic relations and the observation of sufficient experimental data. In a previous work, the phenomenological modeling of planetary gearboxes was primarily based on harmonic analysis of deterministic signals. McFadden et al. [21] explained asymmetric sidebands around the meshing frequency in the spectrum of planetary gear vibration using a phenomenon model. Feng et al. [22] derived the spectral structure of planetary gear sets using a periodic deterministic model, which considered both the amplitude and frequency modu-

4

lation effects of gear damage. Further, the vibration signal models for planet bearing fault diagnosis were deduced by a similar method [23]. Lei et al. [24] also adopted a Fourier series analysis of algebraic equations to describe the phenomena of fault spectrum of planetary gearboxes, which considered all possible vibration transfer paths from meshing points to the transducer. These models are proven to be effective in practice, but provide the fault information only in the frequency domain. In the cyclostationary paradigm, the frequency of both modulated and carrier signals can be represented in cyclic frequency and spectral frequency domains. This representation is ideal for analyzing complex modulation relations of planetary gear sets. The structure of this paper is as follows. In section 2, the basic concepts of cyclostationarity are reviewed. Next, two primary factors are analyzed in detail: the meshing vibration of gears and the impulse response caused by a local fault. It is worth mentioning that the stationary fluctuation of the shaft speed makes original angle periodic signal cyclostationary in time. In section 3, we propose an analysis procedure of planetary gear signals using cyclostationarity. First, removal of the periodic components using the self-adaptive noise cancellation method is recommended. Then, we utilize the second-order cyclostationary descriptors, cyclic spectral correlation, or cyclic spectral coherence to identify fault feature. Sections 4 and 5 provide the simulation analysis and experimental verification, respectively. Finally, section 6 concludes the paper.

5

Figure 1: The two types of vibrations with amplitude modulation because of carrier rotation.

2. Cyclostationary modeling of planetary gear sets Gear meshing is the primary vibration source of planetary gearboxes. Besides, a local defect of tooth changes its meshing stiffness. To all appearances, stiffness change excites impulsive oscillations in the vibration response, which is a distinct characteristic of local faults [25]. Both types of vibrations show amplitude modulation, because the meshing point periodically moves relative to the fixed transducer, as shown in Fig. 1. They have a low correlation in the frequency domain. Therefore, we ignore the cross-term interference and analyze their cyclostationary feature separately in the following analysis for the sake of simplicity. 2.1. Basic concepts Similar to the definition of strict-sense stationary, a stochastic process X(t) is called strict-sense cyclostationary of order k, if arbitrary k time in6

stants t1 , t2 , · · · , tk have a kth-order probability distribution functionwith period T . Fx (x1 , x2 , · · · , xk ; t1 , t2 , · · · , tk ) = Fx (x1 , x2 , · · · , xk ; t1 + T, t2 + T, · · · , tk + T ) (1) Furthermore, the stochastic process is wide-sense cyclostationary of order k only if the moments or cumulants of X(t) periodically vary with time: E{

k Y

X(ti )} = E{

i=1

k Y

X(ti + T )}

(2)

i=1

where E{·} denotes expectation operator, which also means ensemble average. In general, the study of cyclostationary stochastic signals is considered as wide-sense, and specifically, the first-order and second-order statistical analyses are the most commonly used methods in practical engineering. The first-order moment is the expected value function of a signal. The expected value being periodic means the signal shows first-order cyclostationarity (CS1). ∆

mx (t) = E{X(t)} = mx (t + T )

(3)

This feature is used to describe a signal model of the period components along with the additive stationary noise, a standard assumption in rotating machinery. In stationary paradigm, the period components should be carefully explained to satisfy the wide-sense stationary requirement, such as sine waves with a random phase of uniform distribution. From the CS1 perspective, the periodic components would have a broader interpretation, which could be considered as deterministic signals. The second-order cyclostationarity (CS2) is a more useful definition of 7

fault diagnosis, which implies the second-order moment, that is, the autocorrelation function is periodic. ∆

R2x (t1 , t2 ) = E{X ∗ (t1 )X(t2 )} = R2x (t1 + T, t2 + T )

(4)

CS2 is exceptionally useful for expressing modulation characteristics. The amplitude and frequency modulation phenomenon of vibrations can be precisely modeled as CS2, and the multiplicative noise can also be regarded as a type of modulation. It is easy to verify that a periodic signal with stationary noise is not only CS1 but also CS2. But this type of signals can be simply analyzed by a Fourier transform. Moreover, such components would interfere with the demodulation analysis of second-order statistics. To solve this issue, the cumulant functions are used to define purely kth-order cyclostationary (PCSk) process [26], which means that only the kth-order cumulant is periodic. As for the periodic signal x(t), the second-order cumulant (autocovariance function) is zero. Cov2˜x (t1 , t2 ) , E {(˜ x(t1 ) − E{˜ x(t1 )})∗ (˜ x(t2 ) − E{˜ x(t2 )})}

(5)

= R2˜x (t1 , t2 ) − m∗x˜ (t1 )mx˜ (t2 ) = 0 where ∗ denotes the complex conjugation. Since the higher cumulants hold the same results, x˜(t) is the PCS1 signal. To analyze the modulation information of vibrations, it is vital to subtract PCS1 components in advance [16]. It is more convenient to comprehend the CS2 signals in the frequency domain. Using the Wiener-Khinchine relationship, the spectral correlation density function can be gained by the 2D Fourier transform of autocorrelation

8

function:

Z∞ Z∞ SC2x (α, f ) =

R2x (t, τ )e−j2παt e−j2πf τ dtdτ

(6)

−∞ −∞

The practice is accustomed to taking central time t = (t1 + t2 )/2, and τ = t2 − t1 . α is named as the cyclic frequency corresponding to t, which gives the modulated information of the signal, and f is the spectral frequency corresponding to τ , which represents the frequency distribution of carrier wave. 2.2. Gear meshing vibrations Gear vibration signals inherently periodic with the rotating angle θ of the gear shaft. For this reason, gear vibrations are considered as sine waves with amplitude and frequency modulation, as follows: x(θ) =

X

ak (θ) cos(2πkfm θ + bk (θ)),

k ∈ N+

(7)

k

where kfm is the kth-order meshing frequency of gears, and ak (θ) and bk (θ) denote the amplitude and frequency modulation of the kth harmonic respectively. ak (θ) and bk (θ) can be expanded as complex exponential signals by the Fourier series, whose frequency set A is composed of several orders of rotating frequencies of the gear shafts fs or of the planet carrier fc . ak (θ) =

X

αi ej$i θ

i;$i ∈A

bk (θ) =

X

(8) j$i θ

βi e

i;$i ∈A

From the convolution property of the Fourier transform, that is, when the product in the time domain is equivalent to the convolution in the frequency 9

domain, the spectrum of this signal model has the spectral line at kfm ± ifc (k, i ∈ N+ ). Here, the shaft frequency component is omitted because its modulation effect is much weaker than planet carrier. This model is the classic deterministic periodic signal model of gear vibrations, a PCS1 process [13]. In practice, the rotating speed cannot be perfectly constant because of driving power and load fluctuations. This slight fluctuation would alter the statistical property of angular periodic signals [16]. Assume the angular displacement has the following relationship with the rotating speed: Z t V (u)du θ(t) = ω0 t + Φ(t) = ω0 t +

(9)

−∞

where ω0 is the nominal rotating speed and V (t) is the random fluctuation of speed. Therefore, the angular position fluctuation Φ(t) is also a stochastic process. Consider a periodic signal of θ having the Fourier series form as given below: x(θ) =

X

ci ej$i θ

(10)

i

After substituting θ(t) into this equation to replace the constant theta, a stochastic signal relative to time variable t becomes: X(t) =

X

ci ej$i Φ(t) ejωi t

(11)

i

where ωi = ω0 $i and $i is the angular frequency relative to the angle θ. Then, the first-order and second-order moments are calculated to verify the cyclostationary property of X(t). mx (t) =

X i

j$i Φ(t)

E{e

}ci e

jωi t

=

X i

10

jωi t

ci e

Z

ej$i ϕ(t) pϕ (ϕ; t)dϕ

(12)

with pϕ (ϕ; t) being the probability density function of Φ(t). It is obvious that mean function is periodic if pϕ (ϕ; t) is a periodic function relative to t: pϕ (ϕ; t) = pϕ (ϕ; t+T ). In a special case, if the speed fluctuation is stationary, that is, if pϕ (ϕ; t) is constant relative to t, the signal X(t) is CS1. Then, the second-order moment of the signal is: R2x (t1 , t2 ) =

XX i

=

XX i

ci cj e

E{ej$j Φ(t2 )−j$i Φ(t1 ) }ci cj ejωj τ ej(ωj −ωi )t1

j

jωj τ j(ωj −ωi )t1

e

ZZ

(13) ej$j ϕ2 −j$i ϕ1 p2ϕ (ϕ1 , ϕ2 ; t1 , t2 )dϕ1 dϕ2

j

where τ = t2 − t1 ; p2ϕ (ϕ1 , ϕ2 ; t1 , t2 ) is the joint probability density function of Φ(t1 ) and Φ(t2 ). If p2ϕ (ϕ1 , ϕ2 ; t1 , t2 ) is periodic or if the special case of speed fluctuation is stationary, the second-order moment is periodic, that is, the signal is CS2. It can be further proven that the second-order cumulant Cov2x (t1 , t2 ) is periodic, when random variables Φ(t1 ) and Φ(t2 ) are correlated. From the above derivation, the periodic signal of angular displacement becomes a CS2 signal relative to time after considering the influence of the speed fluctuation. Going a step further, the determination term in autocorrelation function (Eq. 13) is the same as the analysis of the periodic signal. The integral term causes some changes in cyclostationary characteristics. It has been summarized that features of the modulated signal model in spectral correlation density SC2x (α, f ) are some discrete intersections at a series of lines: f = ±α/2 ± (kfm ± ifc ), (k, i ∈ N+ ) [12]. The integral term blurs the discrete point, which is explicitly determined by statistical characteristics of Φ(t).

11

2.3. Impulses of local fault The local fault of gear tooth, such as a break, crack, and chip, would cause the meshing stiffness to be discontinuous with gear rotation. The vibration response is an impulsive oscillation sequence with interval timing jitter caused by speed fluctuations. The signal could be modeled as: Y (t) =

X

Ai s(t − iT − τi )

(14)

i

where s(t) is a single attenuated waveform excited by the broken tooth participating in the mesh. T denotes the impulses period determined by the rotating period of the shaft with the fault gear and Ai implies the amplitude of the ith impact. Owing to rotating speed fluctuation, the interval time has a random jitter when compared with the average period T . The stochastic stationary point process {τi }i∈Z is used to describe this time jitter, and it is assumed that {τi }i∈Z has the following conditional probability: P [iT + τi |jT + τj ] = P [iT + τi ], j < i

(15)

That is to say, the following equivalent statistical properties of moments are satisfied: E{τi } = 0

(16)

E{τi τj } = σδ2 δ(i − j) with σδ being the standard deviation of the process and δ(n) being the delta function. Besides, the meshing position of a defective tooth is time-varying along with the carrier rotation. The impulse response of the system between different points on the ring gear to the sensor is different in magnitude. Therefore, the amplitude of impulses {Ai }i∈Z is also modeled as a stochastic point process. However, the statistical properties of {Ai }i∈Z , which are 12

slightly different from stationary process {τi }i∈Z , have the Q period corresponding to the carrier rotation: E{Ai } = E{Ai+Q } = A¯i

(17)

E{Ai Aj } = E{Ai+Q Aj+Q } = A¯2i δ(i − j) Next, the cyclostationary feature of Y (t) is analyzed. The mean function is derived as: my (t) =

X

=

X

A¯i E{s(t − iT − τi )}

i

A¯i

(18)

Z s(t − iT − τ )pτ (τ ; i)dτ

i

where pτ (τ ; i) denotes the probability density function of τi . Since {τi }i∈Z is stationary, it can be concluded that my (t) is a periodic function. The autocorrelation function is calculated as: X X τ τ R2y (t, τ ) = E{[ Ai s(t − − iT − τi )]* [ Aj s(t + − jT − τj )]} 2 2 i j XX τ τ = E{Ai Aj }E{s∗ (t − − iT − τi )s(t + − jT − τj )} 2 2 i j X τ τ = A¯2i [s∗ (t − − iT )s(t + − iT )] ⊗ pτ (t) 2 2 i (19) where the symbol ⊗ means the convolution operation. It is easy to find that the autocorrelation function has a common cycle of T and Q, which proves that Y (t) is a CS2 signal model. To further formulate the frequency feature

13

of the model, the spectral correlation density function is given as: SC2y (α, f ) =

F

2

t→α,τ →f

Z∞ Z∞ {R2y (t, τ )} =

R2y (t, τ )e−j2παt e−j2πf τ dtdτ

−∞ −∞

=

Z∞ Z∞ X −∞ −∞

i

(20)

τ τ A¯2i e−j2πiT α [s* (t − )s(t + )]e−j2παt e−j2πf τ pˆτ (α)dtdτ 2 2

with pˆτ (α) being the Fourier transform of pτ (t). Because of A¯2i is periodic, it could be represented as the following Fourier series: A¯2i =

X

ak ejk(2π/Q)(T i) =

k∈Z

X

ak ej2πT kα2 i

(21)

k∈Z

According to the Poisson formula with α1 = 1/T and α2 = 1/Q: T

X

ej2πT i(kα2 −α) =

i∈Z

X

δ(α − iα1 − kα2 )

(22)

i∈Z

and the Fourier transform relation: τ τ α α F 2 {s∗ (t − )s(t + )} = sˆ∗ (f − )ˆ s(f + ) t→α,τ →f 2 2 2 2

(23)

Eq. 20 can be simplified as:   ∞ ∞ Z Z  τ τ SC2y (α, f ) = [s* (t − )s(t + )]e−j2παt e−j2πf τ dtdτ   2 2 −∞ −∞

×

=

     

1 T

PP k∈Z i∈Z

1 XX ak δ(α − iα1 − kα2 )ˆ pτ (α) T k∈Z i∈Z

sˆ∗ (f − α2 )ˆ s(f + α2 )ak pˆτ (α), α = iα1 + kα2 , i 6= 0 or k 6= 0

SC2y (f ),      0,

i = 0 and k = 0 others (24) 14

It shows that when cyclic frequency equals iα1 ± kα2 , (i, k ∈ N+ ), there are spectral lines on around the frequency band of oscillation in the spectral frequency domain. Moreover, the Fourier transform of the probability density function pˆτ (α) acts as a low-pass filter, making the modulated feature focus on the low-frequency part in the cyclic frequency domain. In the case of α = 0, the spectral correlation is reduced to a power spectral density function of the signal SC2y (f ). To summarize, the superposition of the shaft frequency with the faulty gear and with the carrier rotation frequency is the modulated fault feature. This can help obtain a comprehensive feature expression of both the modulated information and the carrier waves in the cyclostationary model. 3. Analysis procedure As previously mentioned, periodic (PCS1) components of planetary gear vibration signals should be extracted and removed prior to the second-order analysis. The most popular first-order cyclostationary analysis method is the time-domain synchronous average [27], which is widely used in fault diagnosis of gearboxes. However, only a single periodic component can be separated at one time. This is similar to utilizing the following operator to extract the component at frequency α [3].   Z 1 −j2παt Pα {·} = lim (·)e dt ej2παt T →∞ T T

(25)

The frequency of all periodic components should be known beforehand, or one must use statistical test methods to determine periodic components in the noise, which may not yield satisfying results for complex modulation signals. Therefore, another strategy, self-adaptive noise cancellation (SANC) 15

[28], is adopted to separate the deterministic part from other random-like vibrations. The SANC applies an appropriate delay as a reference signal to filter out the correlated component in the primary input. Thus, the periodic components which are correlated between the inputs are persisted with, and the other uncorrelated random part will be separated at the error output. The schematic diagram of the SANC is depicted in Fig. 2. Some adaptive

Figure 2: The block diagram of the SANC.

filter methods, such as those using the least mean square (LMS) or the recursive least square (RLS) algorithms, can be used to design the filter. It is worth noting that a relatively longer time delay and a larger order of the filter will guarantee better results [28]. Then, the second-order cyclostationary descriptors in the frequency domain are utilized to find the existence of fault features. The spectral correlation density function is the most commonly used frequency statistics. The SC2x (α, f ) defined before (Eq. 6) is a 2D Fourier transform of R2x (t, τ ). However, the SC2x function is discrete in the α domain and continuous in the f domain because of a CS2 signal is periodic in variable t. Therefore, the Fourier series expansion can also be applied to t, to avoid using the singular Dirac delta function δ(t) in the spectrum. The resulting descriptor is named the cyclic power spectrum S2x (f ; α), which has the following relationship 16

with the spectral correlation. SC2x (α, f ) =

X

S2x (f ; αi )δ(α − αi )

(26)

αi ∈A

where A is the cyclic frequency set. From another perspective, the cyclic power spectrum measures the linear dependence between a signal with symmetric frequency shift as given in Eq. 23. Similar to the coherence function of two stationary signals, the correlation between different frequency shift signals can be normalized by their respective energies. The descriptor is coined as spectral coherence, which is calculated by S2x (f ; α) γ2x (f ; α) = p S2x (f + α/2; 0)S2x (f − α/2; 0)

(27)

The estimation theory of power spectrum of stationary signals can be generalized into the cyclic spectral estimation. The classic non-parametric estimators, such as the Welch method, are the most popular with a relatively complete theoretical system. The estimation performance can be improved by smoothing frequency bands and averaging between different spectra. In [13], the performance of several cyclic estimators are discussed in detail. On the whole, the choice of window function of the estimator must balance the resolution and estimation bias. Besides, the estimation variance is inversely proportional to recorded cycles of the signal. As derived in section 2, if a local fault occurs in a component of the planetary set, two types of fault features will be generated, in the spectral correlation and/or in the spectral coherence. First, the feature lines will be on the meshing frequency of planetary sets (and its higher orders) in the spectral frequency domain fm , and simultaneously, these features will lie on 17

the carrier rotation frequency fc (and its higher orders) in the cyclic frequency domain. Besides, the impulses of local fault show features on the impulse response frequency in the spectral frequency domain, and on the rotation frequency of the fault gear relative to the carrier, that is, fs , fp , and fr , (and their higher orders) in the cyclic frequency domain. The meshing frequency, carrier rotation frequency, and different gear rotation frequencies relative to the carrier are deduced in [22], as shown in the following formulas: Zs Zr (Zs + Zr )

(28)

fc = fr/c = fm /Zr

(29)

fm = fi

fp/c = k

fm , (k = 1, 2) Zp

(30)

fm Zs

(31)

fs/c = N

where fi is the input shaft frequency (shaft frequency of the sun); fm and, fc represent the meshing frequency and the carrier rotation frequency respectively; fp/c , fr/c , and fs/c are rotation frequencies of the sun, planet, and ring, respectively; Zs , Zr , and Zp represent the tooth number of sun, ring, and planet, respectively; and N denotes planet number; Moreover in Eq. 30, if the impacts of sun-planet and ring-planet pairs are similar, then k = 2, otherwise k = 1. 4. Numerical simulation In this section, two numerical simulations validate the developed model by the suggested analysis procedure. The first simulation analyzes the periodic amplitude modulation signal of rotating angle, which turns out to be CS2 18

signals with rotating speed fluctuation. The angular periodic simulation signal is formulated as: x(θ) = (1 + cos(2πf2 θ)) · cos(2πf1 θ)

(32)

which is a simple amplitude modulation function; besides, no higher harmonics are included in the signal. The carrier frequency f1 = 10 (1/rad), and the modulated frequency f2 = 1 (1/rad). The rotating speed is set as ω0 = 2πf0 = 60π rad/s with sampling frequency Fs = 1000 Hz. The speed fluctuation is stationary with the standard deviation equal to 5%ω0 and 15%ω0 , respectively. If there is no velocity fluctuation, the resulting time signal x(t) theoretically has the carrier frequency of 300 Hz, and the modulated frequency is 30 Hz. Fig. 3 displays the segments of temporal waveforms and the corresponding amplitude spectrum of Fourier transform with no speed variation and the fluctuation of σ = 5%ω0 , and σ = 15%ω0 , respectively. It can be found that there is not a clear distinction in time domain, which is caused by a relatively low sampling frequency. The spectrum, however, clearly shows the impact of speed fluctuation. The frequency becomes increasingly blurred with the increase in fluctuation variance. Then, the signals are processed by the suggested procedure. The resulting cyclic power spectra are shown in Fig. 4. It is worth noting that the signal is deterministically periodic without fluctuation, and thus, does not require the periodic extraction procedure. Hence, Fig. 4 (a) is the cyclic power spectrum of the periodic signal. It can be seen that the angular periodic signal becomes CS2 signal with speed fluctuation. The feature points, which lie in a triangular region, provide the information for both the carrier and modulated frequency. Moreover, the 19

Figure 3: The waveform and spectrum of amplitude modulation signals with different speed fluctuation. (a, d) the periodic signal; (b, e) the standard deviation of fluctuation σ = 5%ω0 ; (c, f) the standard deviation of fluctuation σ = 15%ω0 .

20

Figure 4: The cyclic power spectrum S2x (f ; α) of the simulation signals (unit: g2 /Hz ). (a) the periodic signal; (b) σ = 5%ω0 ; (c) σ = 15%ω0 .

21

fluctuations blur the discrete feature point, which is in agreement with the previous model analysis. The other numerical analysis is the simulation of impulsive oscillation sequence caused by a local fault. The simulation signal y(t) is a sample function of the stochastic signal model in Eq. 14, where s(t) is the free response of a damping system with a single degree of freedom. s(t) = e−at cos(2πf0 t)

(33)

Among them, a denotes the attenuation coefficient, which is set as 0.15; the natural frequency of the oscillation f0 is set as 400 Hz. The nominal period of impulses T is equal to 0.05 s, and the time jitter {τi } is stationary with standard deviation σ = 5%T = 2.5 × 10−3 s. The magnitude of impulses {Ai } is the sampling of periodic Hanning window when the impacts occur. The period of the window function is 0.2 s; and the sampling frequency of the discrete signal is 2000 Hz. A segment of the sampling function y(t) and its corresponding amplitude spectrum is given in Fig. 5. The ambiguity of the spectrum, which differs from

Figure 5: (a) A segment of sampling function of the stochastic signal Y (t) with interval timing jitter. (b) The amplitude spectrum of the sample y(t).

22

Figure 6: (a, c) The periodic components and its amplitude spectrum extracted by the SANC. (b, d) The residual part of the simulation signal and its amplitude spectrum.

periodic signals is caused by the time jitter. First, the periodic components is extracted by the SANC, as plotted in Fig. 6. Fig. 6(c) shows that most of the spectral lines of extracted periodic components are in evenly spaced intervals, although some non-periodic features still can be seen. Then, the cyclic power spectrum of the residual signal is given in Fig. 7. It is clear that the spectral frequency shows the natural frequency of oscillations; besides, in cyclic frequency domain, the feature is the frequency of impulses, and the sideband of interval equals the modulation frequency. It is verified that the impulses caused by a local fault are the CS2 signal, and the fault features are the same as in the model analysis.

23

Figure 7: The cyclic power spectrum of the residual part of the simulation signal (unit: g2 /Hz ).

5. Experimental analysis We conducted an experiment on the Spectra Quest (SQ) fault simulator to further verify the cyclostationary model of planetary gear vibration signals. The test rig is composed of an alternating current motor, a planetary gearbox, a fixed-axis gearbox, and a magnetic powder brake, as shown in Fig. 8. The

Figure 8: (a) The SQI fault simulator. (b) The internal view of the planetary gearbox.

24

Table 1: The gear tooth number and the rotating frequency relative to the carrier rotation when the input shaft rotating frequency is 1 Hz.

Tooth number

Sun

Planet

Ring

Carrier

fm

28

36

100





0.61

0.22

0

21.88

Relative frequency 0.78

planetary gearbox has only one stage, whose tooth number is listed in Table 1. The speed reduction system sets the sun gear shaft as the input shaft and the carrier as the output, and fixes the ring gear. In the meantime, the meshing frequency and the rotating frequency relative to the carrier of each component are also calculated using Eq. 28 — Eq. 31. A tooth-chipped local fault is prefabricated on a planet gear. The vibration signals are sampled by accelerometers arranged on the gearbox housing, and processed through a commercial data acquisition system. The input shaft speed is around 1800 RPM with stationary rotating fluctuation, and the sampling frequency is set at 6000 Hz. Fig. 9 depicts the estimated shaft rotating frequency using a tachometer which acquires one pulse per shaft revolution. We simply calculate the average frequency of each revolution (between two consecutive pulses), so there are jump break points in the plot, which should not exist in an actual situation because of the rotational inertia.

A seg-

ment of the acquired acceleration signal is depicted in Fig. 10. The signal shows a large degree of randomness, accompanied by some distinct impulses. According to the suggested cyclostationary analysis procedure, the signal is processed by the SANC method first. The extracted periodic components and the random-like residuals are given in Fig. 11. From the amplitude differ25

Figure 9: The estimated shaft rotating frequency through a tachometer.

Figure 10: A segment of acquired signal on the planetary gearbox.

Figure 11: The extracted periodic components and the random-like residuals of the experimental signal. (a) The periodic component; (b) the random components.

26

ence of components, the residual random part occupies a larger proportion in the experimental signal. Furthermore, the second-order cyclostationary descriptor is calculated. To highlight the weak energy components, the normalized statistic and the spectral coherence are given in Fig. 12. First, the

Figure 12: The spectral coherence of the residual part of the simulation signal.

carrier rotating frequency fc and its higher orders can be observed in the cyclic frequency domain. Besides, these components are mainly distributed around the fourth-order meshing frequency. These features accord with the cyclostationary modeling of gear meshing vibrations. Here, an unusual phenomenon is that the modulated frequency is the rotating frequency of the carrier rather than the planet, which is different from the common fault feature of fixed-axis gearboxes. We suppose this to be because of the broken tooth on a single planet, which causes an imbalance in the whole carrier system. The rotating frequency is precisely the fault feature of an imbalanced 27

shaft. Moreover, the ifp/c ± kfc (i, k ∈ N+ ) features can be seen in the cyclic frequency domain, where fp/c is the rotating frequency of planet shaft relative to the carrier. However, the spectral frequency of these components lies in a broader frequency band. This is because these fault frequencies are caused by impulses excited through meshing stiffness change. The oscillation frequency of the impulses is different from the meshing frequency. This experimental analysis validates the effectiveness of the proposed cyclostationary model. To further explain the superiority of cyclostationary analysis, a comparison with the classic planetary gear diagnosis method based on the harmonic analysis is carried out. To begin with, we plot the amplitude spectrum of the signal and its detailed view around the meshing frequency (Fig. 13). It can be found that the energy of the signal mainly concentrates around the higher orders of the meshing frequency, which agrees with the spectral coherence analysis. The sideband around the meshing frequency has some confusions that are hard to identify. We further plot the envelope demodulation spectrum of the signal as shown in Fig. 14. The carrier frequency and the planet shaft rotating frequency can be observed here. However, the features of gear mesh and impulses are mixed up. In fact, the envelope demodulation spectrum is equal to integration along the spectral frequency axis in the spectral coherence. The spectral coherence function, in contrast, clearly shows the different features of these components. 6. Conclusion Based on the observation of a kinematic pattern and the vibration phenomenon, a cyclostationary modeling of planetary gear vibration signals is 28

Figure 13: The amplitude spectrum of the signal and its detail view around the meshing frequency.

Figure 14: The envelope demodulation spectrum of the simulation signal.

29

proposed in this paper. The vibrations of planetary gear sets are divided into two main categories: the gear meshing vibrations and the impulses caused by local fault. Considering the modulation of a carrier rotation and the speed fluctuation, both components are modeled in the cyclostationary paradigm. The stationary speed fluctuation, or more generally, the periodic probability density function of fluctuation will cause the angular periodic (PCS1) signals to turn into CS2. The feature of these components in the spectral correlation function is the blurring of the discrete feature point of the corresponding periodic signal. Moreover, the fluctuation makes the equal intervals of original periodic impulses non-equal with timing jitter. In second-order cyclostationary analysis, the oscillation frequency is shown in spectral frequency domain, and the nominal period of impulse and the frequency of amplitude modulation are given in the cyclic frequency domain. As for the process of planetary gear vibrations, we suggest that periodic components should be removed by the SANC method first. Then, the cyclostationary fault feature can be revealed by the cyclic power spectrum or the spectral coherence function. The numerical simulation and experimental analysis verifies the effectiveness and superiority of the proposed model. Acknowledgement We thank the supports given by the National Natural Science Foundation of China (No. 51875433), the Young Talent fund of University Association for Science and Technology in Shaanxi of China (No. 20170502), the Natural Science Foundation of Shaanxi province (No. 2019KJXX-043), the Fundamental Research Funds for the Central Universities, and the China 30

Scholarship Council. [1] Z. Peng, F. Chu, Application of the wavelet transform in machine condition monitoring and fault diagnostics: a review with bibliography, Mechanical systems and signal processing 18 (2) (2004) 199–221. URL https://doi.org/10.1016/S0888-3270(03)00075-X [2] W. A. Gardner, A. Napolitano, L. Paura, Cyclostationarity: Half a century of research, Signal processing 86 (4) (2006) 639–697. URL https://doi.org/10.1016/j.sigpro.2005.06.016 [3] J. Antoni, Cyclostationarity by examples, Mechanical Systems and Signal Processing 23 (4) (2009) 987–1036. URL https://doi.org/10.1016/j.ymssp.2008.10.010 [4] W. Gardner, L. Franks, Characterization of cyclostationary random signal processes, IEEE Transactions on information theory 21 (1) (1975) 4–14. URL https://doi.org/10.1109/TIT.1975.1055338 [5] J. Antoni, Cyclic spectral analysis of rolling-element bearing signals: Facts and fictions, Journal of Sound and vibration 304 (3-5) (2007) 497– 529. URL https://doi.org/10.1016/j.jsv.2007.02.029 [6] D. Abboud, J. Antoni, M. Eltabach, S. Sieg-Zieba, Angle-time cyclostationarity for the analysis of rolling element bearing vibrations, Measurement 75 (2015) 29–39. URL https://doi.org/10.1016/j.measurement.2015.07.017 31

[7] W. Teng, X. Ding, Y. Zhang, Y. Liu, Z. Ma, A. Kusiak, Application of cyclic coherence function to bearing fault detection in a wind turbine generator under electromagnetic vibration, Mechanical Systems and Signal Processing 87 (2017) 279–293. URL https://doi.org/10.1016/j.ymssp.2016.10.026 [8] D. Wang, X. Zhao, L.-L. Kou, Y. Qin, Y. Zhao, K.-L. Tsui, A simple and fast guideline for generating enhanced/squared envelope spectra from spectral coherence for bearing fault diagnosis, Mechanical Systems and Signal Processing 122 (2019) 754–768. URL https://doi.org/10.1016/j.ymssp.2018.12.055 [9] Z. K. Zhu, Z. H. Feng, F. R. Kong, Cyclostationarity analysis for gearbox condition monitoring: Approaches and effectiveness, Mechanical Systems & Signal Processing 19 (3) (2005) 467–482. URL https://doi.org/10.1016/j.ymssp.2004.02.007 [10] S. Baudin, D. R´emond, J. Antoni, O. Sauvage, Non-intrusive rattle noise detection in non-stationary conditions by an angle/time cyclostationary approach, Journal of Sound and Vibration 366 (2016) 501–513. URL https://doi.org/10.1016/j.jsv.2015.11.044 [11] S. Delvecchio, G. D’Elia, G. Dalpiaz, On the use of cyclostationary indicators in ic engine quality control by cold tests, Mechanical Systems & Signal Processing 60-61 (2015) 208–228. URL https://doi.org/10.1016/j.ymssp.2014.09.015 [12] J. Antoni, R. Randall, Differential diagnosis of gear and bearing faults, 32

Journal of Vibration and Acoustics 124 (2) (2002) 165–171. URL https://doi.org/10.1115/1.1456906 [13] J. Antoni, Cyclic spectral analysis in practice, Mechanical Systems & Signal Processing 21 (2) (2007) 597–630. URL https://doi.org/10.1016/j.ymssp.2006.08.007 [14] I. Antoniadis, G. Glossiotis, Cyclostationary analysis of rolling-element bearing vibration signals, Journal of sound and vibration 248 (5) (2001) 829–845. URL https://doi.org/10.1006/jsvi.2001.3815 [15] M. Poirier, M. Gagnon, A. Tahan, A. Coutu, J. Chamberland-Lauzon, Extrapolation of dynamic load behaviour on hydroelectric turbine blades with cyclostationary modelling, Mechanical Systems & Signal Processing 82 (2017) 193–205. URL https://doi.org/10.1016/j.ymssp.2016.05.018 [16] J. Antoni, F. Bonnardot, A. Raad, M. E. Badaoui, Cyclostationary modelling of rotating machine vibration signals, Mechanical Systems & Signal Processing 18 (6) (2004) 1285–1314. URL https://doi.org/10.1016/S0888-3270(03)00088-8 [17] A. Kahraman, Natural modes of planetary gear trains, Journal of Sound Vibration 173 (1994) 125–130. URL https://doi.org/10.1006/jsvi.1994.1222 [18] R. G. Parker, V. Agashe, S. M. Vijayakar, Dynamic response of a planetary gear system using a finite element/contact mechanics model, Jour33

nal of Mechanical Design 122 (3) (2000) 304–310. URL https://doi.org/10.1115/1.1286189 [19] Z. Shen, B. Qiao, L. Yang, W. Luo, X. Chen, Evaluating the influence of tooth surface wear on tvms of planetary gear set, Mechanism and Machine Theory 136 (2019) 206–223. URL https://doi.org/10.1016/j.mechmachtheory.2019.03.014 [20] X. Liang, M. J. Zuo, M. R. Hoseini, Vibration signal modeling of a planetary gear set for tooth crack detection, Engineering Failure Analysis 48 (2015) 185–200. [21] P. McFadden, J. Smith, An explanation for the asymmetry of the modulation sidebands about the tooth meshing frequency in epicyclic gear vibration, Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science 199 (1) (1985) 65–70. URL https://doi.org/10.1243/PIME PROC 1985 199 092 02 [22] Z. Feng, M. J. Zuo, Vibration signal models for fault diagnosis of planetary gearboxes, Journal of Sound & Vibration 331 (22) (2012) 4919– 4939. URL https://doi.org/10.1016/j.jsv.2012.05.039 [23] Z. Feng, H. Ma, M. J. Zuo, Vibration signal models for fault diagnosis of planet bearings, Journal of Sound & Vibration 370 (2016) 372–393. URL https://doi.org/10.1016/j.jsv.2016.01.041 [24] Y. Lei, Z. Liu, L. Jing, F. Lu, Phenomenological models of vibration signals for condition monitoring and fault diagnosis of epicyclic gearboxes, 34

Journal of Sound & Vibration 369 (2016) 266–281. URL https://doi.org/10.1016/j.jsv.2016.01.016 [25] R.-B. Sun, Z.-B. Yang, W. Luo, B.-J. Qiao, X.-F. Chen, Weighted sparse representation based on failure dynamics simulation for planetary gearbox fault diagnosis, Measurement Science and Technology 30 (4) (2019) 045008. URL https://doi.org/10.1088/1361-6501/ab02d8 [26] W. A. Gardner, C. M. Spooner, The cumulant theory of cyclostationary time-series. i. foundation, Signal Processing IEEE Transactions on 42 (12) (1994) 3387–3408. URL https://doi.org/10.1109/78.340775 [27] P. McFadden, Window functions for the calculation of the time domain averages of the vibration of the individual planet gears and sun gear in an epicyclic gearbox, Journal of Vibration and Acoustics 116 (2) (1994) 179–187. URL https://doi.org/10.1115/1.2930410 [28] D. Ho, R. Randall, Effects of time delay, order of fir filter and convergence factor on self-adaptive noise cancellation, in: International Conference on Sound and Vibration (ICSV5), Adelaide, 1997, pp. 11–18.

35

Highlights •

A cyclostationary phenomenological signal model of the planetary set is proposed.



Periodic components are removed by the self-adaptive noise cancellation method.



Modulation features are revealed by the second-order cyclostationary descriptor.



The model is verified by the numerical and experimental analysis.

Author Contribution Statement

Ruo-Bin Sun: Conceptualization, Methodology, Formal analysis, Writing - Original Draft, Writing Review & Editing Zhi-Bo Yang: Conceptualization, Funding acquisition Konstantinos Gryllias: Resources, Writing - Review & Editing, Supervision Xue-Feng Chen: Supervision, Project administration, Supervision

Declaration of interests

☑The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ☐The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: