Wind generating set pivotal bearing fault diagnosis method based on high-order statistic analysis

Abstract

The invention provides a wind generating set pivotal bearing fault diagnosis method based on high-order statistic analysis, which takes skewness and kurtosis as measurement standards for denoising of a statistical filter, inhibits interference of background noise and enhances expression of fault features so as to relieve the problem that a fault signal is covered by the background noise. According to the method, deep features caused by bearing faults are reflected in a high-dimensional space by applying phase-space reconstruction, and feature fault diagnosis is performed on non-Gaussian and non-linear signals in the high-dimensional space by adopting adaptive multi-spectral estimation. Compared with a traditional linear feature diagnosis method, the problems that background noise and fault signals are mixed together and are difficult to separate and non-linear and non-Gaussian features interfere with diagnosis are solved, and a new solution is provided for fault diagnosis of the slewing bearing of the wind generating set.

Description

Technical Field

[0001] This method belongs to the field of wind power generation, and specifically relates to a fault diagnosis method for the slewing bearing of a wind turbine based on high-order statistic analysis. Background Technique

[0002] As a sustainable energy source, wind energy has become one of the important sources of global electricity production. Wind turbines extract wind energy from the air flow in the atmosphere. The slewing bearings in wind turbine generators are divided into two categories: pitch bearings and yaw bearings. The pitch bearing plays a key role in the process of wind turbines extracting wind energy from the air flow in the atmosphere. It is responsible for adjusting the blade pitch angle to obtain the optimal output power of the wind turbine. The yaw bearing is responsible for adjusting the direction of the wind turbine to the wind direction to capture the maximum wind energy. Due to the complex working environment of the slewing bearing and the long-term state of load change, the bearing is prone to failure. The failure of the slewing bearing usually shows as weak vibration signals. When running at low speed, due to sliding, the vibration of other machine components, noise and interference, the signal characteristics related to the failure are blurred and covered up. Therefore, it is very difficult to identify the characteristic frequency of the slewing bearing from the vibration signal, which leads to difficult accurate fault diagnosis; in addition, the bearing fault signal has complex and time-varying dynamic characteristics, which makes the traditional methods based on linear or Gaussian assumptions unable to capture these characteristics.

[0003] As an effective signal processing technology, high-order statistics perform excellently especially in the analysis of non-linear and non-Gaussian signals. It can provide more information than traditional linear analysis methods. Therefore, how to use high-order statistic analysis to conduct fault diagnosis on non-linear and non-Gaussian vibration signals in a complex noise environment is still a technical problem to be solved urgently. Summary of the Invention

[0004] This method proposes a fault diagnosis method for the slewing bearing of a wind turbine based on high-order statistic analysis. First, aiming at the current situation that the weak fault signal of the bearing is covered by noise, skewness and kurtosis are used as the measurement criteria for statistical filter denoising to suppress the interference of background noise and enhance the expression of fault characteristics. The signal after denoising processing will be used as the input of the feature extraction module for subsequent fault diagnosis; secondly, aiming at the non-linear and non-Gaussian characteristics of the fault signal, this invention uses phase space reconstruction to reflect the deep-level characteristics caused by bearing faults in a high-dimensional space, and uses adaptive polyspectrum estimation to analyze the signal in the high-dimensional space, and then uses it for fault diagnosis.

[0005] The technical solution adopted by this method is as follows:

[0006] S1: Obtain the vibration signal of the slewing bearing and set it as the original signal;

[0007] S2: Segment the original signal. The number of segments is selected as an integer multiple of the signal period to ensure that each segment of the signal can contain a complete signal period and avoid the loss of periodic characteristics. Calculate the skewness and kurtosis of each segment of the signal in turn;

[0008] S3: For the problem that the fault signal of the slewing bearing is weak and easily affected by background noise, set the initial thresholds of skewness and kurtosis according to the statistical distribution differences of skewness and kurtosis in different wind turbine data. The skewness value is used to measure the statistical asymmetry of the signal, and the kurtosis value is used to measure the sharpness of the signal distribution. The vibration signal of a normal bearing shows a stable and symmetric distribution, with the skewness value close to 0 and the kurtosis value close to 3; the vibration signal of a faulty bearing has impacts and asymmetric fluctuations due to defects, the skewness value deviates from 0, and the kurtosis value is greater than 3. Design a statistical filter through skewness and kurtosis in high-order statistic analysis to suppress the interference of background noise.

[0009] S4: Compare and calculate the skewness value and kurtosis value of each segment of the signal with the initial thresholds of skewness and kurtosis to obtain the weight coefficient. When the calculated value is greater than 1, it is considered that the signal segment contains fault characteristics, highlight the signal segment, otherwise perform denoising processing to suppress the interference of background noise on fault diagnosis, and then obtain the vibration signal after denoising;

[0010] S5: Take the denoised signal as the input data for the next step of phase space reconstruction. Through the time-delay embedding method, transform the one-dimensional time series into an expression in a high-dimensional phase space to characterize the dynamic characteristics of the signal in the high-dimensional space. In the reconstructed phase space, perform adaptive polyspectrum estimation on the signal to reveal the phase information and high-order non-linear characteristics of the signal, extract the interaction relationship between frequencies in different dimensions, identify potential fault modes in the slewing bearing, and obtain the data diagnosis result;

[0011] S6: Denote the error between the data diagnosis result and the actual state result as ε, set the error requirement between the data diagnosis result and the actual state result as ε0, and dynamically update the high-order statistic parameters and iterative weight parameters through the gradient descent method so that the error between the data diagnosis result and the actual state result reaches the set error requirement between the data diagnosis result and the actual state result, gradually approaching a solution with the minimum error, and improving the accuracy of fault diagnosis. Description of the Drawings

[0012] Figure 1 Flowchart of the Fault Diagnosis Method for the Slewing Bearing of a Wind Turbine Based on High-Order Statistic Analysis Detailed Implementation Manner

[0013] First, install a suitable vibration sensor on the slewing bearing for signal acquisition. The selection of the sensor should be optimized according to the working characteristics of the paddle bearing. The sensor is required to capture the low-frequency signals during the operation of the slewing bearing and the high-frequency signals caused by friction and impact factors, and have sufficient sensitivity to capture small vibration signals. A sensor with a wide dynamic range is used to handle signals with vibration intensities ranging from low to high, so as to ensure that the collected data contains various dynamic information of the bearing operation and can reflect the vibration conditions at all stages during the operation of the slewing bearing. The original signal collected is x(t), and x(t) contains fault signals and noise signals. Further, the specific steps for data denoising based on the statistical characteristics of background noise and fault signals are as follows:

[0014] Step 1: Calculate the skewness and kurtosis values of the vibration signal of the bearing in the normal state:

[0015]

[0016] In the formula

[0017] x represents the vibration signal;

[0018] S represents the signal skewness value;

[0019] E represents the mathematical expectation;

[0020] μ represents the signal mean value;

[0021] σ represents the signal standard deviation.

[0022]

[0023] In the formula

[0024] K represents the signal kurtosis value.

[0025] Set thresholds S0 and K0 according to the statistical distribution diagrams of the skewness and kurtosis values.

[0026] Step 2: Segment the original signal with a length of L. The number of segments is selected as an integer multiple of the signal period number to ensure that each segment of the signal can contain a complete signal period and avoid the loss of periodic characteristics. Calculate the skewness and kurtosis of each segment of the signal, and the denoising weight coefficient is automatically adjusted according to the calculation results of the skewness and kurtosis. The weight parameter is defined as:

[0027]

[0028] In the formula

[0029] ω i is the weight coefficient of the i-th segment of the signal;

[0030] K i is the kurtosis value of the i-th segment of the signal;

[0031] S i is the skewness value of the i-th segment of the signal;

[0032] If the skewness and kurtosis exceed the initial thresholds K0 and S0, it is considered that the signal segment contains fault characteristics, and the weight coefficient is multiplied by the signal segment to highlight the fault characteristics; if the skewness and kurtosis are lower than the initial thresholds K0 and S0, it is considered that the signal segment belongs to noise, and median filtering is performed. The filtered signal of the i-th segment is expressed as:

[0033]

[0034] where

[0035] x i (t) is the original signal of the i-th segment that has not been processed;

[0036] is the vibration signal of the i-th segment after denoising;

[0037] MF represents the median filter.

[0038] Furthermore, the specific steps of phase space reconstruction are as follows:

[0039] The first step: Select the denoised time series data

[0040] The second step: Calculate a suitable time delay τ to determine the element interval of each new state vector. Observe the change trend of the signal autocorrelation function with the time delay, and find the time delay corresponding to the first significant decline as the suitable time delay, so that adjacent state points in the phase space do not coincide as much as possible while maintaining the dynamic characteristics of the system. The calculation formula of the autocorrelation function is:

[0041]

[0042] where

[0043] R(τ) is the autocorrelation function;

[0044] L is the signal length.

[0045] R(τ) reflects the similarity between the signal after a time delay τ and the original signal. Calculate the autocorrelation function of different time delays τ in turn, and find the position where the autocorrelation function decays to the specified threshold (set to 1 / e) for the first time, that is, the time point when the signals lose linear correlation.

[0046] Step 3: Select the embedding dimension. The embedding dimension represents how many independent variables are needed in the phase space to describe the dynamics of the system. When the embedding dimension increases, the local structure of the data points, i.e., the relative distances between data points, should remain consistent. If the embedding dimension is insufficient, some neighboring data points will be wrongly mapped to a farther region, resulting in the appearance of false nearest neighbors. When the number of false nearest neighbors is small, it indicates that the embedding dimension for phase space reconstruction is appropriate. Therefore, the adequacy of the embedding dimension can be judged by monitoring the number of false nearest neighbors to ensure that the local structure of the data is retained, thus determining the embedding dimension. In the m-dimensional space, select a data point x m (t) and its neighboring data point x m (t0), and calculate the Euclidean distance between the two:

[0047] d m (t) = ||x m (t) - x m (t0)||

[0048] In the formula

[0049] d m (t) represents the Euclidean distance between two adjacent data points in the m-dimensional space;

[0050] ||·|| represents the Euclidean norm.

[0051] Embed the data point from the m-dimensional space into the m + 1-dimensional space to obtain a new embedding vector:

[0052] x m+1 (t) = [x(t), x(t + τ), …, x(t + (m - 1)τ), x(t + mτ)]

[0053] In the m + 1-dimensional space, select a data point x m+1 (t) and its neighboring data point x m+1 (t0), and calculate the new Euclidean distance between the two adjacent data points:

[0054] d m+1 (t) = ||x m+1 (t) - x m+1 (t0)||

[0055] In the formula

[0056] d m+1 (t) represents the Euclidean distance between two adjacent data points in the m + 1-dimensional space;

[0057] If the data points are close in the low-dimensional space but become far apart in the high-dimensional space, then this point can be regarded as a false nearest neighbor. Assume that in the m-dimensional space, the data point x m (t) and xm The Euclidean distance d of (t0) m (t) is very small, and in m+1 dimensional space, the Euclidean distance d m+1 (t) becomes larger, it means that the data point may be a false neighbor in the low-dimensional space. For each pair of adjacent points, calculate their distance ratio r(t) in the low-dimensional space and the high-dimensional space:

[0058]

[0059] When r(t) is greater than 20, the point is considered to be a false neighbor point.

[0060] Calculate the proportion of false neighbor points to the total number of data points in the corresponding dimension. As the dimension m of the phase space increases, the proportion of false neighbors gradually decreases until the proportion of false neighbor points tends to be stable, thereby determining the phase space dimension D of the data:

[0061]

[0062] In the formula

[0063] D is the appropriate embedding dimension for phase space reconstruction;

[0064] N f (m) is the proportion of false neighbor points in m dimension, that is, the number of false neighbor data points;

[0065] N(m) is the total number of data points in the m dimension.

[0066] Step 4: The data constructs a new state vector in the phase space according to the delay time and the embedded dimension. Through phase space reconstruction, the multi-dimensional dynamic characteristics of the original system are reconstructed from a single vector, the dynamic information of the original system is restored, and the nonlinear characteristics of the data are analyzed more intuitively, providing a basis for subsequent complex system processing.

[0067] Furthermore, an adaptive multi-spectral estimation method is used to reconstruct the signal x in D-dimensional space. D (t) Perform feature extraction, analyze the frequency components of the signal from different dimensions, and identify the spectrum anomalies caused by the fault, namely:

[0068]

[0069] In the formula

[0070] θ i is the time delay variable;

[0071] ω i is a frequency variable;

[0072] S D (ω1,ω2,…,ω D-1 ) is the signal xD The D - order spectrum of (t);

[0073] C D (θ1, θ2, …, θ D-1 ) is the (D - 1) - order cumulant of the signal x D (t).

[0074] Furthermore, through polyspectrum estimation, high - order information of the signal is obtained, revealing the tiny and rapidly changing characteristics existing during bearing faults, capturing the interaction relationships between various faults, and identifying the diagnostic results.

[0075] Furthermore, denote the error between the data diagnostic result and the actual state result as ε, set the error requirement between the data diagnostic result and the actual state result as ε0, and make the error between the data diagnostic result and the actual state result reach the set error requirement by continuously updating the high - order statistic parameters (skewness and kurtosis thresholds) and iterating the weight parameters. The specific steps are as follows:

[0076] Step 1: Take the mean - square error as the target error function ε for the error between the data diagnostic result and the actual state result;

[0077] Step 2: Calculate the gradients of the error function with respect to the skewness and kurtosis thresholds, and update the skewness and kurtosis values:

[0078]

[0079] Where

[0080] S (k+1) represents the skewness threshold at the (k + 1) - th iteration;

[0081] S (k) represents the skewness threshold at the k - th iteration;

[0082] K (k+1) represents the kurtosis threshold at the (k + 1) - th iteration;

[0083] K (k) represents the kurtosis threshold at the k - th iteration;

[0084] γ represents the learning rate;

[0085] represents the gradient of the error function with respect to the skewness;

[0086] represents the gradient of the error function with respect to the kurtosis.

[0087] Step 3: Calculate the new error ε (k+1) , when ε (k+1)The error requirement between the set data diagnosis result and the actual state result is met, and the iteration process stops.

[0088] Dynamically adjust the skewness threshold and kurtosis threshold according to the change of the error, optimize the denoising weight coefficient, reduce the error between the diagnosis result and the actual component operation state result after multiple iterations, and improve the fault diagnosis accuracy.

Claims

1. This method proposes a fault diagnosis method for the slewing bearing of a wind turbine based on high-order statistic analysis, characterized by Skewness and kurtosis are used as the criteria for statistical filter denoising to suppress the interference of background noise and highlight the expression of fault features; On the other hand, phase space reconstruction is used to reflect the deep features caused by bearing faults in a high-dimensional space, and adaptive polyspectrum estimation is adopted in the high-dimensional space for feature fault diagnosis of non-Gaussian and non-linear signals.

2. According to claim 1, characterized in that Data denoising is carried out according to the statistical characteristics of background noise and fault signals. The specific steps are as follows: The first step: Calculate the skewness and kurtosis values of the bearing vibration signal in the normal state: In the formula x represents the vibration signal; S represents the signal skewness value; E represents the mathematical expectation; μ represents the signal mean; σ represents the signal standard deviation, In the formula K represents the signal kurtosis value, Set thresholds S0 and K0 according to the statistical distribution diagrams of skewness and kurtosis values. The second step: Segment the original signal with a length of L. The number of segments is selected as an integer multiple of the signal period number to ensure that each segment of the signal can contain a complete signal cycle and avoid the loss of periodic characteristics. Calculate the skewness and kurtosis of each segment of the signal. The denoising weight coefficient is automatically adjusted according to the calculation results of skewness and kurtosis. The weight parameter is defined as: In the formula ω i is the weight coefficient of the i-th segment of the signal; K i is the kurtosis value of the i-th segment of the signal; S i is the skewness value of the i-th segment of the signal; If the skewness and kurtosis exceed the initial thresholds K0 and S0, it is considered that this segment of the signal contains fault features, and the weight coefficient is multiplied by the signal segment to highlight the fault features; If the skewness and kurtosis are lower than the initial thresholds K0 and S0, the signal segment is considered to be noise and median filtering is performed. The filtered signal of the i-th segment is expressed as: In the formula x i (t) is the original signal of the i-th segment that has not been processed; is the vibration signal of the i-th segment after denoising processing; MF represents the median filter.

3. According to claim 1, characterized in that Perform phase space reconstruction on the denoised signal. The specific steps are as follows: Step 1: Select the denoised time series data The second step: Calculate a suitable time delay τ to determine the element interval of each new state vector. Observe the change trend of the signal autocorrelation function with the time delay, and find the time delay corresponding to the first significant decrease as the suitable time delay, so that adjacent state points in the phase space do not overlap as much as possible while maintaining the dynamic characteristics of the system. The calculation formula of the autocorrelation function is: In the formula R(τ) is the autocorrelation function; L is the signal length, R(τ) reflects the similarity between the signal after a time delay τ and the original signal. Calculate the autocorrelation function of different time delays τ in turn, and find the position where the autocorrelation function decays to a specified threshold (set as 1 / e) for the first time, that is, the time point when the signals lose linear correlation. Step 3: Select the embedding dimension. The embedding dimension indicates how many independent variables are needed in the phase space to describe the dynamics of the system. When the embedding dimension increases, the local structure of the data points, i.e., the relative distances between the data points, should remain consistent. If the embedding dimension is insufficient, some neighboring data points will be wrongly mapped to farther regions, resulting in the appearance of false nearest neighbors. When the number of false nearest neighbors is small, it indicates that the embedding dimension for phase space reconstruction is appropriate. Therefore, the sufficiency of the embedding dimension can be judged by monitoring the number of false nearest neighbors to ensure that the local structure of the data is retained, thereby determining the embedding dimension. In the m-dimensional space, select a data point x m (t) and its neighboring data point x m (t0), and calculate the Euclidean distance between the two: d m (t) = ||x m (t) - x m (t0)|| In the formula d m (t) represents the Euclidean distance between two adjacent data points in the m-dimensional space; ||·|| represents the Euclidean norm, Embed the data points from the m-dimensional space into the m + 1-dimensional space to obtain a new embedding vector: x m+1 (t) = [x(t), x(t + τ), …, x(t + (m - 1)τ), x(t + mτ)] In the (m + 1)-dimensional space, select a data point x m+1 (t) and its adjacent data point x m+1 (t0), and calculate the new Euclidean distance between two adjacent data points: d m+1 (t) = ||x m+1 (t) - x m+1 (t0)|| In the formula d m+1 (t) represents the Euclidean distance between two adjacent data points in an (m + 1)-dimensional space; If data points are close in a low-dimensional space but become far apart in a high-dimensional space, then this point can be regarded as a false nearest neighbor. Assume that in the m-dimensional space, the data point x m (t) and x m (t0) has a small Euclidean distance d m (t), while in the m+1-dimensional space, the Euclidean distance d m+1 (t) becomes larger, indicating that this data point may be a false nearest neighbor in the low-dimensional space. For each pair of adjacent points, calculate the distance ratio r(t) between them in the low-dimensional space and the high-dimensional space: When r(t) is greater than 20, then this point is considered a false nearest neighbor point. Calculate the proportion of the number of false nearest neighbor points to the total number of data points in the corresponding dimension. As the dimension m of the phase space increases, the false nearest neighbor ratio gradually decreases until the false nearest neighbor point ratio tends to be stable, so as to determine the phase space dimension D of the data: In the formula D is the appropriate embedding dimension for phase space reconstruction; N f (m) is the proportion of false nearest neighbors in the m - dimensional space, that is, the number of data points of false nearest neighbors; N(m) is the total number of data points in the m dimension, The fourth step: The data constructs new state vectors in the phase space according to the delay time and embedding dimension. Through phase space reconstruction, the multi-dimensional dynamic characteristics of the original system are reconstructed from a single vector, the dynamic information of the original system is restored, and the non-linear characteristics of the data are analyzed more intuitively, providing a basis for subsequent complex system processing.

4. According to claim 1, characterized in that The adaptive multi-spectrum estimation method is used to extract the features of the signal x D (t) reconstructed in the D-dimensional space, analyze the frequency components of the signal from different dimensions, and distinguish the spectral anomalies caused by faults, that is: In the formula θ i is the time delay variable; ω i is a frequency variable; S D (ω1, ω2, …, ω D-1 ) is the D-order spectrum of the signal x D (t); C D (θ1, θ2, …, θ D-1 ) is the (D - 1)th order cumulant of the signal x D (t), Furthermore, through multi-spectrum estimation, the high-order information of the signal is obtained to reveal the tiny and rapidly changing features existing in bearing faults, capture the interaction relationships between various faults, and identify the diagnostic results.

5. According to claim 1, characterized in that Furthermore, the error between the data diagnostic result and the actual state result is denoted as ε, and the error requirement between the data diagnostic result and the actual state result is set as ε0. By continuously updating the high-order statistic parameters (skewness and kurtosis thresholds) to iterate the weight parameters, the error between the data diagnostic result and the actual state result reaches the set error requirement between the data diagnostic result and the actual state result. The specific steps are as follows: The first step: The mean square error is used as the target error function ε for the error between the data diagnostic result and the actual state result. The second step: Calculate the gradients of the error function with respect to the skewness and kurtosis thresholds, and update the skewness and kurtosis values: In the formula S (k+1) represents the skewness threshold at the (k + 1)-th iteration; S (k) represents the skewness threshold at the k-th iteration; K (k+1) represents the kurtosis threshold at the (k + 1)-th iteration; K (k) represents the kurtosis threshold at the k-th iteration; γ represents the learning rate. Denote the gradient of the error function with respect to skewness; Denotes the gradient of the error function with respect to kurtosis, Step 3: Calculate the new error ε (k+1) , when ε (k+1) reaches the error requirement between the set data diagnosis result and the actual state result, the iteration process stops. Dynamically adjust the skewness threshold and kurtosis threshold according to the change of the error, optimize the denoising weight coefficient, and reduce the error between the diagnostic result and the actual component operation state result after multiple iterations to improve the accuracy of fault diagnosis.

Images