A distributed electric vehicle mechanical transmission system planetary wheel bearing fault detection method

Abstract

The application discloses a kind of distributed electric vehicle mechanical transmission system planetary wheel bearing fault detection methods, belong to the field of fault diagnosis and signal processing;The method is for the problems existing in the planetary wheel bearing fault detection of distributed electric vehicle mechanical transmission system, such as the limited installation position of sensor, the weak fault impact feature is difficult to extract, with the encoder signal as the signal source, based on the instantaneous angular velocity signal calculation cycle spectrum coherence, the feature saliency of spectral frequency component is evaluated by introducing the morphological outlier index, and the spectral frequency containing rich fault information is fused by adaptive multi-spectrum proportion fusion method, to obtain the fusion spectrum, to extract the fault feature of planetary bearing;The method of the application solves the problem that the weak fault impact feature of planetary wheel bearing is difficult to extract, and can realize the fault detection of planetary wheel bearing fault under variable speed condition.

Description

Technical Field

[0001] This invention relates to a method for detecting faults in planetary gear bearings of a distributed electric vehicle mechanical transmission system, belonging to the field of fault diagnosis and signal processing. Background Technology

[0002] As an important development direction for new energy vehicles, distributed electric vehicles often adopt integrated and lightweight mechanical transmission structures in their drive systems. Among them, planetary gearboxes are widely used in drive axles or wheel-side drive systems due to their advantages such as compact structure, large transmission ratio, and strong load-bearing capacity. However, as key support components, planetary gear bearings are subjected to complex working conditions such as high speed and variable load for a long time, and are simultaneously constrained by the dual motion of planetary gear rotation and revolution. They also bear the combined load from the sun gear and the ring gear, resulting in a relatively high failure rate and posing special challenges to their fault diagnosis. On the one hand, the compact transmission structure of distributed electric vehicles limits the installation location of sensors. Existing methods, such as vibration signal analysis, rely on the arrangement of multiple measurement points, which are not suitable for distributed drive scenarios with limited sensors. On the other hand, vibration signals need to pass through a time-varying transmission path before they can be received by the sensors, which causes severe amplitude and phase modulation, and the fault characteristics are drowned out by noise. On the other hand, although encoder-based instantaneous angular velocity analysis can avoid installation limitations, the impact energy generated by local faults in planetary gear bearings in the instantaneous angular velocity signal is weak, and the signal transmission path is complex. The fault characteristic frequency overlaps with gear meshing frequency, shaft rotation frequency, etc. Traditional spectrum analysis or envelope demodulation methods are easily affected by background noise and interference components, leading to feature extraction failure. Weak fault impacts are difficult to excite significant torsional resonance, which limits the effectiveness of traditional demodulation techniques. Summary of the Invention

[0003] To address the problems existing in the prior art, this invention provides a method for fault detection of planetary gear bearings in the mechanical transmission system of a distributed electric vehicle. This method uses instantaneous angular velocity signals as a basis, and evaluates the saliency of fault information in each frequency component of the cyclic spectrum through cyclic spectrum coherence analysis and the introduction of a morphological outlier index. Based on this, an adaptive proportional multispectral fusion method is constructed to obtain an optimized fused spectrum for extracting fault features of the planetary bearings. This method avoids the limitations of traditional cyclic spectrum analysis methods that rely on resonance demodulation and frequency band selection, and achieves feature extraction of weak planetary bearing faults, effectively solving the problem of fault detection in planetary gear bearings of distributed electric vehicles.

[0004] The fault detection method for planetary gear bearings in the distributed electric vehicle mechanical transmission system of the present invention is as follows: Step 1: Acquire the instantaneous angular displacement sequence θ of the optical encoder using the PicoScope signal acquisition system. m and the corresponding time series t mThe instantaneous angular velocity signal w(θ) is calculated using the forward difference method. m The calculation formula is as follows: ; In the formula: w(θ) m ) represents the time m at angular displacement θ m The instantaneous angular velocity, m=1, 2, 3,…; L represents the number of gratings in the encoder; Step 2: Set the cyclic spectral window width N w =2 k For k=1, 2, ..., calculate w(θ) according to the following steps. m Cyclic spectral coherence of CSCoh(α, f): (1) w(θ) m The cyclic spectral correlation (CSC) function can be defined as the two discrete Fourier transforms of the instantaneous autocorrelation function: ; In the formula, α represents the cycle frequency, f represents the spectral frequency, Δθ=2π / L represents the angular domain sampling interval, L represents the number of gratings in the encoder, and R... w (θ m , ) represents the signal w(θ) m The autocorrelation function of ) This represents the angular delay, where j is the imaginary unit; (2) To eliminate the influence of non-uniform noise, CSC is used for normalization to obtain CSCoh, whose expression is: ; Step 3: Based on the cyclic spectral coherence CSCoh(α, f), the optimized fused spectrum f is obtained using an adaptive proportional multispectral fusion method. co op The specific process is as follows: (1) The spectral frequency component f(α) is obtained from the cyclic spectral coherence CSCoh(α, f), and its expression is: ; In the formula, f i (α) represents the i-th spectral frequency component; (2) The frequency band of interest β in the spectral frequency component f(α) is extracted using the following formula. h (α); ; In the formula, h represents the fault harmonic order, h=1,2,...,H, where H is the maximum fault harmonic order considered (for bearing faults, H is usually taken as 3); α tThis indicates the range of background noise that needs to be evaluated (since the fault characteristics of encoder signals are typically of a low order, α). t It is usually set to 0.5–0.8α. fault α fault The order of theoretical fault characteristics; (3) To eliminate β h The difference in amplitude (α) is normalized using the following formula to obtain the normalized amplitude spectrum β. h-std (α): ; In the formula, min(·) represents the minimum value, and max(·) represents the maximum value; (4) For the normalized amplitude spectrum β h-std (α) Arrange in descending order to obtain an ordered morphological distribution U h (α), expressed as: ; In the formula, sort(·) represents sorting in descending order; (5) To distinguish between normal fluctuations and abnormally high amplitude points caused by faults or other factors, based on the ordered morphological distribution U h (α) Introduce an adaptive normal upper bound d h The expression is: ; ; In the formula, Cov(·) represents the covariance, Var(·) represents the variance, and M[·] represents the median. d represents the mean. h It's U h Simple linear trend prediction A of (α) h (α)+B h The maximum predicted value, A h B represents the predicted slope term. h Represents the intercept term; (6) Further, calculate the upper bound of d beyond the normal state. h Outlier impact integral P h The expression is: ; (7) Calculate each spectral frequency component f using the following formula. i The Morphological Outlier Index (MOI) of (α) evaluates the richness of fault information by the ratio of the amplitude of each fault harmonic order to the integral of the outlier impulse. The expression is: ; (8) For each spectral frequency component f i(α) Sort the sequence {f} in descending order according to the corresponding MOI values. n (α)}, satisfying MOI(f1(α))≥MOI(f2(α))...≥MOI(f Nw / 2 Based on this, the spectral frequency components f(α) containing rich fault information are linearly superimposed according to the fusion ratio r, thereby enhancing the rolling bearing fault components and weakening the non-periodic harmonic components, and obtaining the fusion spectrum f corresponding to the fusion ratio r. co r (α), expressed as: ; In the formula, r = 0.01, 0.02, ... 1, and Q represents the number of spectral frequency components selected by the proportion r; Indicates rounding down; (9) Calculate the fusion spectrum f corresponding to different fusion ratios r. co r For the MOI of (α), select f with the largest MOI. co r (α) as the optimized fusion spectrum f co op (α), expressed as: ; Step 4: Referring to Step 3, calculate the width N of each window. w f below co op (α), denoted as f Nw op (α), choose f Nw op The fusion spectrum with the largest MOI value in (α) is f op (α) is used to reveal the fault characteristics of planetary gear bearings, and its expression is: .

[0005] The beneficial effects of this invention are: (1) The present invention uses an optical encoder as a signal source to realize the fault detection of planetary gear bearings in the distributed electric vehicle mechanical transmission system, which has the advantage of convenient signal acquisition; (2) By selecting a cyclic spectrum slice containing rich fault information instead of a frequency band, this invention solves the problem that the cyclic spectrum demodulation method fails due to the low-pass filtering characteristics of instantaneous angular velocity signals and the difficulty of exciting torsional resonance by weak faults. It avoids the limitations of traditional cyclic spectrum analysis methods that rely on resonance demodulation and frequency band selection, and expands the application of cyclic spectrum analysis in instantaneous angular velocity signals. (3) Based on instantaneous angular velocity signal and cyclic spectrum coherence analysis, this invention introduces morphological outlier index (MOI) to evaluate the characteristic significance of spectral frequency components, and obtains the fused spectrum through adaptive proportional multispectral fusion method, thereby realizing the extraction of weak fault features of planetary gear bearings. Attached Figure Description

[0006] Figure 1 This is a schematic diagram of the method flow of the present invention; Figure 2 Photograph of the experimental platform in Example 1; Figure 3 The instantaneous angular velocity experimental signal w(θ) in Example 1 m ) and its corresponding order spectrum, where Figure (a) is the experimental signal w(θ) of instantaneous angular velocity. m The waveform of θ is shown in Figure (b), which is the experimental signal of instantaneous angular velocity w(θ). m The corresponding order spectrum; Figure 4 This is a schematic diagram of the spectral frequency component f(α) in the present method; Figure 5 The results of the method of the present invention in Example 1 are shown in Figure (a), where the fusion spectrum f is... op The window width, Figure (b) is the fusion spectrum f op The fusion ratio of (α) is shown in Figure (c), which is the fusion spectrum f. op The distribution of the spectral frequency components contained in (α) on the α-f plane and the MOI intensity are shown in Figure (d). op (α); Figure 6 The results of the SAM method analysis in Example 1 are shown in Figure (a), where the optimized weight index MO is shown in Figure (b), and the order spectrum corresponding to the optimized weight index MO is shown in Figure (c). Figure 7 The results of the analysis using the Infogram method in Example 1 are as follows: Figure (a) is the square envelope negative entropy distribution, Figure (b) is the order spectrum determined by Figure (a), Figure (c) is the square envelope spectrum negative entropy distribution, Figure (d) is the order spectrum determined by Figure (c), Figure (e) is the average negative entropy distribution, and Figure (f) is the order spectrum determined by (e). Detailed Implementation

[0007] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative effort are within the scope of protection of the present invention. The methods not specifically described in this embodiment are all conventional methods. Example 1: As Figure 1 As shown, the steps for extracting fault features of the outer ring of the planetary bearing on a distributed electric vehicle test bench under variable speed conditions using the method of the present invention are as follows: Adopting such Figure 2 The distributed electric vehicle planetary reducer experimental platform shown was used for the experiment. The experimental data was acquired through the Picoscope acquisition system, and the sampling frequency f of the high-speed counter was [not specified]. s 5×10 6 Hz, and the instantaneous angular velocity signal w(θ) is calculated using the forward difference method. m The encoder, model K1024G-5000BM-K526, is mounted on the input shaft and has a grating count L=5000. To simulate a fault in the outer ring of the planetary bearing, a small groove approximately 0.8mm wide and 0.8mm deep is machined on the outer ring of the planetary bearing using wire cutting to simulate an outer ring fault. Figure 2 The specific parameters of the planetary gearbox are shown in Table 1. With the input shaft as the reference shaft, the characteristic order of each component of the planetary gearbox is shown in Table 2. The fault characteristic order of the outer ring of the planetary bearing is O. bpfo The calculation formula is as follows: ; ; Table 1 Planetary Gearbox Parameters

[0008] Table 2 Characteristic Order

[0009] Step 1: Set the motor speed to vary within the range of 0-1000 rpm to simulate the typical operating conditions of a distributed electric vehicle from acceleration and smooth driving to deceleration. The instantaneous angular displacement sequence θ of the optical encoder is acquired through the PicoScope signal acquisition system. m and the corresponding time series t m The instantaneous angular velocity signal w(θ) is calculated using the forward difference method. m The calculation formula is as follows: ; In the formula: w(θ) m ) represents the time m at angular displacement θ m The instantaneous angular velocity, m=1, 2, 3,…; L represents the number of gratings in the encoder; the obtained instantaneous angular velocity signal w(θ) m )like Figure 3 As shown in a, the corresponding order spectrum is as follows: Figure 3 As shown in b, it can be seen that in the order spectrum of the original signal, the velocity trend component and the frequency order O sDominant, planetary gear bearing outer ring failure characteristic order O bpfo Spectral lines and higher harmonics are submerged in background noise and cannot be effectively identified.

[0010] Step 2: Set the cyclic spectral window width N w =2 k For k=6, 7, ... 16, calculate w(θ) according to the following steps. m Cyclic spectral coherence of CSCoh(α, f): (1) w(θ) m The cyclic spectral correlation (CSC) function can be defined as the two discrete Fourier transforms of the instantaneous autocorrelation function: ; In the formula, α represents the cycle frequency, f represents the spectral frequency, Δθ=2π / L represents the angular domain sampling interval, and R w (θ m , ) represents the signal w(θ) m The autocorrelation function of ) This represents the angular delay, where j is the imaginary unit; (2) To eliminate the influence of non-uniform noise, CSC is used for normalization to obtain CSCoh, whose expression is: ; Step 3: Set parameter H=3, α fault =O bpfo =4.55, α t =0.5O bpfo Based on CSCoh(α, f), an adaptive proportional multispectral fusion method is used to obtain the optimized fused spectrum f. co op The specific process is as follows: (1) Obtain the spectral frequency component f(α) of the bivariate spectrum, such as Figure 4 As shown, the expression is as follows: ; In the formula, f i (α) represents the i-th spectral frequency component; (2) The frequency band of interest β in the spectral frequency component f(α) is extracted using the following formula. h (α); ; In the formula, h = 1, 2, ..., H represents the order of the fault harmonics, and H is the maximum order of the fault harmonics considered; α t α represents the range of background noise that needs to be evaluated. fault =O bpfo =4.55, α t=0.5O bpfo ; (3) To eliminate β h The difference in amplitude (α) is normalized using the following formula to obtain the normalized amplitude spectrum β. h-std (α): ; In the formula, min(·) represents the minimum value, and max(·) represents the maximum value; (4) For the normalized amplitude spectrum β h-std (α) Arrange in descending order to obtain an ordered morphological distribution U h (α), expressed as: ; In the formula, sort(·) represents sorting in descending order; (5) To distinguish between normal fluctuations and abnormally high amplitude points caused by faults or other factors, based on the ordered morphological distribution U h (α) Introduce an adaptive normal upper bound d h The expression is: ; ; In the formula, Cov(·) represents the covariance, Var(·) represents the variance, and M[·] represents the median. d represents the mean. h It's U h Simple linear trend prediction A of (α) h (α)+B h The maximum predicted value, A h B represents the predicted slope term. h Represents the intercept term; (6) Further, calculate the upper bound of d beyond the normal state. h Outlier impact integral P h The expression is: ; (7) Calculate each spectral frequency component f using the following formula. i The Morphological Outlier Index (MOI) of (α) evaluates the richness of fault information by the ratio of the amplitude of each fault harmonic order to the integral of the outlier impulse. The expression is: ; (8) For each spectral frequency component f i (α) Sort the sequence {f} in descending order according to the corresponding MOI values. n (α)}, satisfying MOI(f1(α))≥MOI(f2(α))...≥MOI(f Nw / 2(α)) is obtained by linearly superimposing the spectral frequency components f(α) containing rich fault information according to the fusion ratio r, so as to obtain the fusion spectrum f corresponding to the ratio r. co r (α), expressed as: ; In the formula, r = 0.01, 0.02...1, and Q represents the number of spectral frequency components selected by the ratio r; Indicates rounding down; (9) Calculate the fusion spectrum f corresponding to different fusion ratios r. co r For the MOI of (α), select f with the largest MOI. co r (α) as the optimized fusion spectrum f co op (α), expressed as: ; Step 4: Calculate the width N of each window w f below co op (α), denoted as f Nw op (α), choose f Nw op The fusion spectrum with the largest MOI value in (α) is f op (α) is used to reveal the fault characteristics of planetary gear bearings, and its expression is: ; Different window widths N w =2 k (k=6,7,..,16) condition, f Nw op (α) MOI index as Figure 5 As shown in (a), it can be seen that in window width N w =2 12 The fusion spectrum f obtained at that time op (α), with an MOI value of 138.4. When the window width N w =2 12 At that time, the fusion spectrum f corresponding to different fusion ratios r co r The MOI value of (α) is as follows Figure 5 As shown in (b), f op The optimized fusion ratio r=0.19, at this ratio, f op The number of spectral frequency components f(α) included in (α) is Q, which is 389. Their distribution on the α-f plane and the MOI intensity are as follows: Figure 5 As shown in (c), the fusion spectrum f obtained by the proposed method is... op(α) can effectively weaken non-periodic harmonic components, enhance the fault characteristics of planetary gear bearings, and fuse the spectrum f op The MOI of (α) increased by 26% compared to the spectral frequency f(α) component with the largest MOI (110) in the cyclic spectrum. Figure 5 The fused spectrum f(α) obtained by linearly superimposing the spectral frequency f(α) components marked in (c) is... op (α) such as Figure 5 As shown in (d), the first-order fault spectrum of the planetary gear bearing (O) is evident. bpfo =4.55×) and its higher-order spectral lines can be effectively identified, verifying the effectiveness of the proposed method under variable speed conditions.

[0011] To further verify the effectiveness of the method of this invention in extracting fault features of planetary bearings, the original signal was first analyzed using the SAM method, with the weighting index MO set to -0.5 ≤ MO ≤ 1.5. The optimized MO was determined as follows: Figure 6 As shown in (a), MO is 0.8, and the order spectrum corresponding to the optimized parameters is as follows. Figure 6 As shown in (b), the frequency order O s Dominant, while the failure characteristic order of the outer ring of the planetary bearing is O. bpfo Unable to effectively identify. On the other hand, the Infogram method was used to analyze the original IAS signal, with the decomposition level set to 7 levels, and the squared envelope negative entropy ( ), square envelope spectrum negative entropy (ΔI) E ) and average negative entropy (ΔI) 1 / 2 The resonance region (center frequency f) determined by the distribution index c Bandwidth B w )like Figure 7 (a) Figure 7 (c) and Figure 7 As shown in (e), the corresponding order spectrum is as follows: Figure 7 (b) Figure 7 (d) and Figure 7 As shown in (f), the fault characteristic order of the outer ring of the planetary gear bearing remains O. bpfo Unable to be effectively identified.

[0012] This document uses specific examples to illustrate the principles and implementation methods of the present invention. The above examples are only for the purpose of helping to understand the present invention and its core ideas. At the same time, for those skilled in the art, based on the ideas of the present invention, there will be changes in the specific implementation methods and application scope. Therefore, the content of this specification should not be construed as a limitation of the present invention.

Claims

1. A method for detecting faults in planetary gear bearings of a distributed electric vehicle mechanical transmission system, characterized in that, The steps are as follows: (1) Acquire the instantaneous angular displacement and corresponding time signal of the encoder, and calculate the instantaneous angular velocity signal w(θ) using the forward difference method. m ); (2) Set the cyclic spectral window width N w =2 k , k=1, 2, ..., calculate the instantaneous angular velocity signal w(θ) m The cyclic spectrum coherence of CSCoh(α, f) is obtained to obtain a bivariate spectrum of spectral frequency f and cyclic frequency α. (3) The spectral frequency component f(α) is obtained from the cyclic spectral coherence CSCoh(α, f), and its expression is: ； In the formula, f i (α) represents the i-th spectral frequency component; The frequency band β in the spectral frequency component f(α) is extracted using the following formula. h (α); ； In the formula, h represents the fault harmonic order, h = 1, 2, ..., H, where H is the maximum fault harmonic order considered; α t α represents the range of background noise that needs to be evaluated. t The value is (0.5~0.8)α. fault α fault The order of theoretical fault characteristics; To eliminate β h The difference in amplitude (α) is normalized using the following formula to obtain the normalized amplitude spectrum β. h-std (α): ； In the formula, min(·) represents the minimum value, and max(·) represents the maximum value; For the normalized amplitude spectrum β h-std (α) Arrange in descending order to obtain an ordered morphological distribution U h (α), expressed as: ； In the formula, sort(·) represents sorting in descending order; To distinguish between normal fluctuations and abnormally high amplitude points caused by faults or other factors, based on the ordered morphological distribution U h (α) Introduce an adaptive normal upper bound d h The expression is: ； ； In the formula, Cov(·) represents the covariance, Var(·) represents the variance, and M[·] represents the median. d represents the mean. h It's U h Simple linear trend prediction A of (α) h (α)+B h The maximum predicted value, A h B represents the predicted slope term. h Represents the intercept term; Furthermore, the calculation of the upper bound d beyond the normal range is performed. h Outlier impact integral P h The expression is: ； The frequency component f of each spectrum is calculated using the following formula. i The Morphological Outlier Index (MOI) of (α) evaluates the richness of fault information by the ratio of the amplitude of each fault harmonic order to the integral of the outlier impulse. The expression is: ； For each spectral frequency component f i (α) Sort the sequence {f} in descending order according to the corresponding MOI values. n (α)}, the sequence satisfies MOI(f1(α))≥MOI(f2(α))...≥MOI(f Nw / 2 (α)) is obtained by linearly superimposing the spectral frequency components f(α) containing rich fault information according to the fusion ratio r, so as to obtain the fused spectrum f corresponding to the fusion ratio r. co r (α), expressed as: ； In the formula, r = 0.01, 0.02, ... 1, and Q represents the number of spectral frequency components selected by the fusion ratio r; Indicates rounding down; Calculate the fusion spectrum f corresponding to different fusion ratios r. co r For the MOI of (α), select the fusion spectrum f with the largest MOI. co r (α) as the optimized fusion spectrum f co op (α), expressed as: ； Calculate the width N of each window w f below co op (α), denoted as f Nw op (α), choose f Nw op The fusion spectrum with the largest MOI value in (α) is f op (α) is used to reveal the fault characteristics of planetary gear bearings, and its expression is: 。 2. The method for detecting planetary gear bearing faults in a distributed electric vehicle mechanical transmission system according to claim 1, characterized in that, In step (2), w(θ) m The cyclic spectral coherence CSCoh(α, f) is calculated as follows: w(θ m The cyclic spectral correlation (CSC) function is defined as the two discrete Fourier transforms of the instantaneous autocorrelation function: ； In the formula, α represents the cycle frequency, f represents the spectral frequency, Δθ=2π / L represents the angular domain sampling interval, L represents the number of gratings in the encoder, and R... w (θ m , ) represents the signal w(θ) m The autocorrelation function of ) This represents the angular delay, where j is the imaginary unit, and w(θ) m ) represents the time m at angular displacement θ m The instantaneous angular velocity, m=1, 2, 3,…; To eliminate the influence of non-uniform noise, CSC is used for normalization to obtain CSCoh, whose expression is: 。

Images