# Fault Diagnosis of Intershaft Bearing Using Variational Mode Decomposition with TAGA Optimization

Authors: Jing Tian; Shu-Guang Wang; Jie Zhou; Yan-Ting Ai; Yu-Wei Zhang; Cheng-Wei Fei
Journal: Shock and Vibration
 (2021)

---

## Abstract

To efficiently extract the features of aeroengine intershaft bearing faults with weak signal, the variational mode decomposition (VMD) method based on the tolerant adaptive genetic algorithm (TAGA) (TAGA-VMD) is proposed by introducing the idea of tolerance into the traditional adaptive genetic algorithm in this paper. In this method, the tolerant genetic algorithm was adopted to find the optimum empirical parametersK and α of VMD. A fault simulation experiment system of intershaft bearings was built for the inner ring fault and outer ring fault of bearings to verify the proposed TAGA-VMD method. The results show that the proposed method can effectively extract the fault feature frequency of intershaft bearings, and the error between the extracted fault feature frequency and the theoretical value of fault frequency is less than 0.1%. The efforts of this study provide one promising fault feature extraction approach for aeroengine intershaft bearing fault diagnosis with weak signal.


---

## Body

## 1. Introduction

Rotor is a key system of aeroengine and directly affects the operating conditions of aeroengine. In rotor systems, rolling bearing works in a severe environment and is prone to failure. Once rolling bearing fault occurs, it directly affects the working state and operation efficiency of a rotor system. In the online monitoring of aeroengine, sensors are mainly installed on the outer surface of casings to monitor the health status of rolling bearing, so that too long signal propagation path induces vibration signals attenuation, low signal to noise ratio (SNR) of fault signals, and difficult extraction of fault features [1]. Therefore, it is urgent to study the weak fault diagnosis of intershaft bearings [2–4].For the signal process of rolling bearing fault diagnosis, numerous signal processing methods, such as deconvolution [5, 6] and wavelet transform [7, 8], emerged for reducing noise and signal filtering [9, 10]. Besides, a new performance degradation prediction (HMEPEM) method was proposed to track the evolution of degradation and predict the performance degradation trend of rolling bearings [11]. Stochastic resonance (SR) [12] can enhance the SNR of weak signal and was widely used in the field of weak signal detection [13, 14]. The fault signal of intershaft bearing is weak and attenuates seriously. If the SR method is adopted to process the fault signal, the feature of fault signals can be effectively extracted further. In the existing fault diagnosis methods of rotating machinery, various signal decomposition algorithms are widely applied to extract the key fault information and filter background noise by processing the original fault signals. As a nonstationary and nonlinear signal processing method, empirical mode decomposition (EMD) [15] was widely used in bearing fault diagnosis [16, 17], because it could adaptively decompose fault signals into different modal components. Ensemble EMD (EEMD) was proposed [18] to solve the problem of mode aliasing and noise sensitivity of EMD in signal decomposition. Then, more improvements have been made based on EMD and EEMD [19, 20]. Wavelet transform was utilized to remove the irrelevant components of fault information, and it is combined with EEMD to improve the operation efficiency of EEMD [21].Variational mode decomposition (VMD) is an adaptive modal decomposition method based on variation ideas [22]. The frequency center and bandwidth of each modal component can be determined by searching the optimal solution of the variational mode, which is promising to effectively overcome the phenomenon of modal aliasing. Therefore, the VMD has a significant effect in processing fault signals of rotating machinery [23, 24]. Mohanty et al. [25] used the VMD and EMD to analyze the vibration and acoustic characteristics of bearings under variable speed conditions by developing a new fault identification method of decomposed signals using correlation coefficient (CC) and Hurst exponent. Xu et al. [26] enhanced the weak impulse component of the original signal by the Teager energy operator and processed the fault signal with the VMD. However, the decomposition effect of VMD on fault signal depends on the selected empirical parameters, especially the number of modal components K and penalty factor α. At present, many researchers have proposed some intelligent optimization methods to solve parameter optimization problems [27, 28]. Wu et al. [29] proposed an improved quantum-inspired differential evolution (MSIQDE) algorithm which can avoid premature convergence and improve the global search ability to optimize the parameters of deep belief network (DBN) and constructed an optimal DBN model and use this model to classify bearing faults. The optimization algorithm is also applied to search for the optimal parameters of VMD [30, 31]. Zhang et al. [32] proposed an adaptive parameter VMD method based on a grasshopper optimization algorithm, which took kurtosis and correlation coefficient as optimization objective function and introduced a weighted kurtosis index, to acquire the optimum parameters of VMD. Cancan et al. used particle swarm optimization (PSO) to optimize the number of modal components K and penalty factor α of VMD, by constructing the correlation coefficient between modal components and the original signal to define fitness function as the evaluation index. However, it has not been clear for the relationship between the correlation coefficient and optimization parameters yet up to now [33].To find the adaptive parameters more precisely, the improved genetic algorithm is utilized to synchronously optimize the number of modal componentsK and penalty factor α of VMD. Firstly, a tolerance genetic algorithm is developed by introducing tolerance into the traditional genetic algorithm to improve search ability and convergence speed. Then, the parameter K and α of VMD are optimized synchronously by the tolerant adaptive genetic algorithm (TAGA). Next, the parameter-adaptive VMD based on TAGA is designed. Finally, the tolerance genetic algorithm-based parameter adaptive VMD is used to process simulation signals and experimental signals, respectively, to extract the fault characteristics of intershaft bearing and realize the monitoring and diagnosis of intershaft bearing conditions. The experimental results show that TAGA can greatly improve the global search ability of genetic algorithm and avoid being trapped in local optimum and TAGA-VMD can extract the feature frequency of intershaft bearing faults more efficiently than general VMD.The outline of this paper is structured as follows. Section2 introduces the basic principle of VMD. The validity of VMD in signal decomposition is discussed in Section 3. In Section 4, the parameter adaptive VMD based on the tolerance genetic algorithm is established. Section 5 adopts the proposed method to process the experimental signal. The conclusions of this study are given in Section 6.


## 2. Variational Mode Decomposition Method

Variational mode decomposition (VMD) is an adaptive mode decomposition method based on variation ideas [22]. Unlike recursive mode decomposition, VMD decomposes complex nonlinear signals into intrinsic mode function (IMF) components with different central frequencies and limited bandwidths by iteratively searching for the optimal solution of the variational mode. As dealing with the frequency domain information of complex signals, the VMD has the advantages of high decomposition accuracy, strong adaptability, and good robustness.


### 2.1. Variational Mode Decomposition Theory

The IMF is defined as a signal with adjustable frequency and amplitude in VMD algorithm. Thekth IMF μk (t) of an original signal can be denoted by(1)μkt=Aktcosϕkt,in which t indicates time; Ak (t) is the instantaneous amplitude of the IMF μk (t); and ϕk (t) is the phase function of time t.Regarding equation (1), the instantaneous frequency ωk (t) of the IMF μk (t) is derived as(2)ωkt=ϕk′t=dϕktdt.AssumingK IMFs in VMD, the original signal may be decomposed into K modal components, in respect of the below-detailed steps:Step 1.
in terms of the Hilbert transform method [22], the analyzed signal and its unilateral spectrum of kth IMF are obtained by(3)δt+jπtμkt,in which j ∈ Z indicates integer multiples of fundamental frequency and δ (t) is the unit pulse function or signal.Step 2.
to estimate the central frequency ofkth IMF, an exponential term is introduced into equation (3) to adjust the frequency spectrum of kth IMF to the corresponding baseband. Equation (3) can be rewritten as(4)δt+jπtμkte−jωkt.Step 3.
compute the square norm of the gradient of the demodulated signal and estimate the bandwidth of each IMF. The constrained variational model is expressed as(5)minμk,ωk,∑k=1K∂tδt+jπtμkte−jωkt22,s.t.,∑k=1Kμk=f,where ωk (k = 1, 2, …, K) is the kth central frequency corresponding to μk and f is the original signal.
To obtain the optimal solution of equation (5), penalty factor α and Lagrange multiplier operator λ (t) are introduced, in which α guarantees the accuracy of signal reconstruction and λ (t) maintains the efficiency of constraint. Therefore, a new solution expression of Lagrange function L (μk, ωk, λ) to the variational issue can be obtained as follows:(6)Lμk,ωk,λ=α∑k∂tδt+j/πtμkte−jωkt22+ft−∑kμkt22+λt,ft−∑kμkt.
To solve the above variational problem, the VMD uses alternate direction method of multipliers [22] to find the saddle point in equation (6) by (n + 1)th iteratively updating μkn+1, λkn+1, and ωkn+1 which are converted to frequency domain through Parseval/Plancherel Fourier equidistant transformation, the equation of the modal components is(7)μ^kn+1ω=f^ω−∑i≠kμ^iω+λ^ω/21+2αω−ωk2.
Similarly, the expression of center frequency in frequency domain can be obtained as follows:(8)ωkn+1=∫0∞ωμ^kω2dω∫0∞μ^kω2dω,where μ^kn+1ω is the Wiener filter of the current residual f^ω−∑i≠kμ^iω and ωkn+1 is the gravity center of the kth IMF power spectrum.


### 2.2. Verification of Variational Mode Decomposition with Simulated Signals

To verify the validity of VMD for signal decomposition, the signal is structured as(9)ft=f1t+f2t+f3t,where f1t, f2t, and f3t represent cosine functions with different frequencies. To simulate the real experimental conditions, Gaussian white noise Г(t) with noise intensity D = 15 and sampling frequency fs = 12 000 Hz is overlapped to the original signal. Then, the original input signal with the Gaussian white noise is expressed as(10)xt=ft+Γt.The time-domain waveforms of the input signalf (t), the components f1 (t), f2 (t), and f3 (t), and the signal x (t) with Gaussian white noise are shown in Figure 1. Their spectrums are shown in Figure 2.Figure 1
Time-domain waveform off1 (t), f2 (t), and f3 (t) as well as x (t) and f (t).Figure 2
The frequency-domain spectrums off1 (t), f2 (t), and f3 (t) as well as x (t) and f (t).As the modal parameterK = 4 and penalty factor α = 2000, the VMD is applied to process the simulation signals. The time-domain waveforms and spectrums of the decomposed output signal are displayed in Figures 3 and 4.Figure 3
Time-domain waveforms of IMF decomposed with the VMD.Figure 4
IMF spectrograms decomposed with the VMD.Comparing Figure2 with Figure 4, it can be seen that the VMD can effectively extract the central frequency of modal components, and the noise energy is mainly concentrated on the third IMF signal (IMF3). It is demonstrated that the VMD has a good effect in extracting feature signals under a noise background.


## 2.1. Variational Mode Decomposition Theory

The IMF is defined as a signal with adjustable frequency and amplitude in VMD algorithm. Thekth IMF μk (t) of an original signal can be denoted by(1)μkt=Aktcosϕkt,in which t indicates time; Ak (t) is the instantaneous amplitude of the IMF μk (t); and ϕk (t) is the phase function of time t.Regarding equation (1), the instantaneous frequency ωk (t) of the IMF μk (t) is derived as(2)ωkt=ϕk′t=dϕktdt.AssumingK IMFs in VMD, the original signal may be decomposed into K modal components, in respect of the below-detailed steps:Step 1.
in terms of the Hilbert transform method [22], the analyzed signal and its unilateral spectrum of kth IMF are obtained by(3)δt+jπtμkt,in which j ∈ Z indicates integer multiples of fundamental frequency and δ (t) is the unit pulse function or signal.Step 2.
to estimate the central frequency ofkth IMF, an exponential term is introduced into equation (3) to adjust the frequency spectrum of kth IMF to the corresponding baseband. Equation (3) can be rewritten as(4)δt+jπtμkte−jωkt.Step 3.
compute the square norm of the gradient of the demodulated signal and estimate the bandwidth of each IMF. The constrained variational model is expressed as(5)minμk,ωk,∑k=1K∂tδt+jπtμkte−jωkt22,s.t.,∑k=1Kμk=f,where ωk (k = 1, 2, …, K) is the kth central frequency corresponding to μk and f is the original signal.
To obtain the optimal solution of equation (5), penalty factor α and Lagrange multiplier operator λ (t) are introduced, in which α guarantees the accuracy of signal reconstruction and λ (t) maintains the efficiency of constraint. Therefore, a new solution expression of Lagrange function L (μk, ωk, λ) to the variational issue can be obtained as follows:(6)Lμk,ωk,λ=α∑k∂tδt+j/πtμkte−jωkt22+ft−∑kμkt22+λt,ft−∑kμkt.
To solve the above variational problem, the VMD uses alternate direction method of multipliers [22] to find the saddle point in equation (6) by (n + 1)th iteratively updating μkn+1, λkn+1, and ωkn+1 which are converted to frequency domain through Parseval/Plancherel Fourier equidistant transformation, the equation of the modal components is(7)μ^kn+1ω=f^ω−∑i≠kμ^iω+λ^ω/21+2αω−ωk2.
Similarly, the expression of center frequency in frequency domain can be obtained as follows:(8)ωkn+1=∫0∞ωμ^kω2dω∫0∞μ^kω2dω,where μ^kn+1ω is the Wiener filter of the current residual f^ω−∑i≠kμ^iω and ωkn+1 is the gravity center of the kth IMF power spectrum.


## 2.2. Verification of Variational Mode Decomposition with Simulated Signals

To verify the validity of VMD for signal decomposition, the signal is structured as(9)ft=f1t+f2t+f3t,where f1t, f2t, and f3t represent cosine functions with different frequencies. To simulate the real experimental conditions, Gaussian white noise Г(t) with noise intensity D = 15 and sampling frequency fs = 12 000 Hz is overlapped to the original signal. Then, the original input signal with the Gaussian white noise is expressed as(10)xt=ft+Γt.The time-domain waveforms of the input signalf (t), the components f1 (t), f2 (t), and f3 (t), and the signal x (t) with Gaussian white noise are shown in Figure 1. Their spectrums are shown in Figure 2.Figure 1
Time-domain waveform off1 (t), f2 (t), and f3 (t) as well as x (t) and f (t).Figure 2
The frequency-domain spectrums off1 (t), f2 (t), and f3 (t) as well as x (t) and f (t).As the modal parameterK = 4 and penalty factor α = 2000, the VMD is applied to process the simulation signals. The time-domain waveforms and spectrums of the decomposed output signal are displayed in Figures 3 and 4.Figure 3
Time-domain waveforms of IMF decomposed with the VMD.Figure 4
IMF spectrograms decomposed with the VMD.Comparing Figure2 with Figure 4, it can be seen that the VMD can effectively extract the central frequency of modal components, and the noise energy is mainly concentrated on the third IMF signal (IMF3). It is demonstrated that the VMD has a good effect in extracting feature signals under a noise background.


## 3. Variational Mode Decomposition Based on Tolerance Adaptive Genetic Algorithms

In the VMD algorithm, the parameters of VMD seriously influence the decomposition. The most influential parameters are the number of modal decompositionK and penalty factor α. Small K will lead to the difficulty of effectively separating the central frequency and the occurrence of mode aliasing phenomenon, while large K will cause the over-decomposition of the signal which reduces the energy of the characteristic signal in the effective IMF component and decrease the amplitude of characteristic frequency. The penalty factor α mainly affects the bandwidth of modal function. Larger penalty factor α reduces the sum of the bandwidth of modal functions and prolongs the running time of the VMD algorithm to cut down the applicability of the VMD algorithm. To find the reasonable K and α, an improved genetic algorithm-TAGA is proposed to search for the optimum parameters of VMD.Genetic algorithms can encode parameters without experience and search globally through multiple routes with strong flexibility. In line with the adaptive genetic algorithm [34, 35], the idea of adaptive tolerance is introduced to the proposed TAGA. The TAGA is to add or subtract the standard deviation of fitness fluctuation from adaptive tolerance fitness limit (TFL), to adjust the points (mutation probability Pm and crossover probability Pc) of piecewise function to improve the performance of genetic algorithm iteration, respectively.The flowchart of optimizing parametersK and α in the VMD by TAGA is shown in Figure 5.Figure 5
The flowchart of parameter optimization in VMD by using TAGA.


### 3.1. Initialization Coding

The classical binary coding method [36] is used to initialize the code of parameter K and α optimization. First, the thresholds of parameter K and α are determined and coded as binary digits according to the precision of optimization parameters. Generally, the number of modal decomposition K does not exceed 7 and α is in the range (1000 and 3000) in respect of experience, to ensure the reconstruction accuracy and validity of the VMD algorithm. In this study, the accuracy of K and α is set as ε = 0.01.


### 3.2. Construction of Tolerant Genetic Operators

Genetic operator is the process of population evolution in genetic algorithm, including selection, crossover, and mutation.It is assumed that the iteration generations of genetics areN and the number of individuals per generation is M. In the first iteration, M individuals are randomly selected from the feasible region composed of K and α for selective inheritance. In the first N/2 generation, the strategy of retaining dominant individuals is directly adopted to obtain a faster convergence rate and improve the efficiency of the algorithm. The latter N/2 generations adopt the traditional roulette strategy [37]. Based on equation (11), the section probability of the individual is obtained to ensure the randomness of the genetic system and prevent falling into the local optimum solution.(11)pi=fi∑i=1Mfi,where pi represents the section probability of ith individual and fi is the fitness value of ith individual.To greatly improve the global search ability of the genetic algorithm and avoid being trapped in local optimum, a tolerance adaptive crossover operator is utilized as the crossover probability is denoted byPc. Then, the individuals of the selected K and α are truncated and crossed by binary bits, to generate new individuals. The crossover probability Pc is(12)Pc=Pcu−Pcu−Pclfc−fmeanfmax−fmean,fc≥TFLc,Pcu,fc<TFLc,where Pcu and Pcl is the upper limit and lower limit of crossover probability; fc indicates the larger fitness value of two crossover individuals; fck is the larger fitness value of two crossover individuals for kth IMF; fmax denotes the maximum fitness value of the current generation; fmean is the average fitness value of the current generation; TFLc is the tolerance fitness limit.(13)TFLc=fmean−1n−1∑k=1nfck−fmean2.When the mutation probability isPm, the binary gene coding of randomly selected population K and α is mutated to maintain individual diversity. To get the optimal approximate solution, a tolerant adaptive mutation operator [38] is proposed to make the system easier jump out of the local optimal solution at the end of iteration to increase the search depth of the algorithm. When Pmu is the upper limit of mutation probability, and Pml is the lower limit of mutation probability, and fm is the fitness value of mutation individual; and fmk is the fitness value of mutation individual for kth IMF and TFLm is the tolerance fitness limit, the mutation probability Pm is(14)Pm=Pmu−Pmu−Pmlfmean−fmfmax−fmean,fm≥TFLm,Pmu,fm<TFLm,(15)TFLm=fmean+1n−1∑k=1nfmk−fmean2.


### 3.3. Evaluation of Fitness Function

The core of the genetic algorithm lies in the establishment of the fitness evaluation system, which is the goal and direction of algorithm optimization. The value of individual fitness reflects the degree of individual fitness. Therefore, the correctness of fitness function design directly affects the performance of the algorithm. The kurtosis is dimensionless parameters and is sensitive to an impulse signal. When the bearing fails, the impulse components of the collected signal become more obvious and the correlated kurtosis becomes higher [39]. The kurtosis of kth IMF component of the fault signal decomposed by VMD is gained by(16)krk=1n∑t=1nxkt−μkσk4,where xk (t) is the kth decomposed IMF of the VMD; μk is the average value of xk (t); σk is the standard deviation of xk (t); and krk is the kurtosis of kth IMF component of the fault signal.The kurtosis of IMFs can be used to construct the fitness functionFn, i.e.,(17)Fn=−MAXkrixit,i=1,2,…,K.The value ofFn depends on the selection of K and α, i.e., the smaller the function value is, the greater obvious the fault components extracted by the VMD from the original signal are.


### 3.4. Termination Conditions

With the continuous optimization of individuals, the evolution terminates an output of the optimal VMD parameters when the deviation between average fitness and the fitness of the optimal individualsK and α satisfies δ < 10−6 or the number of iterations reaches 100 generations.As shown in Figure5, the parameters K and α of VMD are optimized by TAGA. First of all, K and α are coded by the classical binary coding method in the set interval to realize population initialization. In the line with equations (16) and (17), fitness function and calculate individual fitness value of K and α and judge whether the fitness value meets termination conditions. When termination conditions are not satisfied, the genetic operation is carried out on the individuals of K and α with equations (11)–(15) to get new generation individuals. And the process continues until the optimal individuals K and α are output when termination conditions are satisfied. So the optimal parameter K and α of VMD can be obtained by TAGA.


## 3.1. Initialization Coding

The classical binary coding method [36] is used to initialize the code of parameter K and α optimization. First, the thresholds of parameter K and α are determined and coded as binary digits according to the precision of optimization parameters. Generally, the number of modal decomposition K does not exceed 7 and α is in the range (1000 and 3000) in respect of experience, to ensure the reconstruction accuracy and validity of the VMD algorithm. In this study, the accuracy of K and α is set as ε = 0.01.


## 3.2. Construction of Tolerant Genetic Operators

Genetic operator is the process of population evolution in genetic algorithm, including selection, crossover, and mutation.It is assumed that the iteration generations of genetics areN and the number of individuals per generation is M. In the first iteration, M individuals are randomly selected from the feasible region composed of K and α for selective inheritance. In the first N/2 generation, the strategy of retaining dominant individuals is directly adopted to obtain a faster convergence rate and improve the efficiency of the algorithm. The latter N/2 generations adopt the traditional roulette strategy [37]. Based on equation (11), the section probability of the individual is obtained to ensure the randomness of the genetic system and prevent falling into the local optimum solution.(11)pi=fi∑i=1Mfi,where pi represents the section probability of ith individual and fi is the fitness value of ith individual.To greatly improve the global search ability of the genetic algorithm and avoid being trapped in local optimum, a tolerance adaptive crossover operator is utilized as the crossover probability is denoted byPc. Then, the individuals of the selected K and α are truncated and crossed by binary bits, to generate new individuals. The crossover probability Pc is(12)Pc=Pcu−Pcu−Pclfc−fmeanfmax−fmean,fc≥TFLc,Pcu,fc<TFLc,where Pcu and Pcl is the upper limit and lower limit of crossover probability; fc indicates the larger fitness value of two crossover individuals; fck is the larger fitness value of two crossover individuals for kth IMF; fmax denotes the maximum fitness value of the current generation; fmean is the average fitness value of the current generation; TFLc is the tolerance fitness limit.(13)TFLc=fmean−1n−1∑k=1nfck−fmean2.When the mutation probability isPm, the binary gene coding of randomly selected population K and α is mutated to maintain individual diversity. To get the optimal approximate solution, a tolerant adaptive mutation operator [38] is proposed to make the system easier jump out of the local optimal solution at the end of iteration to increase the search depth of the algorithm. When Pmu is the upper limit of mutation probability, and Pml is the lower limit of mutation probability, and fm is the fitness value of mutation individual; and fmk is the fitness value of mutation individual for kth IMF and TFLm is the tolerance fitness limit, the mutation probability Pm is(14)Pm=Pmu−Pmu−Pmlfmean−fmfmax−fmean,fm≥TFLm,Pmu,fm<TFLm,(15)TFLm=fmean+1n−1∑k=1nfmk−fmean2.


## 3.3. Evaluation of Fitness Function

The core of the genetic algorithm lies in the establishment of the fitness evaluation system, which is the goal and direction of algorithm optimization. The value of individual fitness reflects the degree of individual fitness. Therefore, the correctness of fitness function design directly affects the performance of the algorithm. The kurtosis is dimensionless parameters and is sensitive to an impulse signal. When the bearing fails, the impulse components of the collected signal become more obvious and the correlated kurtosis becomes higher [39]. The kurtosis of kth IMF component of the fault signal decomposed by VMD is gained by(16)krk=1n∑t=1nxkt−μkσk4,where xk (t) is the kth decomposed IMF of the VMD; μk is the average value of xk (t); σk is the standard deviation of xk (t); and krk is the kurtosis of kth IMF component of the fault signal.The kurtosis of IMFs can be used to construct the fitness functionFn, i.e.,(17)Fn=−MAXkrixit,i=1,2,…,K.The value ofFn depends on the selection of K and α, i.e., the smaller the function value is, the greater obvious the fault components extracted by the VMD from the original signal are.


## 3.4. Termination Conditions

With the continuous optimization of individuals, the evolution terminates an output of the optimal VMD parameters when the deviation between average fitness and the fitness of the optimal individualsK and α satisfies δ < 10−6 or the number of iterations reaches 100 generations.As shown in Figure5, the parameters K and α of VMD are optimized by TAGA. First of all, K and α are coded by the classical binary coding method in the set interval to realize population initialization. In the line with equations (16) and (17), fitness function and calculate individual fitness value of K and α and judge whether the fitness value meets termination conditions. When termination conditions are not satisfied, the genetic operation is carried out on the individuals of K and α with equations (11)–(15) to get new generation individuals. And the process continues until the optimal individuals K and α are output when termination conditions are satisfied. So the optimal parameter K and α of VMD can be obtained by TAGA.


## 4. Experimental Verification

### 4.1. Test Rig of Intershaft Bearing

A birotor test rig as shown in Figure6 was adopted to perform the simulation experiments of intershaft bearing faults, from which the experimental data were acquired to verify the practicability and effectiveness of the VMD method based on TAGA (TAGA-VMD method) in signal processing. The experimental system of the birotor test rig is composed of a motor drive system, rotor system, support system, and data collection system. This paper adopted a cylindrical roller bearing (mode no. 202EM NSK, Japan). The surface defect of the outer ring and inner ring in the bearing is made by wire cutting, and their width and depth are 0.5 mm as shown in Figure 7. The defect runs through the whole inner ring and outer ring longitudinally. Five integrated circuits piezoelectric (ICP) accelerometers (mode no. 33B30, PCB Company, USA) were used to measure the vibration signals of intershaft bearing faults.Figure 6
Birotor test rig of intershaft bearing.Figure 7
Intershaft bearing with faults: (a) outer ring fault and (b) inner ring fault.
(a)(b)During the experiment, the high-pressure rotor reversely rotates against the low-pressure rotor. When the outer ring fails, the speeds of the high-pressure rotor and low-pressure rotor are 1000r/min and 600 r/min, respectively. When the inner ring fails, the speeds of high-pressure rotor and low-pressure rotor are 1200 r/min and 300 r/min, respectively. The sampling frequency and time are 16384 Hz and the sampling time is 20 s. The parameters of the intershaft bearing are shown in Table 1.Table 1
Parameters of intershaft bearing.
Number of roll rodsInner diameter (mm)Outside diameter (mm)Diameter of the roll rod (mm)1125355In light of bearing failure frequency as shown in equations (18) and (19), the fault frequency of intershaft bearing with outer ring fault and inner ring fault is calculated by equations (18) and (19).Outer ring:(18)f0=n0±ni2×601−dDmcosθZ.Inner ring:(19)fi=n0±ni2×601+dDmcosθZ,where fi and fo indicate the fault frequencies of the inner ring and outer ring; ni and no are the rotational speeds of the inner ring and outer ring, respectively; Z denotes the number of bearings rollers; d stands for the diameter of rollers; Dm is the diameter of bearings; θ explains the pressure angle between rollers and outer ring in the bearings.


### 4.2. Verification of Outer Ring Fault

To validate the adaptive parameter VMD based on TAGA, the vibration signals of outer ring fault with 1000r/min and 600 r/min were collected. In line with equation (18), the fault frequency of outer ring faults is 117.3 Hz. In addition, the frequencies of the high-pressure rotor and low-pressure rotor are 16.6 Hz and 10 Hz, respectively.The fault signal of outer ring collected by the sensor is shown in Figure8. In Figure 8, impulse components are found in the time-domain signal, which are induced by the defect of the outer ring. However, the periodicity of the impulse waveform is not obvious, because the signal contains vibration and noise from other parts of the experimental system. The envelope spectrum of the outer ring fault signal is displayed in Figure 9. As seen from Figure 9, there is a lot of noise in the low-frequency region and the fault signal is relatively weak, so that the characteristic frequency of the bearing fault cannot be extracted effectively.Figure 8
Fault signal of outer ring of intershaft bearing.Figure 9
Envelope spectrum of intershaft bearing outer ring fault signal.To effectively extract the characteristic frequency of bearing fault, the VMD is applied to decompose the collected fault signals. Before using TAGA optimization, we randomly selectK = 3 and α = 2000. After decomposition, the IMF component with the highest kurtosis is selected for envelope analysis. The envelope spectrum is shown in Figure 10.Figure 10
IMF signal envelope spectrum of outer ring.As shown in Figure10, a large amount of noises in the signal are filtered, and the fluctuation of signal peak tends to be stable, and the periodicity of impulse component in the fault signal is more obvious after VMD decomposition. Moreover, both the rotating frequency fr (16.7 Hz) and the fault frequency fo (117.3 Hz) of the outer ring are clear, including two frequency multiplication. This is because the frequency domain information of the vibration signal can be effectively analyzed through the envelope demodulation method, to find the fault characteristic frequency and modulation frequencies in the envelope spectrum.However, the display of other modulation frequencies is incomplete so that no enough frequency characteristics are obtained from the envelope spectrum. For this issue, the parameters of VMD are optimized by the TAGA. We gainedK = 5 and α = 1806. According to the optimized results, similarly, the parameters are set and the IMF component with the highest kurtosis after decomposition is selected for envelope analysis. The envelope spectrum is shown in Figure 11.Figure 11
IMF signal envelope spectrum of outer ring with the VMD by TAGA optimization.As illustrated in Figure11, the frequency components of envelope spectrum obtained by VMD with TAGA optimization are clearer, which includes not only outer loop rotational frequency fr (16.7 Hz) and its triple frequency multiplication but also outer ring fault frequency fo (117.3 Hz) and its quadruplicated frequency multiplication. Besides, it is found that various modulation frequencies (f0 ± fr, f0 ± 2fr, 2f0 ± fr, and 2f0 ± 2fr) center on fault characteristic frequencies and sidebands of rotational frequency. Moreover, the amplitude of the characteristic frequency is obviously higher than that of K = 3. Therefore, the outer ring fault of bearing is clearer. The reason is that when the outer ring of rolling bearing has partial defect fault, the position of the ball passing through the defect always varies due to the working load on the bearing, and the load density at the defect changes periodically with the rotation of the outer ring. In this case, the periodic load causes the amplitude modulation of the fault signal, in which the modulation frequency is equal to rotational frequency of the outer ring.


### 4.3. Verification of Inner Ring Fault

To verify the adaptability of the proposed method against other rotational speeds and fault states, the vibration signals of inner ring fault withno = 1200 r/min and ni = 300 r/min are collected. With equation (18), the characteristic frequency of inner ring faults was determined to be 165 Hz, and the frequencies of the high-pressure rotor and low-pressure rotor were 20 Hz and 5 Hz, respectively. The fault signal of the inner ring is shown in Figure 12. As displayed in Figure 12, impulse components and other interferential information exist in the time-domain signal. Moreover, the peaks of impulse components are unstable and accompanied with high-frequency noise. Like the outer ring fault, most of the energy in the frequency spectrum is concentrated in the high-frequency region.Figure 12
Fault signal of inner ring of intershaft bearing.
(a)(b)By using the VMD, the fault signal of inner ring was decomposed asK = 4 to fin the IMF component with the highest kurtosis for envelope analysis. The envelope spectrum is shown in Figure 13:Figure 13
IMF signal envelope spectrum of inner ring.As illustrated in Figure13, the periodicity of the impulse component in the fault signal is more obvious after VMD. Additionally, the rotating frequency fr and fault frequency fi of the inner ring are 5 Hz and 165 Hz, respectively, including two frequency multiplication. Therefore, fi ± fr is 170 Hz. However, the frequency multiplication of fault bearing and other modulation frequency information are incomplete yet. To gain more frequency characteristics from the envelope spectrum, the parameters of VMD are optimized by TAGA. The optimization result is K = 6 and α = 1865. The IMF component with the highest kurtosis after decomposition is selected for envelope analysis. The envelope spectrum is shown in Figure 14:Figure 14
IMF signal envelope spectrum of inner ring using the VMD with TAGA optimization.As illustrated in Figure14, (1) the frequency components of envelope spectrum of fault signal processed by adaptive parameter VMD based on TAGA (TAGA-VMD) method are clear, which include inner ring rotational frequency fr (5 Hz) and its double frequency multiplication, and inner ring fault frequency fi (165 Hz) and its double frequency multiplication as well. (2) Various modulation frequencies (fi ± fr, fi ± 2fr, fi ± 2fr) center on fault characteristic frequencies and sidebands of rotational frequency and the characteristic of frequency components is more obvious than that of K = 4. Therefore, the inner ring fault of bearing is clearer. (3) When the inner ring of the rolling bearing has partial defect fault, the position of the ball passing through the defect always varies due to the working load on the bearing, and the load density on the defect changes periodically with the rotation of the inner ring. The periodic load causes the amplitude modulation of fault signal, in which the modulation frequency is equal to rotational frequency of inner ring. The inner ring fault of bearing can be clearly seen by the envelope spectrum. Therefore, fault frequency can be more effectively detected in the fault diagnosis of intershaft bearing.


## 4.1. Test Rig of Intershaft Bearing

A birotor test rig as shown in Figure6 was adopted to perform the simulation experiments of intershaft bearing faults, from which the experimental data were acquired to verify the practicability and effectiveness of the VMD method based on TAGA (TAGA-VMD method) in signal processing. The experimental system of the birotor test rig is composed of a motor drive system, rotor system, support system, and data collection system. This paper adopted a cylindrical roller bearing (mode no. 202EM NSK, Japan). The surface defect of the outer ring and inner ring in the bearing is made by wire cutting, and their width and depth are 0.5 mm as shown in Figure 7. The defect runs through the whole inner ring and outer ring longitudinally. Five integrated circuits piezoelectric (ICP) accelerometers (mode no. 33B30, PCB Company, USA) were used to measure the vibration signals of intershaft bearing faults.Figure 6
Birotor test rig of intershaft bearing.Figure 7
Intershaft bearing with faults: (a) outer ring fault and (b) inner ring fault.
(a)(b)During the experiment, the high-pressure rotor reversely rotates against the low-pressure rotor. When the outer ring fails, the speeds of the high-pressure rotor and low-pressure rotor are 1000r/min and 600 r/min, respectively. When the inner ring fails, the speeds of high-pressure rotor and low-pressure rotor are 1200 r/min and 300 r/min, respectively. The sampling frequency and time are 16384 Hz and the sampling time is 20 s. The parameters of the intershaft bearing are shown in Table 1.Table 1
Parameters of intershaft bearing.
Number of roll rodsInner diameter (mm)Outside diameter (mm)Diameter of the roll rod (mm)1125355In light of bearing failure frequency as shown in equations (18) and (19), the fault frequency of intershaft bearing with outer ring fault and inner ring fault is calculated by equations (18) and (19).Outer ring:(18)f0=n0±ni2×601−dDmcosθZ.Inner ring:(19)fi=n0±ni2×601+dDmcosθZ,where fi and fo indicate the fault frequencies of the inner ring and outer ring; ni and no are the rotational speeds of the inner ring and outer ring, respectively; Z denotes the number of bearings rollers; d stands for the diameter of rollers; Dm is the diameter of bearings; θ explains the pressure angle between rollers and outer ring in the bearings.


## 4.2. Verification of Outer Ring Fault

To validate the adaptive parameter VMD based on TAGA, the vibration signals of outer ring fault with 1000r/min and 600 r/min were collected. In line with equation (18), the fault frequency of outer ring faults is 117.3 Hz. In addition, the frequencies of the high-pressure rotor and low-pressure rotor are 16.6 Hz and 10 Hz, respectively.The fault signal of outer ring collected by the sensor is shown in Figure8. In Figure 8, impulse components are found in the time-domain signal, which are induced by the defect of the outer ring. However, the periodicity of the impulse waveform is not obvious, because the signal contains vibration and noise from other parts of the experimental system. The envelope spectrum of the outer ring fault signal is displayed in Figure 9. As seen from Figure 9, there is a lot of noise in the low-frequency region and the fault signal is relatively weak, so that the characteristic frequency of the bearing fault cannot be extracted effectively.Figure 8
Fault signal of outer ring of intershaft bearing.Figure 9
Envelope spectrum of intershaft bearing outer ring fault signal.To effectively extract the characteristic frequency of bearing fault, the VMD is applied to decompose the collected fault signals. Before using TAGA optimization, we randomly selectK = 3 and α = 2000. After decomposition, the IMF component with the highest kurtosis is selected for envelope analysis. The envelope spectrum is shown in Figure 10.Figure 10
IMF signal envelope spectrum of outer ring.As shown in Figure10, a large amount of noises in the signal are filtered, and the fluctuation of signal peak tends to be stable, and the periodicity of impulse component in the fault signal is more obvious after VMD decomposition. Moreover, both the rotating frequency fr (16.7 Hz) and the fault frequency fo (117.3 Hz) of the outer ring are clear, including two frequency multiplication. This is because the frequency domain information of the vibration signal can be effectively analyzed through the envelope demodulation method, to find the fault characteristic frequency and modulation frequencies in the envelope spectrum.However, the display of other modulation frequencies is incomplete so that no enough frequency characteristics are obtained from the envelope spectrum. For this issue, the parameters of VMD are optimized by the TAGA. We gainedK = 5 and α = 1806. According to the optimized results, similarly, the parameters are set and the IMF component with the highest kurtosis after decomposition is selected for envelope analysis. The envelope spectrum is shown in Figure 11.Figure 11
IMF signal envelope spectrum of outer ring with the VMD by TAGA optimization.As illustrated in Figure11, the frequency components of envelope spectrum obtained by VMD with TAGA optimization are clearer, which includes not only outer loop rotational frequency fr (16.7 Hz) and its triple frequency multiplication but also outer ring fault frequency fo (117.3 Hz) and its quadruplicated frequency multiplication. Besides, it is found that various modulation frequencies (f0 ± fr, f0 ± 2fr, 2f0 ± fr, and 2f0 ± 2fr) center on fault characteristic frequencies and sidebands of rotational frequency. Moreover, the amplitude of the characteristic frequency is obviously higher than that of K = 3. Therefore, the outer ring fault of bearing is clearer. The reason is that when the outer ring of rolling bearing has partial defect fault, the position of the ball passing through the defect always varies due to the working load on the bearing, and the load density at the defect changes periodically with the rotation of the outer ring. In this case, the periodic load causes the amplitude modulation of the fault signal, in which the modulation frequency is equal to rotational frequency of the outer ring.


## 4.3. Verification of Inner Ring Fault

To verify the adaptability of the proposed method against other rotational speeds and fault states, the vibration signals of inner ring fault withno = 1200 r/min and ni = 300 r/min are collected. With equation (18), the characteristic frequency of inner ring faults was determined to be 165 Hz, and the frequencies of the high-pressure rotor and low-pressure rotor were 20 Hz and 5 Hz, respectively. The fault signal of the inner ring is shown in Figure 12. As displayed in Figure 12, impulse components and other interferential information exist in the time-domain signal. Moreover, the peaks of impulse components are unstable and accompanied with high-frequency noise. Like the outer ring fault, most of the energy in the frequency spectrum is concentrated in the high-frequency region.Figure 12
Fault signal of inner ring of intershaft bearing.
(a)(b)By using the VMD, the fault signal of inner ring was decomposed asK = 4 to fin the IMF component with the highest kurtosis for envelope analysis. The envelope spectrum is shown in Figure 13:Figure 13
IMF signal envelope spectrum of inner ring.As illustrated in Figure13, the periodicity of the impulse component in the fault signal is more obvious after VMD. Additionally, the rotating frequency fr and fault frequency fi of the inner ring are 5 Hz and 165 Hz, respectively, including two frequency multiplication. Therefore, fi ± fr is 170 Hz. However, the frequency multiplication of fault bearing and other modulation frequency information are incomplete yet. To gain more frequency characteristics from the envelope spectrum, the parameters of VMD are optimized by TAGA. The optimization result is K = 6 and α = 1865. The IMF component with the highest kurtosis after decomposition is selected for envelope analysis. The envelope spectrum is shown in Figure 14:Figure 14
IMF signal envelope spectrum of inner ring using the VMD with TAGA optimization.As illustrated in Figure14, (1) the frequency components of envelope spectrum of fault signal processed by adaptive parameter VMD based on TAGA (TAGA-VMD) method are clear, which include inner ring rotational frequency fr (5 Hz) and its double frequency multiplication, and inner ring fault frequency fi (165 Hz) and its double frequency multiplication as well. (2) Various modulation frequencies (fi ± fr, fi ± 2fr, fi ± 2fr) center on fault characteristic frequencies and sidebands of rotational frequency and the characteristic of frequency components is more obvious than that of K = 4. Therefore, the inner ring fault of bearing is clearer. (3) When the inner ring of the rolling bearing has partial defect fault, the position of the ball passing through the defect always varies due to the working load on the bearing, and the load density on the defect changes periodically with the rotation of the inner ring. The periodic load causes the amplitude modulation of fault signal, in which the modulation frequency is equal to rotational frequency of inner ring. The inner ring fault of bearing can be clearly seen by the envelope spectrum. Therefore, fault frequency can be more effectively detected in the fault diagnosis of intershaft bearing.


## 5. Conclusions

The aim of this study is to present an adaptive parameters-variational mode decomposition (VMD) method, called the VMD method based on the tolerant adaptive genetic algorithm (TAGA) (TAGA-VMD), to accurately extract the features of aeroengine intershaft bearing faults with weak signal. The proposed TAGA-VMD method introduces the idea of tolerance into the traditional adaptive genetic algorithm to find the optimum empirical parametersK and α of VMD. Based on the simulation experiments on intershaft bearing test rig, the fault signals of inner ring fault and outer ring fault are collected and analyzed to verify the proposed TAGA-VMD method. In respect to this study, some conclusions are summarized as follows:(1)
The VMD is illustrated to have a good effect in extracting the fault feature of simulated signals under Gauss white noise background(2)
TAGA can adaptively optimize the empirical parameters of VMD and is promising to be applied to improve the feature extraction precision of aeroengine intershaft bearing faults with weak signal(3)
The proposed TAGA-VMD method can precisely extract the impulse components of fault signals by gaining more obvious characteristics of fault frequency(4)
The proposed TAGA-VMD method of intershaft bearing can effectively extract the feature frequency of intershaft bearing faults with the error of less than 0.1%In this paper, the proposed method is demonstrated to be accurate in identifying fault frequency for the fault diagnosis of intershaft bearing and other rotors. The efforts of this study provide one promising fault feature extraction approach for aeroengine intershaft bearing fault diagnosis with weak signal.



---

*Source: 8828317-2021-01-13.xml*