Gear case coupling fault diagnosis method based on periodic dictionary sparse representation

Abstract

The invention aims to provide a gearbox coupling fault diagnosis method based on periodic dictionary sparse representation, and belongs to the field of gear fault diagnosis. Comprising the following steps: acquiring a gearbox coupling fault signal; calculating the fault impact frequency of each component; calculating an envelope spectrum; selecting a fault component; constructing a periodic dictionary; decoupling the fault signal based on sparse representation; extracting fault features; and determining the fault type of each component. According to the method, whether each component of the gearbox breaks down or not is judged based on the envelope spectrum and the component characteristics, signal decoupling is performed on different fault components in sequence, and the phenomenon of false components or mode aliasing generated in the decoupling process is avoided. According to the invention, based on the periodic dictionary of the periodic characteristics of the fault component, the recognition capability of the sparse representation on the fault pulse is improved, so that the reconstruction signal can accurately restore the fault characteristics of the corresponding fault component. The method has high feature extraction and fault diagnosis precision.

Description

Technical Field

[0001] The present invention relates to a mechanical fault diagnosis method, in particular to a gear fault diagnosis method. Background Art

[0002] Gearbox is a common transmission device, which is widely used in vehicles, ships, generator sets and other environments. Since gearboxes are often affected by high loads and variable working conditions, gearboxes often fail. In order to ensure production efficiency and personnel safety, gearbox fault diagnosis is required.

[0003] The most frequently faulty parts of a gearbox are gears and bearings. Gear faults include broken teeth and wear, and bearing faults include inner ring faults, outer ring faults, rolling element faults, and cage faults. These faults will destroy the original shape of the parts and cause them to impact during operation. The fault diagnosis method based on vibration signals can effectively collect pulse signals generated by impacts. By extracting features from vibration signals containing pulses and using machine learning methods for classification, it can effectively identify the fault types of each component.

[0004] However, when a multi-component coupling fault occurs in a gearbox, the fault pulses of each component overlap with each other, and there are multiple characteristics in the vibration signal, making it difficult to determine the fault type of each component. Summary of the invention

[0005] The object of the present invention is to provide a gearbox coupling fault diagnosis method based on periodic dictionary sparse representation, which can solve the problems of aliasing of fault signals of various components and difficulty in fault diagnosis when a gearbox coupling fault occurs.

[0006] The object of the present invention is achieved in that:

[0007] The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation of the present invention is characterized by comprising the following steps:

[0008] (1) Obtaining gearbox coupling fault signal;

[0009] (2) Calculate the fault impact frequency of each component;

[0010] (3) Calculate the envelope spectrum;

[0011] (4) Select the faulty component;

[0012] (5) Construct a period dictionary;

[0013] (6) Decoupling fault signals based on sparse representation;

[0014] (7) Extract fault features;

[0015] (8) Confirm the fault type of each component.

[0016] The present invention may also include:

[0017] 1. Step (1) is to place the sensor on the surface of the gearbox body to collect the vibration signal generated during the operation of the gearbox. The collected original signal is recorded as x = [x(1), x(2), ..., x(l), ..., x(L)], where x(l) is the amplitude of the lth sampling point, L is the signal length, and the selected sampling frequency is recorded as F s .

[0018] 2. The components in step (2) include gears and bearings, and the fault impact frequency is determined according to the speed of the gearbox and the transmission ratio between the gears:

[0019] First, calculate the fault impact frequency of the gear. In the primary transmission gearbox, the gears include the driving gear and the driven gear. The driving gear is located on the input shaft. Assume that the input shaft speed is n main , then the fault impact frequency f of the driving gear main_g for:

[0020]

[0021] Assuming the transmission ratio between the driving and driven gears is g, the fault impact frequency of the driven gear is f slave_g for:

[0022] f slave_g =f main_g ·g(2)

[0023] Then calculate the bearing fault impact frequency. The bearing fault types include inner ring fault, outer ring fault, rolling element fault, and cage fault. Assume that the diameter of the bearing rolling element is d, the diameter of the circle where the rolling element center is located is D, the number of rolling elements is z, the bearing contact angle is θ, and the speed of the shaft where the bearing is located is n. a , the calculation formula of bearing fault impact frequency is as follows:

[0024] Bearing inner ring fault impact frequency:

[0025]

[0026] Bearing outer ring fault impact frequency:

[0027]

[0028] Bearing rolling element failure impact frequency:

[0029]

[0030] Bearing cage failure impact frequency:

[0031]

[0032] The first-stage transmission gearbox contains input shaft bearings and output shaft bearings. For the input shaft bearings, n a =n main , for output shaft bearing, n a =n main ·g, n a Substitute the above formulas for the impact frequency of the inner ring fault of the bearing, the impact frequency of the outer ring fault of the bearing, the impact frequency of the rolling element fault of the bearing, and the impact frequency of the bearing cage fault into the above formulas, and the impact frequency of the bearing fault includes: the impact frequency of the inner ring fault of the input shaft bearing f main_bi , input shaft bearing outer ring fault impact frequency f main_bo , input shaft bearing rolling element fault impact frequency f main_br , Input shaft bearing retainer fault impact frequency f main_bc , output shaft bearing inner ring fault impact frequency f slave_bi , output shaft bearing outer ring fault impact frequency f slave_bo , output shaft bearing rolling element failure impact frequency f slave_br , output shaft bearing retainer failure impact frequency f slave_bc ;

[0033] The fault impact frequency of each component is calculated as [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ].

[0034] 3. The specific steps of step (3) are:

[0035] Perform envelope spectrum analysis on the original signal x = [x(1), x(2), ..., x(l), ..., x(L)]:

[0036] Perform the Hilbert transform on the signal x:

[0037]

[0038] get

[0039] Construct the parsing signal:

[0040]

[0041] Get the analysis signal

[0042] beg Module length:

[0043]

[0044] Get envelope signal env=[env(1),env(2),…,env(l),…,env(L)];

[0045] Take Fourier transform of env:

[0046]

[0047] The envelope spectrum e = [e(1), e(2),…, e(r),…e(L)] is obtained.

[0048] 4. The process of step (4) is: take the frequency with an amplitude greater than five times the mean of the envelope spectrum as the envelope peak frequency, and record the envelope peak frequency as [f p1 ,f p2 ,f p3 ,…,f pn ,…,f pN ], N is the number of envelope spectrum peak frequencies;

[0049] Select an envelope peak frequency f pn , compared with f pn The impact frequency of each component failure [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ], if f pn In [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ], it indicates that the component corresponding to the fault impulse frequency has failed. If the same value is not found, it indicates that f pn If it is a noise component or fault frequency doubling, the next envelope peak frequency will be reselected for comparison.

[0050] 5. For the components determined to be faulty, a dictionary is constructed for each component: if the faulty component is a bearing, a periodic dictionary is constructed based on the Laplace wavelet; if the faulty component is a gear, a dictionary is constructed based on the Morlet wavelet:

[0051] The process of constructing the Laplace wavelet dictionary is:

[0052] The analytical expression of Laplace wavelet is:

[0053]

[0054] In the formula, τ l is the delay parameter, ζ l is the damping coefficient, f l is the natural frequency, T L is the Laplace wavelet length;

[0055] Let T L =F s / f pn , τ l =1,2,`…,T L , l =0.01,0.02,…,0.1,f l =1kHz,2kHz,…,10kHz; set the parameter [τ l ,ζ l ,f l ] are combined to generate the Laplace wavelet dictionary atom, denoted by d L1 ,d L2 ,...,d Lk ,...,d LK , where d Lk is the kth dictionary atom, K is the number of dictionary atoms, and its value is the total number of variable combinations, that is, K = 100·T L ; Construct the Laplace wavelet dictionary D based on each dictionary atom L :

[0056] D L =[d L1 ,d L2 ,...,d Lk ,...,d LK ];

[0057] The process of constructing the Morlet wavelet dictionary is:

[0058] The analytical expression of Morlet wavelet is:

[0059]

[0060] In the formula, s m 、u m , m ,θ m is the time-frequency parameter of Morlet wavelet, T M is the Morlet wavelet length;

[0061] Let T M =F s / f pn ,u m =1,2,…,T M ,s m =1kHz,2kHz,…,10kHz,ξ m =0.1,0.2,…,1; set the parameter [s m ,u m ,ξ m ] are combined to generate the Morlet wavelet dictionary atom, denoted as d M1 ,d M2 ,...,d Mj ,...,d MJ , where d MJ is the jth dictionary atom, J is the number of dictionary atoms, J = 100T M ; Construct each Moelet dictionary atom into a Morlet dictionary D M :

[0062] D M =[d M1 ,d M2 ,...,d Mj ,...,d MJ ];

[0063] Extend with the impact cycle of the component as the interval:

[0064] The Laplace wavelet dictionary or the Morlet wavelet dictionary is converted to f pn To extend the cycle;

[0065] When the extended dictionary is a Laplace wavelet dictionary, the generated periodic dictionary atom d pk for:

[0066]

[0067] In the formula, m is the number of period extensions, and the period dictionary obtained by extension is:

[0068] D p =[d p1 ,d p2 ,...,d pk,...,d pK ]

[0069] In the formula, L d = m·T L is the length of the dictionary;

[0070] When the Morlet wavelet dictionary is extended, d pj By m d Mj The periodic dictionary obtained by extension is:

[0071] D p =[d p1 ,d p2 ,...,d pj ,...,d pJ ]

[0072] In the formula, L d = m·T M is the length of the dictionary.

[0073] 6. The process of step (6) is: Reconstruct the original signal based on sparse representation segmentation: Let the dictionary D p Size L d ×K, substitute the dictionaries of other sizes into the number of rows and columns to solve;

[0074] First, take the step size as 1 and cut the length L from the original signal d fragments, denoted by R1, R2, …, R i ,…,

[0075] Then the signal fragment R is reconstructed using the sparse representation i , the sparse representation is obtained by solving the sparse coefficient α∈R k×1 , R i Indicated as D p The product of the sparse coefficient α is used to obtain the reconstructed signal fragment

[0076]

[0077] α is solved by the orthogonal matching pursuit (OMP) algorithm, the steps are:

[0078] Set the initial residual ε=R i First, calculate the inner product of each dictionary atom and the residual and find the optimal atom d b :

[0079]

[0080] Setting support sets Db Merge into S p :

[0081] S p =S p ∪d b

[0082] Determine α so that S p Accurately restore R i :

[0083]

[0084] Update residuals:

[0085] ε=R i -S p α

[0086] Set the iteration threshold to δ o , if ||ε|| 2 >δ o , substituting ε into the formula Calculate and update the parameters in the above formula; if ||ε|| 2 <δ o , then the loop ends and α is obtained;

[0087] Finally, the first digits of each reconstructed signal segment are added in sequence and the average is calculated to obtain the reconstructed signal The vibration signal of the faulty component is also obtained by decoupling:

[0088]

[0089] 7. The process of step (7) is: get Then, calculate its RMS:

[0090]

[0091] calculate The kurtosis is:

[0092]

[0093] The fault characteristic of the corresponding component is [rms, kur];

[0094] After completing the feature extraction of the current faulty component, if there is still an envelope peak to be identified, the next envelope peak to be identified is selected, the corresponding faulty component is determined, the periodic dictionary is constructed, and decoupling and feature extraction are performed.

[0095] 8. The process of step (8) is:

[0096] Assume that q groups of features [rms, kur] are collected, and the K-means clustering algorithm is used to determine the fault type of each component. The process is as follows:

[0097] ① Set the maximum number of iterations r max , let the number of clusters be P. For a primary transmission gearbox, there are four components, each of which includes four types of faults: broken teeth, wear, pitting, and other faults of gears, and inner ring, outer ring, rolling element, and cage faults of bearings. Therefore, P = 4 × 4 = 16;

[0098] ② Represent each group of features in the rectangular coordinate plane composed of the rms axis and the kur axis, and randomly initialize P cluster centers c1 (0) ,c2 (0) ,…,c P (0) ;

[0099] ③ Calculate the Euclidean distance between each sample and the cluster center, and assign the sample to the cluster closest to it. The Euclidean distance calculation formula is:

[0100]

[0101] ④ Update the position of each cluster center. The coordinates of the new cluster center on the rms axis and the kur axis are the mean values ​​of all points in the cluster on the corresponding coordinate axes.

[0102] ⑤If the updated cluster center position has changed, go to step ②; if it has not changed, go to step ⑥;

[0103] ⑥When the cluster center position stops changing or reaches the maximum number of iterations, output the clustering results;

[0104] After clustering, the q samples are divided into P classes, 20% of the samples are taken as the training set, the fault type of the training set is input, and the fault type of all samples in each class is defined as the fault type of the majority of training set samples in the class, then the fault type of the remaining 80% of the test set samples is determined.

[0105] The advantages of the present invention are:

[0106] (1) The present invention determines whether each component of the gearbox is faulty based on the envelope spectrum and component characteristics, and sequentially decouples signals of different faulty components, thereby avoiding the generation of false components or modal aliasing during the decoupling process.

[0107] (2) The present invention improves the analytical dictionary of sparse representation and proposes a periodic dictionary based on the periodic characteristics of the faulty component, thereby improving the sparse representation's ability to identify fault pulses, so that the reconstructed signal can accurately restore the fault characteristics of the corresponding faulty component.

[0108] (3) The present invention performs feature extraction based on the decoupled fault signal and determines the fault type according to the K-means algorithm, and has high feature extraction and fault diagnosis accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0109] Figure 1 is a flow chart of the present invention;

[0110] Figure 2a It is the fault signal of the driving gear broken tooth. Figure 2b This is the driven gear wear fault signal. Figure 2c It is the output shaft bearing rolling element fault signal. Figure 2d is the gear meshing signal, Figure 2e is the noise signal, Figure 2f is the coupling fault signal x;

[0111] Figure 3a is the envelope spectrum of the original signal, Figure 3b It is a local enlarged view of the original signal envelope spectrum;

[0112] Figure 4 Partial atoms for periodic dictionaries;

[0113] Figure 5a Reconstruct the signal for the bearing fault (Reconstructed Signal), Figure 5b Reconstruct the signal (residual) for the bearing fault;

[0114] Figure 6a Reconstruct the signal for the active gear failure (Reconstruct the signal), Figure 6b Reconstruct the signal (residual) for the active gear fault;

[0115] Figure 7a Reconstruct the signal for the driven gear fault (Reconstruct the signal), Figure 7b Reconstruct the signal (residual) for the driven gear fault;

[0116] Figure 8 is the fault type of each decoupled sample. DETAILED DESCRIPTION

[0117] The present invention is described in more detail below with reference to the accompanying drawings:

[0118] Combination Figure 1-8 The present invention is a gearbox coupling fault diagnosis method based on sparse representation of periodic dictionary. First, the gearbox coupling fault signal is collected, and then the fault impact frequency of each component is calculated to determine the fault condition of each component. A periodic dictionary is constructed based on the periodic characteristics of the faulty component, and then the fault signal is reconstructed using sparse representation to achieve fault signal decoupling. Then the fault characteristics of each component are extracted, and the fault type of each component is determined through a clustering algorithm to complete the gearbox coupling fault diagnosis. Specifically, the following steps are included: 1. Obtain the gearbox coupling fault signal

[0119] The sensor is placed on the surface of the gearbox body to collect the vibration signal generated during the operation of the gearbox. The collected original signal is recorded as x = [x(1), x(2), ..., x(l), ..., x(L)], x(l) is the amplitude of the lth sampling point, L is the signal length, and the selected sampling frequency is recorded as F s .

[0120] 2 Calculate the failure impact frequency of each component

[0121] The high-frequency failure components of the gearbox are gears and bearings, and their failure impact frequency can be determined according to the speed of the gearbox and the transmission ratio between gears.

[0122] First, calculate the fault impact frequency of the gear. Take the first-stage transmission gearbox as an example. The gears include the driving gear and the driven gear. The driving gear is located on the input shaft. Assume that the input shaft speed is n main , then the fault impact frequency f of the driving gear main_g for:

[0123]

[0124] Assuming the transmission ratio between the driving and driven gears is g, the fault impact frequency of the driven gear is f slave_g for:

[0125] f slave_g =f main_g ·g(2)

[0126] Then calculate the bearing fault impact frequency. Bearing fault types include inner ring fault, outer ring fault, rolling element fault, and cage fault. Assume that the diameter of the bearing rolling element is d, the diameter of the circle where the rolling element center is located is D, the number of rolling elements is z, the bearing contact angle is θ, and the speed of the shaft where the bearing is located is n a The calculation formula for the bearing fault impact frequency is as follows.

[0127] Bearing inner ring fault impact frequency:

[0128]

[0129] Bearing outer ring fault impact frequency:

[0130]

[0131] Bearing rolling element failure impact frequency:

[0132]

[0133] Bearing cage failure impact frequency:

[0134]

[0135] The primary transmission gearbox contains input shaft bearings and output shaft bearings. For the input shaft bearings, n a =n main , for output shaft bearing, n a =n main ·g. a Substituting into equations (3)-(6) respectively, the bearing fault impact frequency includes: input shaft bearing inner ring fault impact frequency f main_bi , input shaft bearing outer ring fault impact frequency f main_bo , input shaft bearing rolling element fault impact frequency f main_br , Input shaft bearing retainer fault impact frequency f main_bc , output shaft bearing inner ring fault impact frequency f slave_bi , output shaft bearing outer ring fault impact frequency f slave_bo , output shaft bearing rolling element failure impact frequency f slave_br , output shaft bearing retainer failure impact frequency f slave_bc There are 8 types in total.

[0136] In summary, the fault impact frequency of each component is calculated as [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ].

[0137] 3 Calculate the envelope spectrum

[0138] The envelope spectrum analysis of the original signal x = [x(1), x(2), ..., x(l), ..., x(L)] is performed as follows:

[0139] ① Perform Hilbert transform on signal x:

[0140]

[0141] get

[0142] ②Construct analysis signal:

[0143]

[0144] Get the analysis signal

[0145] ③Seek Module length:

[0146]

[0147] The envelope signal env=[env(1),env(2),…,env(l),…,env(L)].

[0148] ④ Perform Fourier transform on env:

[0149]

[0150] The envelope spectrum e = [e(1), e(2),…, e(r),…e(L)] is obtained.

[0151] 4 Select the faulty component

[0152] After calculating the envelope spectrum, the frequency with an amplitude greater than five times the mean of the envelope spectrum is taken as the envelope peak frequency, and the envelope peak frequency is recorded as [f p1 ,f p2 ,f p3 ,…,f pn ,…,f pN ], N is the number of peak frequencies of the envelope spectrum.

[0153] Select an envelope peak frequency f pn , compared with f pn The impact frequency of each component failure [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ], if f pn In [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ], it indicates that the component corresponding to the fault impulse frequency has failed. If the same frequency is not found, it indicates that f pn If it is a noise component or fault frequency doubling, the next envelope peak frequency will be reselected for comparison.

[0154] For the components that are determined to be faulty, a dictionary is constructed for each component. If the faulty component is a bearing, a periodic dictionary is constructed based on the Laplace wavelet. If the faulty component is a gear, a dictionary is constructed based on the Morlet wavelet. The specific method is introduced in the next section.

[0155] 5. Build a cycle dictionary

[0156] 5.1 Constructing Laplace wavelet dictionary

[0157] If the component determined to be faulty is a bearing, a dictionary is constructed based on the Laplace wavelet. The analytical expression of the Laplace wavelet is:

[0158]

[0159] In the formula, τ l is the delay parameter, ζ l is the damping coefficient, f l is the natural frequency, T L is the Laplace wavelet length.

[0160] Let T L =F s / f pn , τ l =1,2,`…,T L , l =0.01,0.02,…,0.1,f l =1kHz,2kHz,…,10kHz. l ,ζ l ,f l ] are combined to generate multiple Laplace wavelet dictionary atoms, denoted as d L1 ,d L2 ,...,d Lk ,...,d LK , where d Lk is the kth dictionary atom, K is the number of dictionary atoms, and its value is the total number of variable combinations, that is, K = 100·T L . Construct the Laplace wavelet dictionary D based on each dictionary atom L :

[0161] D L =[d L1 ,d L2 ,...,d Lk ,...,d LK ](12)

[0162] 5.2 Constructing the Morlet Wavelet Dictionary

[0163] If f pn The corresponding component is a gear, so the dictionary is constructed based on Morlet wavelet. The analytical expression of Morlet wavelet is:

[0164]

[0165] In the formula, s m 、u m , m ,θ m is the time-frequency parameter of Morlet wavelet, T M is the Morlet wavelet length.

[0166] Let T M =F s / f pn ,u m =1,2,…,T M ,s m =1kHz,2kHz,…,10kHz,ξ m =0.1,0.2,…,1. m ,u m ,ξ m ] are combined to generate multiple Morlet wavelet dictionary atoms, denoted as d M1 ,d M2 ,...,d Mj ,...,d MJ , where d MJ is the jth dictionary atom, J is the number of dictionary atoms, J = 100T M . Each Moelet dictionary atom is constructed into a Morlet dictionary D M :

[0167] D M =[d M1 ,d M2 ,...,d Mj ,...,d MJ ](14)

[0168] 5.3 Extending the impact cycle of the component

[0169] The Laplace wavelet dictionary or the Morlet wavelet dictionary is converted to f pn When the extended dictionary is a Laplace wavelet dictionary, the generated periodic dictionary atom d pk for:

[0170]

[0171] Where m is the number of period extensions. The period dictionary obtained by extension is:

[0172] D p =[d p1 ,d p2 ,...,dpk ,...,d pK ](16)

[0173] In the formula, L d = m·T L is the length of the dictionary.

[0174] When the Morlet wavelet dictionary is extended, d pj By m d Mj The periodic dictionary obtained by extension is:

[0175] D p =[d p1 ,d p2 ,...,d pj ,...,d pJ ](16)

[0176] In the formula, L d = m·T M is the length of the dictionary.

[0177] 6 Decoupling fault signals based on sparse representation

[0178] The original signal is reconstructed based on sparse representation. Since the periodic dictionary sizes of different components are different, in this section, the dictionary D p Size L d ×K. Dictionaries of other sizes can be substituted into the number of rows and columns to find the solution.

[0179] First, take the step size as 1 and cut the length L from the original signal d fragments, denoted by R1, R2, …, R i ,…,

[0180] Then the signal fragment R is reconstructed using the sparse representation i The sparse representation is obtained by solving the sparse coefficient α∈R k×1 , R i Indicated as D p The product of the sparse coefficient α is used to obtain the reconstructed signal fragment

[0181]

[0182] Among them, α is solved by the orthogonal matching pursuit (OMP) algorithm, the steps are:

[0183] ①Set the initial residual ε=R i First, calculate the inner product of each dictionary atom and the residual and find the optimal atom d b :

[0184]

[0185] ②Set support set D b Merge into S p :

[0186] S p =S p ∪d b (19)

[0187] ③Determine α so that S p Accurately restore R i :

[0188]

[0189] ④Update residual:

[0190] ε=R i -S p α(21)

[0191] ⑤Set the iteration threshold to δ o , if ||ε|| 2 >δ o , then substitute ε into equation (18), calculate and update the parameters in equations (18) to (21); if ||ε|| 2 <δ o , then the loop ends and α is obtained.

[0192] Finally, the first digits of each reconstructed signal segment are added in sequence and the average is calculated to obtain the reconstructed signal The vibration signal of the faulty component is also obtained by decoupling:

[0193]

[0194] 7 Extracting fault features

[0195] get Then, calculate its RMS:

[0196]

[0197] calculate The kurtosis is:

[0198]

[0199] The fault characteristic of the corresponding component is [rms, kur].

[0200] After completing the feature extraction of the current faulty component, if there is still an envelope peak to be identified, the next envelope peak to be identified is selected, the corresponding faulty component is determined, the cycle dictionary is constructed, and decoupling and feature extraction are performed. For a primary transmission gearbox, all the feature types that can be extracted are expressed as follows:

[0201] Driving gear characteristics [rms main_g ,kur main_g ]、Driven gear characteristics [rms slave_g ,kur slave_g ]、Input shaft bearing characteristics [rms main_b ,kur main_b ]、Output shaft bearing characteristics [rms slave_b ,kur slave_b ].

[0202] 8. Determine the fault type of each component

[0203] Assume that q groups of features [rms, kur] are collected, and the K-means clustering algorithm is used to determine the fault type of each component. The process is as follows:

[0204] ① Set the maximum number of iterations r max =100, let the number of clusters be P. For the primary transmission gearbox, there are four components, each of which includes four types of faults (tooth breakage, wear, pitting, other faults of gears, inner ring, outer ring, rolling element, and cage faults of bearings), so P = 4 × 4 = 16.

[0205] ② Represent each group of features in the rectangular coordinate plane composed of the rms axis and the kur axis, and randomly initialize P cluster centers c1 (0) ,c2 (0) ,…,c P (0) .

[0206] ③ Calculate the Euclidean distance between each sample and the cluster center, and assign the sample to the cluster closest to it. The Euclidean distance calculation formula is:

[0207]

[0208] ④ Update the position of each cluster center. The coordinates of the new cluster center on the rms axis and the kur axis are the means of all points in the cluster on the corresponding coordinate axes.

[0209] ⑤If the updated cluster center position has changed, go to step ②; if it has not changed, go to step ⑥.

[0210] ⑥When the cluster center position stops changing or reaches the maximum number of iterations, output the clustering results.

[0211] After clustering, the q samples are divided into P classes, 20% of the samples are taken as the training set, the fault type of the training set is input, and the fault type of all samples in each class is defined as the fault type of the majority of training set samples in the class, then the fault type of the remaining 80% of the test set samples is determined.

[0212] Example:

[0213] Step 1: Get the gear coupling fault signal

[0214] Assume that a certain level of transmission gearbox has three faults: driving gear broken teeth, driven gear broken teeth, and output shaft bearing rolling element failure. Assume that the structural parameters of each component of the gearbox are as shown in Table 1, and the sampling frequency F s =20.48kHz, the signal sample duration is 1s, then the signal sample contains 20 broken tooth fault pulses, 17 wear fault pulses, and 30 bearing fault pulses. The pulses of each fault type have different waveforms, and also contain gear meshing signals and background noise. Based on the above information, the gearbox coupling fault simulation signal x∈R is constructed. 20480×1 ,The waveforms of the coupled fault signal x and its components are shown in Figure 2.

[0215] Table 1 Gearbox structural parameters

[0216] Structural parameters Value Structural parameters Value Input shaft speed 1200rpm Number of teeth on driving gear 17 Driven gear teeth number 20 Driving and driven wheel transmission ratio 17:20 Bearing contact angle 0° Number of bearing rolling elements 9 Rolling element center diameter 60mm Rolling element diameter 27mm

[0217] Step 2: Calculate the fault impulse frequency

[0218] Fault impact frequency f of driving gear main_g for:

[0219]

[0220] Fault impact frequency f of the driven gear slave_g for:

[0221]

[0222] Bearing inner ring fault impact frequency:

[0223]

[0224] Bearing outer ring fault impact frequency:

[0225]

[0226] Bearing rolling element failure impact frequency:

[0227]

[0228] Bearing cage failure impact frequency:

[0229]

[0230] For input shaft bearings, n a =n main , n a Substituting into equations (103)-(106) respectively, the input shaft bearing inner ring fault impact frequency f main_bi =130Hz, input shaft bearing outer ring fault impact frequency f main_bo =49Hz, input shaft bearing rolling element fault impact frequency f main_br =35Hz, input shaft bearing retainer fault impact frequency f main_bc =9Hz.

[0231] For output shaft bearings, n a =n main ·g, n a Substituting into equations (103)-(106) respectively, the output shaft bearing inner ring fault impact frequency f slave_bi =111Hz, output shaft bearing outer ring fault impact frequency f slave_bo =42Hz, output shaft bearing rolling element fault impact frequency f slave_br =30Hz, output shaft bearing retainer fault impact frequency f slave_bc =8Hz.

[0232] The statistics of the fault impact frequency of each component calculated above are shown in Table 2.

[0233] Table 2 Fault impact frequency of each component

[0234] Fault impact frequency Value Fault impact frequency Value <![CDATA[f main_g ]]> 20Hz <![CDATA[f slave_g ]]> 17Hz <![CDATA[f main_bi ]]> 130Hz <![CDATA[f main_bo ]]> 49Hz <![CDATA[f main_br ]]> 36Hz <![CDATA[f main_bc ]]> 9Hz <![CDATA[f slave_bi ]]> 111Hz <![CDATA[f slave_bo ]]> 41Hz <![CDATA[f slave_br ]]> 30Hz <![CDATA[f slave_bc ]]> 8Hz

[0235] Step 3 Calculate the envelope spectrum

[0236] The steps for envelope spectrum analysis of the original signal x are as follows:

[0237] ① Perform Hilbert transform on signal x:

[0238]

[0239] ②Construct analysis signal:

[0240]

[0241] ③Seek Module length:

[0242]

[0243] 4④ Perform Fourier transform on env:

[0244]

[0245] The envelope spectrum e = [e(1), e(2), …, e(r), …e(L)] is obtained, as shown in FIG3 .

[0246] Step 4: Select the faulty component

[0247] Ten envelope peaks can be extracted from the envelope spectrum, namely [f n1 =17Hz,f n2 =20Hz,f n3 =30Hz,f n4 =34Hz,f n5 =40Hz,f n6 =60Hz,f n7 =90Hz,f n8 =120Hz,f n9 =150Hz,f n10 =180Hz]. Select f n1 =17Hz, comparing it with the fault impulse frequency of each component in Table 2, we can know that f n1 =f slave_g = 17Hz indicates that the driven gear is faulty. In addition, f n1 =f main_g =20Hz, indicating that the driving gear is faulty; f n3 =f slave_br =30Hz, f n4 =f slave_bo =40Hz indicates output shaft bearing failure, but the specific fault type cannot be determined. The remaining envelope peaks [f n4 ,f n6 ,f n7 ,f n8 ,f n9 ,f n10 ]There is no corresponding fault characteristic frequency, indicating that it is a multiple of the noise peak or fault frequency.

[0248] Step 5: Build a cycle dictionary

[0249] 5.1 Constructing Laplace wavelet dictionary

[0250] If the faulty component is a bearing, a dictionary is constructed based on Laplace wavelet. The analytical expression of Laplace wavelet is:

[0251]

[0252] In the formula, τ l is the delay parameter, ζ l is the damping coefficient, f l is the natural frequency, T L is the Laplace wavelet length.

[0253] Let T L =F s / f pn , τ l =1,2,`…,T L , l =0.01,0.02,…,0.1,f l =1kHz,2kHz,…,10kHz. l ,ζ l ,f l ] are combined to generate multiple Laplace wavelet dictionary atoms, denoted as d L1 ,d L2 ,...,d Lk ,...,d LK , where d Lk is the kth dictionary atom, K is the number of dictionary atoms, and its value is the total number of variable combinations, that is, K = 100·T L . Construct the Laplace wavelet dictionary D based on each dictionary atom L :

[0254] D L =[d L1 ,d L2 ,...,d Lk ,...,d LK ](112)

[0255] 5.2 Constructing the Morlet Wavelet Dictionary

[0256] If the faulty component is a gear, the dictionary is constructed based on Morlet wavelet. The analytical expression of Morlet wavelet is:

[0257]

[0258] In the formula, s m 、u m , m ,θ m is the time-frequency parameter of Morlet wavelet, T M is the Morlet wavelet length.

[0259] Let T M =F s / f pn ,u m =1,2,…,T M ,s m =1kHz,2kHz,…,10kHz,ξ m =0.1,0.2,…,1.m ,u m ,ξ m ] are combined to generate multiple Morlet wavelet dictionary atoms, denoted as d M1 ,d M2 ,...,d Mj ,...,d MJ , where d MJ is the jth dictionary atom, J is the number of dictionary atoms, J = 100T M . Each Moelet dictionary atom is constructed into a Morlet dictionary D M :

[0260] D M =[d M1 ,d M2 ,...,d Mj ,...,d MJ ](114)

[0261] 5.3 Extending the impact cycle of the component

[0262] Extend the Laplace wavelet dictionary to generate the periodic dictionary atom d pk :

[0263]

[0264] In the formula, m is the number of period extensions, and m=3. The period dictionary obtained by extension is:

[0265] D p =[d p1 ,d p2 ,...,d pk ,...,d pK ](116)

[0266] In the formula, L d = m·T L is the length of the dictionary. Cyclic dictionary atoms such as Figure 4 .

[0267] Take m = 3, extend the Morlet wavelet dictionary, and generate the periodic dictionary atom d pj :

[0268]

[0269] The periodic dictionary obtained by extension is:

[0270] D p =[d p1 ,d p2 ,...,dpj ,...,d pJ ](116)

[0271] In the formula, L d = m·T M is the length of the dictionary. The resulting partial periodic dictionary atoms are as follows Figure 4 Step 6: Decouple fault signals based on sparse representation

[0272] First, take the step size as 1 and cut the length L from the original signal d fragments, denoted by R1, R2, …, R i ,…,

[0273] Then the signal fragment R is reconstructed using the sparse representation i The sparse representation is obtained by solving the sparse coefficient α∈R k×1 , R i Indicated as D p The product of the sparse coefficient α is used to obtain the reconstructed signal fragment

[0274]

[0275] Among them, α is solved by the orthogonal matching pursuit (OMP) algorithm, the steps are:

[0276] ①Set the initial residual ε=R i First, calculate the inner product of each dictionary atom and the residual and find the optimal atom d b :

[0277]

[0278] ②Set support set D b Merge into S p :

[0279] S p =S p ∪d b (119)

[0280] ③Determine α so that S p Accurately restore R i :

[0281]

[0282] ④Update residual:

[0283] ε=R i -S p α(121)

[0284] ⑤Set the iteration threshold to δ o , if ||ε|| 2 >δ o , then substitute ε into equation (118), calculate and update the parameters in equations (118) to (121); if ||ε|| 2 <δ o , then the loop ends and α is obtained.

[0285] Finally, the first digits of each reconstructed signal segment are added in sequence and the average is calculated to obtain the reconstructed signal The vibration signal of the faulty component is also obtained by decoupling:

[0286]

[0287] The reconstructed signal of the output shaft bearing fault is shown in Figure 5, the reconstructed signal of the driving gear is shown in Figure 6, and the reconstructed signal of the driven gear is shown in Figure 7.

[0288] Step 7: Extract fault features

[0289] get Then, calculate its RMS:

[0290]

[0291] calculate The kurtosis is:

[0292]

[0293] The calculated driving gear characteristic of this simulation signal [rms main_g =0.27, kur main_g =6.82], driven gear characteristics [rms slave_g =0.21,kur slave_g =10.22], output shaft bearing characteristics [rms slave_b =0.33, kur slave_b =18.74].

[0294] Step 8: Determine the fault type of each component

[0295] Continue to generate 30 fault signal samples, each of which is a gearbox coupling fault signal. The fault types include broken teeth, wear, pitting, cracks, inner ring, outer ring, rolling element, and cage faults of the bearing. The faulty components include the driving gear, driven gear, input shaft bearing, and output shaft bearing. According to the method proposed in steps 1 to 6, the fault features are extracted, and finally 160 sets of features [rms, kur] are collected. The K-means clustering algorithm is used to determine the fault type of each component. The process is as follows:

[0296] ① Set the maximum number of iterations r max =100, and the number of clusters is 16.

[0297] ② Represent each group of features in the rectangular coordinate plane composed of the rms axis and the kur axis, and randomly initialize P cluster centers c1 (0) ,c2 (0) ,…,c P (0) .

[0298] ③ Calculate the Euclidean distance between each sample and the cluster center, and assign the sample to the cluster closest to it. The Euclidean distance calculation formula is:

[0299]

[0300] ④ Update the position of each cluster center. The coordinates of the new cluster center on the rms axis and the kur axis are the means of all points in the cluster on the corresponding coordinate axes.

[0301] ⑤If the updated cluster center position has changed, go to step ②; if it has not changed, go to step ⑥.

[0302] ⑥When the cluster center position stops changing or reaches the maximum number of iterations, output the clustering results.

[0303] The 160 samples after clustering are divided into 16 classes, 20% of the samples are taken as the training set, the fault type of the training set is input, and the fault type of all samples in each class is defined as the fault type of the majority of training set samples in the class, then the fault type of the remaining 80% of the test set samples is determined. Figure 8 shown.

Claims

1. A gearbox coupling fault diagnosis method based on periodic dictionary sparse representation is characterized by: The steps include: (1) Obtaining gearbox coupling fault signal; (2) Calculate the fault impact frequency of each component; (3) Calculate the envelope spectrum; (4) Select the faulty component; (5) Construct a period dictionary; (6) Decoupling fault signals based on sparse representation; (7) Extract fault features; (8) Confirm the fault type of each component.

2. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: Step (1) specifically comprises placing a sensor on the surface of the gearbox body to collect the vibration signal generated during the operation of the gearbox. The collected original signal is recorded as x = [x(1), x(2), ..., x(l), ..., x(L)], where x(l) is the amplitude of the lth sampling point, L is the signal length, and the selected sampling frequency is recorded as F s .

3. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: The components in step (2) include gears and bearings, and the fault impact frequency is determined according to the rotation speed of the gearbox and the transmission ratio between the gears: First, calculate the fault impact frequency of the gear. In the primary transmission gearbox, the gears include the driving gear and the driven gear. The driving gear is located on the input shaft. Assume that the input shaft speed is n main , then the fault impact frequency f of the driving gear main_g for: Assuming the transmission ratio between the driving and driven gears is g, the fault impact frequency of the driven gear is f slave_g for: f slave_g =f main_g ·g Then calculate the bearing fault impact frequency. The bearing fault types include inner ring fault, outer ring fault, rolling element fault, and cage fault. Assume that the diameter of the bearing rolling element is d, the diameter of the circle where the rolling element center is located is D, the number of rolling elements is z, the bearing contact angle is θ, and the speed of the shaft where the bearing is located is n. a , the calculation formula of bearing fault impact frequency is as follows: Bearing inner ring fault impact frequency: Bearing outer ring fault impact frequency: Bearing rolling element failure impact frequency: Bearing cage failure impact frequency: The first-stage transmission gearbox contains input shaft bearings and output shaft bearings. For the input shaft bearings, n a =n main , for output shaft bearing, n a =n main ·g, n a Substitute the above formulas for the impact frequency of the inner ring fault of the bearing, the impact frequency of the outer ring fault of the bearing, the impact frequency of the rolling element fault of the bearing, and the impact frequency of the bearing cage fault into the above formulas, and the impact frequency of the bearing fault includes: the impact frequency of the inner ring fault of the input shaft bearing f main_bi , input shaft bearing outer ring fault impact frequency f main_bo , input shaft bearing rolling element fault impact frequency f main_br , Input shaft bearing retainer fault impact frequency f main_bc , output shaft bearing inner ring fault impact frequency f slave_bi , output shaft bearing outer ring fault impact frequency f slave_bo , output shaft bearing rolling element failure impact frequency f slave_br , output shaft bearing retainer failure impact frequency f slave_bc ; The fault impact frequency of each component is calculated as [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ].

4. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: The specific steps of step (3) are: Perform envelope spectrum analysis on the original signal x = [x(1), x(2), ..., x(l), ..., x(L)]: Perform the Hilbert transform on the signal x: get Construct the parsing signal: Get the analysis signal beg Module length: Get envelope signal env=[env(1),env(2),…,env(l),…,env(L)]; Take Fourier transform of env: The envelope spectrum e = [e(1), e(2),…, e(r),…e(L)] is obtained.

5. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: The process of step (4) is: take the frequency with an amplitude greater than five times the mean of the envelope spectrum as the envelope peak frequency, and record the envelope peak frequency as [f p1 ,f p2 ,f p3 ,…,f pn ,…,f pN ], N is the number of envelope spectrum peak frequencies; Select an envelope peak frequency f pn , compared with f pn The impact frequency of each component failure [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ], if f pn In [f main_g ,f slave_g ,f main_bi ,f main_bo ,f main_br ,f main_bc ,f slave_bi ,f slave_bo ,f slave_br ,f slave_bc ], it indicates that the component corresponding to the fault impulse frequency has failed. If the same value is not found, it indicates that f pn If it is a noise component or fault frequency doubling, the next envelope peak frequency will be reselected for comparison.

6. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: For the components determined to be faulty, a dictionary is constructed for each component: if the faulty component is a bearing, a periodic dictionary is constructed based on the Laplace wavelet; If the faulty component is a gear, a dictionary is constructed based on the Morlet wavelet: The process of constructing the Laplace wavelet dictionary is: The analytical expression of Laplace wavelet is: In the formula, τ l is the delay parameter, ζ l is the damping coefficient, f l is the natural frequency, T L is the Laplace wavelet length; Let T L =F s / f pn , τ l =1,2,`…,T L , l =0.01,0.02,…,0.1,f l =1kHz,2kHz,…,10kHz; set the parameter [τ l ,ζ l ,f l ] are combined to generate the Laplace wavelet dictionary atom, denoted by d L1 ,d L2 ,...,d Lk ,...,d LK , where d Lk is the kth dictionary atom, K is the number of dictionary atoms, and its value is the total number of variable combinations, that is, K = 100·T L ; Construct the Laplace wavelet dictionary D based on each dictionary atom L : D L =[d L1 ,d L2 ,...,d Lk ,...,d LK ]; The process of constructing the Morlet wavelet dictionary is: The analytical expression of Morlet wavelet is: In the formula, s m 、u m , m ,θ m is the time-frequency parameter of Morlet wavelet, T M is the Morlet wavelet length; Let T M =F s / f pn ,u m =1,2,…,T M ,s m =1kHz,2kHz,…,10kHz,ξ m =0.1,0.2,…,1; set the parameter [s m ,u m ,ξ m ] are combined to generate the Morlet wavelet dictionary atom, denoted as d M1 ,d M2 ,...,d Mj ,...,d MJ , where d MJ is the jth dictionary atom, J is the number of dictionary atoms, J = 100T M ; Construct each Moelet dictionary atom into a Morlet dictionary D M : D M =[d M1 ,d M2 ,...,d Mj ,...,d MJ ]; Extend with the impact cycle of the component as the interval: The Laplace wavelet dictionary or the Morlet wavelet dictionary is converted to f pn To extend the cycle; When the extended dictionary is a Laplace wavelet dictionary, the generated periodic dictionary atom d pk for: In the formula, m is the number of period extensions, and the period dictionary obtained by extension is: D p =[d p1 ,d p2 ,...,d pk ,...,d pK ] In the formula, L d = m·T L is the length of the dictionary; When the Morlet wavelet dictionary is extended, d pj By m d Mj The periodic dictionary obtained by extension is: D p =[d p1 ,d p2 ,...,d pj ,...,d pJ ] In the formula, L d = m·T M is the length of the dictionary.

7. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: The process of step (6) is: reconstruct the original signal based on sparse representation segmentation: let the dictionary D p Size L d ×K, substitute the dictionaries of other sizes into the number of rows and columns to solve; First, take the step size as 1 and cut the length L from the original signal d fragments, denoted as R1, R2, …, Then the signal fragment R is reconstructed using the sparse representation i , the sparse representation is obtained by solving the sparse coefficient α∈R k×1 , R i Indicated as D p The product of the sparse coefficient α is used to obtain the reconstructed signal fragment α is solved by the orthogonal matching pursuit (OMP) algorithm, the steps are: Set the initial residual ε=R i First, calculate the inner product of each dictionary atom and the residual and find the optimal atom d b : Setting support sets D b Merge into S p : S p =S p ∪d b Determine α so that S p Accurately restore R i : Update residuals: e=R i -S p a Set the iteration threshold to δ o , if ||ε|| 2 >δ o , substituting ε into the formula Calculate and update the parameters in the above formula; if ||ε|| 2 <δ o , then the loop ends and α is obtained; Finally, the first digits of each reconstructed signal segment are added in sequence and the average is calculated to obtain the reconstructed signal The vibration signal of the faulty component is also obtained by decoupling:

8. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: The process of step (7) is: get Then, calculate its RMS: calculate The kurtosis is: The fault characteristic of the corresponding component is [rms, kur]; After completing the feature extraction of the current faulty component, if there is still an envelope peak to be identified, the next envelope peak to be identified is selected, the corresponding faulty component is determined, the periodic dictionary is constructed, and decoupling and feature extraction are performed.

9. The gearbox coupling fault diagnosis method based on periodic dictionary sparse representation according to claim 1 is characterized by: The process of step (8) is: Assume that q groups of features [rms, kur] are collected, and the K-means clustering algorithm is used to determine the fault type of each component. The process is as follows: ① Set the maximum number of iterations r max , let the number of clusters be P. For a primary transmission gearbox, there are four components, each of which includes four types of faults: broken teeth, wear, pitting, and other faults of gears, and inner ring, outer ring, rolling element, and cage faults of bearings. Therefore, P = 4 × 4 = 16; ② Represent each group of features in the rectangular coordinate plane composed of the rms axis and the kur axis, and randomly initialize P cluster centers c1 (0) ,c2 (0) ,…,c P (0) ; ③ Calculate the Euclidean distance between each sample and the cluster center, and assign the sample to the cluster closest to it. The Euclidean distance calculation formula is: ④ Update the position of each cluster center. The coordinates of the new cluster center on the rms axis and the kur axis are the mean values ​​of all points in the cluster on the corresponding coordinate axes. ⑤If the updated cluster center position has changed, go to step ②; if it has not changed, go to step ⑥; ⑥When the cluster center position stops changing or reaches the maximum number of iterations, output the clustering results; After clustering, the q samples are divided into P classes, 20% of the samples are taken as the training set, the fault type of the training set is input, and the fault type of all samples in each class is defined as the fault type of the majority of training set samples in the class, then the fault type of the remaining 80% of the test set samples is determined.

Images