CN106323635A - Rolling bearing fault on-line detection and state assessment method - Google Patents

Abstract

A rolling bearing fault on-line detection and state assessment method is disclosed. The method comprises the following steps: twelve dimensional dimensionless parameters are extracted; the twelve dimensional dimensionless parameters comprise six dimensional time domain statistical parameters, three dimensional frequency domain statistical parameters and three dimensional dimensionless parameters in a small wave envelope spectrum; standardized reconstruction characteristic vectors can be obtained; whether a rolling bearing malfunctions is determined, and a state of the rolling bearing is assessed. Via the rolling bearing fault on-line detection and state assessment method, the twelve dimensional dimensionless parameters which can be used for effectively representing the state of the rolling bearing can be automatically extracted, the twelve dimensional dimensionless parameters are subjected to decorrelation and standardization operation, standardized reconstruction characteristic vectors that are distributed to form a hypersphere with an original point being a sphere center, and fault detection and state assessment of the rolling bearing can be realized via 2-norms of the standardized reconstruction characteristic vectors; difficult problems of long on line training time, low efficiency, and hard-to-obtain fault samples and the like of a rolling bearing state assessing model can be solved.

Description

Rolling bearing fault online detection and state evaluation method

Technical Field

The invention belongs to the field of rolling bearing fault intelligent detection and state evaluation methods, and particularly relates to a rolling bearing fault online detection and state evaluation method.

Background

Because rolling bearing fault samples are often difficult to obtain, bearing fault types are complex and various, and several types of faults may exist simultaneously, the state evaluation of the rolling bearing often faces a data field description problem, namely, the adopted characteristic evaluation method should be suitable for the situation that only normal samples exist.

Since on-line monitoring of rolling bearings can be essentially regarded as a description of the normal data domain boundaries of the bearings, it is necessary to study the spatial distribution of multidimensional feature vectors in order to better exploit a priori knowledge to establish a suitable model. Usually, the dimension of the feature is very high (more than 3 dimensions), and the distribution of the feature cannot be displayed visually, but according to the popularization of the low-dimension situation, the distribution of the feature vector is easily imagined to be in a super-ellipsoid shape, and the direction and the length of the main axis of the super-ellipsoid may be changed according to the selected feature. Such boundaries are complex, if the complex boundaries are described, the algorithm needs to have strong nonlinear capability, and methods such as support vector data description, self-organizing neural network and Gaussian mixture model have such capability, and have been used for solving the problems and achieving good effects on experiments.

However, due to the strong non-linear capability of the above method, the training inevitably faces a problem: the computational complexity is large and thus dynamic training is difficult to achieve. However, taking an aircraft engine rolling bearing as an example, the state evaluation model is expected to be embedded in an aircraft engine control system, and the engine control system needs to carefully optimize and allocate computing resources to ensure the safety of important tasks (such as engine control, etc.), and the limited computing resources and the complexity of the model form a pair of contradictions. Therefore, an evaluation model which is effective and low in calculation complexity is constructed, dynamic training and real-time evaluation of the rolling bearing can be achieved, and the method has important engineering significance.

201310015619 discloses a rolling bearing fault detection method based on vibration detection, which comprises the steps of firstly carrying out 3-layer wavelet packet decomposition on rolling bearing data acquired by an acceleration sensor, then solving the energy of a third-layer wavelet packet coefficient reconstruction signal, and then selecting a frequency band with concentrated energy according to the change of the energy value of each frequency band of the third layer to reconstruct the approximate estimation of an original signal; and further analyzing the reconstructed signal by utilizing a cepstrum, and finally comparing the reconstructed signal with the fault characteristic frequency and the side frequency characteristic which are calculated theoretically so as to diagnose the fault. The method needs fault data of the rolling bearing, but no matter actual detection or theoretical calculation, the fault data is often difficult to obtain, the diagnosis process needs human participation, and the method is not beneficial to automatic online training and detection.

Disclosure of Invention

The invention aims to provide a rolling bearing fault detection and state evaluation method for online training, which can solve the difficult problems of long time consumption and low efficiency of online training of a rolling bearing state evaluation model.

In order to realize the purpose, the invention is realized by the following technical scheme:

a rolling bearing fault online detection and state evaluation method comprises the following steps:

s1: extracting 12-dimensional dimensionless parameters, wherein the 12-dimensional dimensionless parameters comprise: the method comprises the following steps of (1) carrying out 6-dimensional time domain statistical parameters, 3-dimensional frequency domain statistical parameters and 3-dimensional dimensionless parameters in a wavelet envelope spectrum;

s2: obtaining a standardized reconstruction feature vector; and

s3: determining whether the rolling bearing is in failure; wherein,

the step S1 specifically includes:

s1-1: collecting vibration acceleration signals, and storing the collected vibration acceleration signals in a segmented manner to obtain a plurality of samples;

s1-2: extracting 6-dimensional time domain statistical parameters from the sample through the time domain statistical parameters, wherein the 6-dimensional time domain statistical parameters comprise: form factor T_SIPeak index T_CIPulse index T_MIMargin index T_CLIKurtosis T_KUAnd degree of distortion T_SK；

S1-3: extracting 3-dimensional frequency domain statistical parameters from the sample through frequency domain statistical parameters, wherein the 3-dimensional frequency domain statistical parameters comprise gravity center frequency F_FCMean square frequency F_MSFSum frequency variance F_VF(ii) a And

s1-4: obtaining a wavelet envelope spectrum of the sample through wavelet transformation, and extracting 3-dimensional dimensionless parameters in the wavelet envelope spectrum, wherein the 3-dimensional dimensionless parameters in the wavelet envelope spectrum comprise: inner ring failure frequency f_ICorresponding characteristic W_BPFIOuter ring fault frequency f_OCorresponding characteristic W_BPFOAnd ball failure frequency f_BCorresponding characteristic W_BSF；

The step S2 specifically includes:

s2-1: based on a minimum reconstruction error criterion, applying a constraint condition to the 12-dimensional dimensionless parameters to obtain 12-dimensional decorrelation reconstruction characteristic vectors;

s2-2: normalizing the 12-dimensional decorrelated reconstructed feature vector to obtain a normalized reconstructed feature vector with a 12-dimensional sample mean value of 0 and a sample standard deviation of 1; and is

The step S3 specifically includes:

s3-1: training to obtain a sample mean value and a sample standard deviation of the 2 norm of the standardized reconstruction characteristic vector under the normal operation state of the rolling bearing, and setting a threshold value of the rolling bearing state abnormity;

s3-2: testing, namely taking the Euclidean distance from the characteristic vector of the sample in the unknown state of the rolling bearing to the characteristic vector of the sample mean value in the normal running state of the rolling bearing as a degradation index, comparing the degradation index with the threshold value, and judging that the rolling bearing is in a fault state when the degradation index is greater than the threshold value; and when the degradation index is less than or equal to the threshold value, determining that the rolling bearing is in a normal operation state.

Further, the greater the deterioration index is, the greater the degree of deterioration of the rolling bearing is determined.

Further, the specific operation of extracting the 6-dimensional time domain statistical parameter in step S1-2 is as follows:

form factor T_SIIs given by the formula (1), peak index T_CIGiven by equation (2), the pulse index T_MIIs given by the formula (3), margin index T_CLIGiven by the formula (4), kurtosis T_KUGiven by formula (5), the skew T_SKIs given by the formula (6),

$T_{SI}=\frac{\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left(y_i^2\right)}}{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left|y_i\right|}---\left(1\right)$

$T_{CI}=\frac{\underset{i=1}{\overset{10}{\mathrm{\Σ}}}{y}_{pi}}{\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\left(y_i\right)}^2}}---\left(2\right)$

$T_{MI}=\frac{\underset{i=1}{\overset{10}{\mathrm{\Σ}}}{y}_{pi}}{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left|y_i\right|}---\left(3\right)$

$T_{CLI}=\frac{\underset{i=1}{\overset{10}{\mathrm{\Σ}}}{y}_{pi}}{{\[\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\sqrt{\left|y_i\right|}\]}^2}---\left(4\right)$

$T_{KU}=\frac{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^4}{{\left(\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^2\right)}^2}---\left(5\right)$

$T_{SK}=\frac{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^3}{{\left(\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^2\right)}^{\frac{3}{2}}}---\left(6\right)$

wherein, y_iIs data in the vibro-acceleration signal; y is_piIs the maximum value of the absolute value of each data segment after the data is divided into 10 data segments.

Further, the specific operation of extracting the 3-dimensional time domain statistical parameter in step S1-3 is as follows:

center of gravity frequency F_FCMean square frequency F given by equation (7)_MSFGiven by equation (8), the frequency variance F_VFIs given by the formula (9),

$F_{FC}=\frac{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}{f}_{i}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}S\left(f_i\right)}---\left(7\right)$

$F_{MSF}=\frac{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}{f}_i^{2}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}S\left(f_i\right)}---\left(8\right)$

$F_{VF}=\frac{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}{(f_i-F_{FC})}^{2}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}S\left(f_i\right)}---\left(9\right)$

wherein, S (f)_i) As a function of the vibration acceleration signal frequency spectrum.

Further, the specific operation of step S1-4 includes the following three steps:

s1-4-1, calculating the inner ring fault frequency f according to the formulas (10) - (12)_IOuter ring fault frequency f_OAnd ball failure frequency f_B，

$f_I=\frac{1}{2}Z(1+\frac{d}{D}cos\mathrm{\α})f_r---\left(10\right)$

$f_O=\frac{1}{2}Z(1-\frac{d}{D}cos\mathrm{\α})f_r---\left(11\right)$

$f_B=\frac{D}{2d}\[1-{\left(\frac{d}{D}\right)}^2{\cos }^2\mathrm{\α}\]f_r---\left(12\right)$

Wherein D represents the diameter of the rolling elements, D represents the pitch diameter of the bearing, Z represents the number of the rolling elements, α represents the contact angle, f_rRepresenting the relative rotation frequency of the inner ring and the outer ring of the rolling bearing;

s1-4-2, carrying out wavelet decomposition on the vibration acceleration signal by adopting a db8 wavelet to obtain 6 signals, wherein the 6 signals comprise 5 detail signals d1, d2, d3, d4, d5 and 1 approximate signal a 5;

s1-4-3, extracting 3-dimensional dimensionless parameters in the wavelet envelope spectrum: setting the frequency spectrum in wavelet packet envelope spectrum, the fault characteristic frequency and the characteristic spectrum peak near each order of frequency multiplication, and setting the envelope spectrum analysis bandwidth as f_eThe envelope spectrum is W (f), and the number of spectral lines is N_eThen S is_eaIs composed of

$S_{ea}=\frac{1}{N_e}\underset{i=0}{\overset{N_e}{\mathrm{\Σ}}}W\left(f_i\right)---\left(13\right)$

Then order S_edSetting the number of spectral lines of fault frequency in the envelope spectrum as n for the average value of spectral lines at each order of multiple frequency of fault characteristic frequency in the envelope spectrum_eThen, then

$S_{ed}=\frac{1}{n_e}\underset{i=0}{\overset{n_e}{\mathrm{\Σ}}}W\left({\mathrm{if}}_d\right)---\left(14\right)$

Constructing a dimensionless parameter:

${\mathrm{\ΔS}}_e=\frac{S_{ed}}{S_{ea}}---\left(15\right)$

respectively carrying out envelope spectrum analysis on the 6 signals in the step S1-4-2 to obtain the fault frequency f of the inner ring_ICorresponding characteristic W_BPFIOuter ring fault frequency f_OCorresponding characteristic W_BPFOFrequency of ball failure f_BCorresponding characteristic W_BSF。

Further, the specific operation of step S2-1 includes the following three steps:

s2-1-1: standardizing the 12-dimensional dimensionless parameters to obtain x_iThen there is ∑_ix_i＝0；

S2-1-2: assuming a new coordinate system { w ] obtained after projective transformation₁,w₂,…,w_dWhere d is a feature dimension; the constraint imposed is that wi is a normal orthogonal basis vector, satisfying | | w_i||₂＝1，w_i ^Tw_j0(i ≠ j), the optimization objective function is

$\begin{array}{cc}\underset{W}{\min }& -tr\left(W^T{\mathrm{XX}}^{T}W\right)\\ s.t.& W^{T}W=I\end{array}---\left(16\right)$

S2-1-3: solving the equation (16) to obtain a projection matrix W, and calculating a decorrelated reconstructed eigenvector z according to the equation (17)_i

z_i＝(z_i1,z_i2,…,z_id)，z_id＝w_j ^Tx_i(17)。

Further, the reconstructed feature vector z is normalized in step S2-2_i ^*Is given by the formula (18)

$z_i^{*}=(z_{i1}^{*},z_{i2}^{*},\mathrm{...},z_{id}^{*}),z_{ij}^{*}=\frac{z_{ij}-{\widehat{\mathrm{\μ}}}_j}{{\widehat{\mathrm{\σ}}}_j}---\left(18\right)$

Wherein

${\widehat{\mathrm{\μ}}}_j=\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_{ij},{\widehat{\mathrm{\σ}}}_j^2=\frac{1}{n-1}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{(z_{ij}-{\widehat{\mathrm{\μ}}}_j)}^2---\left(19\right).$

Further, the threshold value in step S3-1 is given by equation (20)

$D_{max}={\widehat{\mathrm{\μ}}}_D+3{\widehat{\mathrm{\σ}}}_D---\left(20\right)$

Wherein,andare each D_iSample mean and sample standard deviation.

Compared with the prior art, the invention has the following beneficial effects:

the invention provides an online detection method for rolling bearing faults, which can automatically extract multidimensional dimensionless parameters for effectively representing the bearing state, obtain a super-spherically distributed standardized reconstruction characteristic vector by performing decorrelation and standardization processing on original characteristics, and realize the fault detection and state evaluation of a rolling bearing by only calculating 2 norms of the standardized reconstruction characteristic vector, and has the following remarkable advantages different from the traditional method:

1) in the aspect of feature fusion, the invention fully combines the multidimensional information of time domain features, frequency domain features and wavelet envelope spectrum features and expresses the multidimensional information through one-dimensional deterioration indexes, thereby effectively solving the problem of fault diagnosis caused by inconsistent sensitivity of different features to different faults in the prior art;

2) in the aspect of model complexity, the invention improves the spatial distribution of the characteristic vector through linear projection transformation in a new way, and therefore greatly simplifies the complexity of a classifier model, and meanwhile, the invention can well realize the on-line training and dynamic updating of the rolling bearing fault detection and state evaluation model because parameters do not need to be adjusted and iterative calculation is not needed;

3) in the aspect of data sources of model training, the invention only needs normal operation data of the rolling bearing to be used for model training, thereby avoiding the defect that the fault data of the rolling bearing is often difficult to obtain in the traditional method training;

4) in a whole view, the method fully combines the multi-dimensional characteristic information, improves the spatial distribution of characteristic vectors through characteristic transformation, can further improve the fault recognition rate while effectively simplifying the complexity of a model, and is closer to the engineering requirement.

Drawings

FIG. 1 is a flow chart of a method of the present invention;

FIG. 2(a) is a two-dimensional spatial distribution of features that is not relevant;

FIG. 2(b) is a related two-dimensional feature spatial distribution;

FIG. 3 is a schematic diagram illustrating the error caused by the inconsistency between the data distribution and the description;

FIG. 4 is a schematic diagram of distance discrimination analysis;

FIG. 5 is a schematic view of a rolling bearing state evaluation;

FIG. 6(a) is a scatter diagram of the different fault detection results of the bearing according to the method of the present invention;

FIG. 6(b) is a scatter diagram of different fault detection results of a bearing by a support vector data description method; and

fig. 6(c) is a scatter diagram of the different fault detection results of the bearing by the self-organizing neural network method.

Detailed Description

The contents of the rolling bearing fault on-line detection and state evaluation method of the present invention will be further described with reference to the accompanying drawings.

As shown in fig. 1, the invention discloses an online detection and state evaluation method for rolling bearing faults, which is implemented by the following steps:

s1: the 12-dimensional dimensionless parameter extraction method specifically comprises the following four steps:

s1-1: collecting vibration acceleration signals, and storing the collected vibration acceleration signals in a segmented manner to obtain a plurality of samples;

s1-2: extracting 6-dimensional time domain statistical parameters including form factor T by calculating time domain statistical parameters_SIPeak index T_CIPulse index T_MIMargin index T_CLIDegree of kurtosis T_KUDegree of skewness T_SK；

S1-3: extracting 3-dimensional time domain statistical parameters including center-of-gravity frequency F by frequency domain statistical parameter calculation_FCMean square frequency F_MSFFrequency variance F_VF；

S1-4: obtaining a wavelet envelope spectrum through wavelet transformation, and extracting the fault frequency f of the inner ring corresponding to the wavelet envelope spectrum_ICorresponding characteristic W_BPFIOuter ring fault frequency f_OCorresponding characteristic W_BPFOFrequency of ball failure f_BCorresponding characteristic W_BSF；

S2: obtaining a standardized reconstruction feature vector, specifically comprising the following two steps:

s2-1: based on the minimum reconstruction error criterion, applying corresponding constraint conditions to obtain a decorrelated 12-dimensional reconstruction feature vector;

s2-2: carrying out standardization processing on the 12-dimensional reconstruction characteristic vector to obtain a standardized reconstruction characteristic vector with the average value of each characteristic sample being 0 and the standard deviation of the sample being 1;

s3: the rolling bearing fault detection and state evaluation based on distance discrimination analysis specifically comprises the following two steps:

s3-1: training, calculating to obtain a sample mean value and a sample standard deviation of a 2-norm of a sample characteristic vector of the rolling bearing in a normal running state, and defining a threshold value of rolling bearing state abnormity;

s3-2: and testing, namely calculating to obtain an Euclidean distance from the sample characteristic vector to the normal sample mean vector of the rolling bearing in an unknown state as a degradation index, comparing the degradation index with a threshold set in S3-1, judging whether the bearing fails or not, and evaluating the degradation degree of the bearing.

The specific operation of step S1-2 is:

form factor T_SIThe peak index T is given by equation (1)_CIGiven by equation (2), the pulse index T_MIGiven by equation (3), the margin index T_CLIGiven by the formula (4), kurtosis T_KUGiven by formula (5), the skew T_SKIs given by formula (6), wherein y_iThe vibration data in the vibration acceleration signal is original data; y is_piIs the maximum value of the absolute value of each data segment after dividing the original data into 10 data segments.

$T_{SI}=\frac{\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left(y_i^2\right)}}{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left|y_i\right|}---\left(1\right)$

$T_{CI}=\frac{\underset{i=1}{\overset{10}{\mathrm{\Σ}}}{y}_{pi}}{\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\left(y_i\right)}^2}}---\left(2\right)$

$T_{MI}=\frac{\underset{i=1}{\overset{10}{\mathrm{\Σ}}}{y}_{pi}}{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left|y_i\right|}---\left(3\right)$

$T_{CLI}=\frac{\underset{i=1}{\overset{10}{\mathrm{\Σ}}}{y}_{pi}}{{\[\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\sqrt{\left|y_i\right|}\]}^2}---\left(4\right)$

$T_{KU}=\frac{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^4}{{\left(\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^2\right)}^2}---\left(5\right)$

$T_{SK}=\frac{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^3}{{\left(\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{y_i}^2\right)}^{\frac{3}{2}}}---\left(6\right)$

The specific operation of step S1-3 is:

center of gravity frequency F_FCMean square frequency F given by equation (7)_MSFGiven by equation (8), the frequency variance F_VFIs given by formula (9), wherein S (f)_i) As a function of the frequency spectrum of the vibration acceleration signal

$F_{FC}=\frac{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}{f}_{i}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}S\left(f_i\right)}---\left(7\right)$

$F_{MSF}=\frac{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}{f}_i^{2}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}S\left(f_i\right)}---\left(8\right)$

$F_{VF}=\frac{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}{(f_i-F_{FC})}^{2}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\mathrm{\Σ}}}S\left(f_i\right)}---\left(9\right)$

The specific operation of the step S1-3 comprises the following three steps:

s1-3-1) calculating fault characteristic frequency according to the formulas (10) - (12), wherein D represents the diameter of the rolling body, D represents the pitch diameter of the bearing, Z represents the number of the rolling bodies, α represents the contact angle, f_rIndicating the relative rotational frequency of the inner and outer rings of a rolling bearing

$f_I=\frac{1}{2}Z(1+\frac{d}{D}cos\mathrm{\α})f_r---\left(10\right)$

$f_O=\frac{1}{2}Z(1-\frac{d}{D}cos\mathrm{\α})f_r---\left(11\right)$

$f_B=\frac{D}{2d}\[1-{\left(\frac{d}{D}\right)}^2{\cos }^2\mathrm{\α}\]f_r---\left(12\right)$

S1-3-2) carrying out wavelet decomposition on the vibration acceleration signal by adopting a db8 wavelet to obtain 5 detail signals d1, d2, d3, d4, d5 and 1 approximate signal a 5;

s1-3-3) automatic extraction of wavelet envelope spectrum features: setting the frequency spectrum in wavelet packet envelope spectrum, the fault characteristic frequency and the characteristic spectrum peak near each order of frequency multiplication, and setting the envelope spectrum analysis bandwidth as f_eThe envelope spectrum is W (f), and the number of spectral lines is N_eThen S is_eaIs composed of

$S_{ea}=\frac{1}{N_e}\underset{i=0}{\overset{N_e}{\mathrm{\Σ}}}W\left(f_i\right)---\left(13\right)$

Then order S_edSetting the number of spectral lines of fault frequency in the envelope spectrum as n for the average value of spectral lines at each order of multiple frequency of fault characteristic frequency in the envelope spectrum_eThen, then

$S_{ed}=\frac{1}{n_e}\underset{i=0}{\overset{n_e}{\mathrm{\Σ}}}W\left({\mathrm{if}}_d\right)---\left(14\right)$

Constructing a dimensionless parameter:

${\mathrm{\ΔS}}_e=\frac{S_{ed}}{S_{ea}}---\left(15\right)$

respectively carrying out envelope spectrum analysis on the 6 signals in the step S1-3-2), and obtaining the corresponding inner ring fault frequency f through automatic calculation_ICorresponding characteristic W_BPFIOuter ring fault frequency f_OCorresponding characteristic W_BPFOFrequency of ball failure f_BCorresponding characteristic W_BSF。

After the extracted features are normalized, the sample mean μ of each feature is 0, and the sample standard deviation σ is 1, so that the distributions can be preferably compared in the same scale. It can be seen that the boundary of the feature distribution in fig. 2(a) can be described approximately by a circle (the circle in the figure has a radius of 2.5), but the feature distribution in fig. 2(b) appears to be significantly elliptical, with the major axis of the ellipse being at an angle of approximately 45 ° to the major axis of the coordinates. There are 66 cases in which any two of the 12-dimensional features (i.e., 12-dimensional dimensionless parameters) are combined, but neither distribution exceeds that shown in fig. 2. The distribution of fig. 2(a) can be better described by a circle, but when describing a distribution like fig. 2(b), two types of errors are inevitably introduced, as shown in fig. 3. Wherein the first type of error is caused by erroneously identifying the faulty sample as a normal sample, i.e., a false positive example; the second type of error is caused by erroneously determining a normal sample as a faulty sample, i.e., false negative.

Further observing fig. 2(b), it is found that the distribution appears elliptical because there is a large correlation between features, i.e. when the center of gravity frequency is large, the frequency variance is also large, and the frequency variance formula includes the center of gravity frequency, which is the root cause of the existence of the correlation. More generally, inter-feature correlation can be measured by a correlation coefficient, ρ_xyIs as shown in formula (16):

${\mathrm{\ρ}}_{xy}=\frac{n\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x}_i{y}_i-\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x}_i\·\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i}{\sqrt{n\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x}_i^2-{\left(\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x}_i\right)}^2}\·\sqrt{n\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i^2-{\left(\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{y}_i\right)}^2}}---\left(16\right)$

in the formula, x_i,y_iCorresponding to the characteristic values of different characteristics of the ith sample, and n is the sample capacity. The correlation coefficients between corresponding features in fig. 2 have been labeled in the figure.

The conclusions of the above analysis can be generalized to a higher dimensional space: if the described features are all not related, the distribution of the feature vectors in the multi-dimensional space is in a hypersphere shape and can be simply described by a hypersphere; on the contrary, the distribution of the feature vectors is in a super-ellipsoid shape and is not easy to describe.

The specific operation of the step S2-1 comprises the following three steps:

s2-1-1) assuming that the normalization processing is carried out on the 12-dimensional dimensionless parameters, namely the original characteristic vectors, to obtain x_iThen there is ∑_ix_i＝0；

S2-1-2) assuming that a new coordinate system obtained after projective transformation is { w₁,w₂,…,w_d}. Wherein d is a feature dimension; the constraint imposed is that wi is a normal orthogonal basis vector, satisfying | | w_i||₂＝1，w_i ^Tw_jLet 0(i ≠ j), set sample point x_iThe projection in the new coordinate system is z_i＝(z_i1,z_i2,…,z_id) Wherein z is_id＝w_j ^Tx_iIs x_iCoordinates of the j-th dimension in the new coordinate system. Based on z_iTo reconstruct x_iThen, there are:

${\hat{x}}_i=\underset{j=1}{\overset{d}{\mathrm{\Σ}}}{z}_{ij}{w}_j---\left(17\right)$

considering the entire training set, the original feature vector x_iAnd reconstructing the feature vectorThe distance between them is:

$\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\left|\right|\underset{j=1}{\overset{d}{\mathrm{\Σ}}}{z}_{ij}{w}_j-x_i|{|}_2^2=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_i^T{z}_i-2\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_i^T{W}^T{x}_i+const\∝-tr\left(W^T\right(\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x}_i{x}_i^T\left)W\right)---\left(18\right)$

equation (18) should be minimized according to the minimum reconstruction error criterion. Considering w_jIs a standard orthogonal basis, sigma_ix_ix_i ^TIs a covariance matrix, and therefore, the optimization objective function is

$\begin{array}{cc}\underset{W}{\min }& -tr\left(W^T{\mathrm{XX}}^{T}W\right)\\ s.t.& W^{T}W=I\end{array}---\left(19\right)$

S2-1-3) solving the formula (19) to obtain a projection matrix W, and calculating according to the formula (20) to obtain a decorrelated reconstructed feature vector z_i

z_i＝(z_i1,z_i2,…,z_id)，z_id＝w_j ^Tx_i(20)

Reconstructed feature vector z obtained by decorrelation_iThe distribution of (2) is still in a shape of a super ellipsoid, but the direction of the main axis of the super ellipsoid is already approximately parallel to the coordinate axis. Therefore, by further normalizing the reconstructed feature vector, the feature vector distribution can be converted into a hypersphere distribution having the origin of coordinates as the hypersphere center.

The normalized reconstructed feature vector z in step S2-2_iIs given by formula (21)

$z_i^{*}=(z_{i1}^{*},z_{i2}^{*},\mathrm{...},z_{id}^{*}),z_{ij}^{*}=\frac{z_{ij}-{\widehat{\mathrm{\μ}}}_j}{{\widehat{\mathrm{\σ}}}_j}---\left(21\right)$

Wherein,

${\widehat{\mathrm{\μ}}}_j=\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_{ij},{\widehat{\mathrm{\σ}}}_j^2=\frac{1}{n-1}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{(z_{ij}-{\widehat{\mathrm{\μ}}}_j)}^2---\left(22\right)$

the distance discrimination analysis is to adopt a statistical analysis means to discriminate the attribution of the new sample according to the distance between each characteristic value of the new sample and the known class on the premise of determining the classification. The basic idea is to calculate the distance between an individual and each population and consider that the individual belongs to the population closest to the individual.

In a d-dimensional feature space, consider the use of a d-dimensional vector x₀Representing n sample points of a sample set, we want this vector x₀The smaller the sum of squared distances from each sample i (i ═ 1, …, n), the better, the square error criterion function J is defined₀(x₀) The following were used:

$J_0\left(x_0\right)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\left|\right|x_0-x_i|{|}^2---\left(23\right)$

then

$x_0=\widehat{\mathrm{\μ}}=\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x}_i---\left(24\right)$

This conclusion can be demonstrated as follows:

$\begin{array}{c}{J}_0\left(x_0\right)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\left|\right|(x_0-\widehat{\mathrm{\μ}})-(x_i-\widehat{\mathrm{\μ}})|{|}^2\\ =\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\left|\right|x_0-\widehat{\mathrm{\μ}}|{|}^2-2{(x_0-\widehat{\mathrm{\μ}})}^T\underset{i=1}{\overset{n}{\mathrm{\Σ}}}(x_i-\widehat{\mathrm{\μ}})+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}||x_i-\widehat{\mathrm{\μ}}|{|}^2\\ =\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\left|\right|x_0-\widehat{\mathrm{\μ}}|{|}^2+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}||x_i-\widehat{\mathrm{\μ}}|{|}^2\end{array}---\left(25\right)$

wherein the second term on the right of equation (25) is independent of x₀Therefore, the expression takes a minimum value under the condition of expression (24). Thus, the sample mean of the normal samples can be used as a zero-dimensional description of the normal samples.

In the case of only normal samples, the problem translates into: calculating the distance D from the feature vector to the normal sample mean vector according to a certain measurement criterion_iAnd is in contact with a set distance threshold D_maxComparing to judge the bearing state (sample attribution), namely for any sample, if D is satisfied_i≤D_maxThen the sample is classified as normal, otherwise as abnormal. As shown in FIG. 4, the boundary described by the distance discrimination is a radius D_maxIf Euclidean distance is adopted as a measurement criterion, the hypersphere of (1) is

$D_i=\left|\right|z_i^{*}|{|}_2---\left(26\right)$

At the same time, D_iIt is clear and can reflect the degree of deterioration of the rolling bearing, and can be used as an index for evaluating the health state of the rolling bearing, as shown in fig. 5. The more severe the bearing damage, the corresponding D of the sample_iThe larger will be. Further, setting up a plurality of thresholds by respective criteria may give warning limits, abnormal limits, etc. of the bearing.

The formulation of the threshold value in step S3-1 is given by equation (27)

$D_{max}={\widehat{\mathrm{\μ}}}_D+3{\widehat{\mathrm{\σ}}}_D---\left(27\right)$

Wherein,andare each D_iSample mean and sample standard deviation.

The aircraft engine rotor tester with the casing, which is developed by Shenyang engine design research institute, is adopted to carry out a fault simulation test so as to verify the effectiveness of the invention.

The tester is designed in a structure, firstly, the appearance of the tester is consistent with that of a casing of an engine core engine, the size of the tester is reduced to 1/3, and the internal structure is necessarily simplified: the core machine is simplified into a 0-2-0 supporting structure form, the multi-stage compressor is simplified into a single-stage disc structure, and the blades are simplified into inclined planes. Considering that on a real aero-engine, the sensors are difficult to be mounted on the bearing housing, the vibration acceleration sensors were arranged in the vertical and horizontal directions of the turbine casing in the test. In the test, vibration signals are collected through an NI USB9234 data collector, the model of the acceleration sensor is B & K4805, and the sampling frequency is 10.24 kHz.

The test object was a 6206 type ball bearing, and the basic geometrical parameters are shown in table 1. Grooves are respectively machined on an outer ring, an inner ring raceway and a rolling body of the rolling bearing by adopting wire cut electrical discharge machining to simulate the impact generated by bearing damage. The rated speed in the test is 1500 rpm.

TABLE 1 bearing geometry (unit: mm)

Taking a vibration acceleration signal measured at a measuring point above the casing as an example, 12 (dimensional) dimensionless parameters are extracted according to the step 1).

To verify the validity of the proposed method, the method of the invention was compared with a support vector data description and a self-organizing neural network. The evaluation models were trained using half of the normal samples, and the remaining samples were used as unknown samples for testing (110 samples in each state).

In the invention, D in formula (26)_iAs an index of the state deterioration of the rolling bearing, a scatter diagram was obtained as shown in fig. 6 (a); the support vector data description method is realized by an MATLAB toolbox LibSVM-3.18-SVDD (matrix support vector machine) expanded by a LibSVM software package, a kernel function takes a radial basis kernel, corresponding parameters are selected by cross validation, a decision value is taken as a bearing degradation index, and an obtained scatter diagram is shown in fig. 6 (b); the SOM method is implemented by an MATLAB-SOM toolbox, the number of neurons and the number of iterations are selected by cross validation, the minimum matching distance from the feature vector to the neuron is used as a bearing state evaluation index, and the obtained scatter diagram is shown in fig. 6 (c). The classification accuracy of each method for samples in different states of the bearing is shown in table 2.

TABLE 2 accuracy of different method classifications

It can be seen that the method of the present invention ensures a higher level in the classification accuracy and is significantly better than the results of support vector data description and self-organizing neural networks.

The method has the advantages of no need of parameter adjustment and small model complexity, and can better solve the problem of dynamic update of the evaluation model in online monitoring of the rolling bearing; meanwhile, the method fully integrates effective information in 12-dimensional characteristics of a time domain, a frequency domain and a wavelet envelope spectrum, so that the accuracy of fault detection is further improved.

Finally, it should be noted that: the above-mentioned embodiments are only used for illustrating the technical solution of the present invention, and not for limiting the same; although the present invention has been described in detail with reference to the foregoing embodiments, it will be understood by those of ordinary skill in the art that: the technical solutions described in the foregoing embodiments may still be modified, or some or all of the technical features may be equivalently replaced; and the modifications or the substitutions do not make the essence of the corresponding technical solutions depart from the scope of the technical solutions of the embodiments of the present invention.

Claims (8)

1. A rolling bearing fault on-line detection and state evaluation method is characterized in that: the method comprises the following steps:

s1: extracting 12-dimensional dimensionless parameters, wherein the 12-dimensional dimensionless parameters comprise: the method comprises the following steps of (1) carrying out 6-dimensional time domain statistical parameters, 3-dimensional frequency domain statistical parameters and 3-dimensional dimensionless parameters in a wavelet envelope spectrum;

s2: obtaining a standardized reconstruction feature vector; and

s3: determining whether the rolling bearing is in failure; wherein,

the step S1 specifically includes:

s1-1: collecting vibration acceleration signals, and storing the collected vibration acceleration signals in a segmented manner to obtain a plurality of samples;

s1-2: extracting 6-dimensional time domain statistical parameters from the sample through the time domain statistical parameters, wherein the 6-dimensional time domain statistical parameters comprise: form factor T_SIPeak index T_CIPulse index T_MIMargin index T_CLIKurtosis T_KUAnd degree of distortion T_SK；

S1-3: extracting 3-dimensional frequency domain statistical parameters from the sample through frequency domain statistical parameters, wherein the 3-dimensional frequency domain statistical parameters comprise gravity center frequency F_FCMean square frequency F_MSFSum frequency variance F_VF(ii) a And

s1-4: obtaining a wavelet envelope spectrum of the sample through wavelet transformation, and extracting 3-dimensional dimensionless parameters in the wavelet envelope spectrum, wherein the 3-dimensional dimensionless parameters in the wavelet envelope spectrum comprise: inner ring failure frequency f_ICorresponding characteristic W_BPFIOuter ring fault frequency f_OCorresponding characteristic W_BPFOAnd ball failure frequency f_BCorresponding characteristic W_BSF；

The step S2 specifically includes:

s2-1: based on a minimum reconstruction error criterion, applying a constraint condition to the 12-dimensional dimensionless parameters to obtain 12-dimensional decorrelation reconstruction characteristic vectors;

s2-2: normalizing the 12-dimensional decorrelated reconstructed feature vector to obtain a normalized reconstructed feature vector with a 12-dimensional sample mean value of 0 and a sample standard deviation of 1; and is

The step S3 specifically includes:

s3-1: training to obtain a sample mean value and a sample standard deviation of the 2 norm of the standardized reconstruction characteristic vector under the normal operation state of the rolling bearing, and setting a threshold value of the rolling bearing state abnormity;

s3-2: testing, namely taking the Euclidean distance from the characteristic vector of the sample in the unknown state of the rolling bearing to the characteristic vector of the sample mean value in the normal running state of the rolling bearing as a degradation index, comparing the degradation index with the threshold value, and judging that the rolling bearing is in a fault state when the degradation index is greater than the threshold value; and when the degradation index is less than or equal to the threshold value, determining that the rolling bearing is in a normal operation state.

2. The rolling bearing fault online detection and state evaluation method according to claim 1, characterized in that: the greater the degradation index is, the greater the degree of degradation of the rolling bearing is determined to be.

3. The rolling bearing fault online detection and state evaluation method according to claim 1, characterized in that: the specific operation of extracting the 6-dimensional time domain statistical parameter in the step S1-2 is as follows:

form factor T_SIIs given by the formula (1), peak index T_CIGiven by equation (2), the pulse index T_MIIs given by the formula (3), margin index T_CLIGiven by the formula (4), kurtosis T_KUGiven by formula (5), the skew T_SKIs given by the formula (6),

$T_{SI}=\frac{\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}\left(y_i^2\right)}}{\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}\left|y_i\right|}---\left(1\right)$

$T_{CI}=\frac{\underset{i=1}{\overset{10}{\Σ}}{y}_{pi}}{\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}{\left(y_i\right)}^2}}---\left(2\right)$

$T_{MI}=\frac{\underset{i=1}{\overset{10}{\Σ}}{y}_{pi}}{\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}\left|y_i\right|}---\left(3\right)$

$T_{CLI}=\frac{\underset{i=1}{\overset{10}{\Σ}}{y}_{pi}}{{\[\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}\sqrt{\left|y_i\right|}\]}^2}---\left(4\right)$

$T_{KU}=\frac{\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}{y_i}^4}{{\left(\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}{y_i}^2\right)}^2}---\left(5\right)$

$T_{SK}=\frac{\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}{y_i}^3}{{\left(\frac{1}{N}\underset{i=1}{\overset{N}{\Σ}}{y_i}^2\right)}^{\frac{3}{2}}}---\left(6\right)$

wherein, y_iIs data in the vibro-acceleration signal; y is_piIs the maximum value of the absolute value of each data segment after the data is divided into 10 data segments.

4. The rolling bearing fault online detection and state evaluation method according to claim 3, characterized in that: the specific operation of extracting the 3-dimensional time domain statistical parameter in the step S1-3 is as follows:

center of gravity frequency F_FCMean square frequency F given by equation (7)_MSFGiven by equation (8), the frequency variance F_VFIs given by the formula (9),

$F_{FC}=\frac{\underset{i=0}{\overset{n}{\Σ}}{f}_{i}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\Σ}}S\left(f_i\right)}---\left(7\right)$

$F_{MSF}=\frac{\underset{i=0}{\overset{n}{\Σ}}{f}_i^{2}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\Σ}}S\left(f_i\right)}---\left(8\right)$

$F_{VF}=\frac{\underset{i=0}{\overset{n}{\Σ}}{(f_i-F_{FC})}^{2}S\left(f_i\right)}{\underset{i=0}{\overset{n}{\Σ}}S\left(f_i\right)}---\left(9\right)$

wherein, S (f)_i) As a function of the vibration acceleration signal frequency spectrum.

5. The rolling bearing fault online detection and state evaluation method according to claim 3, characterized in that: the specific operation of the step S1-4 comprises the following three steps:

s1-4-1, calculating the inner ring fault frequency f according to the formulas (10) - (12)_IOuter ring fault frequency f_OAnd ball failure frequency f_B，

$f_I=\frac{1}{2}Z(1+\frac{d}{D}cos\mathrm{\α})f_r---\left(10\right)$

$f_O=\frac{1}{2}Z(1-\frac{d}{D}cos\mathrm{\α})f_r---\left(11\right)$

$f_B=\frac{D}{2d}\[1-{\left(\frac{d}{D}\right)}^2{\cos }^2\mathrm{\α}\]f_r---\left(12\right)$

Wherein D represents the diameter of the rolling elements, D represents the pitch diameter of the bearing, Z represents the number of the rolling elements, α represents the contact angle, f_rRepresenting the relative rotation frequency of the inner ring and the outer ring of the rolling bearing;

s1-4-2, carrying out wavelet decomposition on the vibration acceleration signal by adopting a db8 wavelet to obtain 6 signals, wherein the 6 signals comprise 5 detail signals d1, d2, d3, d4, d5 and 1 approximate signal a 5;

s1-4-3, extracting 3-dimensional dimensionless parameters in the wavelet envelope spectrum: setting the frequency spectrum in wavelet packet envelope spectrum, the fault characteristic frequency and the characteristic spectrum peak near each order of frequency multiplication, and setting the envelope spectrum analysis bandwidth as f_eThe envelope spectrum is W (f), and the number of spectral lines is N_eThen S is_eaIs composed of

$S_{ea}=\frac{1}{N_e}\underset{i=0}{\overset{N_e}{\Σ}}W\left(f_i\right)---\left(13\right)$

Then order S_edSetting the number of spectral lines of fault frequency in the envelope spectrum as n for the average value of spectral lines at each order of multiple frequency of fault characteristic frequency in the envelope spectrum_eThen, then

$S_{ed}=\frac{1}{n_e}\underset{i=0}{\overset{n_e}{\Σ}}W\left({\mathrm{if}}_d\right)---\left(14\right)$

Constructing a dimensionless parameter:

${\mathrm{\ΔS}}_e=\frac{S_{ed}}{S_{ea}}---\left(15\right)$

respectively carrying out envelope spectrum analysis on the 6 signals in the step S1-4-2 to obtain the fault frequency f of the inner ring_ICorresponding characteristic W_BPFIOuter ring fault frequency f_OCorresponding characteristic W_BPFOFrequency of ball failure f_BCorresponding characteristic W_BSF。

6. The rolling bearing fault online detection and state evaluation method according to claim 1, characterized in that: the specific operation of the step S2-1 comprises the following three steps:

s2-1-1: standardizing the 12-dimensional dimensionless parameters to obtain x_iThen there is ∑_ix_i＝0；

S2-1-2: assuming a new coordinate system { w ] obtained after projective transformation₁,w₂,…,w_dWhere d is a feature dimension; imposed constraintsIs wi is a normal orthogonal basis vector satisfying | | w_i||₂＝1，w_i ^Tw_j0(i ≠ j), the optimization objective function is

$\begin{array}{cc}\underset{W}{min}& -tr\left(W^T{\mathrm{XX}}^{T}W\right)\\ s.t.& W^{T}W=I\end{array}---\left(16\right)$

S2-1-3: solving the equation (16) to obtain a projection matrix W, and calculating a decorrelated reconstructed eigenvector z according to the equation (17)_i

z_i＝(z_i1,z_i2,…,z_id)，z_id＝w_j ^Tx_i(17)。

7. The rolling bearing fault online detection and state evaluation method according to claim 6, characterized in that: the normalized reconstructed feature vector z in step S2-2_i ^*Is given by the formula (18)

$z_i^{*}=(z_{i1}^{*},z_{i2}^{*},...,z_{id}^{*}),z_{ij}^{*}=\frac{z_{ij}-{\widehat{\mathrm{\μ}}}_j}{{\widehat{\mathrm{\σ}}}_j}---\left(18\right)$

Wherein,

${\widehat{\mathrm{\μ}}}_j=\frac{1}{n}\underset{i=1}{\overset{n}{\Σ}}{z}_{ij},{\widehat{\mathrm{\σ}}}_j^2=\frac{1}{n-1}\underset{i=1}{\overset{n}{\Σ}}{(z_{ij}-{\widehat{\mathrm{\μ}}}_j)}^2---\left(19\right).$

8. the rolling bearing fault online detection and state evaluation method according to claim 1, characterized in that: the threshold value in step S3-1 is given by equation (20)

$D_{max}={\widehat{\mathrm{\μ}}}_D+3{\widehat{\mathrm{\σ}}}_D---\left(20\right)$

Wherein,andare each D_iSample mean and sample standard deviation.