Rotor Fault Diagnosis Method Based on VMD Symmetrical Polar Image and Fuzzy Neural Network

Abstract

Rotor fault diagnosis has attracted much attention due to its difficulties such as non-stationarity of fault signals, difficulty in fault feature extraction and low diagnostic accuracy of small samples. In order to extract fault feature information of rotors more effectively and to improve fault diagnosis precision, this paper proposed a fault diagnosis method based on variational mode decomposition (VMD) symmetrical polar image and fuzzy neural network. Firstly, the original rotor vibration signal is decomposed by using the VMD method and the relevant parameter selection algorithm of the VMD method is also proposed. Secondly, the intrinsic mode functions (IMF), which are sensitive to the signal characteristics, are selected for signal reconstruction based on a comprehensive evaluation factor method. As well, the reconstructed signal is transformed into a two-dimensional snowflake image through using the symmetrical polar coordinate method. Finally, the image features are extracted by the gray level co-occurrence matrix to form the state feature vector, which is input into the fuzzy neural network to realize the rotor fault diagnosis. Through the analysis of measured signals, the experimental results show that the proposed method can reach a higher recognition rate of 98% and the k-cross-validation experiment is used to demonstrate the robustness of the fuzzy neural network, and the average recognition accuracy of this experiment is 99.2%. Compared with some similar methods, the proposed method still has the highest fault recognition precision 98.4%, and the smallest standard deviation 0.5477.

1. Introduction

Rotating machinery is widely used in many industries such as aerospace and railway transportation. Rotor is an important component of rotating machinery equipment. Due to the increase in equipment complexity and precision, and the continuity of work tasks, local faults will inevitably occur on the rotor. According to statistical results, rotor local faults account for more than 50% of the causes leading to the failure of rotating machinery [1]. If the rotor faults cannot be detected in time, the working performance of the whole equipment will be affected, and huge economic losses will be caused in serious cases. Therefore, rotor fault diagnosis has always been a hot spot in the field of fault diagnosis [2].

Generally, the extraction of rotor fault features is mainly based on the time domain, frequency domain or time-frequency domain of the signal. In recent years, image processing has been widely used to diagnose rotor faults. Method of amplitude standardization based on signal sampling length can be used to transform the rotor vibration signal into a two-dimensional gray image, in which the texture feature, the feature frequency, and the brightness of the gray image are used to establish the characteristic matrix of the image. Thereby different degrees of unbalanced faults can be recognized [3]. Jiang et al. [4] used the nonlinear chirp mode decomposition method to decompose the gas turbine rotor vibration signal into components, then the instantaneous frequency of each frequency component was obtained, and the fault diagnosis of the rotor was carried out by the convolutional neural network (CNN). Zhu et al. [5] proposed a rotor fault diagnosis approach that transformed multiple rotor vibration signals into symmetrized dot pattern (SDP) image, and identified the SDP graphical feature characteristic of different states using a CNN. Hui et al. [6] expressed the time-frequency image of the rotor vibration signal as a texture feature tensor for the pattern recognition of rotor fault states with the linear support higher-tensor machine (SHTM). The experimental results showed that the method of classifying time-frequency texture feature tensor can achieve ideal diagnosis precision. The core of the fault diagnosis method based on image processing technology is to generate a simple and intuitive image, which can take abundant characteristic information of the rotor fault signals. Among image generation methods, the symmetrical polar coordinate method transforms the signal into a polar coordinate and expresses it in the form of a two-dimensional snowflake image, which is directly obtained from the time-domain signal and without the need for frequency analysis. This method is simple to calculate and it is an ideal image generation method.

However, due to the diversity of working conditions and the complexity of the mechanical system, the rotor vibration signal shows obvious signs of non-linear and non-stationary [7]. In addition, due to the influence of environmental noise, acquisition devices and other factors, the collected rotor vibration signals contain various disturbing signals [8]. It is difficult to ensure the precision of image feature acquisition by directly obtaining the symmetric polar image of this signal. Therefore, in order to strengthen the fault characteristic components, it is necessary to eliminate the disturbing signals before generating snowflake images.

Variational modal decomposition (VMD) [9] method is an adaptive and nonrecursive signal decomposition method. VMD method decomposes the non-linear and non-stationary signal into a series of intrinsic modal functions (IMFs) with limited bandwidth by solving the optimal solution of the variational problem.

Table 1. The characteristics of common signal decomposition methods.

Method	Shortcomings
WT	➀ Non-adaptive; ➁ Lack of selection criteria for wavelet basis functions; ➂ Different wavelet basis analysis results are different.
EMD	➀ End effect; ➁ Mode mixing; ➂ Envelopes fitting; ➃ Illusive mode component; ➄ Poor anti-noise performance.
LMD	➀ End effect; ➁ Mode mixing; ➂ Illusive mode component ➃ Criteria for iteration termination conditions; ➄ Poor anti-noise performance; ➅ computing efficiency.
EEMD	➀ End effect; ➁ Noise in IMF component; ➂ Envelopes fitting; ➃ Illusive mode component; ➄ Poor anti-noise performance.
VMD	➀ Lack of criteria for the selection of key parameters; ➁ Illusive mode component.

lists the shortcomings of common signal decomposition methods. From Table 1, we can know that compared with empirical mode decomposition (EMD) [10], local mean decomposition (LMD) [11], wavelet transform (WT) [12] and ensemble empirical mode decomposition (EEMD), VMD method could effectively avoid the mode mixing problem in the decomposition process, suppress disturbing components, and ensure the precision of subsequent signal feature extraction. At present, this method is widely applied to lots of fields, such as rotating machinery fault diagnosis and signal denoising. In view of the fact that vibration signals of rolling bearings are much contaminated by noise in the early failure period, Wang et al. [13] proposed a denoising method based on VMD and singular value decomposition (SVD), and the experimental result shows that the proposed method can effectively retain the useful signals and denoise the bearing signals in extremely noisy backgrounds. Fu et al. [14] combined the advantages of VMD, wavelet threshold denoising (WTD), and singular spectrum analysis (SSA), proposed a measured load signal denoising method and experimental results of synthetic signals demonstrate that the presented method outperforms the conventional digital signal denoising methods. Meng et al. [15] applied VMD and soft wavelet-thresholding to denoise rotor vibration signals. Through noise reduction experiments of rotor vibration signals under different working conditions, the results show that the signal-to-noise ratio of the VMD method basically doubles the EMD one. Liu et al. [16] proposed a novel denoising method that combines VMD and detrended fluctuation analysis (DFA). Through simulated and real signals analysis, the results show that the superior performance of this proposed filtering over EMD-based denoising and discrete wavelet threshold filtering. Unfortunately, from Table 1 we can also see that the selection of key parameters and the elimination of disturbing components in the process of VMD decomposition have always been the main factor that affects the decomposition precision and the efficiency of feature extraction.

The pattern recognition method based on deep learning could automatically extract image features. But in order to obtain appropriate parameters, the deep network requires a large number of sample data to train the model. Insufficient, samples lead to overfitting of the model training, resulting in low generalization ability of the model and poor fault diagnosis precision. Even though, the present research on fault diagnosis methods of rotating machinery are still based on small samples. As a common intelligent information processing method, fuzzy neural network has strong self-learning and data direct processing ability, which can express the structure of the result clearly [17], making itself an effective method for small sample fault diagnosis. Si et al. [18] established a fault signal acquisition engine test stand and simulated six kinds of engine working conditions. WT method was used to extract signal features and compose state eigenvectors and the fuzzy neural network was used to realize the effective diagnosis of engine faults with small samples. Hu et al. [19] proposed a fault diagnosis method for mine hoist braking system based on fuzzy neural networks, the experimental results show that the fuzzy neural network can diagnose faults accurately. Zhang et al. [20] applied fuzzy neural network to diagnose ventilator fault, through the example comparison of the fuzzy neural network and the BP neural network, the results illustrate that the fault diagnosis method of the fuzzy neural network can recognize the faults rapidly, accurately and steadily, which provides an efficient way for the diagnosis.

Aiming at the non-stationary rotor fault vibration signal and it is hard to accurately extract image features by directly obtaining the symmetric polar image of the signal with environmental noise. In this paper, a novel method for fault feature extraction and pattern recognition is proposed. The signal feature information is extracted by VMD symmetrical polar image, and fault recognition is achieved using fuzzy neural network. The major contributions of this paper compared to the existing approaches are that:

(1)

Combining the advantages of the VMD method in non-stationary signal processing with the effectiveness of symmetric polar image in feature information representation, it provides a new idea for introducing image processing technology into the field of rotating machinery fault diagnosis.

(2)

Aiming at the problem of selecting key parameters in the VMD decomposition process, a parameter selection method based on the maximum central frequency criterion and the minimum multi-scale fuzzy entropy criterion is proposed. The proposed method is simple, efficient and accurate, which can effectively improve the decomposition accuracy of the VMD method.

(3)

In order to solve the problem that false IMF components affect the accuracy of signal feature extraction, a comprehensive evaluation factor method is proposed. This method fully considers the characteristics of each IMF component in the time domain and frequency domain and the degree of correlation with the original signal, which can effectively filter out noise and other interference components, and select IMF components sensitive to signal feature information.

(4)

In view of the fact that the fault diagnosis is a usually small sample, the fuzzy neural network is used to identify the rotor fault features obtained from the VMD symmetrical polar image. Compared with BP neural network, SVM and Random Forest, this method has higher average recognition accuracy of 98.4%, which provides a solution to the fault diagnosis problem in the case of a small sample.

2. Rotor Fault Diagnosis Based on VMD Symmetric Polar Image and Fuzzy Neural Network

A novel method for fault feature extraction and pattern recognition is proposed to eliminate the influence of environmental noise on the precision of extracting feature information of symmetrical mirror image of rotor vibration signals in this paper.

First of all, VMD is applied to decompose the rotor vibration signal and comprehensive evaluation factor (CEF) method is used to reconstruct the signal. Secondly, two-dimensional feature parameters of the snowflake image of the reconstructed signal are calculated to construct multidimensional feature vectors. In the end, the fuzzy neural network is applied for rotor fault diagnosis. Flow chart of the proposed method is shown in Figure 1.

• Collection of rotor vibration signal. The vibration signals of various fault states of rotor are collected by the ZT-3 rotor fault-simulating experimental device.
• Decomposition of the signal. VMD is applied to decompose rotor vibration signals into several IMF components. During VMD decomposition, the parameters K and α are selected according to the principle of the maximum value of center frequency and the principle of minimum value of multi-scale fuzzy entropy (MFE).
• Signal reconstruction. In order to eliminate environmental noise, the background signal and other disturbing components, and retain useful IMF components, the original rotor vibration signal is reconstructed by the CEF method.
• Feature extraction. The reconstructed signal is transformed into a two-dimensional snowflake image through the symmetrical polar coordinate method. Two-dimensional feature parameters of the snowflake image are used to form a state feature vector, which can be obtained by the gray level co-occurrence matrix.
• Fuzzy neural network recognition. The extracted state feature vectors are randomly allocated into two sets, which are the training sample set and the test sample set, respectively. The training sample set is used for training the fuzzy neural network. The testing data is imported into the well-trained fuzzy neural network. Then different working conditions of the rotor can be recognized and the recognition results can be automatically output.

Figure 1. Flow chart of rotor fault diagnosis based on VMD symmetric polar image and fuzzy neural network.

Figure 1. Flow chart of rotor fault diagnosis based on VMD symmetric polar image and fuzzy neural network.

3. Model Establishment

The mathematical model for this novel fault diagnosis method includes two parts: fault feature extraction and fault recognition. In the part of fault feature extraction, the original rotor vibration signal is decomposed by VMD and the signal is reconstructed by the CEF method, which can eliminate the disturbing IMF component. Then, a method is adopted to extract the fault feature of the snowflake image of the reconstructed signal and compose state feature vectors. In the part of fault recognition, the fuzzy neural network is used to perform fault recognition so that the accurate fault classification result is obtained.

3.1. Rotor Fault Feature Extraction Based on VMD Symmetrical Polar Image

3.1.1. Variational Mode Decomposition

VMD method is an adaptive signal decomposition method, which decomposes the original rotor vibration signal f(t) into a series of modal functions u₁, u₂, ..., u_k, with different center frequencies. The sum of modal functions is equal to the input signal under the constraint, and the constrained variational problem is constructed as follows:

$$\begin{array}{l}\underset{u_k,{\omega }_k}{\min }\left({{\displaystyle \sum _{k=1}^K\Vert {\partial }_t\left[\left(\delta (t)+\frac{j}{\pi t}\right)u_k(t)\right]e^{-j{\omega }_{k}t}\Vert }}_2^2\right)\\ {\displaystyle \sum _{k=1}^K{u}_k=f(t)}\end{array}\}$$

(1)

where, {u_k} = {u₁, u₂,‧‧‧, u_K} and {ωk} = {ω₁, ω₂,‧‧‧, ω_K} are modal functions, and each mode function’s central frequencies.

To solve the problem, the quadratic penalty factor α and Lagrange multiplication operator λ are used to ensure the reconstruction precision of the signal and the constraint strictness in the presence of Gaussian noise, respectively.

$$L(u_k,{\omega }_k,\lambda )=\alpha {{\displaystyle \sum _{k=1}^K\Vert {\partial }_t\left[\left(\delta (t)+\frac{j}{\pi t}\right)u_k(t)\right]e^{-j{\omega }_{k}t}\Vert }}_2^2+{\Vert f(t)-{\displaystyle \sum _{k=1}^K{u}_k(t)}\Vert }_2^2+\langle \lambda (t),f(t)-{\displaystyle \sum _{k=1}^K{u}_k(t)}\rangle$$

(2)

Initialize u_k, ω_k and λ. The solution to the constrained variational problem is using the alternate direction method of multipliers to find the saddle point of the augmented Lagrangian in a sequence of iterative sub-optimization. Then modal component u_k and the center frequency ω_k can be expressed as follows:

$${\widehat{u}}_k^{n+1}(\omega )=\frac{\widehat{f}(\omega )-{\displaystyle \sum _{i\ne k}{\widehat{u}}_i(\omega )}+\frac{\widehat{\lambda }(\omega )}{2}}{1+2\alpha {(\omega -{\omega }_k)}^2}$$

(3)

$${\omega }_k^{n+1}=\frac{{\displaystyle {\int }_0^{\infty }\omega {\left|{\widehat{u}}_k(\omega )\right|}^{2}d\omega }}{{\displaystyle {\int }_0^{\infty }{\left|{\widehat{u}}_k(\omega )\right|}^{2}d\omega }}$$

(4)

In the VMD algorithm, modal decomposition parameter K and quadratic penalty factor α directly affect the precision of signal decomposition. At present, the parameter selection method based on the default parameter setting method is difficult to highlight the advantages of this signal decomposition method and limits its performance. Although intelligent search algorithms can optimize the parameter values given above, such methods take a long time and are difficult to achieve the purpose of actual detection. Therefore, an efficient method of the VMD decomposition key parameter selection is necessary.

3.1.2. Selection Method of Modal Decomposition Parameter K

The modal decomposition parameter K determines the number of IMF components obtained by VMD. If the value of K is not selected properly, it is easy to produce modal mixing, which affects the precision of signal analysis. The center frequency of each IMF component obtained by the VMD decomposition of the signal is distributed from low to high. If the best-preset ratio parameter K is obtained, it means that the center frequency distribution between adjacent IMF components is reasonable [21]. Moreover, the final result value K will not be similar or mixed. If the optimal value of K is selected, it means that the center frequency value of the highest order IMF component is the largest and the value of the maximum center frequency still tends to be stable during subsequent decomposition. In this paper, the optimal value of K is selected based on the principle that the center frequency of the IMF component is the maximum value for the first time.

3.1.3. Selection Method of Quadratic Penalty Factor α

The choice of the quadratic penalty factor α determines the bandwidth of each IMF component obtained by the VMD. The smaller the α value, the larger the bandwidth of the IMF component.

An important parameter for evaluating the complexity of non-stationary signals is fuzzy entropy. Compared to other methods such as sample entropy and approximate entropy, the value of fuzzy entropy is less affected by related parameters and has stronger stability. The rotating machinery vibration is random vibration, and the signal is highly irregular and complex. When a fault occurs, the impact component in the vibration signal increases, and the regularity and self-similarity of the signal also increase [22], which results in a small fuzzy entropy value of the signal. In order to describe the complexity of the signal more accurately, MFE is selected as the evaluation factor of the penalty factor α. When the optimal value of α is selected, the MFE value of the reconstructed signal based on the VMD should be the minimum value.

For signal {x(i), i = 1, 2,…, N}, given the embedding dimension of the fuzzy entropy and the similarity tolerance r, a new coarse-grained signal y_j(τ) can be established according to the original signal x(i):

$$y_j(\tau )=\frac{1}{\tau }{\displaystyle \sum _{i=(j-1)\tau +1}^{j\tau }{x}_i}, 1\le j\le \frac{N}{\tau }$$

(5)

where, τ = 1, 2,…, N, is the scale factor.

3.1.4. Comprehensive Evaluation Factor Method

When the rotating machinery has a fault, the vibration signal energy will change. Compared with the fault information component, the noise component has less energy [23]. The energy fluctuation coefficient η_k can effectively represent the degree of energy change in the signal [24]. The Time-domain correlation coefficient selection method is the most commonly used method to distinguish false disturbing components in the signal [25]. The noise disturbing component in the signal not only affects the signal’s pattern but also affects the calculation precision of the time-domain correlation coefficient. But in the frequency domain, the correlation coefficient ρ_k of the signal power spectrum is hardly affected by noise. Therefore, in order to eliminate the influence of noise and other disturbing components on the precision of fault feature extraction, the above parameters are comprehensively characterized for signal patterns, a CEF P_k is constructed here:

$$\{\begin{array}{l}{P}_k=\left|\frac{{\eta }_k+{\rho }_k}{2}\right|\\ {\eta }_k=\frac{{{\displaystyle \sum _{t=1}^n\left|u_k(t)-{\overline{u}}_k(t)\right|}}^2}{{\displaystyle \sum _{t=1}^n{\left|x(t)-\overline{x}(t)\right|}^2}}\\ {\rho }_k=\frac{{\displaystyle \sum _{t=1}^n(u_k(t)-{\overline{u}}_k(t))\cdot (x(t)-\overline{x}(t))}}{\sqrt{{\displaystyle \sum _{t=1}^n{(u_k(t)-{\overline{u}}_k(t))}^2\cdot {(x(t)-\overline{x}(t))}^2}}}\end{array}$$

(6)

where, u_k(t) is the k-th IMF component, x(t) is the original signal. P_k is dimensionless, and its numerical range is between [0, 2]. Based on [24,25], the threshold of P_k is set to be 0.15.

3.2. The Basic Principle of Symmetric Polar Image

The symmetrical polar coordinate method transforms the signal under polar coordination and expresses it in the form of a two-dimensional snowflake pattern. The snowflake image is directly obtained from the time-domain signal, without the need for frequency analysis. In addition, because of the symmetry of the snowflake image, it can show the differences between each other more clearly than the time-domain waveform of the signal, and it is handy to compute. Hence, this method is an ideal image generation method.

The amplitude of signal x(t) at time i is x_i, at time i + L it is x_i + L. Figure 2 shows the basic schematic symmetrical polar coordinate method [26,27].

In Figure 2, r(i) is the polar coordinate radius, δ(i) is the rotation angle in the counterclockwise direction from the initial line, and θ(i) is the clockwise rotation angle from the initial line. By changing the rotation angle of the initial line, a set of signals (x_i, x_i + L) can be converted to a hexagonal snowflake mirror symmetric image in polar coordinates. The specific conversion equation is as follows:

$$r(i)=\frac{x_i-x_{\min }}{x_{\max }-x_{\min }}$$

(7)

$$\delta (i)=\varphi +\frac{x_{(i+L)}-x_{\min }}{x_{\max }-x_{\min }}g$$

(8)

$$\theta (i)=\varphi -\frac{x_{(i+L)}-x_{\min }}{x_{\max }-x_{\min }}g$$

(9)

where x_max is the maximum of the signal amplitude; x_min is the minimum of the signal amplitude; L is the time interval; φ is the rotation angle of the initial line; g is the amplification factor of angle. Generally, L = 3~10, φ = 60°, g = 20°~60° is applicable [26,27].

3.3. Fault Feature Extraction

The gray level co-occurrence matrix with texture information in the spatial distribution relationship between pixels, can accurately describe the roughness, vertical direction and complexity of the image. Therefore, it is often used as a feature parameter in image analysis [27]. The gray level co-occurrence matrix is generated as follows: take one pixel A (x, y) and another pixel B (x + Dx, y + Dy) in the gray image, then make statistics on the probability of their simultaneous occurrence P(i, j, d, θ). The equation of P(i, j, d, θ) is as follows:

$$\begin{array}{l}P(i,j,d,\theta )=\{[(x,y),(x+Dx,y+Dy)|\\ f(x,y)=i,f(x+Dx,y+Dy)=j]\}\end{array}$$

(10)

where i, j = 0, 1, 2, ‧‧‧, N−1, and N is the gray value of the image; x, y are the horizontal and vertical values of coordinate of pixel A in the image; the deviation between pixel A and pixel B in x direction and y direction are Dx and Dy, respectively; d is the distance between A and B; θ is the angle between the connecting line of two pixels and the horizontal line [28]. Generally, θ = 0°, 45°, 90° and 135° [29], as shown in Figure 3.

Due to the large amount of data in the gray level co-occurrence matrix, it is generally not used as a feature to distinguish image texture directly, but based on some statistics obtained from its calculation as a feature of image texture classification. Therefore, the researchers proposed 14 statistical feature quantities obtained from the calculation of gray level co-occurrence matrix [26]. In the study, six commonly and effective statistical features were used as the characteristic parameter: max probability (mp), entropy (ent), contrast (con), correlation (cor), energy (ene) and inverse difference moment (idm).

Before extracting the feature statistics of the gray level co-occurrence matrix, the normalization process is generally performed according to Equation (11).

$$g(i,j)=\frac{P(i,j)}{{\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}P(i,j)}}}$$

(11)

3.3.1. Max Probability

It reflects the maximum probability of the simultaneous occurrence of pixel A and pixel B in the gray level co-occurrence matrix. The equation is

$$mp=\underset{i,j}{\max }(g(i,j))$$

(12)

3.3.2. Entropy

It can represent the information amount of the image and the non-uniformity of the image texture. If the image is full of texture, the entropy value can reach the top. If there is no texture, the entropy value is equal to 0. The equation is

$$ent=-{\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}g(i,j)\log (g(i,j))}}$$

(13)

3.3.3. Contrast

It reflects the clarity of the image texture. The deeper the image texture groove, the higher the contrast value, and vice versa. The equation is

$$con={\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}{(i-j)}^{2}g(i,j)}}$$

(14)

3.3.4. Correlation

It reflects the degree of similarity of gray level co-occurrence matrix pixels in row and column directions. It is also a measure of the linear relationship of images. The equations are

$$\{\begin{array}{l}cor=\frac{{\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{i=0}^{N-1}(i-{\mu }_i)(j-{\mu }_j)g(i,j)}}}{{\sigma }_i{\sigma }_j}\\ {\mu }_i={\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}i\cdot (i,j)}}\\ {\mu }_j={\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}j\cdot (i,j)}}\\ {\sigma }_i^2={\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}{(i-{\mu }_i)}^2\cdot g(i,j)}}\\ {\sigma }_j^2={\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}{(j-{\mu }_i)}^2\cdot g(i,j)}}\end{array}$$

(15)

3.3.5. Energy

Its value is the sum of the squares of each gray level co-occurrence matrix element. It reflects the uniformity of the gray level change of image texture and the thickness of the texture. The equation is

$$ene={\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}{g}^2(i,j)}}$$

(16)

3.3.6. Inverse Difference Moment

It reflects the regularity and local variation of image texture. The more regular the image texture, the higher the inverse difference moment value, and vice versa. The equation is

$$idm={\displaystyle \sum _{i=0}^{N-1}{\displaystyle \sum _{j=0}^{N-1}\frac{1}{1+{(i-j)}^2}}}g(i,j)$$

(17)

3.4. Fuzzy Neural Network

Fuzzy neural network combines the advantages of fuzzy logic and neural networks well, and this method is widely used in fault diagnosis [30,31,32]. The fuzzy set represented by the membership function value in the fuzzy theory transforms the information with certain similarities. Meanwhile, a neural network has a strong ability for self-learning and information storage. Therefore, the fuzzy set obtained by fuzzy theory can be trained by neural networks to realize fault diagnosis. Furthermore, the fuzzy theory changes the fuzzy fault information of the recognition into the deterministic information, which can efficiently obtain more accurate diagnosis results.

The learning of fuzzy neural networks is a process of updating the network connection parameters constantly, and finally makes the network performance reach optimal effect. The learning error is calculated by the actual output value and the target value, then the backpropagation error is used to adjust the model condition parameters. The calculation steps of the fuzzy neural network used in this paper are as follows.

Step 1. For the k-dimensional input vector x = [x₁, x₂, ‧‧‧, x_k], the membership degree of each input variable x_j is calculated according to the fuzzy rules, and the membership function is Gaussian type.

$$\mu A_j^i=\exp [-(x_i-c_j^i)2/b_j^i]$$

(18)

where,$j=1,2,\cdot \cdot \cdot ,k$, $i=1,2,\cdot \cdot \cdot ,n$, $c_j^i$, $b_j^i$ are the center and width of the membership function respectively; k is the dimension of the input parameter (the number of eigenvectors); n is the number of fuzzy subsets.

Step 2. Each membership degree is fuzzy calculated, and the fuzzy operator is continuous multiplication operator.

$${\omega }^i=\mu A_j^1(x_1)\mu A_j^2(x_2)\cdot \cdot \cdot \mu A_j^k(x_k), i=1,2,\cdot \cdot \cdot ,n$$

(19)

Step 3. According to the result of fuzzy calculation, the output value of the fuzzy model is calculated.

$$y_i={\displaystyle \sum _{i=1}^n{\omega }^i(p_0^i+p_1^i{x}_1+\cdot \cdot \cdot +p_k^i{x}_k)}/{\displaystyle \sum _{i=1}^n{\omega }^i}$$

(20)

Step 4. Calculate the error between the expected output and the actual output of the network.

$$e=\frac{1}{2}{(y_d-y_c)}^2$$

(21)

where, y_d is the expected output of the network; y_c is the actual output of the network.

Step 5. Adjust the network connection parameters by Equation (22).

$$\{\begin{array}{l}{p}_j^i(k)=p_j^i(k-1)-\alpha \frac{\partial e}{\partial p_j^i},\frac{\partial e}{\partial p_j^i}=(y_d-y_c){\omega }^i/{\displaystyle \sum _{i=1}^n{\omega }^i{x}_j}\\ c_j^i(k)=c_j^i(k-1)-\beta \frac{\partial e}{\partial c_j^i},b_j^i(k)=b_j^i(k-1)-\beta \frac{\partial e}{\partial b_j^i} \end{array}$$

(22)

4. Experimental Analysis

4.1. Signal Acquisition

For the sake of verifying the effectiveness of this proposed method, the rotor vibration signals used in this section are generated by ZT-3 rotor fault simulating experimental device, which could be seen in Figure 4. It is mainly consisted of one drive motor, onephotoelectric sensor, an acceleration sensor and a signal acquisition system. The motor speed is controlled by a speed controller, whose rated current of the motor is 2.5 A and the output power is 250 W. The rotor vibration signal is measured by using an AI005 acceleration sensor. The signal is processed by USB-5936 dynamic signal collector and it is stored by the computer.

Table 2. Equipment list of the experiment.

Order Number	Equipment Name	Main Parameters (or Components)	Quantity
1	ZT-3 rotor fault simulating experimental device	➀ Speed regulation range of governor: 0–10,000 RPM; ➁ Shaft length: 320 mm, diameter: 9.5 mm; ➂ Rotor width is 19 mm, diameter is 76 mm.	1
2	AI005 acceleration sensor	➀ Sensitivity: (20 ± 5 ℃) 50 mv/mm/s2; ➁ Frequency range: 4–2000 Hz.	1
3	USB-5936 dynamic signal collector	➀ Number of channels: 4; ➁ Resolution ratio: 24 bit; ➂ Dynamic range: 102 dB.	1
4	Desktop computer	➀ CPU: i5-12490F; ➁ Graphic card: GeForce RTX 3060; ➂ Memory:512 G + 1 T; ➃ Hard disk: SSD + HHD.	1
5	Auxiliary device	➀ Metal gasket; ➁ M2 screws; ➂ Rub-impact screw; ➃ Data cable.	Several

shows the equipment list of the experiment.

During the experiment, four fault states were simulated: bearing block looseness fault (BBLF), misalignment fault (MF), imbalance fault (IF) and local rub-impact fault (LRIF).

Table 3. The fault simulation settings.

Order Number	The Fault Types	Fault Simulation Settings
1	BBLF	The screw used for fixing on the bearing block is loosened to simulate the looseness fault of the bearing block.
2	MF	Insert a metal gasket under the bearing block to raise the rotating shaft, and then tighten the screw to simulate the misalignment fault.
3	IF	Put two standard M2 screws (each screw weighs 0.3 g) into the circumference-shaped groove of the rotor disc to simulate the imbalance fault.
4	LRIF	A support is installed on the base of the test bench, and a threaded hole is reserved above the support, a long screw made of plastic material is screwed into the threaded hole to make it in contact with the rotating shaft of the rotor, as the rotation speed increases, friction will occur between the rotating shaft and the screw, so as to simulate the local rub-impact fault.

shows the fault simulation settings. Figure 5 shows the four faults in this experiment.

During signal acquisition, the signal sampling frequency was set to 2000 Hz and the motor speed was set to 2700 r/min. For each working condition of the rotor, 50 groups of samples were collected as a sample set in which 30 groups of signal samples were randomly chosen as the training sample set, the remaining samples were test sample set. The sampling length of each signal was 2000 points.

Table 4. Sample information of rotor in the experiment.

Fault Types	The Quantity of Training Samples	The Quantity of Test Samples	Category Label
Normal(Nor)	30	20	1
BBLF	30	20	2
MF	30	20	3
IF	30	20	4
LRIF	30	20	5

shows the relevant information on the sample signals used in this experiment.

4.2. Feature Extraction of the Snowflake Image Based on VMD Reconstruction Signal

Figure 6 and Figure 7 show the time-domain waveform and spectrum of a group of rotor signals randomly selected under different states. From Figure 6 and Figure 7, noise-disturbing components are clearly shown in the signal and the time-domain waveform of the rotor signal in some states have similar patterns, such as BBLF and MF. Therefore, the rotor fault state could not be effectually distinguished by the time-frequency waveform or signal spectrum. For removing noise and extracting the characteristics of the signal accurately, it is very necessary to preprocess the signal.

The center frequency of each IMF component of an MF signal obtained by VMD under different K values is shown in

Table 5. The center frequency of each IMF component corresponds to different K values.

The Modal Decomposition Parameter K	Center Frequency/Hz
	IMF1	IMF2	IMF3	IMF4	IMF5	IMF6
K = 2	22	302	—	—	—	—
K = 3	22	301	566	—	—	—
K = 4	21	95	451	687	—	—
K = 5	21	95	450	686	688	—
K = 6	21	95	313	458	687	703

.

From Table 5, when K = 4, the IMF component center frequency obtains the maximum for the first time. When K = 5, the center frequency 686 Hz and 688 Hz appear, indicating that over-decomposition occurs during the VMD decomposition process, and the maximum growth rate Δ is defined here to evaluate the fluctuation of the center frequency:

$$\Delta =\frac{{\omega }_k-{\omega }_{k-1}}{{\omega }_{k-1}}\times 100\%$$

(23)

where, when K = k, ω_k and ω_k-1 are the k-th and (k−1)th IMF components of the original rotor vibration signal obtained by VMD. Figure 8 shows the change of maximum growth rate Δ under different values of K.

From Figure 8 we can know that when K = 4, Δ = 52.33%. With the increase of K, the value of Δ does not increase significantly. For example, when K = 5, Δ is only 0.29%, indicating that the maximum center frequency does not increase obviously. Therefore, it can be considered that the optimal value of K is 4.

For each scale factor τ, the fuzzy entropy of the signal y_j(τ) at each scale is calculated first, and then the mean value of the fuzzy entropy is obtained. When K = 4 and with different values of α, the mean value of MFE of the reconstructed MF signal is shown in Figure 9, where, m = 2, r = 0.15‧SD (standard deviation of the signal to be decomposed), scale factor s = 10.

From Figure 9, when α = 2500, the MFE mean value of the reconstructed signal after the VMD decomposition is the smallest, which indicates that the impact component related to the fault feature in the reconstructed signal contains the most, and also shows strong regularity and self-similarity. So the penalty factor α = 2500 is taken to decompose the MF signal.

When K = 4 and α = 2500, the VMD decomposition results of the MF signal are shown in Figure 10.

It could be seen that the VMD method overcomes the modal mixing problem from Figure 10. The frequency of each IMF component, with the pattern of concentrating near its center frequency, shows that this method has effectively reduced information leakage between the modal components, which proved its effectiveness.

Calculate P_k for the determination of the false IMF components. The calculation results are displayed in Figure 11.

As shown in Figure 11, the sensitive components obtained from the discrimination results of comprehensive evaluation factor Pk of the IMF components are the IMF1 and IMF2. Following the previous analysis, IMF3 and IMF4 are the noise disturbing component. Therefore, IMF1 and IMF2 are selected as sensitive components to reconstruct the MF signal. According to the above methods, the rotor vibration signals in other states in Figure 6 can be reconstructed. Figure 12 and Figure 13 show the time-domain waveform and spectrum of the rotor reconstructed signal. It can be seen from Figure 12 and Figure 13 that the noise components in the MF reconstructed signal are effectively filtered out, and the impact characteristics in the signal are highlighted, which ensures the precise composition of two-dimensional snowflake images.

60 groups of vibration signals under each state collected before were constructed by the VMD method, and 300 snowflake images were obtained by substituting these signals into the polar coordinate algorithm. In this experiment, L was set 3, g was set 30°. snowflake images corresponding to different states were illustrated in Figure 14.

It can be seen from Figure 14 that the two-dimensional snowflake image of the reconstructed signal in different states is more intuitive than the one-dimensional vibration signal. Rotor faults can be diagnosed by the feature difference of the snowflake image.

In order to distinguish the difference between the snowflake image of the reconstructed rotor signal under different fault conditions, the texture feature of the image is extracted as the key parameter after it is transformed into a gray image. In image texture analysis, extract six feature parameters (mp, ent, con, cor, ene, idm) in four directions (φ = 0°, 45°, 90°, 135° and d = 1) from gray level co-occurrence matrix, altogether 24 parameters obtained. Due to different directions, the co-occurrence matrix will contain different texture information, but when the difference between the snowflake images of the two signals is small, the texture information will have similar values, affecting the accuracy of feature extraction. Therefore, the mean value of the co-occurrence matrix infourdirections is taken as the final gray level co-occurrence matrix, and features are extracted.

Table 6. Characteristic parameters of the rotor under different states.

Working State	Group of Signals	Characteristic Parameter
		mp	ent	con	cor	ene	idm
Nor	1	0.8510	1.5074	0.5883	0.9395	0.7259	0.9330
	2	0.8497	1.5081	0.5892	0.9417	0.7254	0.9288
	3	0.8524	1.5076	0.5887	0.9386	0.7308	0.9347
BBLF	1	0.8826	1.1955	0.4209	0.9486	0.7803	0.9521
	2	0.8829	1.1964	0.4212	0.9497	0.7792	0.9515
	3	0.8831	1.2018	0.4303	0.9510	0.7824	0.9526
MF	1	0.8440	1.5431	0.5718	0.9447	0.7146	0.9338
	2	0.8446	1.5523	0.5825	0.9465	0.7183	0.9346
	3	0.8448	1.5427	0.5730	0.9471	0.7148	0.9327
IF	1	0.8137	1.7647	0.6321	0.9482	0.6653	0.9249
	2	0.8149	1.7708	0.6415	0.9457	0.6687	0.9283
	3	0.8131	1.7652	0.6426	0.9473	0.6719	0.9242
LRIF	1	0.9073	1.0933	0.5399	0.9005	0.8236	0.9421
	2	0.9087	1.1021	0.5412	0.9109	0.8248	0.9468
	3	0.9130	1.0985	0.5381	0.9114	0.8239	0.9416

shows the feature parameters extracted under different conditions.

4.3. Fault Diagnosis Results and Comparative Analysis

The feature parameters extracted from the gray level co-occurrence matrix are used to form the state feature vector which is input to the fuzzy neural network for fault diagnosis. In this experiment, the dimension of input data is 6 and the dimension of output data is 1. It is determined that the number of input nodes in the network is 6 and the number of output nodes is 1. According to this situation, the number of selected membership functions is better to be 12 after the optimization. Therefore, the construction of the fuzzy neural network structure is 6-12-1. The center c and the width b of the membership function are randomly initialized, the maximum number of iterations is set to 1000, the initial network learning rate is 0.01, and the system error is 0.001. Figure 15 shows the fault recognition results of the proposed method. It could be seen that the fault recognition precision of the proposed method was 98%. The reliability of the proposed method in rotor fault feature extraction and fault recognition is verified. In order to further verify the advantages of this method, the following aspects are analyzed.

(1) For examining VMD effects in the proposed method, LMD, EMD and EEMD are applied to decompose the same rotor vibration signal. The rotor signal is reconstructed by the above method according to the steps in this paper, then the fault feature sample set is constructed based on the characteristic parameters of the reconstructed signal symmetric polar image. For eliminating the influence of accidental error, each method is tested five times, and the average recognition precision and standard deviation of the five test results are taken as the evaluation indexes of the method. The comparison of the five test results is shown in Figure 16. The average recognition precision and standard deviation are shown in

Table 7. The average recognition precision and standard deviation of the compared method.

Method	Average Recognition Precision/%	Standard Deviation
VMD	98.4	0.5477
LMD	92.8	0.8367
EMD	92.6	0.8944
EEMD	93.8	1.3038

. It can be seen from Figure 16 and Table 7 that the proposed method has the highest recognition precision and the average recognition precision is 98.4%. The standard deviation of the proposed method is 0.5477, which is far lower than that of the other three methods. It shows that the proposed method can not only effectively reduce the influence of disturbing components on the precision of signal feature extraction, but also improve the precision and stability of fault diagnosis.

(2) In order to study the influence of the partition ratio of the training sample set and test sample set on the recognition precision, the proportion of the training sample set and test sample set in each state is divided into five cases: 2/8, 4/6, 5/5, 6/4, 8/2 and five recognition times for each case. Figure 17 shows the average recognition results. Figure 16 shows that the proposed method still has the highest fault recognition precision under different proportions. The recognition precision of all methods increases with the increase of the proportion of training sample sets. If there are enough training samples, other methods can also obtain higher precision, but the computational efficiency will decrease.

(3) In addition, a k-cross-validation experiment (k = 5) is used to demonstrate the robustness of the proposed model. During the experiment, 50 groups of rotor signal samples in each state are divided into five subsets. Each subset (10 samples) data is used as a verification set and the remaining four groups of subsets (40 samples) data are used as a training set and input into the fuzzy neural network for training. The recognition results of five tests obtained through training are shown in Figure 18. It can be seen from Figure 18 that the recognition accuracy of the fuzzy neural network is ideal, and the average recognition accuracy is 99.2%, which proves that the fuzzy neural network can identify rotor faults in the case of small samples.

(4) To study the influence of different classifiers on recognition precision and calculation time, and to verify the advantages, the sample set extracted from the symmetrical polar coordinate image of the reconstructed signal based on VMD is identified five times by four classifiers: support vector machine (SVM), BP neural network, Random Forest and fuzzy neural network, and record the average recognition accuracy precision and calculation times (this time includes the time used to construct the feature sample set). The results are shown in

Table 8. Average recognition precision and calculation times of the four classifiers.

Classifier	Average Recognition Precision/%	Training Time/s	Test Time/s
SVM	96.4	25.83	10.63
BP neural network	92.2	31.46	11.12
Random Forest	93.6	24.72	10.77
Fuzzy neural network	98.4	26.24	10.85

. Where, the penalty factor C of SVM is 1, and the width coefficient σ = 1.58; for the BP neural network, the neural network structure is 6-14-5, the maximum number of iterations is set to 1000, the initial network learning rate is 0.01, and the system error is 0.001; for Random Forest, the decision tree is 500. This proposed method shows the highest recognition precision. Besides, compared with BP neural network, its computing efficiency has been slightly improved.

5. Conclusions

This paper proposes a rotor fault diagnosis method based on VMD symmetrical polar image and fuzzy neural network. The method performs well in fault diagnosis of rotor vibration signal under different conditions. The conclusions are as follows:

(1)

The two-dimensional snowflake image of the rotor vibration signal generated by the symmetric polar coordinate method is directly obtained from the time-domain signal and without the need for frequency analysis. Based on the parameter selection algorithm and CEF method, VMD can effectively eliminate noise disturbing components, strengthen fault feature components, and ensure the precision of image feature acquisition. Through five times recognition tests, the results show that the average recognition precision and standard deviation of the symmetrical polar coordinate image of the reconstructed rotor signal based on VMD are 98.4% and 0.5477. Compared with some similar signal decomposition methods (e.g., LMD are 92.8% and 0.8367, EMD are 92.6% and 0.8944 and EEMD are 93.8% and 1.3038), VMD could improve precision and stability of fault diagnosis.

(2)

For the rotor vibration signals under different working conditions, even in the case of small samples, the recognition rate of the fault diagnosis method based on VMD symmetric polar image and fuzzy neural network can reach 98%. Through the k-cross-validation experiment, the average recognition accuracy of this experiment is 99.2%, which demonstrates the robustness of the fuzzy neural network.

(3)

Compared with BP neural network, the fuzzy neural network has obvious superiority in identifying precision (e.g., the average recognition precision of the BP neural network is 92.2% and the fuzzy neural network is 98.4%) and saving calculation time (e.g., the training time and test time of BP neural network are 26.24 s and 10.85 s, fuzzy neural network are 26.24 s and 10.85 s). Compared with SVM and Random Forest, the fuzzy neural network also has a higher average recognition precision (SVM is 96.4% and Random Forest is 93.6%), and has the same computing efficiency.

Limited by the laboratory conditions, this experimental device is mainly used to simulate the above four types of rotor faults. Through the study of relevant literature [1,2,3,4,5], these four types of rotor faults are common faults and have certain universality. If necessary requirements can be obtained in the future, the authors will continue to study fault diagnosis methods for other rotor faults (e.g., rotor part falling, oil film whirl, and rotor broken bar fault) on the basis of this study and the generalization ability of the model proposed in this paper.