Rolling bearing fault signal processing method based on GWO parameter searching optimization VMD

Abstract

The invention discloses a rolling bearing fault signal processing method based on GWO parameter searching optimization VMD, and the method comprises the following steps: 1, employing a GWO-VMD combined wavelet threshold function de-noising method, carrying out the reconstruction of a useful noise component and a de-noising component, and obtaining a de-noised signal; 2, in each iteration process, the multi-scale permutation entropy is used as a fitness function, a GWO algorithm is used for automatic optimization to obtain parameters alpha and k, the accuracy of poor parameters caused by artificial experience is effectively improved, an optimal parameter combination is found, and signals are decomposed more accurately; and step 3, setting the re-constrained VMD according to the results of the parameters alpha and k obtained in the above step to obtain each IMF component, calculating the envelope spectrum of each component, and performing envelope analysis on the envelope spectrum to find the fault frequency.

Description

Technical Field

[0001] The invention relates to a rolling bearing fault signal processing method, in particular to a rolling bearing fault signal processing method based on GWO (grey wolf optimization) parameter-seeking optimization VMD. Background Art

[0002] Rolling bearings are widely used in aviation engines, railway transportation equipment and large industrial machinery and equipment due to their low friction resistance, low power consumption and high mechanical efficiency. They are responsible for power transmission and power conversion and are one of the main sources of mechanical failure. Since the health status of rolling bearings is directly related to the safe and stable operation of the equipment, timely and effective discovery of faults in key rolling bearings is of great significance to ensure the safe and stable operation of the machinery. Generally, the changes in the operating status of mechanical equipment can be effectively reflected through its dynamic signals (such as vibration signals), and the hidden fault feature information can be separated by signal processing technology. Therefore, fault signal processing is very important to ensure the healthy operation of rotating machinery. Rotating machinery and equipment will generate vibration signals different from normal operation when they fail. For example, cracks will appear on the outer ring of the rolling bearing, and the cracks on the surface will generate a series of periodic pulses when rotating. How to extract the periodic pulses in the vibration signal and identify the fault type is a hot issue for the fault signal processing of rotating machinery.

[0003] In order to extract fault features and identify fault types in vibration signals, people have designed a variety of methods, such as Fourier transform, wavelet transform, power spectral density estimation, empirical mode decomposition (EMD), local mean decomposition (LMD), etc. However, Fourier transform requires a fixed window length to be specified for frequency domain analysis, and the processing of non-periodic signals requires artificial extension or truncation, which cannot fully display their spectral characteristics. The calculation of power spectral density usually depends on the length of the signal and the selected window function. Improper selection of window function or too short signal may lead to deviation or distortion of spectrum estimation, especially when processing non-steady-state signals or When the signal length is limited, this sensitivity may affect the accuracy of spectral analysis. The effect of wavelet transform is extremely dependent on the wavelet basis function selected empirically, while EMD and LMD are very sensitive to noise and are prone to problems such as poor decomposition effect. Compared with EMD and LMD, VMD can flexibly select the number of modes according to needs, and has a more stable mathematical foundation and less mode overlap. Previous studies have shown that VMD has achieved good results in the field of rotating machinery fault signal processing. However, the performance of VMD is affected by its parameters, especially the penalty factor α and the modal decomposition number k. Choosing the best parameter combination is very important for the accuracy of the algorithm. Summary of the invention

[0004] The purpose of the present invention is to improve the accuracy of the penalty factor α and the modal decomposition number k of manually setting the rolling bearing fault signal type, and to provide a new solution for rolling bearing fault signal processing. A rolling bearing fault signal processing method based on GWO parameter optimization VMD is proposed.

[0005] A rolling bearing fault signal processing method based on GWO parameter optimization VMD, the method comprises the following steps:

[0006] Step 1: A denoising method based on GWO-VMD combined with wavelet threshold function is used. This method uses the correlation coefficient method to distinguish the IMF components decomposed by VMD into useful components and noise components, and then uses the wavelet threshold function to further denoise the noise components, and then reconstructs the useful noise components and denoised components to obtain a denoised signal;

[0007] Step 2: Use the GWO algorithm to automatically optimize the parameters. In each iteration, multi-scale permutation entropy is used as the fitness function, and the GWO algorithm is used to automatically optimize the parameters α and k.

[0008] Step 3: Set the parameters α and k obtained in the above steps to reconstruct VMD, obtain each component, calculate the envelope spectrum of each component, and perform envelope analysis on it to find the fault frequency.

[0009] Preferably, in step 1, a denoising method based on GWO-VMD combined with wavelet threshold function is adopted to obtain a denoised signal, and the steps are:

[0010] The calculation of wavelet coefficients needs to be determined by a threshold, and the selected threshold is shown in formula (6).

[0011]

[0012] Where N is the signal length, σ is the root mean square error of the noise;

[0013] The threshold function is as follows:

[0014]

[0015] Its j,k is the wavelet coefficient, λ is the threshold, is the threshold wavelet coefficient, sgn is the sign function, a and β are the adjustment parameters of the threshold function, where 0≤a<1, and β is a positive constant;

[0016] The optimized characteristic mode function is obtained by using the GWO-VMD algorithm, which is further processed. The high-frequency components containing noise are maintained by wavelet threshold, and the correlation coefficient method is introduced (the larger the absolute value of the correlation coefficient, the stronger the correlation, the closer the correlation coefficient is to 1 or -1, the stronger the correlation, and the closer the correlation coefficient is to 0, the weaker the correlation) to distinguish. As shown below, the correlation coefficient is calculated to distinguish between effective signals and noise signals. The correlation coefficient threshold less than or equal to is high-frequency noise, greater than the correlation coefficient threshold is a useful signal, and the IMF correlation coefficient value less than or equal to the correlation coefficient is processed.

[0017]

[0018] Where IMF i (t) is the i-th component, and N is the number of sampling points.

[0019] Preferably, the initial parameters a, A, and C of GWO are set to As the position of the i-th gray wolf, the range is 100 to 3500, randomly initialize the positions of N gray wolves in the gray wolf population, and then convert the position X of each gray wolf into i As penalty factor [α1, α2, …, α k ] Perform VMD on the signal, calculate the fitness value according to the multi-scale permutation entropy fitness function, and record the optimal three wolf leader individual positions X α , X β , X γ The positions of other gray wolf individuals are updated through formula (5), and the values ​​of parameters a, A, and C are updated through formulas 4 and 5. The above steps are repeated until the iteration is completed to obtain the optimal parameters α and k. The IMF components are obtained based on the VMD decomposition of the optimal parameter combination.

[0020] The gray wolf population hunting formula is:

[0021] D=|C·X p (t)-X(t)| (6)

[0022] X(t+1)=X p (t)-A·D (7)

[0023] Where: C and A are coefficients; X p t is the prey position; X(t) is the wolf position, and the update formulas of C and A are:

[0024] C=2r2 (8)

[0025] A=2ar1-a (9) Where r1 and r2 are random coefficients, whose modulus range is [0, 1], and a is the convergence factor, which decreases linearly from 2 to 0 with the number of iterations. Then, the position X of the three alpha wolves can be obtained by equations (1) and (2): α (t+1), X β (t+1), X γ (t+1), and then determine the positions of other gray wolves in the population through the positions of the three alpha wolves:

[0026]

[0027] The GWO algorithm searches for the global optimal solution based on the first three optimal solutions in the current solution space. At the same time, it changes the weight of an individual with the help of auxiliary coefficients during the search process to reduce the probability of the algorithm falling into the local optimum.

[0028] The rolling bearing fault signal processing method of the present invention, which uses the parameters based on the Grey Wolf Optimization Algorithm (GWO) optimization to adaptively obtain the penalty factor α and the modal decomposition number k, has the following beneficial effects compared with other adaptive signal decomposition methods for faulty bearing signals:

[0029] 1. The wavelet threshold function is used to construct a new filter. Compared with the traditional method, it has a faster convergence speed and stronger fault information separation ability. It not only improves the spectrum aliasing problem of the original VMD, but also can filter out stronger noise signals. Therefore, this method has a stronger filtering ability.

[0030] 2. A rolling bearing fault signal processing method is adopted, which involves optimization based on the GWO (Grey Wolf Optimization) algorithm, adaptively finding the penalty factor α and the modal decomposition number k. Compared with the traditional method, it does not need to rely on the operator's experience to manually set parameters, which effectively solves the problem that VMD needs to manually set parameters based on experience, resulting in poor decomposition effect. BRIEF DESCRIPTION OF THE DRAWINGS

[0031] Figure 1 It is a framework diagram of the rolling bearing fault signal processing method based on GWO parameter optimization VMD of the present invention;

[0032] Figure 2 The time domain and frequency domain decomposition results of the experimental bearing are shown;

[0033] Figure 3 Comparison of SNR graphs of different methods for experimental bearings;

[0034] Figure 4 The SNR diagram of the highest component of signal-to-noise ratio of different methods for experimental bearings;

[0035] Figure 5This is a structural diagram of the fault simulation test bench;

[0036] Figure numerals: 1 - tachometer, 2 - vertical acceleration sensor, 3 - horizontal acceleration sensor, 4 - bearing test seat, 5 - mass disk, 6 - left bearing seat, 7 - driving motor. DETAILED DESCRIPTION

[0037] Specific implementation method 1: Combination Figure 1 The embodiment of the present invention is described, and a rolling bearing fault signal processing method based on GWO parameter optimization VMD is characterized in that the method comprises the following steps:

[0038] Step 1: Use the GWO algorithm to automatically optimize the parameters. In each iteration, multi-scale permutation entropy is used as the fitness function, and the GWO algorithm is used to automatically optimize the parameters α and k.

[0039] Step 2: A denoising method based on GWO-VMD combined with wavelet threshold function is used. This method uses the correlation coefficient method to distinguish the IMF components decomposed by VMD into useful components and noise components, and then uses the wavelet threshold function to further denoise the noise components, reconstruct the useful noise components and denoised components, and obtain the denoised signal.

[0040] Step 3: Set the parameters α and k to construct the reconstructed VMD, obtain each component, calculate the envelope spectrum of each component, and find the fault frequency.

[0041] Specific implementation method 2: Setting the values ​​of the initial parameters a, A, and C of GWO to As the position of the i-th gray wolf, the range is 100 to 3500, randomly initialize the positions of N gray wolves in the gray wolf population, and then convert the position X of each gray wolf into i As penalty factor [α1, α2, …, α k ] Perform VMD on the signal, calculate the fitness value according to the multi-scale permutation entropy fitness function, and record the optimal three wolf leader individual positions X α , X β , X γ The positions of other gray wolf individuals are updated by formula (5), and the values ​​of parameters a, A, and C are updated by formulas 4 and 5. The above steps are repeated until the iteration is completed, and the optimal parameters α and k are obtained. The IMF components are obtained by VMD decomposition based on the optimal parameter combination.

[0042] The gray wolf population hunting formula is:

[0043] D=|C·X p (t)-X(t)| (1)

[0044] X(t+1)=Xp (t)-A·D (2)

[0045] Where: C and A are coefficients; X p t is the prey position; X(t) is the wolf position. The update formulas for C and A are:

[0046] C=2r2 (3)

[0047] A=2ar1-a (4) In the formula, r1 and r2 are random coefficients, and their modulus range is [0, 1]; a is the convergence factor, which decreases linearly from 2 to 0 with the number of iterations. Then, the position X of the three alpha wolves can be obtained by using formula (1) and formula (2): α (t+1), X β (t+1), X γ (t+1), and then determine the positions of other gray wolves in the population through the positions of the three alpha wolves:

[0048]

[0049] The GWO algorithm searches for the global optimal solution based on the first three optimal solutions in the current solution space. At the same time, it makes full use of the auxiliary coefficient to change the weight of an individual in the search process, so the probability of the algorithm falling into the local optimum is low.

[0050] Specific implementation method three: This implementation method is different from the specific implementation method two in that: in step one, a denoising method based on GWO-VMD combined with wavelet threshold function is used to obtain a denoised signal, and the specific steps are:

[0051] The calculation of wavelet coefficients needs to be determined by the threshold. If the threshold is too large, some useful signals will be mistakenly filtered out as noise. On the contrary, there will be residual noise, and the signal will be distorted due to incomplete denoising. The threshold selected in this paper is shown in formula (6).

[0052]

[0053] Where N is the signal length and σ is the root mean square error of the noise.

[0054] The hard threshold function can retain more of the peak features of the real signal, but due to poor continuity, oscillations will occur during the reconstruction process. There is a continuous difference between the estimated value and the true value, which will cause some errors in the reconstructed signal even if the soft threshold function has greater continuity and the reconstructed signal is smoother. The hard threshold and soft threshold expressions are shown in formulas (7) and (8).

[0055] Hard thresholding method:

[0056]

[0057] Soft thresholding method:

[0058]

[0059] Where ω j,k is the wavelet coefficient, is the threshold wavelet coefficient, λ is the threshold, and sgn is the sign function.

[0060] This paper improves the threshold function while ensuring that mathematical reasoning is satisfied. The mathematical formula of the function is shown in formula (9).

[0061]

[0062] Its j,k is the wavelet coefficient, λ is the threshold, is the threshold wavelet coefficient, sgn is the sign function. a and β are the adjustment parameters of the threshold function. Where 0≤a<1, β is a positive constant. The theoretical proof of the continuity, asymptotics and bias of the threshold function proposed in this step is as follows:

[0063] a) Proof of continuity:

[0064] When j,k →λ + hour,

[0065]

[0066] b) Similarly, when hour,

[0067]

[0068] c) From the above formulas (10) and (11), it is not difficult to see that Mathematically prove the continuity of the threshold function.

[0069] d) Prove the asymptotic property of the threshold function:

[0070] When j,k →+∞,

[0071]

[0072] e) When ω j,k →-∞,

[0073]

[0074] f) From the above formulas (12) and (13), we can derive The asymptote of

[0075] g) Prove the deviation of the threshold function: When the odd function property is satisfied, there is ω j,k →±∞,

[0076]

[0077] From the above formula (14), we can see that the new threshold function ω j,k and There is no constant deviation between them.

[0078] The intermittent problem of the hard threshold function and the deviation problem of the soft threshold function can be effectively solved by parameter tuning.

[0079] The optimized characteristic mode function is further processed, and the high-frequency components containing more noise are maintained by wavelet threshold. The correlation coefficient method is introduced for distinction. As shown in formula (15), by calculating the correlation coefficient, the effective signal and the noise signal can be distinguished. Smaller correlation coefficient values ​​are usually considered to be high-frequency noise, and larger correlation coefficient values ​​are useful signals. Smaller IMF correlation coefficient values ​​are processed.

[0080]

[0081] Where IMF i (t) is the i-th component, and N is the number of sampling points.

[0082] The following examples are used to verify the beneficial results of the present invention.

[0083] Embodiment: The steps of the rolling bearing fault signal processing method based on GWO parameter optimization VMD in this embodiment are as follows:

[0084] Step 1: Use the GWO algorithm to automatically optimize the parameters. In each iteration, the maximum multi-scale permutation entropy is used as the fitness function. GWO is used to automatically optimize the parameters. Starting from the mathematical formula, the initial parameters a, A, and C of GWO are set. As the position of the i-th gray wolf, the range is 100 to 8000, randomly initialize the positions of N gray wolves in the gray wolf population, and then convert the position X of each gray wolf into i As penalty factor [α1, α2, …, α k ] Perform VMD on the signal, calculate the fitness value according to the multi-scale permutation entropy formula, and record the optimal three wolf leader individual positions X α , X β , X γ , through the relevant GWO formula, update the positions of other gray wolf individuals and the values ​​of parameters a, A, and C. Repeat the above steps until the iteration is completed and the optimal parameters α and k are obtained.

[0085] Step 2: Further process the optimized characteristic mode function and perform wavelet threshold maintenance on the high-frequency components containing more noise. The correlation coefficient method is used to distinguish the IMF components, and the correlation coefficient threshold is 0.5. Those above the threshold are considered to be valid signals, and those below the threshold are considered to be noise components. These noise components are classified into the improved wavelet threshold processing. The correlation coefficient method is introduced to calculate the correlation coefficient using the wavelet threshold formula, which can distinguish between valid signals and noise signals. Smaller correlation coefficient values ​​are generally considered to be high-frequency noise, and larger correlation coefficient values ​​are useful signals. Smaller IMF correlation coefficient values ​​are processed, and finally the valid signal and the processed IMF components are reconstructed to obtain the denoised signal.

[0086] Step 3, using the rolling bearing fault signal processing method based on GWO-VMD combined with wavelet threshold function, decompose the fault bearing vibration signal and extract the fault characteristics. The fault bearing data collected from the SpectraQuest rotating machinery fault simulation test bench is used in this example. The fault simulation test bench includes a speed 1, a vertical acceleration sensor 2, a horizontal acceleration sensor 3, a bearing test seat 4, a mass disk 5, a left bearing seat 6, and a drive motor 7. The output shaft of the drive motor 7 is coaxially fixedly connected to the left end of the rotating shaft, the left side of the rotating shaft is rotatably installed on the left bearing seat 6, the left and right sides of the middle of the rotating shaft are respectively fixed with clothes moth mass disks, and the right side of the rotating shaft is rotatably installed on the bearing test seat 4. The vertical acceleration sensor 2 and the horizontal acceleration sensor 3 are fixed on the bearing test seat 4, and a tachometer is provided on the drive motor 7 to display the speed of the drive motor 7. The left bearing seat 6, the drive motor 7, and the bearing test seat 4 are all fixed on the base. A test bearing is installed in the bearing test seat 4. The test bearing model is M-BFK-3 / 4, the speed is 1800r / min, the sampling frequency of the vertical acceleration sensor 2 and the horizontal acceleration sensor 3 is set to f=25600Hz, and the outer ring fault frequency of the test bearing is f o =107.9Hz, frequency f r =30Hz.

[0087] The collected bearing signals were analyzed using the proposed method. The number of iterations of GWO was set to 10, the population was set to 20, the optimization range of the penalty factor was [100, 3500], and the iteration range of the decomposition modulus was [2, 12]. The multi-scale permutation entropy was used as the fitness function. When the number of iterations was 2, the fitness function reached the minimum value of 0.6873. The optimal parameter combination was finally obtained as the penalty factor α = 50, the modal function k = 6, Figure 2 The time domain decomposition results of the bearing and the corresponding envelope spectrum are shown in Figure 2. Figure 2 It can be seen that the fault characteristic frequency f of the outer ring can be clearly observed in IMF1 oAnd its multiple harmonic frequencies, other components are also reflected in the fault frequency.

[0088] In order to further verify the effect of this method on rolling bearing fault signal processing, this method is compared with several commonly used adaptive signal decomposition methods, including the original VMD, empirical mode decomposition (EMD), local mean decomposition (LMD), and intrinsic time scale decomposition (ITD). In order to quantify the effect of re-constraining VMD to extract fault feature frequencies, the signal-to-noise ratio (SNR) is used for quantification. Figure 3 , Four The signal-to-noise ratio of each component in the five methods is shown. Among them, the horizontal axis τ = n × m, when τ = 1 to 6, represents the first component fault frequency f o When τ=7~12, it means the second component fault frequency f o From the figure, we can see that except for the second component, the signal-to-noise ratios of the other components of the re-constrained VMD are all greater than those of the EMD. The signal-to-noise ratio of the fifth component of the ITD is higher than that of the LMD, but the overall signal-to-noise ratio is the lowest. In contrast, the re-constrained VMD proposed in this paper has the highest signal-to-noise ratio overall. Among the above four methods, the components with the highest signal-to-noise ratio are Figure 4 As shown, the signal-to-noise ratio of reconstrained VMD is higher than all other methods.

Claims

1. A rolling bearing fault signal processing method based on GWO parameter optimization VMD, characterized by: The method comprises the following steps: Step 1: A denoising method based on GWO-VMD combined with wavelet threshold function is used. This method uses the correlation coefficient method to distinguish the IMF components decomposed by VMD into useful components and noise components, and then uses the wavelet threshold function to further denoise the noise components, and then reconstructs the useful noise components and denoised components to obtain a denoised signal; Step 2: Use the GWO algorithm to automatically optimize the parameters. In each iteration, multi-scale permutation entropy is used as the fitness function, and the GWO algorithm is used to automatically optimize the parameters α and k. Step 3: Set the parameters α and k obtained in the above steps to reconstruct VMD, obtain each component, calculate the envelope spectrum of each component, and perform envelope analysis on it to find the fault frequency.

2. The rolling bearing fault signal processing method based on GWO parameter optimization VMD according to claim 1 is characterized in that: In the step 1, a denoising method based on GWO-VMD combined with wavelet threshold function is adopted to obtain a denoised signal, and the steps are as follows: The calculation of wavelet coefficients needs to be determined by a threshold, and the selected threshold is shown in formula (6). Where N is the signal length, σ is the root mean square error of the noise; The threshold function is as follows: Its j,k is the wavelet coefficient, λ is the threshold, is the threshold wavelet coefficient, sgn is the sign function, a and β are the adjustment parameters of the threshold function, where 0≤a<1, and β is a positive constant; The optimized characteristic mode function is obtained by using the GWO-VMD algorithm, which is further processed. The high-frequency components containing noise are maintained by wavelet threshold, and the correlation coefficient method is introduced for distinction. The following shows that the correlation coefficient is calculated to distinguish the effective signal from the noise signal. The correlation coefficient threshold less than or equal to is high-frequency noise, and the correlation coefficient threshold greater than is a useful signal. The IMF correlation coefficient value less than or equal to the correlation coefficient is processed. Where IMF i (t) is the i-th component, and N is the number of sampling points.

3. The rolling bearing fault signal processing method based on GWO parameter optimization VMD according to claim 1 is characterized in that: Set the initial parameters a, A, and C of GWO to As the position of the i-th gray wolf, the range is 100 to 3500, randomly initialize the positions of N gray wolves in the gray wolf population, and then convert the position X of each gray wolf into i The signal is VMDed as the penalty factor [α1, α2, …, αk], and the fitness value is calculated according to the multi-scale permutation entropy fitness function, and the positions X of the three best wolf leaders are recorded. α , X β , X γ , update the positions of other gray wolf individuals through formula (5), and update the values ​​of parameters a, A, and C through formulas 4 and 5. Repeat the above steps until the iteration is completed, and the optimal parameters α and k are obtained. The IMF components are obtained based on the VMD decomposition of the optimal parameter combination; The gray wolf population hunting formula is: D=|C·X p (t)-X(t)| (6) X(t+1)=X p (t)-A·D (7) Where: C and A are coefficients; X p t is the prey position; X(t) is the wolf position, and the update formulas of C and A are: C＝2r2 (8) A=2ar1-a (9) Where r1 and r2 are random coefficients, whose modulus range is [0, 1], and a is the convergence factor, which decreases linearly from 2 to 0 with the number of iterations. Then, the position X of the three alpha wolves can be obtained by equations (1) and (2): α (t+1), X β (t+1), X γ (t+1), and then determine the positions of other gray wolves in the population through the positions of the three alpha wolves: The GWO algorithm searches for the global optimal solution based on the first three optimal solutions in the current solution space. At the same time, it changes the weight of an individual with the help of auxiliary coefficients during the search process to reduce the probability of the algorithm falling into the local optimum.

Images