CN110686768A - Improved rotating machinery nonstationary vibration signal calculation order ratio analysis method - Google Patents

Abstract

The invention discloses an improved analysis method for calculating the order ratio of a non-stationary vibration signal of a rotating machine, which comprises the following steps: the method comprises the steps of preprocessing a key phase pulse signal to obtain a time point sequence, a time difference sequence and a time gain sequence, performing similarity operation on the time gain sequence obtained through preprocessing, performing acceleration classification on a similarity coefficient sequence obtained through the similarity operation to obtain a gain classification sequence, reversely deducing the gain classification sequence according to the time difference sequence and the time point sequence in sequence to obtain a classification time point sequence, performing angular domain resampling on a non-stationary vibration signal by using the classification time point sequence to obtain an angular domain resampling signal, and performing fast Fourier transform on the angular domain resampling signal to obtain an order spectrum. The method can be effectively used for calculating the order ratio analysis of the non-stationary vibration signal of the rotary machine, and has higher operation efficiency than the traditional method.

Description

Improved rotating machinery nonstationary vibration signal calculation order ratio analysis method

Technical Field

The invention relates to an improved method for analyzing the calculation order ratio of a non-stationary vibration signal of a rotary machine, belonging to the technical field of state monitoring and fault diagnosis of mechanical equipment.

Background

In practical application, the rotating mechanical equipment is in a non-steady motion state due to objective factors such as working environment, time-varying working conditions and the like and even human factors, and the rotating mechanical equipment which can be stably in a steady motion state is fresh in practical application. For example, a planetary gear box in the wind power generation industry may cause time-varying working conditions due to the characteristics of the working environment in practical application, so that the planetary gear box is often in a non-steady motion state with time-varying rotation speed. Therefore, the method has important significance for fault diagnosis and state monitoring of the rotating mechanical equipment in the non-stationary state.

Most fault diagnosis analysis methods today can only be used for vibration signals that are stationary, while for vibration signals that are not stationary, they lose their intended effectiveness. The calculation order ratio analysis method is one of common analysis methods for non-stationary vibration signals, converts non-stationary time domain vibration signals into stationary angular domain signals through angular domain resampling, and can provide support for most characteristic analysis methods for stationary vibration signals. However, the arrival of the big data era and the development of the fields of online monitoring, real-time monitoring and the like put higher demands on the operational efficiency of the applied analysis method. The traditional calculation order ratio analysis method has the problem of low operation efficiency due to excessive iteration times, so that the calculation order ratio analysis method has the problem of insufficient operation efficiency when being practically applied to the fields such as online monitoring, real-time monitoring and the like which need higher operation efficiency and have higher real-time requirements in the current big data era, and is not beneficial to the practical application in the field.

Disclosure of Invention

The invention provides an improved method for analyzing the calculated order ratio of a non-stationary vibration signal of a rotating machine, which can be effectively used for completing the analysis of the calculated order ratio of the non-stationary vibration signal.

The technical scheme of the invention is as follows: an improved method for calculating order ratio analysis of a non-stationary vibration signal of a rotating machine, the method comprising: the method comprises the steps of preprocessing a key phase pulse signal to obtain a time point sequence, a time difference sequence and a time gain sequence, performing similarity operation on the time gain sequence obtained through preprocessing, performing acceleration classification on a similarity coefficient sequence obtained through the similarity operation to obtain a gain classification sequence, performing reverse derivation on the gain classification sequence in sequence according to the direction from the matching time difference sequence to the time point sequence to obtain a classification time point sequence, performing angular domain resampling on a non-stationary vibration signal by using the classification time point sequence to obtain an angular domain resampling signal, and performing fast Fourier transform on the angular domain resampling signal to obtain an order spectrum.

Preprocessing the key-phase pulse signal to obtain a time point sequence, a time difference sequence and a time gain sequence as step S1:

s1, preprocessing the key phase pulse signal S (T) to obtain a time point sequence T from each rotation of the reference axis to the whole circle_n(i)＝{t₁，t₂，t₃，...，t_N}; followed by a sequence of time points T_n(i) Preprocessing is carried out to obtain the time used for constructing each rotation of the reference shaft as a time difference sequence T_{n_1nd}(j) (ii) a Time difference sequence T_{n_1nd}(j) Preprocessing is carried out to obtain the difference value of the time used by the reference shaft for every two continuous rotations as a time gain sequence T_{n_2nd}(k) (ii) a Wherein i is more than or equal to 1 and less than or equal to N, the total number of whole revolution of the reference axis N-1, t₁Is the starting time point of the 1 st rotation of the reference axis, t_NJ is more than or equal to 1 and less than or equal to i-1, and k is more than or equal to 1 and less than or equal to j-1 at the time point when the Nth extracted reference shaft rotates to the whole circle.

The step S1 specifically includes:

s1.1, selecting half of an amplitude pole of a key phase pulse signal S (T) as a mark value, and circularly obtaining a time value of the closest point of the amplitude of each pulse rising edge in the key phase pulse signal S (T) to the mark value to obtain a time point sequence T_n(i)；

S1.2, sequence of time points T_n(i) The pretreatment of (1): sequence of time points T_n(i) Performing difference operation to obtain the time used by the reference shaft for each rotation as a time difference sequence T_{n_1nd}(j)＝{(t₂-t₁)、(t₃-t₂)、...、(t_N-t_N-1) }; wherein (t)_N-t_N-1) The time taken for the last revolution of the reference shaft;

s1.3, sequence of time differences T_{n_1nd}(j) The pretreatment of (1): sequence of time differences T_{n_1nd}(j) Performing difference operation to obtain the difference value of the time used by the reference shaft for each continuous rotation for two circles as a time gain sequence T_{n_2nd}(k)＝{(t₃-t₂)-(t₂-t₁)、...、(t_N-t_N-1)-(t_N-1-t_N-2) }; wherein (t)_N-t_N-1)-(t_N-1-t_N-2) Is the difference in the time taken for the last two revolutions of the reference shaft.

Performing similar operation on the time gain sequence obtained by the preprocessing as step S2, wherein S2 specifically includes:

s2.1, introducing a similar operation formula:

for time gain sequence T_{n_2nd}(k) Carrying out similar operation; where Sc is the coefficient of similarity, dt_sFor a time gain sequence T_{n_2nd}(k) S is more than or equal to 2 and less than or equal to k in the s-th data value;

s2.2, calculating a formula for the time gain sequence T through the similarity coefficient Sc obtained in the step S2.1_{n_2nd}(k) Performing cyclic operation to obtain a time gain sequence T_{n_2nd}(k) The sequence of similarity coefficients Sc.

Performing acceleration classification on the similarity coefficient sequence obtained by the similarity operation to obtain a gain classification sequence as step S3, where S3 specifically includes:

s3.1, establishing an acceleration type mark for the 1 st value in the similar coefficient sequence, and setting a comparison value as 1; as a first comparison;

second to last comparison:

firstly, calculating the magnitude relation between the similarity coefficient Sc of the current value and a threshold value L: if Sc is larger than L, setting the comparison value of the current time as 0, and if Sc is less than or equal to L, setting the comparison value of the current time as 1;

then, the comparison value of the current time is judged to be the same as the comparison value of the previous time: if the two values are the same, judging that the current time and the previous time have the same acceleration type, recording the acceleration type mark of the previous time, and keeping the comparison value of the current time as 1 for the next time of same judgment; if the two values are different, the acceleration type of the current time is judged to be different from that of the previous time, a new acceleration type mark is generated at the moment, and the comparison value of the current time is set to be 1 for the next same judgment;

s3.2, constructing M gain classification sequences T with the same acceleration according to the acceleration type marks obtained in S3.1_{n_2nd}' (k); wherein M is more than or equal to 1 and less than or equal to k.

The step S4 of reversely deriving the gain classification sequence in sequence from the direction of the matching time difference sequence and the time point sequence to obtain a classification time point sequence, where S4 specifically includes:

s4.1, classifying the gain classified sequence T into M types_{n_2nd}' (k) is divided into two types of push-back: the 1 st time and the M th time, wherein M is more than or equal to 2 and less than or equal to M; wherein, the 1 st gain classification sequence refers to a sequence with the addition degree type marked as the 1 st type, and the mth gain classification sequence has the same principle;

s4.2, converting T obtained in S4.1_{n_2nd}' (k) increasing the maximum index value belonging to the index range of the 1 st acceleration type by 1 in the 1 st gain classification sequence as the 1 st acceleration type at T_{n_1nd}' (j) and the 1 st acceleration type at T_{n_1nd}The starting index value in' (j) is¹；

S4.3, converting T obtained in S4.2_{n_1nd}' (j) increasing the maximum index value belonging to the index range of the 1 st acceleration type by 1 in the 1 st time difference classification sequence as the 1 st acceleration type at T_n' (i) last index value, 1 st acceleration type at T_nThe starting index value in' (i) is 1;

s4.4, converting T obtained in S4.1_{n_2nd}' (k) in the m-th gain classification sequence, the initial and final index values belonging to the m-th acceleration type range are both increased by 1, and the m-th acceleration type is represented by T_{n_1nd}The start and end index values in (j);

s4.5, converting T obtained in S4.4_{n_1nd}' (j) largest of the m-th time difference class sequences belonging to the m-th acceleration type index rangeIncreasing the index value by 1 as the mth acceleration type at T_n' (i) the last index value, the initial index value being its value at T_{n_1nd}' (j) initial index values of the index ranges; sorting out time point classification sequence T_n' (i) index range point values of each acceleration degree type, removing overlapped index point values, arranging in order, and obtaining time point sequence T by point values on the arrangement_n(i) The corresponding index value is the time point of the sorting, thereby constructing a sorting time point sequence T which is divided into M types_n″(i)。

The step S5 of performing angular domain resampling on the non-stationary vibration signal using the sorted time point sequence to obtain an angular domain resampled signal, where the step S5 specifically includes:

s5.1, sorting time point sequence T divided into M types_n"(i), the time range for obtaining M-class equivalent accelerations is t₁"time to t₂"time, the angle through which the reference axis rotates in this time frame is 0 to b x 2 pi, the formula exists:

carry in t₁″、t₂"time can be calculated to obtain polynomial coefficient a₀、a₁(ii) a Wherein b is t₁"time to t₂"number of full revolutions of the reference axis during time, t₁To sort time point sequences T_nThe starting time point, t, of a certain type of acceleration time range in (i)₂To sort time point sequences T_nThe end time point of the acceleration time range of certain type in (i);

s5.2, a obtained in S5.1₀、a₁Substitution into the formula:

the time sequence t corresponding to the equal angle sampling in the range of 0 to b x 2 pi can be obtained by calculation_{n_1nd}(ii) a Wherein p is an interpolation coefficient, and delta theta is an equal angle interval;

s5.3, with the time series t obtained in S5.2_{n_1nd}Carrying out interpolation solution on the non-stationary vibration signal x (t) to obtain an angular domain resampling signal w (t)。

The invention has the beneficial effects that: according to the method, the key phase signals are preprocessed for multiple times, the time gain sequence of each rotation process of the reference shaft can be obtained, the sequence is used for subsequent similar operation to obtain a similar sequence, the similar sequence is compared by a threshold value to finally obtain the acceleration type distribution range of all rotation processes of the reference shaft, the time range of the reference shaft with the same acceleration type is deduced through the distribution range, the time range is applied to an improved calculation order ratio analysis method to play a role in reducing the iterative operation times of an improved first-order linear polynomial equation, and further the solving times of the polynomial coefficient are reduced. Finally, the result is reflected in that angular domain resampling time points are obtained more quickly to carry out interpolation resampling, and therefore the improvement of the operation efficiency is achieved. Under the modern big data era, with the development of the fields of on-line monitoring, real-time monitoring and the like, higher requirements are put forward on the operational efficiency of various analysis methods. The improved calculation order ratio analysis method effectively improves the operation efficiency on the basis of keeping the analysis precision of the original method, solves the problem that the calculation order ratio analysis method is insufficient in operation efficiency when being actually applied to the fields such as online monitoring, real-time monitoring and the like which need higher operation efficiency and have higher real-time requirement under the current big data era, and is beneficial to the actual application of the calculation order ratio analysis method in the field.

Drawings

Fig. 1 is a plan structure view of a ZJS50 integrated design type mechanical design experiment table, an acceleration sensor, and an eddy current sensor in embodiment 1 of the present invention;

fig. 2 is a vibration time domain waveform (local part) received by an acceleration sensor when the ZJS50 comprehensive design type mechanical design experiment table in embodiment 1 of the present invention is in operation;

fig. 3 shows that when the ZJS50 integrated design type mechanical design experiment table in embodiment 1 of the present invention is in operation, a key phase pulse signal (local) is received by an eddy current sensor;

fig. 4 is a scale spectrum obtained by analyzing the envelope-preprocessed signal by using a conventional computational scale analysis algorithm in embodiment 1 of the present invention;

fig. 5 is a scale spectrum (threshold L is 0.1) obtained by analyzing an envelope-preprocessed signal by using an improved computational scale analysis algorithm in embodiment 1 of the present invention;

fig. 6 is a scale spectrum (threshold L is 0.5) obtained by analyzing the envelope-preprocessed signal by using an improved computational scale analysis algorithm in embodiment 1 of the present invention;

FIG. 7 is a block diagram of an implementation of the present invention;

in fig. 1, the reference numerals are: the system comprises a motor 1, a gear box 2, a load 3, an eddy current sensor 4, an acceleration sensor I5, an acceleration sensor II 6 and an acceleration sensor III 7.

Detailed Description

The invention will be further described with reference to the following figures and examples, without however restricting the scope of the invention thereto.

Example 1: an improved method for calculating order ratio analysis of a non-stationary vibration signal of a rotating machine, the method comprising: the method comprises the steps of preprocessing a key phase pulse signal to obtain a time point sequence, a time difference sequence and a time gain sequence, performing similarity operation on the time gain sequence obtained through preprocessing, performing acceleration classification on a similarity coefficient sequence obtained through the similarity operation to obtain a gain classification sequence, performing reverse derivation on the gain classification sequence in sequence according to the direction from the matching time difference sequence to the time point sequence to obtain a classification time point sequence, performing angular domain resampling on a non-stationary vibration signal by using the classification time point sequence to obtain an angular domain resampling signal, and performing fast Fourier transform on the angular domain resampling signal to obtain an order spectrum.

The key phase pulse signal s (t) and the vibration signal x (t) sampled at equal time intervals are obtained as shown in fig. 3 and fig. 2, respectively.

The signal is then further processed as follows:

further, a time point sequence, a time difference sequence, and a time gain sequence obtained by preprocessing the key phase pulse signal may be set as step S1:

s1, preprocessing the key phase pulse signal S (t) to obtain the time point sequence of the reference shaft rotating to the whole circle each timeColumn T_n(i)＝{t₁，t₂，t₃，...，t_N}; followed by a sequence of time points T_n(i) Preprocessing is carried out to obtain the time used for constructing each rotation of the reference shaft as a time difference sequence T_{n_1nd}(j) (ii) a Time difference sequence T_{n_1nd}(j) Preprocessing is carried out to obtain the difference value of the time used by the reference shaft for every two continuous rotations as a time gain sequence T_{n_2nd}(k) (ii) a Wherein i is not less than 1 and not more than N, N-1 is the total number of whole revolution of the reference axis, t₁Is the starting time point of the 1 st rotation of the reference axis, t_NJ is more than or equal to 1 and less than or equal to i-1, and k is more than or equal to 1 and less than or equal to j-1 at the time point when the Nth extracted reference shaft rotates to the whole circle.

Further, the step S1 may specifically be:

s1.1, selecting half of an amplitude pole of a key phase pulse signal S (T) as a mark value, and circularly obtaining a time value of the closest point of the amplitude of each pulse rising edge in the key phase pulse signal S (T) to the mark value to obtain a time point sequence T_n(i)；

S1.2, sequence of time points T_n(i) The pretreatment of (1): sequence of time points T_n(i) Performing difference operation to obtain the time used by the reference shaft for each rotation as a time difference sequence T_{n_1nd}(j)＝{(t₂-t₁)、(t₃-t₂)、...、(t_N-t_N-1) }; wherein (t)_N-t_N-1) The time taken for the last revolution of the reference shaft;

s1.3, sequence of time differences T_{n_1nd}(j) The pretreatment of (1): sequence of time differences T_{n_1nd}(j) Performing difference operation to obtain the difference value of the time used by the reference shaft for each continuous rotation for two circles as a time gain sequence T_{n_2nd}(k)＝{(t₃-t₂)-(t₂-t₁)、...、(t_N-t_N-1)-(t_N-1-t_N-2) }; wherein (t)_N-t_N-1)-(t_N-1-t_N-2) Is the difference in the time taken for the last two revolutions of the reference shaft.

Further, the step of performing similar operation on the time gain sequence obtained by the preprocessing may be set as step S2, where the step S2 specifically is:

s2.1, introducing a similar operation formula:

for time gain sequence T_{n_2nd}(k) Carrying out similar operation; where Sc is the coefficient of similarity, dt_sFor a time gain sequence T_{n_2nd}(k) S is more than or equal to 2 and less than or equal to k in the s-th data value; such as dt_s-1＝dt₁＝(t₃-t₂)-(t₂-t₁)；

S2.2, calculating a formula for the time gain sequence T through the similarity coefficient Sc obtained in the step S2.1_{n_2nd}(k) Performing cyclic operation to obtain a time gain sequence T_{n_2nd}(k) The sequence of similarity coefficients Sc.

Further, the step S3 of performing acceleration classification on the similarity coefficient sequence obtained by the similarity operation to obtain a gain classification sequence may be set, where the step S3 specifically is:

s3.1, establishing an acceleration type mark for the 1 st value in the similar coefficient sequence, and setting a comparison value as 1; as a first comparison;

second to last comparison:

firstly, calculating the magnitude relation between the similarity coefficient Sc of the current value and a threshold value L: if Sc is larger than L, setting the comparison value of the current time as 0, and if Sc is less than or equal to L, setting the comparison value of the current time as 1;

then, the comparison value of the current time is judged to be the same as the comparison value of the previous time: if the two values are the same, judging that the current time and the previous time have the same acceleration type, recording the acceleration type mark of the previous time, and keeping the comparison value of the current time as 1 for the next time of same judgment; if the two values are different, the acceleration type of the current time is judged to be different from that of the previous time, a new acceleration type mark is generated at the moment, and the comparison value of the current time is set to be 1 for the next same judgment;

s3.2, constructing M gain classification sequences T with the same acceleration according to the acceleration type marks obtained in S3.1_{n_2nd}' (k); wherein M is more than or equal to 1 and less than or equal tok。

Examples are: t is_{n_2nd}(k) With { a, b, c, d, e }, a similar operation requires two T' s_{n_2nd}(k) The data in (1) gives a Sc value corresponding to the similar Sc sequence as { ab, bc, cd, de }, and the threshold L is set to be f. Comparing the Sc sequence with f:

creating a new acceleration type A for the first time, and setting a comparison value for the first time as 1;

setting the comparison value of the second time as 1 when the ab is less than or equal to f for the second time, performing the same judgment with the comparison value of the first time, and recording the acceleration type A with the same result, wherein the comparison value is kept as 1;

setting the comparison value of the second time as 1 when bc is less than or equal to f for the third time, performing the same judgment on the comparison value of the first time, recording the acceleration type A with the same result, and keeping the comparison value as 1;

setting the comparison value of the third time to be 0 when the cd is larger than the f for the fourth time, carrying out the same judgment on the comparison value of the third time and the comparison value of the second time, and creating a new acceleration type B with different results, wherein the comparison value is reset to be 1;

setting the fourth comparison value to be 1 when de is less than or equal to f for the fifth time, performing the same judgment on the fourth comparison value and the third comparison value, and recording the acceleration type B when the results are the same, wherein the comparison value is kept to be 1;

since the first comparison creates a new acceleration type A, the gain classification sequence T_{n_2nd}' (k) starting position adds the acceleration type A, thus T_{n_2nd}(k) Corresponding to { a, b, c, d, e }, the final gain classification sequence T_{n_2nd}′(k)＝{A，A，A，B，B}。

The setting method of the threshold value L is as follows: setting a larger threshold value L (for example, 1) and a smaller threshold value L (for example, 0.05), and checking whether the precision of the analysis result corresponding to the larger threshold value L is obviously changed compared with the precision of the analysis result corresponding to the smaller threshold value L after the program is run:

a larger threshold L (e.g. 2) may be set if there is no significant change, and a new threshold (e.g. 0.8) lower than the initial larger threshold (1) may need to be set if there is.

Repeating the above steps can determine a more suitable threshold L.

Further, the step S4 may be to sequentially obtain the sorted time point sequence by reversely deriving the gain sorted sequence from the direction of the matching time difference sequence and the time point sequence, where the step S4 specifically is:

s4.1, classifying the gain classified sequence T into M types_{n_2nd}' (k) is divided into two types of push-back: the 1 st time and the M th time, wherein M is more than or equal to 2 and less than or equal to M; wherein, the 1 st gain classification sequence refers to a sequence with the addition degree type marked as the 1 st type, and the mth gain classification sequence has the same principle;

s4.2, converting T obtained in S4.1_{n_2nd}' (k) increasing the maximum index value belonging to the index range of the 1 st acceleration type by 1 in the 1 st gain classification sequence as the 1 st acceleration type at T_{n_1nd}' (j) and the 1 st acceleration type at T_{n_1nd}The starting index value in' (j) is 1;

s4.3, converting T obtained in S4.2_{n_1nd}' (j) increasing the maximum index value belonging to the index range of the 1 st acceleration type by 1 in the 1 st time difference classification sequence as the 1 st acceleration type at T_n' (i) last index value, 1 st acceleration type at T_nThe starting index value in' (i) is 1;

s4.4, converting T obtained in S4.1_{n_2nd}' (k) in the m-th gain classification sequence, the initial and final index values belonging to the m-th acceleration type range are both increased by 1, and the m-th acceleration type is represented by T_{n_1nd}The start and end index values in (j);

s4.5, converting T obtained in S4.4_{n_1nd}' (j) increasing the maximum index value belonging to the index range of the mth acceleration type in the mth time difference class sequence by 1 as the mth acceleration type in T_n' (i) the last index value, the initial index value being its value at T_{n_1nd}' (j) initial index values of the index ranges; sorting out time point classification sequence T_n' (i) index range point values of each acceleration degree type, removing overlapped index point values, arranging in order, and obtaining time point sequence T by point values on the arrangement_n(i) The corresponding index value is the time point of the sorting, thereby constructing the classification into M classesClass time point sequence T_n″(i)。

For T_{n_2nd}' (k) sequence carries out the word description operation mode in the S4.2 step and the S4.4 step, and T is used_{n_2nd}' (k) { a, B } illustrates (a stands for a first acceleration type and B stands for a second acceleration type);

s4.2 step which corresponds to the acceleration type 1 (in this case the letter A), corresponds to T_{n_2nd}' (k) denotes three letters A (representing the first acceleration type) having indices of 1 to 3, the maximum index value of which is 3, the index value obtained by adding 1 to the maximum index value is 4, and 4 is used as the last index value, and the initial index value is held at 1. The result is to change the initial three letters A to four letters A (index range 1-3 to 1-4).

Step S4.4, which is directed to acceleration type 2 and its subsequent acceleration types (in this case, the letter B), the initial and final index values are increased by 1, which means that the index range is not changed, and only the shift process is performed. Therefore, the two letters B are not changed, only the position slides backwards by 1 bit, and the initial index range is changed from 4-5 to 5-6.

Thus obtaining T_{n_1nd}' (j) { a, B }, which sequence is used to deduce T_n′(i)；

For T_{n_1nd}' (j) sequence carries out the character description operation mode in the S4.3 step and the S4.5 step, and T is used_{n_1nd}' (j) { a, B } illustrates (a stands for a first acceleration type, B stands for a second acceleration type);

s4.3 step which corresponds to the acceleration type 1 (in this case the letter A), to T_{n_1nd}' (j) is an index of 1 to 4, and has a maximum index value of 4, an index value obtained by adding 1 to the maximum index value of 4, 5, and an initial index value of 1 with 5 as a final index value. The result is that the initial four letters A are changed to five letters A (the index ranges from 1-4 to 1-5).

Step S4.5, which corresponds to the acceleration type T2 and its subsequent acceleration type (in this case the letter B)_{n_1nd}In the case of' (j) are two with indices of 5 to 6A letter B (representing a second acceleration type) having a maximum index value of 6, an index value of 7 obtained by adding 1 to the maximum index value, and an end index value of 7, the initial index value being the value at T_{n_1nd}' (j) the initial index value of the index range is 5. The result is to change the initial two letters B to three letters B (index range 5-6 to 5-7). T thus obtained_n' (i) the sequence is as follows:

T_n′(i)＝{A，A，A，A，A

B，B，B}

T_n' (i) storing acceleration type classification information, T_n(i) The time points of the entire revolution of the reference axis are stored in the sequence, one for one.

T_n' (i) overlap in the sequence, in this case the index 5, which represents the demarcation between acceleration type A and acceleration type B, corresponding to T_n(i) The sequence is that indexes 1-5 are acceleration types A, and indexes 5-7 are acceleration types B. In each full revolution corresponding to the reference axis, the 1 st to 4 th revolutions are acceleration type A, and the 5 th to 6 th revolutions are acceleration type B.

Through T_n' (i) the sequence can be represented by T_n(i) Sequences are classified (at most into M classes, M is less than or equal to k, i.e., T_{n_2nd}' (k) { M in a, B } is 2).

Through T_nThe sequence (i) can acquire the start and end index values of each acceleration type, and the index values acquired in this example are {1, 5, 7 }. Where 1 and 5 represent the index ranges for acceleration type a and 5 and 7 represent the index ranges for acceleration type a. (if acceleration type C is present, it is arranged in order, in this case not)

T can be acquired through {1, 5, 7}_n(i) The time point corresponding to the index value 1, 5, 7 in the sequence is t₁，t₅，t₇. Thereby constructing T_n″(i)＝{t₁，t₅，t₇}。

Further, the angular resampling may be performed on the non-stationary vibration signal by using the sorted time point sequence to obtain an angular resampled signal, which is set as step S5, where the step S5 specifically is:

s5.1, sorting time point sequence T divided into M types_n"(i), the time range for obtaining M-class equivalent accelerations is t₁"time to t₂"time of day, where t₁To sort time point sequences T_nThe starting time point, t, of a certain type of acceleration time range in (i)₂To sort time point sequences T_nThe end time point of the acceleration time range of certain type in (i). The angle through which the reference axis rotates in this time frame is 0 to b x 2 pi, and the formula exists:

carry in t₁″、t₂"time can be calculated to obtain polynomial coefficient a₀、a₁(ii) a Wherein b is t₁"time to t₂"number of whole revolutions of the reference axis during time;

as explained in the previous example, the acceleration types are only two, a and B (for ease of distinction, the letters are replaced by 1, and 2), and the sequence B here has only two components, i.e., B ═ B₁，b₂}。

Wherein b is₁Representing the number of full-period rotations of the corresponding reference axis in acceleration type 1, b in the above example₁4 (representing the same acceleration of 1 st to 4 th rotation of the reference axis), b₂Then represents the number of full-period rotations of the corresponding reference axis in acceleration type 2, in the above example b₂2 (representing the same acceleration of 5 th to 6 th rotation of the reference shaft).

b₁，b₂The data source of (1) is through T_n' (i) the start and end index values of each acceleration type acquired in the sequence, the index values acquired in this example are {1, 5, 7 }. Here b₁＝5-1＝4，b₂＝7-5＝2。

If there is a third acceleration type 3, then b must be present₃And b is₃Representative is the number of full-period rotations of the reference axis belonging to acceleration type 3.

Thus, the sum of the sequences b in this example is b₁Value of (a) plus b₂The sum of the summed sequences b is 4+2 to 6. I.e. the reference shaft rotates 6 times in total, so corresponds to T_n(i) The sequence has 7 time points. So Σ b ═ N-1(N ═ 7).

S5.2, a obtained in S5.1₀、a₁Substitution into the formula:

the time sequence t corresponding to the equal angle sampling in the range of 0 to b x 2 pi can be obtained by calculation_{n_1nd}(ii) a Wherein p is an interpolation coefficient, and delta theta is an equal angle interval;

s5.3, with the time series t obtained in S5.2_{n_1nd}And (5) carrying out interpolation solution on the non-stationary vibration signal x (t) to obtain an angular domain resampling signal w (t).

And performing fast Fourier transform on the obtained angular domain resampling signals w (t) to obtain an order spectrum, completing the analysis of the calculated order, and performing fast Fourier transform on the angular domain resampling signals w (t) obtained by different thresholds L to obtain the analysis results of fig. 5 and 6.

As shown in fig. 1, the ZJS50 comprehensive design type mechanical design experiment table is used for simulating a speed change working condition occurring in practical application, and a gear tooth breaking fault type is preset for experiment. Fig. 1 shows the mounting positions of three acceleration sensors on the ZJS50 integrated design type mechanical design experiment table and the mounting position of one eddy current sensor on the input shaft of the ZJS50 integrated design type mechanical design experiment table.

Fig. 2 and 3 respectively show the time domain vibration waveform and the key phase pulse signal of the ZJS50 integrated design type mechanical design experiment table in operation, which are acquired by the acceleration sensor and the eddy current sensor. It can be seen from fig. 3 that there is a gradual increase in the key phase pulse interval, indicating that the ZJS50 integrated design type mechanical design laboratory bench is in the deceleration process.

Fig. 4 shows a ratio spectrum obtained by the conventional calculation ratio analysis algorithm after the vibration signal of fig. 2 is subjected to envelope preprocessing and is subjected to expansion analysis on the basis of the key-phase pulse time sequence of fig. 3. The signal is shown to have major components of order 2, 4, 6, 8, etc. with amplitudes of 3.41, 1.33, 1.96, etc., respectively, with sideband components.

Fig. 5 is a ratio spectrum analyzed by a modified computational ratio analysis algorithm at a threshold L of 0.1 during the steps of the method. The signal is shown to have major components of order 2, 4, 6, 8, etc. with amplitudes of 2.97, 0.94, 1.03, 1.66, etc., respectively, with sideband components. Compared with the conventional calculation order ratio analysis algorithm shown in fig. 4, the improved calculation order ratio analysis algorithm shown in fig. 5 has no reduction in analysis accuracy and only has partial difference in amplitude.

Fig. 6 is a ratio spectrum analyzed by a modified computational ratio analysis algorithm at a threshold L of 0.5 during the steps of the method. The signal is shown to have major components of order 2, 4, 6, 8, etc. with amplitudes of 1.47, 0.68, 0,61, 1.01, etc., respectively, with sideband components. Compared with the conventional calculation order ratio analysis algorithm shown in fig. 4, the improved calculation order ratio analysis algorithm shown in fig. 6 has no reduction in analysis accuracy and only has partial difference in amplitude.

Table 1 shows the comparison of the operation efficiency between the conventional calculation order analysis algorithm of fig. 4 and the improved calculation order analysis algorithm of fig. 5 with the threshold L of 0.1, and the improved calculation order analysis algorithm of fig. 6 with the threshold L of 0.5. The analysis of the traditional calculation order ratio analysis algorithm on the experimental data takes 0.282299 seconds, the analysis of the improved calculation order ratio analysis algorithm with the threshold L being 0.1 takes 0.171788 seconds on the same data, the analysis efficiency is improved by 39.15%, the analysis of the improved calculation order ratio analysis algorithm with the threshold L being 0.5 takes 0.129053 seconds on the same data, and the analysis efficiency is improved by 54.28%.

TABLE 1

As can be seen from comparison of fig. 4 to fig. 6 and table 1, compared with the conventional calculation order analysis algorithm, the improved calculation order analysis algorithm shortens the operation time and effectively improves the operation efficiency on the basis of maintaining the data analysis accuracy. And by setting different threshold values L, the operation time consumption of the improved algorithm is shortened again to a certain extent, and the operation efficiency is further improved.

While the present invention has been described in detail with reference to the embodiments shown in the drawings, the present invention is not limited to the embodiments, and various changes can be made without departing from the spirit of the present invention within the knowledge of those skilled in the art.

Claims (7)

1. An improved analysis method for calculating the order ratio of a non-stationary vibration signal of a rotating machine is characterized by comprising the following steps: the method comprises the following steps: the method comprises the steps of preprocessing a key phase pulse signal to obtain a time point sequence, a time difference sequence and a time gain sequence, performing similarity operation on the time gain sequence obtained through preprocessing, performing acceleration classification on a similarity coefficient sequence obtained through the similarity operation to obtain a gain classification sequence, performing reverse derivation on the gain classification sequence in sequence according to the direction from the matching time difference sequence to the time point sequence to obtain a classification time point sequence, performing angular domain resampling on a non-stationary vibration signal by using the classification time point sequence to obtain an angular domain resampling signal, and performing fast Fourier transform on the angular domain resampling signal to obtain an order spectrum.

2. The improved rotating machinery non-stationary vibration signal computational order analysis method of claim 1, wherein: preprocessing the key-phase pulse signal to obtain a time point sequence, a time difference sequence and a time gain sequence as step S1:

s1, preprocessing the key phase pulse signal S (T) to obtain a time point sequence T from each rotation of the reference axis to the whole circle_n(i)＝{t₁，t₂，t₃，...，t_N}; followed by a sequence of time points T_n(i) Preprocessing is carried out to obtain the time used for constructing each rotation of the reference shaft as a time difference sequence T_{n_1nd}(j) (ii) a Time difference sequence T_{n_1nd}(j) Preprocessing is carried out to obtain the difference value of the time used by the reference shaft for every two continuous rotations as a time gain sequence T_{n_2nd}(k) (ii) a Wherein i is not less than 1 and not more than N, N-1 is the total number of whole revolution of the reference axis, t₁Is the starting time point of the 1 st rotation of the reference axis, t_NJ is more than or equal to 1 and less than or equal to i-1, and k is more than or equal to 1 and less than or equal to j-1 at the time point when the Nth extracted reference shaft rotates to the whole circle.

3. The improved rotating machinery non-stationary vibration signal calculation order analysis method as claimed in claim 2, wherein: the step S1 specifically includes:

s1.1, selecting half of an amplitude pole of a key phase pulse signal S (T) as a mark value, and circularly obtaining a time value of the closest point of the amplitude of each pulse rising edge in the key phase pulse signal S (T) to the mark value to obtain a time point sequence T_n(i)；

S1.2, sequence of time points T_n(i) The pretreatment of (1): sequence of time points T_n(i) Performing difference operation to obtain the time used by the reference shaft for each rotation as a time difference sequence T_{n_1nd}(j)＝{(t₂-t₁)、(t₃-t₂)、...、(t_N-t_N-1) }; wherein (t)_N-t_N-1) The time taken for the last revolution of the reference shaft;

s1.3, sequence of time differences T_{n_1nd}(j) The pretreatment of (1): sequence of time differences T_{n_1nd}(j) Performing difference operation to obtain the difference value of the time used by the reference shaft for each continuous rotation for two circles as a time gain sequence T_{n_2nd}(k)＝{(t₃-t₂)-(t₂-t₁)、...、(t_N-t_N-1)-(t_N-1-t_N-2) }; wherein (t)_N-t_N-1)-(t_N-1-t_N-2) Is the difference in the time taken for the last two revolutions of the reference shaft.

4. The improved rotating machinery non-stationary vibration signal computational order analysis method of claim 1, wherein: performing similar operation on the time gain sequence obtained by the preprocessing as step S2, wherein S2 specifically includes:

s2.1, introducing a similar operation formula:for time gain sequence T_{n_2nd}(k) Carrying out similar operation; where Sc is the coefficient of similarity, dt_sFor a time gain sequence T_{n_2nd}(k) S is more than or equal to 2 and less than or equal to k in the s-th data value;

s2.2, calculating a formula for the time gain sequence T through the similarity coefficient Sc obtained in the step S2.1_{n_2nd}(k) Performing cyclic operation to obtain a time gain sequence T_{n_2nd}(k) The sequence of similarity coefficients Sc.

5. The improved rotating machinery non-stationary vibration signal computational order analysis method of claim 1, wherein: performing acceleration classification on the similarity coefficient sequence obtained by the similarity operation to obtain a gain classification sequence as step S3, where S3 specifically includes:

s3.1, establishing an acceleration type mark for the 1 st value in the similar coefficient sequence, and setting a comparison value as 1; as a first comparison;

second to last comparison:

firstly, calculating the magnitude relation between the similarity coefficient Sc of the current value and a threshold value L: if Sc is larger than L, setting the comparison value of the current time as 0, and if Sc is less than or equal to L, setting the comparison value of the current time as 1;

then, the comparison value of the current time is judged to be the same as the comparison value of the previous time: if the two values are the same, judging that the current time and the previous time have the same acceleration type, recording the acceleration type mark of the previous time, and keeping the comparison value of the current time as 1 for the next time of same judgment; if the two values are different, the acceleration type of the current time is judged to be different from that of the previous time, a new acceleration type mark is generated at the moment, and the comparison value of the current time is set to be 1 for the next same judgment;

s3.2, constructing M gain classification sequences T with the same acceleration according to the acceleration type marks obtained in S3.1_{n_2nd}' (k); wherein M is more than or equal to 1 and less than or equal to k.

6. The improved rotating machinery non-stationary vibration signal computational order analysis method of claim 1, wherein: the step S4 of reversely deriving the gain classification sequence in sequence from the direction of the matching time difference sequence and the time point sequence to obtain a classification time point sequence, where S4 specifically includes:

s4.1, classifying the gain classified sequence T into M types_{n_2nd}' (k) is divided into two types of push-back: the 1 st time and the M th time, wherein M is more than or equal to 2 and less than or equal to M; wherein, the 1 st gain classification sequence refers to a sequence with the addition degree type marked as the 1 st type, and the mth gain classification sequence has the same principle;

s4.2, converting T obtained in S4.1_{n_2nd}' (k) increasing the maximum index value belonging to the index range of the 1 st acceleration type by 1 in the 1 st gain classification sequence as the 1 st acceleration type at T_{n_1nd}' (j) and the 1 st acceleration type at T_{n_1nd}The starting index value in' (j) is 1;

s4.3, converting T obtained in S4.2_{n_1nd}' (j) increasing the maximum index value belonging to the index range of the 1 st acceleration type by 1 in the 1 st time difference classification sequence as the 1 st acceleration type at T_n' (i) last index value, 1 st acceleration type at T_nThe starting index value in' (i) is 1;

s4.4, converting T obtained in S4.1_{n_2nd}' (k) in the m-th gain classification sequence, the initial and final index values belonging to the m-th acceleration type range are both increased by 1, and the m-th acceleration type is represented by T_{n_1nd}The start and end index values in (j);

s4.5, converting T obtained in S4.4_{n_1nd}' (j) increasing the maximum index value belonging to the index range of the mth acceleration type in the mth time difference class sequence by 1 as the mth acceleration type in T_n' (i) the last index value, the initial index value being its value at T_{n_1nd}' (j) initial index values of the index ranges; sorting out time point classification sequence T_n' (i) index range point values of each acceleration degree type, removing overlapped index point values, arranging in order, and obtaining time point sequence T by point values on the arrangement_n(i) The corresponding index value is the time point of the sorting, thereby constructing a sorting time point sequence T which is divided into M types_n″(i)。

7. The improved rotating machinery non-stationary vibration signal computational order analysis method of claim 1, wherein: the step S5 of performing angular domain resampling on the non-stationary vibration signal using the sorted time point sequence to obtain an angular domain resampled signal, where the step S5 specifically includes:

s5.1, sorting time point sequence T divided into M types_n"(i), the time range for obtaining M-class equivalent accelerations is t₁"time to t₂"time, the angle through which the reference axis rotates in this time frame is 0 to b x 2 pi, the formula exists:

carry in t₁″、t₂"time can be calculated to obtain polynomial coefficient a₀、a₁(ii) a Wherein b is t₁"time to t₂"number of full revolutions of the reference axis during time, t₁To sort time point sequences T_nThe starting time point, t, of a certain type of acceleration time range in (i)₂To sort time point sequences T_nThe end time point of the acceleration time range of certain type in (i);

s5.2, a obtained in S5.1₀、a₁Substitution into the formula:the time sequence t corresponding to the equal angle sampling in the range of 0 to b x 2 pi can be obtained by calculation_{n_1nd}(ii) a Wherein p is an interpolation coefficient, and delta theta is an equal angle interval;

s5.3, with the time series t obtained in S5.2_{n_1nd}And (5) carrying out interpolation solution on the non-stationary vibration signal x (t) to obtain an angular domain resampling signal w (t).

Images