On-load tap-changer fault diagnosis method based on multi-feature fusion

Abstract

The invention discloses an on-load tap-changer fault diagnosis method based on multi-feature fusion, relates to the technical field of on-load tap-changer fault diagnosis, and is used for solving the problem that complex information in a signal cannot be comprehensively captured due to the fact that an existing vibration signal analysis method generally focuses on extraction of a single feature. Comprising data acquisition, feature extraction and fault diagnosis. According to the method, advanced feature extraction technologies of complementary ensemble empirical mode decomposition (CEEMD), quality factor adjustable wavelet transform (QFAT) and short-time Fourier transform (STFT) are combined, and fault information in a vibration signal is comprehensively extracted. By optimizing the core parameters of the OOSVM classifier and combining the predatory optimization algorithm, the global search capability and local fine search capability of the algorithm are improved, the accuracy, robustness and timeliness of fault diagnosis are further improved, and the requirements of power equipment fault diagnosis for high precision, real-time performance and adaptability are met.

Description

Technical Field

[0001] The invention relates to the technical field of on-load tap changer fault diagnosis, and in particular to an on-load tap changer fault diagnosis method based on multi-feature fusion. Background Art

[0002] The on-load tap changer (OLTC) is an important component of the power transformer and is responsible for the key task of regulating the output voltage. With the increase in the service life of power equipment and the influence of complex operating environment, OLTC is prone to mechanical wear, poor contact and other faults. If these faults are not discovered in time, they will seriously affect the normal operation of the transformer and the stability of the power system. Traditional fault diagnosis methods rely on manual inspection and regular maintenance, cannot achieve real-time monitoring and have low diagnostic accuracy. Especially in the face of complex fault modes and noise interference, traditional methods have poor adaptability and timeliness, and are difficult to meet the needs of modern power equipment fault diagnosis.

[0003] Traditional feature extraction methods, such as time domain analysis, frequency domain analysis, and empirical mode decomposition (EMD), can usually only extract signal features for specific types of faults, lack comprehensive analysis of multi-dimensional and dynamic changes in signals, and are difficult to deal with diverse fault modes in complex electrical equipment operating environments. In addition, existing vibration signal analysis methods usually focus on the extraction of a single feature, are unable to fully capture the complex information in the signal, and lack the ability to comprehensively analyze vibration signals from multiple dimensions, such as time-frequency characteristics, frequency components, and nonlinear characteristics, resulting in low accuracy in fault diagnosis. Summary of the invention

[0004] The purpose of the present invention is to provide an on-load tap changer fault diagnosis method based on multi-feature fusion, which is used to solve the problem that the existing vibration signal analysis method usually focuses on the extraction of a single feature, cannot fully capture the complex information in the signal, and lacks the ability to comprehensively analyze the vibration signal from multiple dimensions.

[0005] The technical solution adopted by the present invention to solve the technical problem is: a method for diagnosing faults of on-load tap changers based on multi-feature fusion, comprising the following steps.

[0006] S1 Data acquisition.

[0007] Collect the original vibration signals of the on-load tap changer under four working conditions: normal state, gear jamming, contact loosening, and fastening screw loosening. i (t).

[0008] S2 feature extraction.

[0009] S2.1 Feature extraction of vibration signal of on-load tapchanger.

[0010] S2.1.1 The complementary ensemble empirical mode decomposition (CEEMD) is used to transform the original vibration signal X i (t) is decomposed into several intrinsic mode components IMF, which are added to the original vibration signal X i (t) adds a positive and a negative white noise m respectively i (t), we get two different new signals P i and N i .

[0011] P i =X i (t)+ m i (t).

[0012] N i =X i (t)- m i (t).

[0013] Perform EMD decomposition on the two different new signals mentioned above to obtain j IMF components and calculate their average values ​​to obtain the final IMF component C j (t).

[0014] .

[0015] C ij + (t) represents the signal P i The jth intrinsic mode component IMF obtained after EMD decomposition; C ij - (t) represents the signal N i The jth intrinsic mode component IMF obtained after EMD decomposition; m is the amount of white noise added during EMD decomposition.

[0016] Repeat the above steps, adding a new normally distributed white noise sequence each time, and taking the IMF with the largest correlation coefficient each time as the final result; calculate the fifth-order component C with the largest correlation coefficient i (t), and find the energy E of each component i , .

[0017] Construct the total energy E of the vibration signal based on the energy of each component ∑ : .

[0018] Find the IMF component energy entropy H EN : .

[0019] S2.1.2 Using quality factor adjustable wavelet transform TQWT to extract the subsequence energy E of the on-load tap changer vibration signal j, the specific steps are as follows: first determine the low-pass scale α and high-pass scale β of the filter to obtain the wavelet basis that meets the signal characteristics. The calculation formulas are: ; Where Q is the quality factor; r is the oversampling rate of the TQWT transformed signal, which is calculated by dividing the sum of the decomposed wavelet coefficients by the signal length; the low-pass output of each filter group of TQWT is used as the input of the continuous filter group, and the maximum number of decomposition layers Jmax of the multi-layer filter group is: .

[0020] Where N is the length of the analytical signal, To round down.

[0021] Then, the TQWT is decomposed; the low-pass filter frequency response function of each layer of wavelet filter is set to H0(ω), and the high-pass filter frequency response function of each layer of wavelet filter is set to H1(ω), w1, w2, …, w J+1 are the wavelet coefficients from the 1st to the J+1th layer, and the calculation formula is as follows.

[0022] .

[0023] .

[0024] Where θ(ω) is the performance function, which is used to construct the frequency response functions H0(ω) and H1(ω). .

[0025] The spectrum X(ω) of the signal is obtained by Fourier transform and is input into the filter bank of TQWT as the input signal; the input signal of the (j-1)th layer is decomposed into two components, denoted as w (j) (ω) and v (j) (ω), w (j) (ω) is the high quality factor component corresponding to the j-th filter bank, v (j) (ω) is the low quality factor component corresponding to the j-th filter bank; v (j) (ω) is used as the input of the next layer of filters, and w (j) (ω) is transformed by inverse Fourier transform, and the time domain component w is obtained (j) (t) is used as the j-th wavelet decomposition coefficient of TQWT; the wavelet decomposition coefficients corresponding to layers 1 to J are stored in the corresponding rows of the decomposition coefficient matrix W to obtain the J-th layer decomposition coefficient matrix W, and the (J+1)-th layer coefficients are the decomposition residues.

[0026] The oscillation energy characteristics of each subsequence of the decomposed OLTC original vibration signal correspond to the filter of this layer, which can characterize the vibration energy characteristics of the original vibration signal. J} is counted as the time-frequency domain feature sequence of the original vibration signal; among which, the energy of the j-th layer subsequence is: .

[0027] S2.1.3 uses the short-time Fourier transform-based vibration frequency component amplitude entropy VFCAE to describe the irregularity of the frequency components of the on-load tapchanger vibration signal in different time periods.

[0028] First, the multi-channel signal is converted into an STFT matrix, and then the VFCAE is extracted to construct the feature vector.

[0029] For a discrete time series s = {s1, s2, .., s N}, its STFT is calculated as follows: .

[0030] Among them, h(nk) is the Gaussian window function, the value of S(f, k) is a complex number; s(n) is the discrete time series of the original signal, e -j2πfn / N is the complex exponential part of the Fourier transform.

[0031] The Euclidean norm of the elements in S is calculated as follows: .

[0032] Where: ‖ ‖ is the Euclidean norm.

[0033] The amplitude of the STFT matrix is ​​divided into L intervals based on the maximum value, and the distribution of the amplitude intervals D(p) is counted, where p is an integer between 0 and L; then, a weight is added to the calculated entropy value. ; Among them, ω(p) is the weight of the pth interval; Finally, the average entropy value in a time period is calculated as the feature value, and the expected feature vector dimension is taken as y, then the feature value is: ; Where k means that the average entropy of the kth feature dimension is being calculated; ; Indicates rounding down.

[0034] S2.2 Feature fusion.

[0035] For the above three feature quantities F(k), E j and H EN Perform standardization.

[0036] .

[0037] Among them, μ F、 μ E、 μ H is the mean of the characteristic quantity; σ F、 σ E、 σ His the standard deviation of the feature quantity.

[0038] The eigenvectors are combined using the weighted summation method.

[0039] ; Among them, ω1, ω2, ω3 are the set weights.

[0040] S3 fault diagnosis.

[0041] The optimization operator support vector machine IVOA-OOSVM based on the improved Viper optimization algorithm is used to perform fault diagnosis on the features.

[0042] The IOWA operator is defined as follows: .

[0043] where b j is IOWA for (u i , a i ) has the jth largest u i a i Value, u i is the order-inducing variable, a i is the independent variable; ω j Represents the weight of the IOWA operator.

[0044] The expression of Shannon entropy is as follows: .

[0045] First, the slack variables ξ are sorted according to the vector u: .

[0046] In the formula, B u is the ordered parameter vector of the slack variable ξ, i.e. b j is the slack variable corresponding to the jth largest IOWA in u; different weights are assigned according to the sorted positions, and the sorted slack variables ξ and the corresponding weight vector ω are defined; the weight vector ω is predetermined according to the definition of the IOWA operator, satisfying and , where n is the number of data points; ω T is the transpose of the weight vector ω.

[0047] Calculate the weighted slack variables and L(ξ), expressed mathematically as: .

[0048] Among them, v(i) is the sorted index, indicating the position of the i-th largest slack variable ξ in the original vector.

[0049] The positive sample optimization equation of the optimization operator support vector machine is: .

[0050] The negative sample optimization equation of the optimization operator support vector machine is: .

[0051] In the formula, matrix A represents positive samples, B represents negative samples; K is the kernel function to be determined; b1 and b2 are biases; ξ and η are slack variables; W1 and W2 are normal vectors; e1 and e2 are unit column vectors, and the number of rows of e1 is the same as that of the kernel function. Same, the number of rows in e2 and the kernel function Same; c1 and c2 are penalty factors.

[0052] The improved Viper optimization algorithm IVOA is used to optimize the core parameters of the support vector machine, the penalty factor C and the Gaussian kernel function parameter g.

[0053] Furthermore, by arranging sensors on the surface of the oil tank to collect the original vibration signal X of the on-load tap changer i (t).

[0054] Furthermore, the mathematical model of the improved Viper optimization algorithm IVOA is as follows: Phase I: prey identification process.

[0055] The prey is freely distributed in the environment, and its behavior simulation expression is: .

[0056] Where P i is the location of the target object of the i-th viper; X k is the state of the viper; k is a natural number [1, N]; N is the number of viper populations.

[0057] The viper randomly selects a prey to attack. The viper behavior simulation expression is as follows: .

[0058] In the formula, F ri is the ideal fitness value; is the new state of the ith viper in the vth dimension; r is an arbitrary value in [1, N], which is used to generate irregular behavior of the viper; the value of I is 1 or 2; x i，v is the original state of the ith viper in the vth dimension; p i，v is the position of the target of the ith viper in the vth dimension, F i is the fitness value of the i-th viper at its current position.

[0059] .

[0060] In the formula, is the new state of the i-th viper; is the fitness value corresponding to the new state of the viper in the vth dimension; calculate the current optimal solution : .

[0061] For each viper The corresponding fitness value.

[0062] Initialize the population using backpropagation initialization to speed up the search process: .

[0063] in, is the new position of the i-th individual in the initialized population; is the adjustment factor; is the distance between the current solution and the optimal solution, indicating the initialization direction of the population individuals: .

[0064] The wild horse algorithm search mechanism is introduced to improve the parameter r and calculate the individual angle coordinate θ1(i): .

[0065] x i is the current iteration individual; x best is the optimal individual in the population; o is the origin; the value range of i is (1, N).

[0066] The iteration area is divided into 4 parts, s is the area number, and the number of individuals in each area is calculated; the probability P(s) is assigned to the area according to the number of individuals, and the probability calculation expression is as follows: .

[0067] Where num(s) is the number of individuals in the sth region.

[0068] Use fitness weight to search the area and calculate the regional fitness weight A(s): .

[0069] Select the search area s by area selection, and select the angle θ in the sth area k , and its generating expression is as follows: .

[0070] Among them, α is the adjustment factor that controls the intensity of the adjustment, α∈[0,1]; d i is the Euclidean distance between the current individual and the optimal solution: ;d max is the maximum possible distance in the entire search space: , that is, the distance between the farthest individual and the optimal solution.

[0071] After the generationk Substitute into the following formula and iterate: .

[0072] In the formula, is the new state of the i-th viper in the v-th dimension after iteration.

[0073] Phase 2: Chasing and escaping.

[0074] After the viper catches its prey, the prey tries to escape; during the chase, assuming that the hunting is close to an attack range with a radius of R, the mathematical expression of this process is as follows: , .

[0075] In the formula, is the new state of the ith viper in the vth dimension in the second stage; t is the current number of iterations; T is the maximum number of iterations.

[0076] .

[0077] In the formula, is the new state of the i-th viper in the second stage; is the fitness value in the new state.

[0078] The particle swarm optimization algorithm PSO is introduced to update the position of the viper. The update expression is as follows: .

[0079] In the formula, is the new state of the ith viper in the vth dimension after the introduction of the PSO algorithm in the second stage; c1, c2 are learning factors that control the degree of guidance of individuals and groups; r1, r2 are random numbers in the range of [0, 1]; X best,v is the individual's historical optimal position, X global,v is the global optimal position of an individual.

[0080] The beneficial effects of the present invention are as follows: the present invention combines the advanced feature extraction technologies of complementary set empirical mode decomposition (CEEMD), quality factor adjustable wavelet transform (QFAT) and short-time Fourier transform (STFT) to comprehensively extract fault information from vibration signals. By optimizing the core parameters of the OOSVM classifier and combining the predator optimization algorithm to enhance the global search capability and local fine search capability of the algorithm, the accuracy, robustness and timeliness of fault diagnosis are further improved, and the requirements of high precision, real-time and adaptability for power equipment fault diagnosis are met. BRIEF DESCRIPTION OF THE DRAWINGS

[0081] Figure 1 The figure is a flow chart of the diagnostic method of the present invention.

[0082] Figure 2 This is a diagram of gear jamming caused by loose transmission shaft.

[0083] Figure 3 This is a diagram of loose contacts.

[0084] Figure 4 Figure 2. Loosening of tightening bolts.

[0085] Figure 5 This is the flow chart of the optimization operator support vector machine based on the improved Viper optimization algorithm. DETAILED DESCRIPTION

[0086] A method for diagnosing faults of an on-load tap changer based on multi-feature fusion according to the present invention is described in detail below with reference to the accompanying drawings.

[0087] like Figure 1 As shown, a method for on-load tap changer fault diagnosis based on multi-feature fusion includes the following steps.

[0088] S1 Data acquisition.

[0089] Sensors are arranged on the surface of the oil tank to collect vibration signals of the on-load tap changer. The vibration signals of the on-load tap changer under four working conditions, namely, normal state, gear jamming, loose contacts, and loose fastening screws, are simulated and collected. Multiple groups of tests are performed under the normal and mechanical failure states of the on-load tap changer to collect multiple groups of data. Figure 2 As shown in the figure, it is a diagram of gear jam; Figure 3 As shown in the figure, it is a contact loosening diagram; Figure 4 The figure shows the loosening of the fastening bolts.

[0090] S2 feature extraction.

[0091] S2.1 Feature extraction of vibration signal of on-load tapchanger.

[0092] S2.1.1 The complementary set empirical mode decomposition (CEEMD) is used to decompose the original vibration signal into several intrinsic mode components (IMFs). i (t) adds a positive and a negative white noise m respectively i (t), we get two different new signals P i and N i .

[0093] P i =X i (t)+ m i (t).

[0094] N i =X i (t)- m i (t).

[0095] Perform EMD decomposition on these two different new signals respectively, obtain j IMF components and calculate their average value to obtain the final IMF component C j (t).

[0096] .

[0097] C ij + (t) represents the signal P i The jth intrinsic mode component IMF obtained after EMD decomposition; C ij - (t) represents the signal N i The jth intrinsic mode component IMF obtained after EMD decomposition; m is the amount of white noise added during EMD decomposition.

[0098] Repeat the above steps, adding a new normally distributed white noise sequence each time, and taking the IMF with the largest correlation coefficient each time as the final result. Find the fifth-order component C with the largest correlation coefficient i (t), and find the energy E of each component i : .

[0099] Construct the total energy E of the vibration signal based on the energy of each component ∑ : .

[0100] Find the IMF component energy entropy H EN : .

[0101] S2.1.2 Using quality factor adjustable wavelet transform TQWT to extract the subsequence energy E of the on-load tap changer vibration signal j The quality factor Q, redundancy parameter r and decomposition level J are the main parameters of TQWT. The quality factor Q is a technical indicator for evaluating the energy concentration of the signal, which is the ratio of the center frequency corresponding to the maximum value in the signal frequency domain to the 3 dB bandwidth.

[0102] The specific steps are as follows: First, determine the low-pass scale α and high-pass scale β of the filter to obtain the wavelet basis that meets the signal characteristics. The calculation formulas are: .

[0103] Where r is the oversampling rate of the TQWT transformed signal, which is calculated by dividing the sum of the decomposed wavelet coefficients by the signal length. The low-pass output of each filter bank of TQWT is used as the input of the continuous filter bank. The maximum number of decomposition layers Jmax of the multi-layer filter bank is: .

[0104] Where N is the length of the analytical signal, To round down.

[0105] Next, we decompose TQWT. The low-pass filter frequency response function of each layer of wavelet filter is set to H0(ω), and the high-pass filter frequency response function of each layer of wavelet filter is set to H1(ω), w1, w2, …, w J+1 are the wavelet coefficients from the 1st to the J+1th layer, and the calculation formula is as follows.

[0106] .

[0107] .

[0108] Where θ(ω) is the performance function used to construct the frequency response function H0 (1) (ω) and H1(ω), . ω is the angular frequency.

[0109] The spectrum X(ω) of the signal is obtained by Fourier transform and is input into the filter bank of TQWT as the input signal. The input signal of the (j-1)th layer is decomposed into two components, denoted as w (j) (ω) and v (j) (ω), respectively, the high and low quality factor components corresponding to the j-th filter bank. (j) (ω) is used as the input of the next layer of filters, and w (j) (ω) is transformed by inverse Fourier transform, and the time domain component w is obtained (j) (t) is used as the j-th wavelet decomposition coefficient of TQWT. The wavelet decomposition coefficients corresponding to layers 1 to J are stored in the corresponding rows of the decomposition coefficient matrix W to obtain the J-th layer decomposition coefficient matrix W, and the (J+1)-th layer coefficients are the decomposition residues.

[0110] The oscillation energy characteristics of each subsequence of the decomposed OLTC vibration signal correspond to the filter of this layer, which can characterize the vibration energy characteristics of the original signal. Therefore, the energy {E1, E2, …, E J} is counted as the time-frequency domain feature sequence of the original vibration signal. The energy of the j-th layer subsequence is: .

[0111] S2.1.3 combines the advantages of time-frequency domain analysis and information entropy technology in strong anti-noise feature extraction, and adopts the short-time Fourier transform vibration frequency component amplitude entropy VFCAE to describe the irregularity of the frequency component of the on-load tap changer vibration signal in different time periods. Among various signal processing technologies, short-time Fourier transform STFT is suitable for long vibration signals to display the signal frequency components that change over time. First, the multi-channel signal is converted into an STFT matrix, and then VFCAE is extracted to construct the feature vector.

[0112] For a discrete time series s = {s1, s2, .., s N}, its STFT is calculated as follows: .

[0113] Among them, h(nk) is the Gaussian window function, the value of S(f, k) is a complex number. s(n) is the discrete time series of the original signal, e -j2πfn / N is the complex exponential part of the Fourier transform.

[0114] To convert it into real matrix form, the Euclidean norm of the elements in S is calculated as follows: .

[0115] Where: ‖ ‖ is the Euclidean norm.

[0116] First, the amplitude of the STFT matrix is ​​divided into L intervals based on the maximum value, and the distribution of the amplitude intervals is counted D(p), where p is an integer between 0 and L. Then, in order to highlight the influence of high-amplitude frequency components, weights are added to the calculated entropy values, and the higher the amplitude, the greater the weight.

[0117] .

[0118] Among them, ω(p) is the weight of the pth interval. Finally, in order to reduce the size of the entropy feature, the average entropy value in a time period is calculated as the feature quantity, and the expected feature vector dimension is taken as y, then the feature quantity is: .

[0119] Among them, k represents the index, indicating that the average entropy of the kth feature dimension is being calculated; ; Indicates rounding down.

[0120] S2.2 Feature fusion.

[0121] Due to the different dimensions of the feature quantities, the above three feature quantities F(k), E j and H EN To standardize: .

[0122] Among them, μ and σ are the mean and standard deviation of the feature value respectively.

[0123] The eigenvectors are combined using the weighted summation method: .

[0124] Among them, ω1, ω2, ω3 are the set weights.

[0125] S3 fault diagnosis.

[0126] like Figure 5 As shown in the figure, the optimization operator support vector machine IVOA-OOSVM based on the improved Viper optimization algorithm is used to perform fault diagnosis on the features.

[0127] The optimization operator support vector machine is improved on the basis of the twin vector machine. Its basic idea is to find two non-parallel hyperplanes by solving two quadratic programming problems and use the IOWA operator to replace the hinge loss function generalization in order to make different trade-offs for the slack variable ξ. The IOWA operator is an average aggregation operator that generalizes the OWA operator by using a reordering process based on induced variables. The n-dimensional IOWA operator is a mapping IOWA: R n →R, which has a and The associated n-dimensional weight vector ω, the IOWA operator is defined as follows: .

[0128] where b j is IOWA for (u i , a i ) has the jth largest u i a i Value, u i is the order-inducing variable, a i is the independent variable. j Represents the weight of the IOWA operator.

[0129] The expression of Shannon entropy is as follows: .

[0130] First, sort the slack variables ξ. According to the previous definition, sort the slack variables ζ according to the vector u: .

[0131] In the formula, B u is the ordered parameter vector of the slack variable ξ, i.e. b j is the slack variable corresponding to the j-th largest IOWA in u. Different weights are assigned according to the sorted positions. These weights will be used to calculate the weighted slack variable sum, thereby replacing ξ in the original hinge loss function. Define the sorted slack variable ξ and the corresponding weight vector ω. ω T is the transpose of the weight vector ω. The weight vector ω is predetermined according to the definition of the IOWA operator and satisfies and , where n is the number of data points.

[0132] Calculate the weighted slack variable sum, expressed mathematically as: .

[0133] Among them, v(i) is the sorted index, indicating the position of the i-th largest slack variable ξ in the original vector.

[0134] From this, the optimization equation for the positive sample of the optimization operator support vector machine is as follows: .

[0135] The negative sample optimization form of the optimization operator support vector machine is similar: .

[0136] In the formula, matrix A represents positive samples, B represents negative samples; K is the kernel function to be determined; b1 and b2 are biases; ξ and η are slack variables; W1 and W2 are normal vectors; e1 and e2 are unit column vectors, and the number of rows of e1 is the same as that of the kernel function. Same, the number of rows in e2 and the kernel function Same; c1 and c2 are penalty factors.

[0137] From the above analysis, we can see that the penalty factors c1, c2 and the Gaussian kernel function parameter g are two very important hyperparameters, which have a significant impact on the performance, generalization ability and computational complexity of the model. The improved Viper optimization algorithm IVOA is used to optimize the core parameters of the support vector machine, the penalty factor C and the Gaussian kernel function parameter g.

[0138] The algorithm establishes a mathematical model based on the different hunting stages of the hypothetical viper as follows: Stage 1: prey identification process.

[0139] The prey is freely distributed in the environment, and its behavior simulation expression is: .

[0140] Where P i is the location of the target object of the i-th viper; X k is the state of the viper; k is a natural number in [1, N]; N is the number of viper populations. The viper randomly selects a prey to attack, and the viper behavior simulation expression is as follows: .

[0141] In the formula, F ri is the ideal fitness value; is the new state of the i-th viper; is the new state of the ith viper in the vth dimension; r is an arbitrary value in [1, N], which is used to generate irregular behavior of the viper; the value of I is 1 or 2. i，vis the original state of the ith viper in the vth dimension, p i，v is the target position of the i-th viper in the v-th dimension, i.e., the position of the prey; F i is the fitness value of the i-th viper at its current position.

[0142] .

[0143] In the formula, is the fitness value corresponding to the new state of the viper in the vth dimension. Calculate the current optimal solution : .

[0144] For each viper In each iteration of the algorithm, back propagation is used to dynamically adjust the population distribution according to the current optimal solution. The optimal solution updated each time is used as a reference for the next round of population initialization. The position of the population is adjusted through error propagation so that the population moves closer to the current optimal solution. , then update the current optimal solution and its corresponding fitness.

[0145] Initialize the population using backpropagation initialization to speed up the search process: .

[0146] in, is the new position of the i-th individual in the initialized population, is the adjustment factor, is the distance between the current solution and the optimal solution, indicating the initialization direction of the population individuals: .

[0147] The wild horse algorithm search mechanism is introduced to improve the parameter r. Calculate the individual angle coordinate θ1(i): .

[0148] x i is the current iteration individual; x best is the optimal individual in the population; o is the origin; the value range of i is (1, N).

[0149] The iteration area is divided into 4 parts, s is the area number, and the number of individuals in each area is calculated. The probability P(s) is assigned to the area according to the number of individuals. The more individuals there are, the greater the probability. The probability calculation expression is as follows: .

[0150] Where num(s) is the number of individuals in the sth region. Use the fitness weight to search the region and calculate the regional fitness weight A(s): .

[0151] Select the search area s by area selection, and select the angle θ in the sth area k , and its generating expression is as follows: .

[0152] Among them, α is the adjustment factor that controls the intensity of the adjustment, α∈[0,1]; d i is the Euclidean distance between the current individual and the optimal solution: ;d max is the maximum possible distance in the entire search space: , that is, the distance between the farthest individual and the optimal solution.

[0153] After the generation k Substitute into the following formula and iterate: .

[0154] In the formula, is the new state of the i-th viper in the v-th dimension after iteration.

[0155] Phase 2: Chasing and escaping.

[0156] After the viper catches its prey, the prey tries to escape. During the chase, assuming that the hunting is close to an attack range with a radius of R, the mathematical expression of this process is as follows: .

[0157] .

[0158] In the formula, is the new state of the ith viper in the vth dimension in the second stage; t is the current number of iterations; T is the maximum number of iterations.

[0159] .

[0160] In the formula, is the new state of the i-th viper in the second stage; is the fitness value in the new state.

[0161] In the second stage of the algorithm, the prey tries to escape. Since the iteration parameter R is limited, the particle swarm optimization algorithm PSO is introduced to update the position of the viper to avoid the algorithm falling into the local optimum. The update expression is as follows: .

[0162] In the formula, is the new state of the ith viper in the vth dimension after the introduction of the PSO algorithm in the second stage; c1, c2 are learning factors that control the degree of guidance of individuals and groups; r1, r2 are random numbers in the range of [0, 1]; X best,v and X global,v are the individual's historical optimal position and the global optimal position respectively.

[0163] The model is trained using known fault samples, and the extracted fault features are input into the trained fault diagnosis model to obtain the fault diagnosis results. The optimization operator support vector machine optimized based on the improved Viper optimization algorithm has better global search capabilities and can perform global exploration in complex environments. By imitating the search behavior of predators in nature, the algorithm can find the global optimal solution in a wider solution space. Compared with other optimization algorithms, this algorithm has better iterative convergence capabilities and faster optimization speed for both multi-modal functions and unimodal functions, demonstrating the superiority of the optimization operator support vector machine based on the improved Viper optimization algorithm.

Claims

1. A fault diagnosis method for on-load tap changer based on multi-feature fusion, characterized in that: The following steps are involved: S1 Data Collection Collect the original vibration signals of the on-load tap changer under four working conditions: normal state, gear jamming, contact loosening, and fastening screw loosening. i (t); S2 feature extraction S2.1 On-load tapchanger vibration signal feature extraction S2.1.1 The complementary ensemble empirical mode decomposition (CEEMD) is used to transform the original vibration signal X i (t) is decomposed into several intrinsic mode components IMF, and the energy entropy H of the intrinsic mode components IMF is calculated. EN ; S2.1.2 Using quality factor adjustable wavelet transform TQWT to extract the subsequence energy E of the on-load tap changer vibration signal j ; S2.1.3 uses the amplitude entropy of the vibration frequency component VFCAE based on short-time Fourier transform to describe the irregularity of the frequency component of the vibration signal of the on-load tapchanger in different time periods; calculates the average entropy value F (k) within a time period; S2.2 Feature Fusion For the above three feature quantities F(k), E j and H EN To standardize: ； Among them, μ F、 μ E、 μ H is the mean of the characteristic quantity; σ F、 σ E、 σ H is the standard deviation of the characteristic quantity; The eigenvectors are combined using the weighted summation method: ; Among them, ω1, ω2, ω3 are the set weights; S3 Troubleshooting The optimization operator support vector machine IVOA-OOSVM based on the improved Viper optimization algorithm is used to perform fault diagnosis on the features.

2. The on-load tap changer fault diagnosis method based on multi-feature fusion according to claim 1, characterized in that: The original vibration signal X of the on-load tap changer is collected by placing sensors on the surface of the oil tank. i (t).

3. The on-load tap changer fault diagnosis method based on multi-feature fusion according to claim 2 is characterized in that: The specific steps of step S2.1.1 are: To the original vibration signal X i (t) adds a positive and a negative white noise m respectively i (t), we get two different new signals P i and N i : P i = X i (t)+ m i (t); N i = X i (t)- m i (t); Perform EMD decomposition on the two different new signals mentioned above to obtain j IMF components and calculate their average values ​​to obtain the final IMF component C j (t); ； C ij + (t) represents the signal P i The jth intrinsic mode component IMF obtained after EMD decomposition; C ij - (t) represents the signal N i The jth intrinsic mode component IMF obtained after EMD decomposition; m is the amount of white noise added during EMD decomposition; Repeat the above steps, adding a new normally distributed white noise sequence each time, and taking the IMF with the largest correlation coefficient each time as the final result; calculate the fifth-order component C with the largest correlation coefficient i (t), and find the energy E of each component i : ； Construct the total energy E of the vibration signal based on the energy of each component ∑ : ; Find the IMF component energy entropy H EN : .

4. The on-load tap changer fault diagnosis method based on multi-feature fusion according to claim 3 is characterized in that: The mathematical model of the improved Viper optimization algorithm IVOA is: Stage 1: Prey Identification Process The prey is freely distributed in the environment, and its behavior simulation expression is: ; Where P i is the location of the target object of the i-th viper; X k is the state of the viper; k is a natural number [1, N]; N is the number of viper populations; The viper randomly selects a prey to attack. The viper behavior simulation expression is as follows: ； In the formula, F ri is the ideal fitness value; is the new state of the ith viper in the vth dimension; r is an arbitrary value in [1, N], which is used to generate irregular behavior of the viper; the value of I is 1 or 2; x i，v is the original state of the ith viper in the vth dimension; p i，v is the position of the target of the ith viper in the vth dimension, F i is the fitness value of the current position of the i-th viper; ； In the formula, is the new state of the i-th viper; is the fitness value corresponding to the new state of the viper in the vth dimension; calculate the current optimal solution : ; For each viper The corresponding fitness value; Initialize the population using backpropagation initialization to speed up the search process: ; in, is the new position of the i-th individual in the initialized population; is the adjustment factor; is the distance between the current solution and the optimal solution, indicating the initialization direction of the population individuals: ; The wild horse algorithm search mechanism is introduced to improve the parameter r and calculate the individual angle coordinate θ1(i): ； x i is the current iteration individual; x best is the optimal individual in the population; o is the origin; the value range of i is (1, N); The iteration area is divided into 4 parts, s is the area number, and the number of individuals in each area is calculated; the probability P(s) is assigned to the area according to the number of individuals, and the probability calculation expression is as follows: ; Where num(s) is the number of individuals in the sth region; Use fitness weight to search the area and calculate the regional fitness weight A(s): ; Select the search area s by area selection, and select the angle θ in the sth area k , and its generating expression is as follows: ; Among them, α is the adjustment factor that controls the intensity of the adjustment, α∈[0,1]; d i is the Euclidean distance between the current individual and the optimal solution: ;d max is the maximum possible distance in the entire search space: , that is, the distance between the farthest individual and the optimal solution; After the generation k Substitute into the following formula and iterate: ; In the formula, is the new state of the i-th viper in the v-th dimension after iteration; Stage 2: Chase and Escape After the viper catches its prey, the prey tries to escape; during the chase, assuming that the hunting is close to an attack range with a radius of R, the mathematical expression of this process is as follows: ； ； In the formula, is the new state of the ith viper in the vth dimension in the second stage; t is the current number of iterations; T is the maximum number of iterations; ； In the formula, is the new state of the i-th viper in the second stage; is the fitness value in the new state; The particle swarm optimization algorithm PSO is introduced to update the position of the viper. The update expression is as follows: ； In the formula, is the new state of the ith viper in the vth dimension after the introduction of the PSO algorithm in the second stage; c1, c2 are learning factors that control the degree of guidance of individuals and groups; r1, r2 are random numbers in the range of [0, 1]; X best,v is the individual's historical optimal position, X global,v is the global optimal position of an individual.

5. The on-load tap changer fault diagnosis method based on multi-feature fusion according to claim 4, characterized in that: The specific steps of step S2.1.2 are as follows: First, determine the low-pass scale α and high-pass scale β of the filter to obtain the wavelet basis that meets the signal characteristics. The calculation formulas are: ; Where Q is the quality factor; r is the oversampling rate of the TQWT transformed signal, which is calculated by dividing the sum of the decomposed wavelet coefficients by the signal length; the low-pass output of each filter group of TQWT is used as the input of the continuous filter group, and the maximum number of decomposition layers Jmax of the multi-layer filter group is: ; Where N is the length of the analytical signal, To round down; Then, the TQWT is decomposed; the low-pass filter frequency response function of each layer of wavelet filter is set to H0(ω), and the high-pass filter frequency response function of each layer of wavelet filter is set to H1(ω), w1, w2, …, w J+1 is the wavelet coefficient from the 1st to the J+1th layer, and the calculation formula is: ； ； Where θ(ω) is the performance function, which is used to construct the frequency response functions H0(ω) and H1(ω). ; The spectrum X(ω) of the signal is obtained by Fourier transform and is input into the filter bank of TQWT as the input signal; the input signal of the (j-1)th layer is decomposed into two components, denoted as w (j) (ω) and v (j) (ω), w (j) (ω) is the high quality factor component corresponding to the j-th filter bank, v (j) (ω) is the low quality factor component corresponding to the j-th filter bank; v (j) (ω) is used as the input of the next layer of filters, and w (j) (ω) is transformed by inverse Fourier transform, and the time domain component w is obtained (j) (t) is used as the j-th wavelet decomposition coefficient of TQWT; the wavelet decomposition coefficients corresponding to layers 1 to J are stored in the corresponding rows of the decomposition coefficient matrix W to obtain the J-th layer decomposition coefficient matrix W, and the (J+1)-th layer coefficient is the decomposition residue; The oscillation energy characteristics of each subsequence of the decomposed OLTC original vibration signal correspond to the filter of this layer, which can characterize the vibration energy characteristics of the original vibration signal. J } is counted as the time-frequency domain feature sequence of the original vibration signal; among which, the energy of the j-th layer subsequence is: .

6. The on-load tap changer fault diagnosis method based on multi-feature fusion according to claim 5, characterized in that: The specific content of step S2.1.3 is: first, convert the multi-channel signal into an STFT matrix, and then extract the VFCAE to construct a feature vector; For a discrete time series s = {s1, s2, .., s N }, its STFT is calculated as follows: ； Among them, h(nk) is the Gaussian window function, the value of S(f, k) is a complex number; s(n) is the discrete time series of the original signal, e -j2πfn / N is the complex exponential part of the Fourier transform; The Euclidean norm of the elements in S is calculated as follows: ; Where: ‖ ‖ is the Euclidean norm; The amplitude of the STFT matrix is ​​divided into L intervals based on the maximum value, and the distribution of the amplitude intervals is counted D(p), where p is an integer between 0 and L; then, a weight is added to the calculated entropy value; ; Among them, ω(p) is the weight of the pth interval; Finally, the average entropy value in a time period is calculated as the feature value, and the expected feature vector dimension is taken as y, then the feature value is: ; Where k means that the average entropy of the kth feature dimension is being calculated; ; Indicates rounding down.

7. The on-load tap changer fault diagnosis method based on multi-feature fusion according to claim 6, characterized in that: The specific content of step S3 is: The IOWA operator is defined as follows: ; where b j is IOWA for (u i , a i ) has the jth largest u i a i Value, u i is the order-inducing variable, a i is the independent variable; ω j Represents the weight of the IOWA operator; The expression of Shannon entropy is as follows: ; First, sort the slack variables ξ according to the vector u: ; In the formula, B u is the ordered parameter vector of the slack variable ξ, i.e. b j is the slack variable corresponding to the jth largest IOWA in u; different weights are assigned according to the sorted positions, and the sorted slack variables ξ and the corresponding weight vector ω are defined; the weight vector ω is predetermined according to the definition of the IOWA operator, satisfying and , where n is the number of data points; ω T is the transpose of the weight vector ω; Calculate the weighted slack variables and L(ξ), expressed mathematically as: ; Among them, v(i) is the sorted index, indicating the position of the i-th largest slack variable ξ in the original vector; The positive sample optimization equation of the optimization operator support vector machine is: ； The negative sample optimization equation of the optimization operator support vector machine is: ； In the formula, matrix A represents positive samples, B represents negative samples; K is the kernel function to be determined; b1 and b2 are biases; ξ and η are slack variables; W1 and W2 are normal vectors; e1 and e2 are unit column vectors, and the number of rows of e1 is the same as that of the kernel function. Same, the number of rows in e2 and the kernel function Same; c1 and c2 are penalty factors; The improved Viper optimization algorithm IVOA is used to optimize the core parameter penalty factor C and Gaussian kernel function parameter g of the support vector machine.

Images