CN113776837A - Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM - Google Patents

Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM Download PDF

Info

Publication number
CN113776837A
CN113776837A CN202111255452.0A CN202111255452A CN113776837A CN 113776837 A CN113776837 A CN 113776837A CN 202111255452 A CN202111255452 A CN 202111255452A CN 113776837 A CN113776837 A CN 113776837A
Authority
CN
China
Prior art keywords
fault
bearing
imf
nlm
rolling bearing
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202111255452.0A
Other languages
Chinese (zh)
Other versions
CN113776837B (en
Inventor
栾孝驰
徐石
沙云东
柳贡民
赵奉同
赵宇
张席
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shenyang Aerospace University
AECC Shenyang Liming Aero Engine Co Ltd
Original Assignee
Shenyang Aerospace University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shenyang Aerospace University filed Critical Shenyang Aerospace University
Priority to CN202111255452.0A priority Critical patent/CN113776837B/en
Publication of CN113776837A publication Critical patent/CN113776837A/en
Application granted granted Critical
Publication of CN113776837B publication Critical patent/CN113776837B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
    • G01M13/00—Testing of machine parts
    • G01M13/04—Bearings
    • G01M13/045—Acoustic or vibration analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/15—Correlation function computation including computation of convolution operations
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Abstract

The invention discloses a rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM, which comprises the steps of firstly decomposing fault signals acquired by a sensor by using CEEMDAN, screening decomposed IMF components by using a correlation coefficient-energy ratio-kurtosis criterion, and reconstructing the screened signals into a new signal; then, a non-local mean filter optimized by a wolf algorithm is used for carrying out parameter optimization selection to realize the optimal denoising effect, SG smooth filtering is carried out on the denoised signals to carry out secondary denoising, and finally obtained signals are subjected to feature extraction to obtain the fault features of the bearing; and finally, judging the fault type of the rolling bearing according to the fault characteristics of the bearing. The method can effectively inhibit noise interference in the vibration signal and improve the accuracy of fault diagnosis of the bearing.

Description

Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM
Technical Field
The invention relates to the technical field of bearing fault diagnosis, in particular to a rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM.
Background
Rolling bearings are widely used in various rotary mechanical systems, and the motion state of the rolling bearings has great influence on the accuracy, reliability and service life of the whole mechanical system. Because the bearing works at high speed and high load, the bearing is easy to break down, and the working state of the bearing can influence the running performance and safety of the whole machine. Therefore, fault diagnosis of the rolling bearing is essential for health monitoring of the rotating machine system. This will help to reduce costs associated with emergency maintenance and production. In the working process of the bearing, the local defects of the bearing are easy to damage due to reasons of improper installation, overload, poor lubrication and the like, and when the bearing fails early, the failure is diagnosed in time, so that the aggravation of the damage can be effectively avoided. Signal analysis is a key issue in mechanical failure diagnosis research and application. It is a tool that extracts fault features and then identifies fault patterns. In this field, there are two main types of bearing fault detection methods widely used at present, namely, acoustic signal analysis and vibration signal analysis. Among them, vibration signal based analysis is the most popular monitoring technique because of its ease of measurement. However, the main obstacle in analyzing the bearing fault from the vibration data is that the signal is easily submerged by background noise and is mixed with the vibration signal of the shaft, gear and other machines, so it is important to design an effective signal enhancement method to reduce the noise and highlight the fault characteristics.
Aiming at the characteristics of pulse property and nonlinearity of a rolling bearing signal, Empirical Mode Decomposition (EMD) can effectively perform self-adaptive processing on the bearing signal to improve the signal-to-noise ratio, but the problem of modal aliasing exists. In order to solve the problem of modal aliasing in EMD, Ensemble Empirical Mode Decomposition (EEMD) is used, and unlike EMD, EEMD solves the problem of modal aliasing by adding white gaussian noise for many times during the decomposition process. However, due to the introduction of the white gaussian noise, the final reconstructed signal is affected by the residual white noise, and the original signal cannot be accurately reconstructed. At present, wavelet transformation and blind source separation have proved to have good effect on separating noise and fault signals, but the problems of local distortion, partial loss of effective information and the like exist. Therefore, there are many defects in extracting the fault signal of the rolling bearing, and a more efficient and accurate bearing fault diagnosis technology is needed.
Disclosure of Invention
Aiming at the defects of the prior art, the invention provides a rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM. In the method, adaptive noise complete integrated empirical mode decomposition (CEEMDAN) is adopted, and the method calculates the only signal margin by adding adaptive white Gaussian noise in each decomposition process. Non-Local mean filtering (NLM) is an improved filtering to the conventional neighborhood filtering method, and uses the self-similarity of images, and the filtering effect depends on three parameters: m (search box radius), λ (bandwidth parameter), P (similar box radius), the choice of parameters is very important.
In order to solve the technical problems, the technical scheme adopted by the invention is as follows: a rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM comprises the following steps:
step 1: collecting vibration signals of a rolling bearing by adopting a vibration sensor, and simultaneously measuring relevant parameters of the bearing;
relevant parameters of the bearing include: contact angle alpha, ball number Z, ball diameter D and pitch diameter D.
Step 2: decomposing the acquired vibration signal of the rolling bearing by adopting a CEEMDAN method to decompose different IMF components;
and step 3: screening the decomposed IMF components by combining the kurtosis value, the correlation coefficient and the energy ratio, and performing linear reconstruction on the screened components;
the process of screening the decomposed IMF components by combining the kurtosis value, the correlation coefficient and the energy ratio is as follows:
step 3.1: the kurtosis value K for each IMF component is calculated as follows:
Figure BDA0003323704000000021
wherein x is the amplitude of the IMF component, and u is the average value of the IMF component amplitudes; σ is the standard deviation of the IMF component amplitude;
step 3.2: calculating a Pearson correlation coefficient r of each IMF component and the original signal;
step 3.3: the energy ratio coefficient of each IMF component is calculated as follows:
Figure BDA0003323704000000022
wherein E isxTotal energy that is the IMF component; eIMF(i) Is the energy possessed by the ith IMF component; epsilon is an energy ratio coefficient;
step 3.4: and carrying out weighted summation on the calculated kurtosis value K, the correlation coefficient r and the energy ratio coefficient epsilon to obtain a comprehensive screening index Kr epsilon value as follows:
Krε＝a1K+a2r+a3ε
wherein, a1、a2、a3Weighted values of kurtosis values, correlation coefficients and energy ratios of the IMF components are respectively;
step 3.5: and sorting the comprehensive screening index Kr epsilon values of all IMF components from large to small, and screening the largest first n IMF components.
And 4, step 4: the bandwidth parameter and the similar box radius of the NLM are optimized by adopting a wolf algorithm GWO, and the process is as follows:
step 4.1: initializing an NLM algorithm, and setting a bandwidth parameter lambda, a similar frame radius P and a search frame radius M; selecting a fixed search box radius M for calculating the speed;
step 4.2: initializing a gray wolf population, a convergence factor a and a coefficient vector A, C in the gray wolf algorithm;
Figure BDA0003323704000000031
A＝2a×r2-α
C＝2r1
where t is the number of iterations, tmax is the maximum number of iterations, r1And r2Is taken as [0, 1]]Random numbers within the interval;
step 4.3: taking the kurtosis value calculation function as a fitness function of a grey wolf algorithm, calculating the fitness value of each grey wolf, and simultaneously storing the parameters of the first three grey wolfs with higher fitness;
step 4.4: updating the position of each gray wolf;
step 4.5: updating a, A and C;
step 4.6: calculating the fitness values of all the gray wolves, and updating the fitness and the position of the first three gray wolves with the best fitness;
step 4.7: judging whether the maximum iteration times is reached, and if so, ending the optimization iteration process; and if not, returning to execute the step 4.4 to the step 4.7 to continue the optimization iteration process until the iteration is finished and outputting the optimal bandwidth parameter lambda and the similar box radius P.
And 5: carrying out noise reduction processing on the signal reconstructed in the step 3 by adopting an optimized GWO-NLM method, and outputting an optimal noise reduction signal;
step 6: carrying out secondary noise reduction on the output optimal noise-reduced signal through SG smooth filtering, carrying out envelope spectrum analysis on the signal subjected to secondary noise reduction, and extracting fault characteristic frequency;
and 7: and (4) calculating theoretical fault characteristic frequency of the bearing by combining related parameters of the bearing, comparing the fault characteristic frequency extracted in the step (6) with the theoretical fault characteristic frequency, and judging the fault type of the rolling bearing.
Adopt the produced beneficial effect of above-mentioned technical scheme to lie in:
1. the self-adaptive noise complete integration empirical mode decomposition is adopted, self-adaptive white Gaussian noise is added in each decomposition process, and the only signal margin is calculated, so that the reconstruction error is almost 0, and the calculation efficiency is high.
2. The kurtosis-correlation coefficient-energy ratio criterion selected by the method provided by the invention screens out the decomposed modal component, can effectively remove the irrelevant component, avoids the one-sidedness of the IMF component selected by a single index, and achieves the purpose of inhibiting the background noise.
3. The kurtosis index adopted by the method provided by the invention is used as the target function of the Hui wolf optimization algorithm, and the kurtosis is sensitive to periodic impact, so that the problem of optimizing parameters of the NLM algorithm can be solved, and the denoising performance of the NLM algorithm is optimized.
4. The method provided by the invention is applied to simulation experiment signals and actual rolling bearing vibration signals, and the result shows that the signals obtained by denoising by the method effectively inhibit the interference of noise, the signal-to-noise ratio is improved, and the method contains abundant and obvious fault information and has a certain engineering application value.
Drawings
FIG. 1 is a flowchart of a rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM in the embodiment of the present invention;
FIG. 2 is a denoising envelope spectrogram for simulating a bearing signal in the embodiment of the present invention;
FIG. 3 is an overall diagram of a denoising envelope spectrogram of an experimental bearing inner ring fault in the embodiment of the invention;
FIG. 4 is a partial detail diagram of a denoising envelope spectrogram of an experimental bearing inner ring fault in the embodiment of the invention;
FIG. 5 is an overall diagram of a denoising envelope spectrogram of an experimental bearing outer ring fault in the embodiment of the invention;
FIG. 6 is a partial detail view of a denoising envelope spectrogram of an experimental bearing outer ring fault in the embodiment of the invention.
Detailed Description
The following detailed description of embodiments of the present invention is provided in connection with the accompanying drawings and examples. The following examples are intended to illustrate the invention but are not intended to limit the scope of the invention.
The embodiment adopts real experimental data for analysis, and is taken from the bearing data center of western university of storage. And selecting a fault bearing for analysis as a 6205-2RJEM SKF type deep groove ball bearing, and carrying out single-part damage processing on faults of the inner ring and the outer ring of the bearing by utilizing an electric spark technology. The sampling frequency of the vibration data was 12000 HZ.
As shown in fig. 1, the rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM in the present embodiment is as follows.
Step 1: collecting vibration signals of a rolling bearing by adopting a vibration sensor, and simultaneously measuring relevant parameters of the bearing;
relevant parameters of the bearing include: contact angle alpha, ball number Z, ball diameter D and pitch diameter D.
Step 2: decomposing the acquired vibration signal of the rolling bearing by adopting a CEEMDAN method to decompose different IMF components;
in this embodiment, white noise is added to the original signal, EMD decomposition is performed on the original signal to which the noise is added, and the decomposed IMF is subjected to EMD decomposition1iThe components are averaged to obtain IMF1：
X(t)＝y(t)+z0ni(t)
Figure BDA0003323704000000041
In the formula: y is the original signal; z is a radical of0The amplitude of the added white noise; i is the number of times the signal is constructed, i belongs to Z+(ii) a N is the average total times of the set, N belongs to Z+；niTo satisfy the white noise of N (0,1) distribution.
The residual component of the first stage can be derived:
x1(t)＝X(t)-IMF1
in the formula: x is the number of1(t) is the residual component of the first stage.
Define aNovel operation EMDk(. cndot.) means the kth IMF component obtained after EMD decomposition. To the residual signal x1(t) plus white noise x1(t)+z1EMD1(ni(t)), and then averaging the first IMF by EMD decomposition to obtain the IMF2:
Figure BDA0003323704000000051
In the formula: IMF2Is the second phase modal component.
For K1, 2.., K, the K-th residual component can be derived:
xk(t)＝xk-1(t)-IMFk
repeating the above steps until the residual signal can not be decomposed by EMD, and obtaining a plurality of modal components and residual signals after the decomposition is completed, and the result is that:
Figure BDA0003323704000000052
in the formula: k is the total number of modes in the mode decomposition process.
And step 3: screening the decomposed IMF components by combining the kurtosis value, the correlation coefficient and the energy ratio, and performing linear reconstruction on the screened components;
the process of screening the decomposed IMF components by combining the kurtosis value, the correlation coefficient and the energy ratio is as follows:
step 3.1: the kurtosis value K for each IMF component is calculated as follows:
Figure BDA0003323704000000053
wherein x is the amplitude of the IMF component, and u is the average value of the IMF component amplitudes; σ is the standard deviation of the IMF component amplitude;
step 3.2: calculating a Pearson correlation coefficient r of each IMF component and the original signal;
step 3.3: the energy ratio coefficient of each IMF component is calculated as follows:
Figure BDA0003323704000000054
wherein E isxTotal energy that is the IMF component; eIMF(i) Is the energy possessed by the ith IMF component; epsilon is an energy ratio coefficient;
step 3.4: and carrying out weighted summation on the calculated kurtosis value K, the correlation coefficient r and the energy ratio coefficient epsilon to obtain a comprehensive screening index Kr epsilon value as follows:
Krε＝a1K+a2r+a3ε
wherein, a1、a2、a3Weighted values of kurtosis values, correlation coefficients and energy ratios of the IMF components are respectively;
step 3.5: and sorting the comprehensive screening index Kr epsilon values of all IMF components from large to small, and screening the largest first n IMF components.
In this embodiment, the pearson correlation coefficient is used to solve the correlation coefficient. The correlation coefficient is a numerical value between [ -1,1], and quantitatively describes the degree of correlation of X and Y, i.e., the larger the correlation coefficient is, and the lowest the degree of correspondence is when the correlation coefficient is 0.
Let two samples be X and Y, respectively, and the correlation coefficient be:
Figure BDA0003323704000000061
in the formula: r is the correlation coefficient of the two samples; cov (X, Y) is the covariance of the two samples;
Figure BDA0003323704000000062
is the variance of X;
Figure BDA0003323704000000063
is the variance of Y.
In this embodiment, when the Kr epsilon value of the comprehensive screening index is calculated, the same weight is taken for the three different indexes to be analyzed and calculated, normalization processing is performed on the indexes, the Kr epsilon value of each IMF component after CEEMDAN decomposition is calculated and sequenced, the first six IMF components of the maximum value are taken for linear reconstruction, and a reconstructed signal y is obtained.
And 4, step 4: bandwidth parameters and similar box radii of the NLM are optimized using the grey wolf algorithm GWO, and the GWO algorithm is inspirational from division of labor and collaborative foraging of wolfs. The method is a new group intelligent algorithm, and simulates the hierarchical structure of the wolf and the predation behavior of the wolf. The highest ranking wolf is the a wolf, located at the top of the food chain, responsible for leadership, decision making and other actions, followed by the B, C and E wolfs. Although the wolf of type B and the wolf of type C are not the highest ranking wolf species, they can take over the wolf of type a as a new leader when they lose leadership. The E kind of wolf is the lowest rank of wolf group and is responsible for balancing the relationship in the group.
The GWO algorithm treats each wolf as a potential solution, where A wolf is the first best solution, and wolfs B and C are the second best solution and the second best solution, respectively. The GWO algorithm is an iterative optimization process that continually updates the position of the wolf seed A, B, C. The wolf pack updates the distance and the position through the following formula to complete the search of the prey.
D＝|C×Xp(t)-X(t)|
X(t+1)＝Xp(t)+D
Wherein D is the distance between the gray wolf and the prey, and t is the iteration number; xpIs the location of the prey, and X is the location of the gray wolf, with its initial location coordinates defined as (c, g). The expression for coefficient vectors a and C is a ═ 2a × r2-α，C＝2r1. When | A |>1 denotes global search, i.e. the wolf population expands the search scope to more easily find the prey. In contrast, | a | < 1, which represents local search, the wolf will contract the enclosure to search for prey.
Figure BDA0003323704000000071
As the number of iterations increases, the convergence factor a decreases linearly from 2 to 0, and tmax is the maximum number of iterations, r1And r2Respectively ofIn [0, 1]]And randomly taking values in the interval.
When the wolf group judges the position of the prey, the head wolf A, the leading wolf B and the wolf C surround the prey, because the wolf A, the wolf B and the wolf C are the closest to the prey, the positions of the three wolfs are gradually close to the prey, and they are described as follows:
Da＝|C1×Xa(t)-X(t)|
Db＝|C2×Xb(t)-X(t)|
Dc＝|C3×Xc(t)-X(t)|
wherein, XaRepresents the current position of A wolf, XbRepresenting the current position of B wolf, XcRepresenting the current location of the C wolf. C1,C2,C3Is a random variable. X (t) is the current position of the wolf pack. The step sizes of wolf E to wolf A, B, C are marked as X respectively1X2 and X3As follows:
X1＝|C3×Xa-A1Da|
X2＝|C2×Xb-A2Db|
X3＝|C3×Xc-A2Dc|
the final position of wolf E is defined as follows:
Figure BDA0003323704000000072
in the hunting of the wolf group, the wolf A, wolf B and wolf C have different degrees of adaptation to the prey. And calculating different fitness values to obtain a first optimal solution, a second optimal solution and a suboptimal solution, and keeping the current position information of the solutions. Meanwhile, the wolf pack judges the moving direction according to the three sets of position information to approach to the prey to finish hunting. The location of the gray wolf is then updated again until the optimal solution is obtained. The position coordinate value corresponding to the optimal solution is defined as (bestc, bestg).
The specific process of step 4 is as follows:
step 4.1: initializing an NLM algorithm, and setting a bandwidth parameter lambda, a similar frame radius P and a search frame radius M; selecting a fixed search box radius M for calculating the speed;
step 4.2: initializing a gray wolf population, a convergence factor a and a coefficient vector A, C in the gray wolf algorithm;
Figure BDA0003323704000000081
A＝2a×r2-α
C＝2r1
where t is the number of iterations, tmax is the maximum number of iterations, r1And r2Is taken as [0, 1]]Random numbers within the interval;
step 4.3: taking the kurtosis value calculation function as a fitness function of a grey wolf algorithm, calculating the fitness value of each grey wolf, and simultaneously storing the parameters of the first three grey wolfs with higher fitness;
step 4.4: updating the position of each gray wolf;
step 4.5: updating a, A and C;
step 4.6: calculating the fitness values of all the gray wolves, and updating the fitness and the position of the first three gray wolves with the best fitness;
step 4.7: judging whether the maximum iteration times is reached, and if so, ending the optimization iteration process; and if not, returning to execute the step 4.4 to the step 4.7 to continue the optimization iteration process until the iteration is finished and outputting the optimal bandwidth parameter lambda and the similar box radius P.
And 5: carrying out noise reduction processing on the signal reconstructed in the step 3 by adopting an optimized GWO-NLM method, and outputting an optimal noise reduction signal;
assume that the noisy signal y consists of the true signal u plus the noisy signal n: y is u + n;
NLM denoised value at sample i
Figure BDA0003323704000000082
Is derived from a weighted average of all samples of its neighborhood ω.
Figure BDA0003323704000000083
Wherein the content of the first and second substances,
Figure BDA0003323704000000084
is a normalization constant. Weight ω (i, j) satisfying 0<ω(i,j)<1 and ∑ andjω (i, j) ═ 1, which compares the neighborhood around samples i and j. If the neighborhood of sample i is similar to the neighborhood of sample j, then take the large value, otherwise take the small value. The weight ω (i, j) is calculated by:
Figure BDA0003323704000000085
where δ is a bandwidth parameter and Δ represents a local sample patch.
Figure BDA0003323704000000086
Represents the sum of the squares of the point-to-point differences for the samples in the patch centered at j.
In the above formula, the special case of i ═ j is not considered, since the weight ω (i, j) can be much larger than the weight of each other sample. Since each neighborhood is similar to itself (ω (i, j) ═ 1), to prevent sample i from weighting itself too high, let ω (i, j) equal the maximum weight of the other samples.
Figure BDA0003323704000000091
The larger the radius of the search box is, the better the denoising effect is, but the longer the algorithm running time can be caused due to the overlarge radius of the search box, so the algorithm running time and the optimization algorithm time are balanced, the radius of the fixed search box is selected, the denoising effect of the filter is improved by optimizing bandwidth parameters and the radius of a similar frame, and the bandwidth parameters play a decisive role in the denoising effect.
Step 6: carrying out secondary noise reduction on the output optimal noise-reduced signal through SG smooth filtering, carrying out envelope spectrum analysis on the signal subjected to secondary noise reduction, and extracting fault characteristic frequency;
and 7: and (4) calculating theoretical fault characteristic frequency of the bearing by combining related parameters of the bearing, comparing the fault characteristic frequency extracted in the step (6) with the theoretical fault characteristic frequency, and judging the fault type of the rolling bearing.
In this embodiment, the denoising envelope spectrum for simulating the bearing signal is shown in fig. 2, and it can be seen that the envelope spectrum of the denoised signal highlights the fault characteristic frequency and the frequency doubling component.
The whole of the denoising envelope spectrogram of the experimental bearing inner ring fault is shown in fig. 3, and it can be seen that the signal envelope spectrogram after the noise reduction by GWO-NLM-SG can extract quite abundant information, and the fault frequency fi161.5Hz, double frequency 2fi323.4Hz triple frequency 3fi484.9Hz quadruple frequency 4fi646.7Hz, frequency quintupling 5fi808.2Hz, six times frequency 6fi970.1Hz, seven times frequency 7fi1132.0 Hz. Frequency conversion FrIs 30.0Hz and double frequency 2Fr59.7Hz, triple frequency 3Fr89.7Hz quadruple frequency 4Fr119.8Hz, frequency quintupling 5Fr149.4Hz, six times frequency 6Fr179.4 Hz. In the interval of 0-1000Hz, the peak value of the fault frequency is obviously seen to be the largest, the amplitude of the corresponding frequency multiplication is firstly reduced, the quadruple frequency is lowest, then the frequency multiplication starts to be increased, and the frequency multiplication starts to be reduced after the frequency multiplication is five times.
The local details of the denoising envelope spectrogram of the experimental bearing inner ring fault are shown in fig. 4, so that not only the rotating frequency, the fault frequency and the corresponding frequency multiplication can be extracted, but also the modulation frequency rich in the center of the frequency multiplication and the fault frequency can be extracted. In which significant frequency conversion and frequency multiplication nF occursr(where n is 1,2, 3, 4, 5, 6, 7), fault frequency and its multiple mfi(where m is 1, 2), and a modulation frequency f centered on the fault frequency and the frequency multiplicationi-aFr(wherein a is 1,2, 3, 4, 5), fi+bFr(wherein b is 1,2, 3, 4, 5, 6), 2fi-cFr(wherein c is 1,2, 3, 4, 5, 6, 7, 8), 2fi+dFr(whereind＝1，2)、3fi-eFr(wherein e is 3, 4, 5, 6, 7, 8, 9). Through comprehensive analysis, the fault can be clearly analyzed as the fault of the inner ring of the rolling bearing.
The whole of the denoising envelope spectrogram of the outer ring fault of the experimental bearing is shown in fig. 5, and it can be seen that the method effectively extracts that the frequency conversion, the frequency doubling of 2, the frequency doubling of 3, the frequency doubling of 4, the frequency doubling of 5, the inner ring fault frequency and the frequency doubling of 2 to 11 are all clearly extracted, the trend of the amplitude of the outer ring fault frequency is that the amplitude is firstly reduced, the local minimum value is reached at the frequency tripling position of the fault frequency, then the frequency doubling amplitude is gradually increased corresponding to the local minimum value, and the local maximum value is reached at the frequency quintupling position of the fault frequency and is gradually reduced. In summary, it can be determined that the outer ring of the rolling bearing is malfunctioning.
The local details of the denoising envelope spectrogram of the experimental bearing outer ring fault are shown in fig. 6, and it can be obviously seen that the denoised signal contains abundant fault information which is extracted.

Claims (4)

1. A rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM is characterized by comprising the following steps:
step 1: collecting vibration signals of a rolling bearing by adopting a vibration sensor, and simultaneously measuring relevant parameters of the bearing;
step 2: decomposing the acquired vibration signal of the rolling bearing by adopting a CEEMDAN method to decompose different IMF components;
and step 3: screening the decomposed IMF components by combining the kurtosis value, the correlation coefficient and the energy ratio, and performing linear reconstruction on the screened components;
and 4, step 4: optimizing the bandwidth parameter and the similar frame radius of the NLM by adopting a wolf algorithm GWO;
and 5: carrying out noise reduction processing on the signal reconstructed in the step 3 by adopting an optimized GWO-NLM method, and outputting an optimal noise reduction signal;
step 6: carrying out secondary noise reduction on the output optimal noise-reduced signal through SG smooth filtering, carrying out envelope spectrum analysis on the signal subjected to secondary noise reduction, and extracting fault characteristic frequency;
and 7: and (4) calculating theoretical fault characteristic frequency of the bearing by combining related parameters of the bearing, comparing the fault characteristic frequency extracted in the step (6) with the theoretical fault characteristic frequency, and judging the fault type of the rolling bearing.
2. The CEEMDAN and GWO-NLM-based rolling bearing fault diagnosis method according to claim 1, wherein: relevant parameters of the bearing include: contact angle alpha, ball number Z, ball diameter D and pitch diameter D.
3. The CEEMDAN and GWO-NLM-based rolling bearing fault diagnosis method according to claim 1, wherein: the process of screening the decomposed IMF components by combining the kurtosis value, the correlation coefficient and the energy ratio is as follows:
step 3.1: the kurtosis value K for each IMF component is calculated as follows:
Figure FDA0003323703990000011
wherein x is the amplitude of the IMF component, and u is the average value of the IMF component amplitudes; σ is the standard deviation of the IMF component amplitude;
step 3.2: calculating a Pearson correlation coefficient r of each IMF component and the original signal;
step 3.3: the energy ratio coefficient of each IMF component is calculated as follows:
Figure FDA0003323703990000012
wherein E isxTotal energy that is the IMF component; eIMF(i) Is the energy possessed by the ith IMF component; epsilon is an energy ratio coefficient;
step 3.4: and carrying out weighted summation on the calculated kurtosis value K, the correlation coefficient r and the energy ratio coefficient epsilon to obtain a comprehensive screening index Kr epsilon value as follows:
Krε＝a1K+a2r+a3ε
wherein, a1、a2、a3Weighted values of kurtosis values, correlation coefficients and energy ratios of the IMF components are respectively;
step 3.5: and sorting the comprehensive screening index Kr epsilon values of all IMF components from large to small, and screening the largest first n IMF components.
4. The CEEMDAN and GWO-NLM-based rolling bearing fault diagnosis method according to claim 1, wherein: the process of the step 4 is as follows:
step 4.1: initializing an NLM algorithm, and setting a bandwidth parameter lambda, a similar frame radius P and a search frame radius M; selecting a fixed search box radius M for calculating the speed;
step 4.2: initializing a gray wolf population, a convergence factor a and a coefficient vector A, C in the gray wolf algorithm;
Figure FDA0003323703990000021
A＝2a×r2-α
C＝2r1
where t is the number of iterations, tmax is the maximum number of iterations, r1And r2Is taken as [0, 1]]Random numbers within the interval;
step 4.3: taking the kurtosis value calculation function as a fitness function of a grey wolf algorithm, calculating the fitness value of each grey wolf, and simultaneously storing the parameters of the first three grey wolfs with higher fitness;
step 4.4: updating the position of each gray wolf;
step 4.5: updating a, A and C;
step 4.6: calculating the fitness values of all the gray wolves, and updating the fitness and the position of the first three gray wolves with the best fitness;
step 4.7: judging whether the maximum iteration times is reached, and if so, ending the optimization iteration process; and if not, returning to execute the step 4.4 to the step 4.7 to continue the optimization iteration process until the iteration is finished and outputting the optimal bandwidth parameter lambda and the similar box radius P.
CN202111255452.0A 2021-10-27 2021-10-27 Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM Active CN113776837B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202111255452.0A CN113776837B (en) 2021-10-27 2021-10-27 Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202111255452.0A CN113776837B (en) 2021-10-27 2021-10-27 Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM

Publications (2)

Publication Number Publication Date
CN113776837A true CN113776837A (en) 2021-12-10
CN113776837B CN113776837B (en) 2023-08-22

Family

ID=78873462

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111255452.0A Active CN113776837B (en) 2021-10-27 2021-10-27 Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM

Country Status (1)

Country Link
CN (1) CN113776837B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116337446A (en) * 2023-04-13 2023-06-27 中国航发沈阳发动机研究所 Rolling bearing fault diagnosis method based on LMD decomposition and GWO-PNN

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120220486A1 (en) * 2009-11-01 2012-08-30 Caerus Molecular Diagnostics Incorporated Methods and apparatus for binding assays
CN109632311A (en) * 2019-01-21 2019-04-16 北京化工大学 A kind of adaptive acoustical signal Method for Bearing Fault Diagnosis

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120220486A1 (en) * 2009-11-01 2012-08-30 Caerus Molecular Diagnostics Incorporated Methods and apparatus for binding assays
CN109632311A (en) * 2019-01-21 2019-04-16 北京化工大学 A kind of adaptive acoustical signal Method for Bearing Fault Diagnosis

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
吴涛;熊新;吴建德;马军;: "基于QH-ITD和AMCKD的滚动轴承故障诊断研究", 电子测量与仪器学报, no. 04, pages 84 - 94 *
姚芳;姜涛;刘明宇;董超群;郑帅;: "基于GWO-ELM的逆变器开路故障诊断", 电源学报, no. 01, pages 49 - 57 *
黄海松;范青松;魏建安;黄东;: "基于CEEMDAN-IGWO-SVM的轴承故障诊断研究", 组合机床与自动化加工技术, no. 03, pages 27 - 30 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116337446A (en) * 2023-04-13 2023-06-27 中国航发沈阳发动机研究所 Rolling bearing fault diagnosis method based on LMD decomposition and GWO-PNN

Also Published As

Publication number Publication date
CN113776837B (en) 2023-08-22

Similar Documents

Publication Publication Date Title
CN109841226B (en) Single-channel real-time noise reduction method based on convolution recurrent neural network
CN110595780B (en) Bearing fault identification method based on vibration gray level image and convolution neural network
CN109827777B (en) Rolling bearing fault prediction method based on partial least square method extreme learning machine
CN113176092B (en) Motor bearing fault diagnosis method based on data fusion and improved experience wavelet transform
CN109404285B (en) Improved shuffled frog leaping algorithm enhanced self-adaptive band-pass filtering method for screw compressor fault diagnosis
CN111504635B (en) Planetary gear fault diagnosis method based on differential evolution probability neural network
CN109946080B (en) Mechanical equipment health state identification method based on embedded circulation network
CN113569990B (en) Strong noise interference environment-oriented performance equipment fault diagnosis model construction method
CN113607415A (en) Bearing fault diagnosis method based on short-time stochastic resonance under variable rotating speed
CN114169377A (en) G-MSCNN-based fault diagnosis method for rolling bearing in noisy environment
CN113642508A (en) Bearing fault diagnosis method based on parameter self-adaptive VMD and optimized SVM
CN112345252A (en) Rolling bearing fault diagnosis method based on EEMD and improved GSA-SOM neural network
CN113776837A (en) Rolling bearing fault diagnosis method based on CEEMDAN and GWO-NLM
CN111881848A (en) Motor fault signal extraction method based on variational modal decomposition and improved particle swarm
CN111811819A (en) Bearing fault diagnosis method and device based on machine learning
CN115062665A (en) Rolling bearing early fault diagnosis method based on self-adaptive variational modal decomposition
CN117030268B (en) Rolling bearing fault diagnosis method
CN114964783B (en) Gearbox fault detection model based on VMD-SSA-LSSVM
CN116973102A (en) Gear box noise fault diagnosis method
CN114705431A (en) Rolling bearing fault diagnosis method based on multi-parameter screening criterion and GWO-PNN
CN115470630A (en) VMD-SSA-LSTM-based rolling bearing residual service life prediction method
CN115146687A (en) Fault feature extraction method based on second-order variable-scale parameter self-matching stochastic resonance
CN114897016A (en) Fan bearing fault intelligent diagnosis method based on multi-source frequency spectrum characteristics
CN113657268B (en) Signal automatic decomposition method applied to wind turbine generator gearbox fault diagnosis
CN113869358A (en) Bearing fault diagnosis method based on cyclic correlation entropy and one-dimensional shallow convolution neural network

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
TA01 Transfer of patent application right

Effective date of registration: 20220815

Address after: 110136, Liaoning, Shenyang moral and Economic Development Zone, No. 37 South Avenue moral

Applicant after: SHENYANG AEROSPACE University

Applicant after: AECC SHENYANG LIMING AERO-ENGINE Co.,Ltd.

Address before: 110136, Liaoning, Shenyang moral and Economic Development Zone, No. 37 South Avenue moral

Applicant before: SHENYANG AEROSPACE University

TA01 Transfer of patent application right
GR01 Patent grant
GR01 Patent grant