Blade vibration parameter testing method based on blade tip timing signal

Abstract

The application provides a kind of blade vibration parameter test method based on blade tip timing signal, it is related to mechanical blade vibration technical field, comprising the following steps: real-time acquisition of the rotational speed function dynamically constructed by rotational speed sensor signal, calculate the theoretical arrival time of each pulse, obtain the arrival time difference caused by blade vibration, calculate the vibration displacement of blade under the current rotational speed fluctuation condition;According to the rotational speed function, the theoretical instantaneous speed of blade passing through the sensor is calculated, the corresponding theoretical pulse width is obtained, the pulse width deviation caused by blade vibration is obtained, and the vibration speed of blade under the current rotational speed fluctuation condition is calculated;Vibration displacement and vibration speed are integrated into multi-measurement vector sparse reconstruction model and jointly solved, and the sparse frequency support set is recovered, to obtain the frequency parameters of blade vibration;Self-adaptive adjustment of the regularization parameter in block sparse bayesian learning method, and the frequency and amplitude of blade vibration are output using the adjusted parameter.

Description

Technical Field

[0001] This invention belongs to the field of mechanical blade vibration technology, specifically relating to a method for testing blade vibration parameters based on blade tip timing signals. Background Technology

[0002] In the field of rotating machinery blade vibration testing, measurement methods based on blade tip timing signals are widely used due to their advantages such as non-contact operation and full blade monitoring. In existing technologies, blade tip timing sensors mounted on the casing are typically used to record the blade's arrival time. The arrival time difference is obtained by comparing the measured arrival time with the theoretical arrival time, and then the blade's vibration displacement is calculated. Finally, spectral analysis or sparse reconstruction of the vibration displacement is performed to obtain the vibration frequency and amplitude.

[0003] However, rotating machinery often operates under variable speed conditions, and speed fluctuations are unavoidable. Existing methods suffer from the following technical problems: First, speed fluctuations lead to significant errors in the calculation of theoretical arrival time and theoretical pulse width, causing a sharp drop in the signal-to-noise ratio of vibration displacement signals extracted solely based on the time difference of arrival. Second, the pulse waveform output by the blade tip timing sensor contains not only arrival time information but also pulse width information closely related to the blade's instantaneous velocity. However, most existing technologies only utilize the time difference of arrival, neglecting the blade vibration velocity information reflected by the pulse width, resulting in insufficient utilization of measurement information. When speed fluctuations are large or the signal-to-noise ratio is low, it is difficult to accurately recover the blade's vibration frequency based solely on vibration displacement, leading to insufficient measurement accuracy and robustness. Summary of the Invention

[0004] In view of the above-mentioned defects or deficiencies in the prior art, a method for testing blade vibration parameters based on blade tip timing signals is provided, including the following steps: During the variable speed operation of rotating machinery, a blade tip timing sensor installed on the casing is used to collect the pulse waveform signal when the blade sweeps across, and the measured arrival time and measured pulse width of each pulse are extracted from the pulse waveform signal. The speed sensor signal is acquired in real time, and a speed function describing the change of rotor speed with time is dynamically constructed based on the speed sensor signal; the theoretical arrival time of each pulse is calculated based on the speed function; the measured arrival time is compared with the theoretical arrival time to obtain the arrival time difference caused by blade vibration; the vibration displacement of the blade under the current speed fluctuation condition is calculated based on the arrival time difference. Obtain the pre-calibrated mapping relationship between pulse width and instantaneous velocity; calculate the theoretical instantaneous velocity of the blade when it passes the sensor based on the rotational speed function, and obtain the corresponding theoretical pulse width using the mapping relationship; compare the measured pulse width with the theoretical pulse width to obtain the pulse width deviation caused by blade vibration; calculate the vibration velocity of the blade under the current rotational speed fluctuation condition based on the pulse width deviation; The vibration displacement and the vibration velocity are integrated into a multi-measurement vector sparse reconstruction model, in which the vibration displacement and vibration velocity share the same sparse frequency support set in the frequency domain. The block sparse Bayesian learning method is used to jointly solve the multi-measurement vector sparse reconstruction model to recover the sparse frequency support set and obtain the frequency parameters of the blade vibration. Based on the current speed fluctuation amplitude and signal-to-noise ratio level, the regularization parameters in the block sparse Bayesian learning method are adaptively adjusted, and the frequency and amplitude of blade vibration are output using the adjusted parameters.

[0005] According to the technical solution provided in this application, before comparing the measured pulse width with the theoretical pulse width, a gap compensation step is further included: Obtain the pre-calibrated mapping relationship between the tip clearance and the pulse width compensation amount; Based on the current rotational speed and the pre-established rotor thermal deformation model, estimate the change in blade tip clearance under the current operating conditions; Using the mapping relationship, calculate the pulse width compensation amount corresponding to the change in blade tip clearance; Subtract the pulse width compensation amount from the measured pulse width to obtain the compensated measured pulse width.

[0006] According to the technical solution provided in this application, before integrating the vibration displacement and the vibration velocity into a multi-measurement vector sparse reconstruction model, the following steps are also included: Vibration displacement data and vibration velocity data are used as non-uniform sampling sequences to be processed; Based on the circumferential installation angle of each blade tip timing sensor and the rotational speed function, the actual sampling time corresponding to each pulse is calculated to form a non-uniform sampling time grid. A uniform sampling time grid is constructed, and the sampling interval of the uniform sampling time grid is adaptively determined according to the current speed fluctuation amplitude; The vibration displacement data and vibration velocity data on the non-uniform sampling sequence are mapped onto the uniform sampling time grid using the distance-weighted interpolation method to obtain the resampled uniform displacement sequence and uniform velocity sequence. The process of integrating the vibration displacement and the vibration velocity into a multi-measurement vector sparse reconstruction model includes the following steps: The resampled uniform displacement sequence and uniform velocity sequence are used as inputs to the multi-measurement vector sparse reconstruction model.

[0007] According to the technical solution provided in this application, the method of jointly solving the multi-measurement vector sparse reconstruction model using a block sparse Bayesian learning method includes the following steps: Initialize the noise variance parameters for the displacement and velocity channels; During the expectation maximization iteration of block sparse Bayesian learning, the noise variance estimates of each channel are updated using displacement residuals and velocity residuals respectively. Based on the currently updated noise variance estimate, adaptive weights are calculated for the displacement and velocity channels, with the channel with smaller noise variance receiving a larger weight. When updating the sparse frequency support set in each iteration, the likelihood terms of the displacement measurement vector and the velocity measurement vector are weighted using the adaptive weights so that the channel with smaller noise variance dominates the joint solution. Iterate until convergence, and output the final sparse frequency support set.

[0008] According to the technical solution provided in this application, the step of calculating the vibration displacement of the blade under the current rotational speed fluctuation condition based on the arrival time difference includes the following steps: The system acquires multiple measured pulse widths of the same blade within the same rotor revolution by using multiple blade tip timing sensors installed at different circumferential angles on the casing. Gap compensation is performed on each measured pulse width to obtain multiple compensated pulse widths; Obtain the pre-calibrated mapping relationship between pulse width and instantaneous velocity, and convert each compensated pulse width into the instantaneous linear velocity of the blade at the corresponding sensor. Based on the instantaneous linear velocity and current rotational speed function at each sensor, calculate the first tip radius at each sensor, where the first tip radius is equal to the instantaneous linear velocity divided by the tip angular velocity at the current rotational speed; The first blade tip radius obtained from multiple sensors is fused to obtain the dynamic blade tip radius at the current rotation speed; The blade tip linear velocity is calculated using the dynamic blade tip radius, and then the arrival time difference is converted into vibration displacement using the blade tip linear velocity.

[0009] According to the technical solution provided in this application, the speed sensor is a sensor installed on a toothed disk, and the toothed disk has multiple teeth evenly distributed around its circumference. Before dynamically constructing the speed function describing the rotor speed change over time based on the speed sensor signal, the method further includes the following steps: The arrival time of each pulse is extracted from the speed sensor signal, and the pulse interval value corresponding to each tooth is calculated based on the difference in arrival times of adjacent pulses. The calculated pulse interval values ​​are sequentially filled into the circular buffer in the circumferential order of the toothed disk. The length of the circular buffer is equal to the total number of teeth on the toothed disk. Each time a new pulse interval value is entered, the standard deviation of all pulse interval values ​​in the circular buffer is calculated, and the sum of the absolute values ​​of the differences between pulse interval values ​​at adjacent positions in the circular buffer is also calculated. The step of dynamically constructing a speed function describing the change of rotor speed over time based on the speed sensor signal includes the following steps: If the standard deviation is greater than or equal to the first threshold, or the sum of the absolute values ​​of the differences is less than or equal to the second threshold, then the original pulse interval value is directly used to construct the rotational speed function through interpolation.

[0010] According to the technical solution provided in this application, after calculating the standard deviation of all pulse interval values ​​within the circular buffer and simultaneously calculating the sum of the absolute values ​​of the differences between pulse interval values ​​at adjacent positions within the circular buffer, the method further includes the following steps: If the standard deviation is less than the first threshold and the sum of the absolute values ​​of the differences is greater than the second threshold, it is determined that there is a periodic measurement error in the speed sensor signal caused by the eccentricity of the gear plate. The pulse interval value in the cyclic buffer is corrected position by position, and the speed function is constructed by interpolation using the corrected pulse interval value.

[0011] According to the technical solution provided in this application, the step of correcting the pulse interval value in the circular buffer position by position includes the following steps: During variable speed operation, whenever a new pulse interval value is filled into the cyclic buffer, the pre-built eccentricity error table is retrieved according to the tooth number of the corresponding tooth, and the eccentricity error compensation value corresponding to the tooth number is read. The pulse interval value is then subtracted from the read eccentricity error compensation value to obtain the corrected pulse interval value. The eccentricity error table is constructed through the following steps: During the startup phase or calibration condition after maintenance of rotating machinery, the machine is run at a constant speed at one-third lower than the first-order bending vibration natural frequency of the blade. The original pulse interval value of each tooth is collected during one complete rotation of the toothed disc. The original pulse interval value of each tooth is subtracted from the arithmetic mean of the pulse interval values ​​of all teeth in that rotation to obtain the eccentricity error compensation value of each tooth. The eccentricity error compensation value is then stored in an eccentricity error table according to the circumferential order of the teeth.

[0012] According to the technical solution provided in this application, after dynamically constructing a speed function describing the rotor speed change over time based on the speed sensor signal, and before calculating the theoretical instantaneous speed of the blade when it passes the sensor based on the speed function, the following steps are further included: Obtain the vibration displacement data of all blades within the same rotor rotation circle, and calculate the statistical average value of the vibration displacement of all blades; If the statistical average value is not within the static balance threshold range, it is determined that the current rotational speed function has a systematic low-frequency drift error. Based on the sign and magnitude of the statistical average, the rotational speed function is subjected to overall translation correction so that the statistical average of all blade vibration displacements after correction is within the range of the static balance threshold.

[0013] According to the technical solution provided in this application, after outputting the frequency and amplitude of blade vibration using the adjusted parameters, the method further includes the following steps: At the same rotational speed, the steps are repeated for all blades on the same stage rotor to obtain the vibration frequency and vibration amplitude of each blade; Calculate the average frequency and standard deviation of all blade vibration frequencies, and the average amplitude and standard deviation of all blade vibration amplitudes; For each blade, the absolute value of the frequency deviation between its vibration frequency and the average frequency is calculated. When the absolute value of the frequency deviation is greater than a preset first multiple multiplied by the frequency standard deviation, the blade is determined to be a frequency abnormal blade. For each blade, calculate the absolute value of the amplitude deviation between its vibration amplitude and the average amplitude. When the absolute value of the amplitude deviation is greater than a preset second multiple multiplied by the amplitude standard deviation, the blade is determined to be an amplitude abnormal blade. Output the number and type of all abnormal blades.

[0014] Compared with existing technologies, the advantages of this application are as follows: This application simultaneously extracts the measured arrival time and measured pulse width of the blade tip timing pulse during variable speed operation, calculates the vibration displacement and vibration velocity of the blade under the current speed fluctuation conditions, and integrates them into a multi-measurement vector sparse reconstruction model, which is then jointly solved using a block sparse Bayesian learning method. The vibration displacement and vibration velocity share the same sparse frequency support set in the frequency domain, which is equivalent to using two measurement information with different physical meanings and noise characteristics to mutually verify and complement each other, significantly improving the accuracy and noise resistance of frequency parameter recovery. Simultaneously, the regularization parameter is adaptively adjusted according to the current speed fluctuation amplitude and signal-to-noise ratio level, enabling the sparse reconstruction process to dynamically adapt to changes in operating conditions, avoiding the over-smoothing or under-sparseness problems caused by fixed parameters under variable speed conditions. This method is particularly suitable for blade vibration monitoring during the start-up, shutdown, or speed change processes of rotating machinery such as aero-engines and gas turbines, effectively suppressing measurement errors caused by speed fluctuations and achieving high-precision blade vibration frequency and amplitude testing. Attached Figure Description

[0015] Other features, objects, and advantages of this application will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings: Figure 1 A flowchart illustrating the steps of the blade vibration parameter testing method based on the blade tip timing signal provided in this application. Detailed Implementation

[0016] The present application will now be described in further detail with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the invention and not intended to limit it. Furthermore, it should be noted that, for ease of description, only the parts relevant to the invention are shown in the accompanying drawings.

[0017] It should be noted that, unless otherwise specified, the embodiments and features described in this application can be combined with each other. This application will now be described in detail with reference to the accompanying drawings and embodiments.

[0018] As mentioned in the background section, this application proposes a method for testing blade vibration parameters based on blade tip timing signals, such as... Figure 1 As shown, it includes the following steps: S1. During the variable speed operation of rotating machinery, the blade tip timing sensor installed on the casing is used to collect the pulse waveform signal when the blade sweeps by, and the measured arrival time and measured pulse width of each pulse are extracted from the pulse waveform signal. S2. Acquire the speed sensor signal in real time, and dynamically construct a speed function describing the change of rotor speed with time based on the speed sensor signal; calculate the theoretical arrival time of each pulse based on the speed function; compare the measured arrival time with the theoretical arrival time to obtain the arrival time difference caused by blade vibration; calculate the vibration displacement of the blade under the current speed fluctuation condition based on the arrival time difference. S3. Obtain the pre-calibrated mapping relationship between pulse width and instantaneous velocity; calculate the theoretical instantaneous velocity of the blade when it passes the sensor according to the rotational speed function, and obtain the corresponding theoretical pulse width using the mapping relationship; compare the measured pulse width with the theoretical pulse width to obtain the pulse width deviation caused by blade vibration; calculate the vibration velocity of the blade under the current rotational speed fluctuation condition based on the pulse width deviation. S4. Integrate the vibration displacement and the vibration velocity into a multi-measurement vector sparse reconstruction model, wherein the vibration displacement and vibration velocity in the multi-measurement vector sparse reconstruction model share the same sparse frequency support set in the frequency domain. S5. The block sparse Bayesian learning method is used to jointly solve the multi-measurement vector sparse reconstruction model to recover the sparse frequency support set and obtain the frequency parameters of the blade vibration. S6. Based on the current speed fluctuation amplitude and signal-to-noise ratio level, adaptively adjust the regularization parameters in the block sparse Bayesian learning method, and use the adjusted parameters to output the frequency and amplitude of blade vibration.

[0019] Specifically, firstly, at least one blade tip timing sensor is installed on the casing of rotating machinery (such as aircraft engines, gas turbines, compressors, etc.). This sensor, typically a fiber optic, capacitive, or eddy current sensor, is used to sense the changes in the physical field generated when the blade sweeps past. When the blade rotates past directly below the sensor probe, the sensor outputs an analog pulse waveform. The rising and falling edges of this pulse correspond to the times when the leading and trailing edges of the blade pass through a specific reference position (such as the center of the sensor spot), respectively. The data acquisition system digitizes this analog signal at a fixed high sampling rate (e.g., 100 MHz) and, through threshold comparison or peak detection algorithms, accurately extracts the measured arrival time of each pulse (usually defined as the time corresponding to the 50% amplitude point of the pulse leading edge) and the measured pulse width (defined as the time difference between the 50% amplitude point of the pulse leading edge and the 50% amplitude point of the pulse trailing edge).

[0020] Simultaneously, the system acquires speed sensor signals in real time. The speed sensor, typically installed near the rotor's toothed disc, employs a magnetoelectric or Hall effect sensor to detect the passage of uniformly distributed teeth (e.g., 60 teeth) on the disc. Each detected tooth pulse is recorded at its arrival time. Since the rotor's speed continuously changes with time during variable-speed operation, the theoretical arrival time cannot be simply calculated using a constant speed. This method dynamically constructs a speed function, i.e., an angular velocity function ω(t), describing the rotor's speed change over time. Specifically, the instantaneous angular velocity ω_i = 2π / (N·Δt_i) corresponding to the time interval Δt_i = t_{i+1} - t_i between two adjacent tooth pulses (where N is the total number of teeth on the toothed disc) is used to obtain a series of discrete angular velocity-time points (t_i, ω_i). Then, cubic spline interpolation or piecewise linear interpolation is used to construct a continuously differentiable speed function ω(t) = dθ / dt, where θ(t) is a function of the rotor angle changing with time. This speed function accurately reflects acceleration, deceleration, or speed fluctuations.

[0021] Based on the rotational speed function, the theoretical arrival time of each blade tip timing pulse can be calculated. Let the circumferential installation angle of the j-th sensor be α_j (starting from a certain reference zero position), and the blade number be k. The theoretical arrival time is defined as: when the leading edge of blade k just reaches the measurement position of sensor j under vibration-free conditions and when the rotor rotates strictly according to the rotational speed function ω(t). The specific calculation method is as follows: given the measured arrival time t0 of blade k passing the reference sensor (such as the sensor installed at 0° position) in the previous revolution, the theoretically required angle increment for blade k to rotate to angle α_j again is α_j (if α_j is less than the cumulative angle of the previous revolution, 2π needs to be added). By numerically integrating the rotational speed function ω(t) starting from t0: θ(t) = Solving the equation θ(t) = α_j + 2π·m (where m is the number of revolutions) yields the theoretical arrival time t_theory. Subtracting the measured arrival time t_actual from t_theory gives the arrival time difference Δt_tof = t_actual - t_theory. This arrival time difference is mainly caused by the circumferential or axial displacement of the blade leading edge due to blade vibration. Given the known tip linear velocity U(t) = ω(t)·R (where R is the tip radius), the vibration displacement d = U(t)·Δt_tof. However, since rotational speed fluctuations affect the instantaneous value of U(t), this method uses a rotational speed function to calculate the instantaneous angular velocity ω(t_actual) at the current moment, then multiplies it by the nominal tip radius R0 to obtain the instantaneous linear velocity, and subsequently calculates the vibration displacement.

[0022] Meanwhile, this method utilizes pulse width information to extract vibration velocity. The pulse width (i.e., the time it takes for the blade to sweep across the sensor's measured spot) is inversely proportional to the blade's instantaneous linear velocity: when blade vibration causes a change in tangential velocity, the pulse width changes accordingly. A mapping relationship f between pulse width and instantaneous velocity, v = f(PW), is obtained beforehand through calibration experiments. This calibration is performed while the rotating machinery is rotating at a low, uniform speed. High-precision equipment such as a laser vibrometer is used to simultaneously measure the actual linear velocity of the blade, record the corresponding pulse width, and fit the mapping curve (usually an inverse proportional relationship or a polynomial relationship). Based on the rotational speed function, the theoretical instantaneous velocity v_theory = ω(t_theory)·R0 when the blade passes the sensor is calculated, and then the corresponding theoretical pulse width PW_theory = f is obtained using the mapping relationship. -1 (v_theory). Comparing the measured pulse width PW_actual with PW_theory, we obtain the pulse width deviation ΔPW = PW_actual - PW_theory. Based on physical relationships, the relationship between the blade tangential vibration velocity v_vib and the pulse width deviation is: v_vib = -(v_theory) 2 The formula is obtained by differentiating PW = K / v (where K is a calibration constant). This allows for the calculation of the vibration velocity under the current rotational speed fluctuation conditions.

[0023] After obtaining the vibration displacement and vibration velocity sequences, they need to be integrated into a multi-measurement vector sparse reconstruction model. Blade vibration typically consists of several dominant frequency components (such as harmonics, blade natural frequencies, etc.), exhibiting sparsity in the frequency domain. Let the vibration displacement signal be x_d(t) and the vibration velocity signal be x_v(t). Both are composed of the same frequency components, differing only in amplitude and phase. Therefore, in the frequency domain, they share the same sparse frequency support set (i.e., the set of non-zero frequency positions). Construct the multi-measurement vector model: Y = AX + E, where Y is the measurement matrix (dimensional M×2, where M is the number of sampling points, and each column represents the displacement and velocity measurements), A is the overcomplete dictionary matrix (composed of sine and cosine basis functions or Fourier basis), X is the sparse coefficient matrix (each column corresponds to the frequency domain coefficients of one channel, but the non-zero rows are in the same position), and E is noise. The sparse reconstruction objective of this model is to recover rows of sparse X (i.e., each row either has two zero columns or two non-zero columns) given Y and A.

[0024] For this model, a block sparse Bayesian learning method is used for joint solution. Block sparse Bayesian learning treats the displacement and velocity coefficients corresponding to each frequency as a block, with coefficients within the block sharing the same hyperparameters (controlling whether the block is activated). Specifically, a sparse prior, such as one determined by autocorrelation, is assigned to each frequency block, and then the hyperparameters and noise variance are iteratively updated using an expectation-maximization algorithm. In each iteration, the reconstructed signal is calculated based on the current posterior mean, and the block activation probability is updated. After iteration convergence, the frequencies corresponding to the non-zero blocks are the frequency parameters of the blade vibration.

[0025] Since the speed fluctuation amplitude and signal-to-noise ratio level vary in actual testing, this method introduces an adaptive regularization parameter adjustment mechanism. The regularization parameter λ is defined as follows: 2 Related to the hyperparameter γ of the sparse prior: λ ∝ σ 2 / γ. In the block sparse Bayesian learning framework, the initial values ​​of the block prior and the update strategy of the noise variance are dynamically adjusted by monitoring the current speed fluctuation amplitude (e.g., the peak-to-peak speed change rate) and the signal-to-noise ratio (SNR) (estimated by the ratio of the pulse waveform amplitude to the root mean square of the background noise). When the speed fluctuation is large or the SNR is low, the regularization parameter is increased to emphasize sparsity and prevent overfitting; conversely, the regularization parameter is decreased to retain more details. Finally, using the adjusted parameters, the frequency (i.e., the frequency value corresponding to the support set) and amplitude (given by the coefficient amplitude of the displacement channel, or obtained by weighting the two channels) of the blade vibration are output.

[0026] The technical principle of this solution lies in the fact that blade vibration simultaneously affects both pulse arrival time and pulse width. The former mainly reflects displacement information, while the latter mainly reflects velocity information, and both have the same frequency components in the frequency domain. Through sparse reconstruction of multiple measurement vectors, the complementarity of displacement and velocity can be utilized to improve the robustness and accuracy of frequency identification. Adaptive regularization parameter adjustment ensures the stability of the algorithm under different operating conditions. The technical benefits include: accurate extraction of blade vibration frequency and amplitude even under conditions of large fluctuations in rotational speed; strong noise resistance; and avoidance of spurious frequencies caused by rotational speed measurement errors or insufficient information from a single physical quantity, as in traditional methods.

[0027] In a preferred embodiment, before comparing the measured pulse width with the theoretical pulse width, a gap compensation step is further included: Obtain the pre-calibrated mapping relationship between the tip clearance and the pulse width compensation amount; Based on the current rotational speed and the pre-established rotor thermal deformation model, estimate the change in blade tip clearance under the current operating conditions; Using the mapping relationship, calculate the pulse width compensation amount corresponding to the change in blade tip clearance; Subtract the pulse width compensation amount from the measured pulse width to obtain the compensated measured pulse width.

[0028] Specifically, a gap compensation step is further introduced to eliminate pulse width measurement errors caused by changes in blade tip clearance, thereby improving the accuracy of vibration velocity calculation. Blade tip clearance refers to the radial distance between the tip of the rotating blade and the inner wall of the casing. In actual operation, the blade tip clearance is not constant but is affected by various factors such as blade elongation caused by centrifugal force, changes in casing and blade dimensions caused by thermal expansion, and rotor radial movement. Changes in blade tip clearance alter the effective spot size or electromagnetic field coupling strength when the blade sweeps across the sensor, directly affecting the measured pulse width: when the gap increases, the sensor signal strength weakens, and the pulse width often increases due to the shift in the threshold crossing point; when the gap decreases, the pulse width decreases accordingly. This pulse width variation caused by gap changes is unrelated to blade vibration; without compensation, it will lead to spurious components in the vibration velocity calculation.

[0029] The gap compensation steps in this embodiment are as follows: First, a mapping relationship between the blade tip clearance and the pulse width compensation is established beforehand. This mapping relationship is obtained through calibration experiments. On a dedicated calibration fixture or test bench of rotating machinery, the blade is mounted on a mechanism capable of precisely controlling radial displacement. Simultaneously, a laser displacement sensor or a capacitive clearance gauge is used to measure the actual blade tip clearance value *g*, and the pulse width *PW* output by the blade tip timing sensor is recorded. With the blade stationary or rotating at extremely low speed (to avoid vibration interference), the actual blade tip clearance value *g* is changed, resulting in a series of (g, PW) data pairs. Since rotational speed also affects the pulse width (velocity effect), calibration should be performed at multiple constant rotational speeds. Finally, the clearance compensation function ΔPW_comp = h(g, ω) is fitted, which represents the required pulse width compensation ΔPW_comp = h(Δg, ω) given an instantaneous angular velocity ω and the deviation Δg of the actual blade tip clearance value *g* from the nominal clearance *g0*. This function can be a lookup table stored in tabular form or a polynomial fitting formula.

[0030] Secondly, based on the current rotational speed and the pre-established rotor thermal deformation model, the change in blade tip clearance under the current operating conditions is estimated. The rotor thermal deformation model includes: the thermal expansion model of the casing (calculating radial deformation based on the thermal expansion coefficient and temperature field distribution of the casing material), the thermal elongation model of the blades (radial growth of the blade material caused by heating and centrifugal force), and the change in bearing clearance with temperature. In actual operation, the characteristic temperature under the current operating conditions is obtained by temperature sensors installed on the casing and rotor (or by obtaining the temperature-speed mapping relationship through thermodynamic simulation), and substituted into the thermal deformation model to calculate the change in blade tip clearance Δg_thermal. In addition, the radial elongation of the blades caused by centrifugal force Δg_centrifugal = (ρ×ω) should also be considered. 2 ×R 2 ) / (2E) (approximate formula), where ρ is density, E is Young's modulus, and R is tip radius. The total clearance change Δg = Δg_thermal + Δg_centrifugal + Δg_other (e.g., wear).

[0031] Then, using the pre-calibrated mapping relationship h, the pulse width compensation amount ΔPW_comp = h(Δg, ω(t)) corresponding to the current blade tip clearance change Δg is calculated, where ω(t) is the instantaneous angular velocity given by the rotational speed function at the current moment. Since the pulse width change caused by the clearance change is gradual (relative to the blade vibration frequency), this compensation amount can be considered to be basically constant within one rotation cycle, or updated once per revolution.

[0032] Finally, the measured pulse width PW_actual is subtracted from the pulse width compensation ΔPW_comp to obtain the compensated measured pulse width PW_corrected = PW_actual - ΔPW_comp. This compensated pulse width eliminates the offset introduced by the gap variation and reflects the pulse width change caused only by the blade tangential vibration velocity. PW_corrected is used for subsequent comparison with the theoretical pulse width and for vibration velocity calculation.

[0033] The technical principle of this solution lies in the fact that the influence of blade tip clearance variation on pulse width and the influence of blade vibration on pulse width are separable in the frequency domain—clearance variation typically manifests as a low-frequency or quasi-static component, while blade vibration frequency is often several times higher than the rotational frequency. By pre-calibrating and estimating the clearance variation using a thermodynamic model, and then performing feedforward compensation, this systematic error can be effectively filtered out. The technical benefits include: significantly improving the accuracy of vibration velocity measurement during speed changes, avoiding misinterpretation of thermal expansion or centrifugal deformation as blade vibration, especially under start-up acceleration or rapid operating condition changes. Clearance compensation can bring the baseline of the vibration velocity signal to zero, thereby ensuring the purity of the frequency support set in subsequent sparse reconstruction. Furthermore, this compensation step is also applicable to scenarios where instantaneous linear velocity is calculated using pulse width to obtain the dynamic blade tip radius, improving the accuracy of vibration displacement calculation.

[0034] In a preferred embodiment, before integrating the vibration displacement and the vibration velocity into a multi-measurement vector sparse reconstruction model, the following steps are further included: Vibration displacement data and vibration velocity data are used as non-uniform sampling sequences to be processed; Based on the circumferential installation angle of each blade tip timing sensor and the rotational speed function, the actual sampling time corresponding to each pulse is calculated to form a non-uniform sampling time grid. A uniform sampling time grid is constructed, and the sampling interval of the uniform sampling time grid is adaptively determined according to the current speed fluctuation amplitude; The vibration displacement data and vibration velocity data on the non-uniform sampling sequence are mapped onto the uniform sampling time grid using the distance-weighted interpolation method to obtain the resampled uniform displacement sequence and uniform velocity sequence. The process of integrating the vibration displacement and the vibration velocity into a multi-measurement vector sparse reconstruction model includes the following steps: The resampled uniform displacement sequence and uniform velocity sequence are used as inputs to the multi-measurement vector sparse reconstruction model.

[0035] Specifically, the blade tip timing sensors are mounted in a circumferentially fixed position (e.g., evenly distributed around the circumference of the casing, or arranged at a specific angle), while the blades rotate at a non-constant speed. Therefore, the time intervals between the blades sweeping past each sensor are not equal. Furthermore, variable speed operation causes the sampling interval to vary over time. This results in a typical non-uniform sampling sequence: the sampling time corresponding to each pulse is determined by the actual physical time when the blade reaches the sensor, rather than by a manually set equal interval. Traditional Fourier transforms or conventional sparse reconstruction methods require the input signal to be sampled on a uniform time grid; otherwise, spectral leakage and spurious peaks will occur.

[0036] The specific steps of this implementation method are as follows: First, the calculated vibration displacement data (each pulse corresponds to a vibration displacement value) and vibration velocity data (each pulse corresponds to a vibration velocity value) are used as non-uniform sampling sequences to be processed. There are N pulses in total, each pulse corresponding to a measured arrival time t_i (extracted from the original pulse, not the theoretical arrival time), as well as vibration displacement d_i and vibration velocity v_i. These (t_i, d_i) and (t_i, v_i) constitute two non-uniformly sampled discrete signals.

[0037] Secondly, based on the circumferential installation angle and rotational speed function of each blade tip timing sensor, the actual sampling time corresponding to each pulse is calculated. It is important to note here that due to rotational speed fluctuations, the time interval for sampling the same blade by the same sensor in different revolutions is not uniform. In fact, the measured arrival time t_i of each pulse has already been accurately recorded by the data acquisition system, so there is no need to recalculate it; t_i can be directly used as the actual sampling time. However, to construct a uniform grid later, the precise time position of each sampling point needs to be determined. The core of this step is to form a non-uniform sampling time grid T_nonuniform = {t_i, i=1..N}.

[0038] Then, a uniform sampling time grid is constructed. The sampling interval ΔT of this grid needs to be adaptively determined according to the current speed fluctuation amplitude. The specific method is as follows: analyze the range of change of the speed function ω(t) within the measurement time interval, let the maximum speed be ω_max, the minimum speed be ω_min, and the speed fluctuation amplitude Δω = ω_max - ω_min. To avoid spectral aliasing, the uniform sampling frequency should be at least twice the highest expected frequency of blade vibration (e.g., 3 times the first-order bending natural frequency of the blade). However, if the speed fluctuation is severe, too small a ΔT will result in too many uniform grid points, increasing the computational burden; too large a ΔT may lose details. The adaptive strategy is: ΔT = min(1 / (2·f_max),0.1·(2π / ω_avg) ), where f_max is the pre-set upper limit of the blade vibration frequency (e.g., 5 kHz), and ω_avg is the average angular velocity. A more refined approach is to dynamically adjust according to the speed fluctuation amplitude: when the fluctuation is large, appropriately reduce ΔT to capture the vibration modulation caused by the speed change; when the fluctuation is small, increase ΔT to reduce the computational load. The uniform sampling time grid is T_uniform = {t_start + q·ΔT, q=0,1,...,K-1}, where t_start and t_end cover the range of all non-uniform sampling times, q is the index of the uniform sampling point (an integer from 0, 1, 2...K-1), and K is the total number of uniform sampling points.

[0039] Next, distance-weighted interpolation is used to map the non-uniform sampling sequence onto a uniform sampling time grid. Commonly used distance-weighted interpolation methods include inverse distance-weighted interpolation and linear interpolation. This implementation recommends piecewise linear interpolation: for each target time τ on the uniform grid, find two adjacent times t_i and t_{i+1} in the non-uniform time sequence such that t_i≤τ≤t_{i+1}, then the interpolated displacement value d(τ) = d_i + (d_{i+1}-d_i)·(τ-t_i) / (t_{i+1}-t_i). The velocity value is handled similarly. This linear interpolation method is simple and efficient, and the error introduced is controllable when the sampling rate is sufficient. For points near the boundary, extrapolation or boundary filling can be used. To reduce interpolation error, the average density of the non-uniform sampling should be more than twice that of the uniform sampling density; otherwise, it is necessary to increase the number of original sensors or increase the sampling rate of the rotational speed signal.

[0040] After the above resampling, a uniform displacement sequence {d(τ_k)} and a uniform velocity sequence {v(τ_k)} are obtained, which have the same time grid and the same number of points K. These two sequences are the resampled uniform displacement sequence and uniform velocity sequence. Subsequently, in the step of "integrating vibration displacement and vibration velocity into a multi-measurement vector sparse reconstruction model", these two uniform sequences are used as two columns of the input Y of the multi-measurement vector sparse reconstruction model.

[0041] The technical advantages of this implementation include: eliminating spectral leakage caused by non-uniform sampling, enabling sparse reconstruction to accurately identify vibration frequencies; adaptive meshing ensures a balance between computational efficiency and accuracy under different rotational speed fluctuations; and providing a standardized data format for subsequent joint solutions of multiple measurement vectors.

[0042] In a preferred embodiment, the step of jointly solving the multi-measurement vector sparse reconstruction model using a block sparse Bayesian learning method includes the following steps: Initialize the noise variance parameters for the displacement and velocity channels; During the expectation maximization iteration of block sparse Bayesian learning, the noise variance estimates of each channel are updated using displacement residuals and velocity residuals respectively. Based on the currently updated noise variance estimate, adaptive weights for the displacement and velocity channels are calculated. If the noise variance of the displacement channel is less than that of the velocity channel, the displacement channel receives a weight greater than that of the velocity channel; otherwise, the velocity channel receives a weight greater than that of the displacement channel. When updating the sparse frequency support set in each iteration, the likelihood terms of the displacement measurement vector and the velocity measurement vector are weighted using the adaptive weights so that the channels with small noise variance dominate in the joint solution. Iterate until convergence, and output the final sparse frequency support set.

[0043] Specifically, the steps are as follows: First, initialize the noise variance parameters of the displacement and velocity channels. Let the measurement matrix Y = [y_d, y_v], where y_d is a uniform displacement sequence vector of length K (the data length of the total number K uniform sampling points after uniform resampling in the same group), and y_v is a uniform velocity sequence vector. Construct a complete dictionary A with a size of K×L, where L is the number of candidate frequencies (e.g., the frequency range is from 0 to f_max, divided according to the frequency resolution Δf). A Fourier basis is usually used: A(k,l) = exp(j·2πf_l·τ_r), which is in complex form. In actual modeling, it is usually split into two columns, a cosine basis and a sine basis, to represent the real and imaginary parts respectively, where j is the imaginary unit, f_l is the l-th candidate frequency, and τ_r is the r-th uniform sampling time. The multi-measurement vector model is: y_d = A x_d + n_d, y_v = A x_v + n_v, where x_d is the L-dimensional complex coefficient vector of the displacement channel, x_v is the L-dimensional complex coefficient vector of the velocity channel, n_d is Gaussian white noise of the displacement channel, n_v is Gaussian white noise of the velocity channel, and the noise variances of the displacement channel and the velocity channel are σ_d and σ_v, respectively. 2 and σ_v 2 During initialization, σ_d is set. 2 = 1, σ_v 2 = 1, or initialized based on an estimate of the signal energy (e.g., 10% of the variance of the measurement vector).

[0044] The second step involves the expectation-maximization iterative process of block sparse Bayesian learning. In each iteration, the posterior distribution is first calculated using the current noise variance and hyperparameter estimates. The block sparse prior is set as follows: a hyperparameter γ_l is introduced for each frequency index l (corresponding to a block), which controls the variance of two coefficients (x_d(l) and x_v(l)) within the block. Specifically, the prior distribution is p(x_d(l), x_v(l)|γ_l) = CN(0,γ_l Σ_0), where Σ_0 is the intra-block correlation matrix, typically set to a 2×2 identity matrix (i.e., assuming the displacement and velocity coefficients are independent and identically distributed). All γ_l values ​​form the vector γ. Furthermore, an inverse gamma prior is introduced for the noise variance. In the expectation step, the posterior mean and covariance of the coefficients are calculated; in the maximization step, γ_l and the noise variance parameter are updated.

[0045] In the maximization step, besides updating the hyperparameter γ_l of each block, the noise variance of the displacement and velocity channels also needs to be updated separately. The noise variance is updated as follows: using the posterior mean of the coefficients calculated in the current iteration, the displacement and velocity signals are reconstructed respectively. That is, the reconstructed displacement is obtained by multiplying the complete dictionary A by the posterior mean of the displacement coefficients, and the reconstructed velocity is obtained similarly. Then, the reconstructed displacement is subtracted from the measured displacement vector to obtain the displacement residual vector; the reconstructed velocity is subtracted from the measured velocity vector to obtain the velocity residual vector. The sum of the squares of each element in the displacement residual vector divided by the total number of uniform sampling points K is the updated displacement noise variance. The velocity noise variance is obtained similarly by dividing the sum of the squares of the velocity residuals by K. This update step allows the noise variance estimate to gradually approach the actual measured noise level with each iteration, and the noise variances of the displacement and velocity channels can change independently without interference.

[0046] After each iteration updates the noise variance, adaptive weights for the displacement and velocity channels are calculated. The principle for weight calculation is that channels with smaller noise variances have more reliable measurement data and should have a larger weight in the joint solution. Specifically, the weight for the displacement channel can be the reciprocal of the displacement noise variance, and the weight for the velocity channel can be the reciprocal of the velocity noise variance. To prevent the values ​​from being too large or too small, the two weights can be normalized so that the sum of the normalized weights equals 1. After normalization, the channel with smaller noise variance receives a weight close to 1, and the channel with larger noise variance receives a weight close to 0. For example, if the noise variance of the displacement channel is much smaller than that of the velocity channel, then the displacement weight will be close to 1, and the velocity weight close to 0, and the algorithm will mainly rely on the displacement data to determine the frequency support set.

[0047] After obtaining the adaptive weights, when updating the hyperparameter γ_l of each block in the subsequent maximization step, the likelihood terms of the displacement and velocity channels need to be weighted. Specifically, when constructing the posterior covariance matrix of each block, the observed data of the displacement channel is multiplied by the square root of the displacement weight, and the observed data of the velocity channel is multiplied by the square root of the velocity weight, effectively scaling the data for each channel separately. Alternatively, when calculating the expected log-likelihood of the complete data for each block, the displacement and velocity weights are multiplied by their respective likelihood terms and then summed. The effect of this is that when the noise variance of a channel is large (i.e., the weight is small), the contribution of the likelihood term of that channel to the hyperparameter update is significantly weakened, while channels with low noise dominate the update process. This weighting mechanism allows the algorithm to automatically adapt to the imbalance of the signal-to-noise ratio between the two channels, avoiding contamination of the sparse reconstruction results by low-quality channels.

[0048] In each iteration, the following steps are performed sequentially: calculate the posterior distribution, update the noise variance, calculate the adaptive weights, and update the hyperparameters of each block. The iterative process continues until a convergence condition is met. The convergence condition can be set as follows: the maximum relative change of all block hyperparameters γ_l between two consecutive iterations is less than a preset threshold (e.g., one-thousandth); or the changes in displacement and velocity noise variances are both less than preset thresholds (e.g., one ten-thousandth); or the number of iterations reaches a preset upper limit (e.g., 500 times). Typically, the expectation-maximization algorithm converges within tens to two hundred iterations.

[0049] After the iteration converges, the final block hyperparameter vector γ is obtained. For each candidate frequency index l, if the corresponding γ_l is greater than a certain threshold (e.g., 1e-3, or significantly greater than the initial value), then the frequency is considered to be a real frequency component in the blade vibration. The frequency values ​​corresponding to all γ_l that satisfy the condition constitute a sparse frequency support set, i.e., the frequency parameters of the blade vibration. For each selected frequency, its vibration amplitude can be obtained by the magnitude of the posterior mean of the displacement coefficients within that block. If a more accurate amplitude estimation is required, the posterior mean of the velocity coefficients can also be combined: since there is a differential relationship between velocity and displacement, i.e., the velocity amplitude equals the displacement amplitude multiplied by the angular frequency, two amplitude estimates can be calculated separately from the displacement and velocity channels, and then weighted and averaged according to the weights of their respective channels to obtain the final amplitude output.

[0050] Through the above steps, the joint solution of the multi-measurement vector sparse reconstruction model is completed, outputting the frequency and amplitude of blade vibration. The advantages of this method are: no manual adjustment of the sparse regularization parameter is required, as it is automatically determined by the Bayesian framework; the adaptive weighting mechanism effectively addresses the actual working conditions where there are significant differences in noise levels between displacement and velocity measurements; and the block sparse structure ensures that the two channels share the same frequency support set, conforming to the physical nature of blade vibration, thereby improving the robustness and accuracy of frequency identification.

[0051] In a preferred embodiment, calculating the vibration displacement of the blade under the current rotational speed fluctuation condition based on the arrival time difference includes the following steps: The system acquires multiple measured pulse widths of the same blade within the same rotor revolution by using multiple blade tip timing sensors installed at different circumferential angles on the casing. Gap compensation is performed on each measured pulse width to obtain multiple compensated pulse widths; Obtain the pre-calibrated mapping relationship between pulse width and instantaneous velocity, and convert each compensated pulse width into the instantaneous linear velocity of the blade at the corresponding sensor. Based on the instantaneous linear velocity and current rotational speed function at each sensor, calculate the first tip radius at each sensor, where the first tip radius is equal to the instantaneous linear velocity divided by the tip angular velocity at the current rotational speed; The first blade tip radius obtained from multiple sensors is fused to obtain the dynamic blade tip radius at the current rotation speed; The blade tip linear velocity is calculated using the dynamic blade tip radius, and then the arrival time difference is converted into vibration displacement using the blade tip linear velocity.

[0052] Specifically, this embodiment further illustrates how to utilize multiple blade tip timing sensors with different circumferential angles to obtain the dynamic blade tip radius through pulse width information, thereby more accurately converting the arrival time difference into vibration displacement. This method is particularly suitable for scenarios where the blade tip radius changes significantly due to centrifugal force and thermal deformation under variable speed operation conditions.

[0053] First, multiple blade tip timing sensors need to be installed circumferentially on the casing of the rotating machinery, for example, three or four sensors, spaced at certain angles (e.g., uniformly distributed at 90 degrees or non-uniformly distributed). When the same blade sweeps past these sensors sequentially within the same revolution of the rotor, each sensor outputs a pulse waveform. After extracting the measured pulse width from each pulse, it needs to be processed according to the previously described clearance compensation steps: estimate the change in blade tip clearance based on the current rotational speed and thermal deformation model, and then subtract the corresponding compensation amount from the measured pulse width according to the pre-calibrated clearance-pulse width compensation mapping relationship to obtain the compensated pulse width. This step eliminates the interference of clearance changes on the pulse width.

[0054] Next, a pre-calibrated mapping relationship between pulse width and instantaneous velocity is needed. This mapping relationship is obtained in the laboratory or under calibration conditions by simultaneously measuring the actual linear velocity of the blade and the corresponding pulse width, and usually exhibits an inverse proportional characteristic. Substituting the compensated pulse width at each sensor location into this mapping relationship, the instantaneous linear velocity of the blade sweeping across the sensor location can be calculated. Note that this instantaneous linear velocity is the tangential velocity of the blade relative to the casing, which is equal to the superposition of the blade tip linear velocity and the blade tangential vibration velocity at the current rotational speed. However, during variable-speed operation, the rotational speed itself also changes, so the linear velocity cannot be obtained directly by multiplying the nominal radius by the angular velocity.

[0055] To determine the true blade tip radius, the angular velocity provided by the current rotational speed function must be used simultaneously. For each sensor, the instantaneous linear velocity obtained in the previous step is divided by the blade tip angular velocity at the current moment (i.e., the moment the pulse occurs). The result is called the first blade tip radius at that sensor. Here, the blade tip angular velocity is directly provided by the rotational speed function, and it is equal to the rotor's instantaneous angular velocity. The physical meaning of the first blade tip radius is: assuming the blade has no tangential vibration (i.e., the linear velocity is entirely contributed by rotation), then what should the blade tip radius be to generate the measured instantaneous linear velocity? In reality, blade tangential vibration causes the measured linear velocity to deviate from the pure rotational velocity; therefore, the first blade tip radius calculated by a single sensor includes the influence of vibration.

[0056] To obtain the dynamic tip radius at the current rotational speed (i.e., the actual elongation of the blade, excluding vibration components), it is necessary to fuse the first tip radii obtained from multiple sensors. The basis for this fusion is that the tangential vibration of the blade is a dynamic quantity that varies with time, while the elongation of the tip radius due to centrifugal force and thermal deformation is essentially constant within one rotational cycle. The average or median of the first tip radii measured by multiple sensors at different angles can effectively cancel out the vibration components, because the mean of the vibration within one cycle is zero (assuming the vibration is a periodic motion around the equilibrium position). Therefore, the arithmetic mean of the first tip radii obtained from all sensors for the same blade within the same rotation can be taken as the dynamic tip radius for that blade in that rotation. Alternatively, a more robust statistical method can be used, such as averaging after removing the maximum and minimum values.

[0057] Once the dynamic blade tip radius is obtained, it can be used to convert the arrival time difference into a more accurate vibration displacement. The arrival time difference is the difference between the actual arrival time of the blade leading edge caused by blade vibration and the theoretical arrival time. Physically, it is the distance the blade deviates from its theoretical position in the circumferential direction divided by the instantaneous linear velocity. However, the instantaneous linear velocity here should be calculated by multiplying the current dynamic blade tip radius by the angular velocity, not the nominal radius. Specifically, for each pulse, calculate the dynamic blade tip radius at that moment (which can be obtained by fusing the values ​​from the cycle or by interpolation), multiply it by the angular velocity at that moment to obtain the true instantaneous linear velocity, and then multiply this linear velocity by the arrival time difference to obtain the circumferential vibration displacement. If radial vibration displacement is required, it also needs to be converted according to the sensor installation angle and the blade vibration mode.

[0058] The technical principle of this solution lies in utilizing redundant information from multiple sensors to independently acquire the instantaneous linear velocity at each sensor via pulse width modulation. This instantaneous velocity is then combined with the rotational speed function to separate the true dynamic blade tip radius, thereby eliminating errors in the time-of-arrival conversion caused by changes in blade tip radius. The technical benefits include: accurate measurement of blade vibration displacement under conditions of variable speed operation and significant thermal deformation, avoiding spurious vibration signals introduced by using a fixed nominal radius. This solution is particularly suitable for applications such as aero-engines where blade tip radius variations can reach several millimeters.

[0059] In a preferred embodiment, the speed sensor is a sensor mounted on a toothed disc, and the toothed disc has multiple teeth evenly distributed around its circumference; Before dynamically constructing the speed function describing the rotor speed change over time based on the speed sensor signal, the method further includes the following steps: The arrival time of each pulse is extracted from the speed sensor signal, and the pulse interval value corresponding to each tooth is calculated based on the difference in arrival times of adjacent pulses. The calculated pulse interval values ​​are sequentially filled into the circular buffer in the circumferential order of the toothed disk. The length of the circular buffer is equal to the total number of teeth on the toothed disk. Each time a new pulse interval value is entered, the standard deviation of all pulse interval values ​​in the circular buffer is calculated, and the sum of the absolute values ​​of the differences between pulse interval values ​​at adjacent positions in the circular buffer is also calculated. The step of dynamically constructing a speed function describing the change of rotor speed over time based on the speed sensor signal includes the following steps: If the standard deviation is greater than or equal to the first threshold, or the sum of the absolute values ​​of the differences is less than or equal to the second threshold, then the original pulse interval value is directly used to construct the rotational speed function through interpolation.

[0060] Specifically, this embodiment addresses the case where the speed sensor is a toothed disc type (with multiple teeth evenly distributed circumferentially), and describes in detail how to dynamically construct a speed function describing the rotor speed change over time. The core innovation lies in utilizing a circular buffer to monitor the consistency of the tooth pulse interval in real time, and based on this, determining whether to directly use the original pulse interval to construct the speed function.

[0061] First, a geared disc is mounted on the rotor shaft of the rotating machinery. Multiple teeth, for example, 60 teeth, are evenly machined around the disc. A fixed speed sensor (such as a magnetoelectric or Hall effect sensor) detects the pulse generated as each tooth passes. A data acquisition system records the arrival time of each pulse. From these times, the time interval between two adjacent pulses can be calculated, called the pulse interval value. Because the geared disc is circumferentially uniform, if the rotor rotates at a constant speed and without eccentricity, all pulse interval values ​​should be equal. However, in variable speed operation, the pulse interval value will change with the rotational speed.

[0062] In this embodiment, these pulse interval values ​​are sequentially filled into a circular buffer in circumferential order of the toothed disc. The circular buffer is a first-in, first-out array, and its length is exactly equal to the total number of teeth on the toothed disc. Whenever a new tooth's pulse interval value is calculated, it is stored in the corresponding position in the buffer, overwriting the old value stored in the previous circumference. In this way, the buffer always stores the pulse interval values ​​of all teeth in the most recent circumference.

[0063] After each new pulse interval value is entered, two statistics need to be calculated simultaneously: the first is the standard deviation of all pulse interval values ​​within the buffer, which reflects the dispersion between each tooth interval. The second is the sum of the absolute values ​​of the differences between adjacent pulse interval values ​​within the buffer, that is, the absolute values ​​of the differences between the first and second, the second and third, and so on, until the absolute values ​​of the differences between the last and the first. This sum reflects the drasticness and periodicity of the pulse interval value along the circumferential direction.

[0064] Next, the calculation of the rotational speed function is determined based on the values ​​of these two statistics. If the calculated standard deviation is greater than or equal to a preset first threshold, or the sum of the absolute values ​​of adjacent differences is less than or equal to a preset second threshold, the original pulse interval value is directly used to construct the rotational speed function through interpolation. The first threshold is set based on the gearbox machining tolerance and the allowable rotational speed fluctuation range; for example, it can be taken as 5% of the average pulse interval value. The second threshold is used to determine whether the pulse interval value exhibits obvious periodic fluctuations. If the sum of the absolute values ​​of adjacent differences is very small, it indicates that the differences between the pulse interval values ​​are small and the changes are gradual, and it can be used directly.

[0065] The principle behind this method is as follows: when the standard deviation is large, it indicates significant differences between pulse interval values, which may be due to actual speed fluctuations rather than measurement errors. Therefore, these raw data should be used directly to reflect speed changes. When the sum of the absolute values ​​of adjacent differences is small, it means that the interval change between adjacent teeth is very smooth, and even if there are some fluctuations, they are in accordance with physical laws and can also be used directly. The specific method of constructing the speed function using the raw pulse interval values ​​is as follows: each pulse interval value is regarded as the reciprocal of the instantaneous angular velocity at the midpoint of that interval, resulting in a series of discrete velocity-time points. Then, linear interpolation or cubic spline interpolation is used to obtain a continuous time-angular velocity function. The technical effects include: the ability to quickly and adaptively determine the quality of the speed sensor signal, avoiding unnecessary corrections that introduce additional errors, while ensuring the real-time performance and accuracy of the speed function during variable speed operation, providing a reliable basis for subsequent calculations of arrival time difference and theoretical pulse width.

[0066] In a preferred embodiment, after calculating the standard deviation of all pulse interval values ​​within the circular buffer and simultaneously calculating the sum of the absolute values ​​of the differences between pulse interval values ​​at adjacent positions within the circular buffer, the method further includes the following step: If the standard deviation is less than the first threshold and the sum of the absolute values ​​of the differences is greater than the second threshold, it is determined that there is a periodic measurement error in the speed sensor signal caused by the eccentricity of the gear plate. The pulse interval value in the cyclic buffer is corrected position by position, and the speed function is constructed by interpolation using the corrected pulse interval value.

[0067] Specifically, if the standard deviation is less than the first threshold and the sum of the absolute values ​​of adjacent differences is greater than the second threshold, this embodiment determines that there is a periodic measurement error in the speed sensor signal caused by the eccentricity of the gear disk, performs position-by-position correction, and then uses the corrected pulse interval value to construct the speed function.

[0068] First, let's explain why this judgment result occurs. A standard deviation less than the first threshold means that the overall dispersion of the pulse interval values ​​of each tooth is not large, that is, their magnitudes are all roughly the same. However, the sum of the absolute values ​​of adjacent differences is greater than the second threshold, which means that although the overall dispersion is small, the interval values ​​between adjacent teeth exhibit a fluctuating, alternating pattern. The typical physical reason for this phenomenon is the eccentricity of the toothed disk: when the mounting axis of the toothed disk does not coincide with the rotation axis, the distance between the teeth on the disk and the sensor changes periodically during rotation, causing the pulse interval value to exhibit a sinusoidal change. Specifically, in the eccentric direction, the linear velocity of the teeth is not constant when passing the sensor, and the interval between adjacent teeth will show an alternating deviation. The period of this deviation is exactly one revolution, and the sum of adjacent differences will be relatively large, while the standard deviation may be insignificant due to the cancellation of positive and negative values.

[0069] Once an eccentricity error is determined, the pulse interval values ​​within the circulating buffer need to be corrected position by position. The correction is based on a pre-constructed eccentricity error table. The table is constructed as follows: During the startup phase of the rotating machinery or under calibration conditions after maintenance, the machine is run at a constant speed one-third lower than the first-order bending vibration natural frequency of the blade. The low speed is chosen to avoid the influence of blade vibration and dynamic stress; uniform speed operation ensures that if there is no eccentricity, all pulse interval values ​​should be identical. At this time, the original pulse interval value of each tooth is collected during one complete rotation of the gear disk, and then the arithmetic mean of the pulse interval values ​​of all teeth in this rotation is calculated. For each tooth, its original pulse interval value is subtracted from this average value; the difference is the eccentricity error compensation value for that tooth. This difference reflects the deviation of the tooth's pulse interval from the average value due to eccentricity. Storing the compensation values ​​of all teeth in circumferential order forms the eccentricity error table.

[0070] During actual variable speed operation, whenever a new pulse interval value is added to the circulating buffer, the eccentricity error compensation value corresponding to that tooth number is first read from the eccentricity error table based on the tooth number (each tooth on the gear sprocket has a unique number, which can be identified through the ring synchronization signal or counting method). Then, this compensation value is subtracted from the original pulse interval value to obtain the corrected pulse interval value. It is important to note the direction of subtraction: if eccentricity causes the pulse interval of that tooth to be too large, then the compensation value is positive, and the result after subtraction is the normal value; conversely, the result is negative. The corrected pulse interval value eliminates the periodic error caused by eccentricity and reflects the true speed change.

[0071] Then, using the corrected pulse interval value, the rotational speed function is constructed through interpolation. The specific interpolation method is the same as described above. In addition, since the eccentricity error may change during long-term operation (such as changes in eccentricity due to thermal deformation), the eccentricity error table can be recalibrated periodically under suitable operating conditions, for example, after each shutdown maintenance.

[0072] The technical effects of this implementation include: significantly improving the accuracy of the rotational speed function, thereby enhancing the accuracy of blade vibration parameter testing; avoiding misjudging eccentricity as rotational speed fluctuations and eliminating the false vibration frequencies caused by it; and being suitable for scenarios where the eccentricity of the gear disc may change during long-term operation, maintaining the correction effect through periodic calibration.

[0073] In a preferred embodiment, the step of position-by-position correction of the pulse interval value within the circular buffer includes the following steps: During variable speed operation, whenever a new pulse interval value is filled into the cyclic buffer, the pre-built eccentricity error table is retrieved according to the tooth number of the corresponding tooth, and the eccentricity error compensation value corresponding to the tooth number is read. The pulse interval value is then subtracted from the read eccentricity error compensation value to obtain the corrected pulse interval value. The eccentricity error table is constructed through the following steps: During the startup phase or calibration condition after maintenance of rotating machinery, the machine is run at a constant speed at one-third lower than the first-order bending vibration natural frequency of the blade. The original pulse interval value of each tooth is collected during one complete rotation of the toothed disc. The original pulse interval value of each tooth is subtracted from the arithmetic mean of the pulse interval values ​​of all teeth in that rotation to obtain the eccentricity error compensation value of each tooth. The eccentricity error compensation value is then stored in an eccentricity error table according to the circumferential order of the teeth.

[0074] Specifically, this embodiment further refines how to perform position-by-position correction on the pulse interval value within the cyclic buffer and describes in detail the construction process of the eccentricity error table. This correction method can effectively eliminate the periodic measurement error caused by gear disc eccentricity and improve the accuracy of the rotational speed function.

[0075] In the actual operation of rotating machinery, eccentricity of the gear disk installation is difficult to completely avoid. When the rotation center of the gear disk does not coincide with its geometric center, the linear velocity of the teeth on the disk as they pass the sensor is no longer constant: on one side of the eccentricity, the linear velocity of the teeth is slightly higher than the average velocity, resulting in a smaller pulse interval; on the opposite side, the linear velocity is slightly lower than the average velocity, resulting in a larger pulse interval. This error exhibits a strict periodicity, with a period of exactly one full revolution, and the error value of each tooth is relatively fixed. To compensate for this error, an eccentricity error table needs to be established in advance, and the table needs to be consulted and corrected in real time during variable speed operation.

[0076] The eccentricity error table is constructed under calibration conditions during the startup phase or after maintenance of rotating machinery. Calibration conditions require the rotating machinery to operate at a constant speed, which must be lower than one-third of the natural frequency of the first-order bending vibration of the blades. This speed condition is chosen because: when the speed is sufficiently low, the centrifugal force on the blades is small, and the amplitude of the blades' own vibrations (such as flutter and forced vibration) can be ignored; simultaneously, at low speeds, aerodynamic excitation is weak, and the blades are essentially in static equilibrium. At this point, the change in the toothed disc pulse interval is almost entirely caused by toothed disc eccentricity, while the influence of blade vibration, speed fluctuations, and other factors can be ignored. In practice, the rotor is allowed to rotate steadily at the selected constant speed for several revolutions. After the speed stabilizes, the original pulse interval value of each tooth within one complete rotation of the toothed disc is collected. Since it is running at a constant speed, all pulse interval values ​​should be identical if there is no eccentricity. However, in reality, due to the presence of eccentricity, the pulse interval values ​​of each tooth will exhibit a sinusoidal variation. Then, the arithmetic mean of the pulse interval values ​​of all teeth in this revolution is calculated. This average represents the pulse interval corresponding to this speed under ideal, eccentric-free conditions. For each tooth, its original pulse interval value is subtracted from the average value; the difference is the eccentricity error compensation value for that tooth. If the pulse interval of a tooth is larger than the average value, it indicates that the linear velocity of that tooth is slower when passing the sensor, and the compensation value is positive; otherwise, it is negative. The compensation values ​​of all teeth are stored sequentially in a table according to the circumferential order of the tooth disk (e.g., numbered clockwise starting from the reference tooth), thus forming the eccentricity error table. The length of this table is equal to the total number of teeth on the tooth disk, and each tooth number corresponds to a compensation value.

[0077] During variable speed operation, the system calculates the pulse interval value for each tooth in real time and fills it into the circular buffer in circumferential order. Each time a new pulse interval value is filled, the tooth number corresponding to the current pulse must first be identified. Tooth number identification is achieved through a circumferential synchronization signal: typically, there is a missing tooth or a specially marked tooth on the tooth sprocket; the system detects this marked tooth and begins counting, thus determining the tooth number corresponding to each pulse. Then, based on the tooth number, the eccentricity error compensation value corresponding to that tooth number is read from a pre-built eccentricity error table. Finally, the original pulse interval value is subtracted from the compensation value to obtain the corrected pulse interval value. This corrected value eliminates the system error caused by eccentricity and reflects the true speed change. It is particularly important to note that the direction of subtraction must be consistent with that during calibration: during calibration, the compensation value is obtained by subtracting the average value from the original value; therefore, during calibration, subtracting the compensation value from the original value restores the ideal value. The corrected pulse interval value is stored in the circular buffer and used for subsequent speed function construction.

[0078] Since the eccentricity of the gear disc may change during long-term operation (e.g., thermal expansion causing shaft bending, loosening of installation, etc.), it is recommended to recalibrate the eccentricity error table periodically. Recalibration can be performed automatically after each shutdown for maintenance or when the rotating machinery is in a stable low-speed warm-up state. Furthermore, if the rotating machinery is equipped with multiple speed sensors, they can be cross-calibrated to update the eccentricity error table online.

[0079] The technical principle of this solution lies in utilizing the periodicity and fixed nature of eccentricity error and the non-periodicity and randomness of rotational speed fluctuations. Through low-speed uniform calibration, the eccentricity error is isolated, forming a lookup table compensation. The correction process is simple and efficient, requiring no complex online identification algorithms. The technical benefits include: significantly improving the measurement accuracy of the rotational speed function, thereby enhancing the accuracy of calculating the blade arrival time difference and theoretical pulse width; eliminating spurious vibration frequencies caused by eccentricity, and avoiding misjudging eccentricity as blade vibration or rotational speed fluctuations.

[0080] In a preferred embodiment, after dynamically constructing a speed function describing the rotor speed change over time based on the speed sensor signal, and before calculating the theoretical instantaneous speed of the blade as it passes the sensor based on the speed function, the following steps are further included: Obtain the vibration displacement data of all blades within the same rotor rotation circle, and calculate the statistical average value of the vibration displacement of all blades; If the statistical average value is not within the static balance threshold range, it is determined that the current rotational speed function has a systematic low-frequency drift error. Based on the sign and magnitude of the statistical average, the rotational speed function is subjected to overall translation correction so that the statistical average of all blade vibration displacements after correction is within the range of the static balance threshold.

[0081] Specifically, this embodiment adds a step of performing an overall translational correction on the rotational speed function. This correction step is performed after the rotational speed function is constructed and before the theoretical instantaneous velocity is calculated, in order to eliminate the overall shift in blade vibration displacement caused by the systematic low-frequency drift error present in the rotational speed sensor signal.

[0082] In actual testing, the construction of the rotational speed function depends on the arrival time of the rotational speed sensor pulse. However, due to factors such as gear disc machining errors, sensor installation angle deviations, signal transmission delays, and clock drift in the data acquisition system, the rotational speed function may contain a slowly changing low-frequency drift component. This drift causes a systematic deviation in the calculated vibration displacement values ​​of all blades, manifested as the statistical average of the vibration displacement of all blades deviating from zero over a long time scale. Physically, when the blades are running stably in rotating machinery, if there is no permanent deformation or overall displacement, the vibration displacement of all blades should fluctuate around zero, and their statistical average should theoretically be zero (or within a very small static equilibrium threshold range). If the calculated average deviates significantly from zero, it indicates that there is a systematic error in the rotational speed function, which needs to be corrected.

[0083] The specific implementation steps are as follows. First, acquire the vibration displacement data of all blades within the same rotor rotation cycle. "Same rotation cycle" refers to one full rotation of the rotor, during which each blade passes the sensor exactly once (if there is only one sensor, multiple rotations of data are needed; if multiple sensors are installed circumferentially, data from multiple blades can be obtained within the same rotation). Collect the vibration displacement values ​​of all blades and calculate their statistical average. The statistical average can be the arithmetic mean or a robust statistical measure such as the median to avoid the influence of outliers from individual blades.

[0084] Next, it is determined whether the statistical average value falls within the pre-set static balance threshold range. The static balance threshold range is set based on the allowable error of blade vibration measurement and common sense physics. For example, it can be set from -1 / 1000th of a millimeter to +1 / 1000th of a millimeter, or based on 1 / 1000th of the blade tip radius. If the statistical average value falls within this range, it indicates that the accuracy of the rotational speed function meets the requirements and no correction is needed. If the statistical average value is not within this range, for example, if it is a positive value and large, it is determined that the current rotational speed function has a systematic low-frequency drift error. This drift may originate from the overall high or low rotational speed function, causing the calculated theoretical arrival time to be systematically earlier or later than the measured arrival time, thus causing the vibration displacement of all blades to shift in the same direction.

[0085] Next, based on the sign and magnitude of the statistical average, the rotational speed function is subjected to an overall translational correction. The specific correction method is to multiply the entire rotational speed function ω(t) by a correction coefficient, or more directly, to add a linear correction term to the rotor angle function θ(t) corresponding to the rotational speed function. Since the drift mainly manifests as a low-frequency component and can be approximated as a constant offset within a relatively short time window, the simplest approach is to calculate the correction amount and then modify the formula for calculating the theoretical arrival time. Let the statistical average of all blade vibration displacements be d_mean, and the nominal radius of the blade tip be R0. Then the corresponding arrival time difference offset is Δt_bias = d_mean / (ω_avg·R0), where ω_avg is the average angular velocity of the current revolution. Using this Δt_bias as the overall correction amount, when calculating the theoretical arrival time of each pulse, this correction amount is uniformly added or subtracted to bring the statistical average of all blade vibration displacements to zero after correction. Another equivalent approach is to directly shift the speed function: multiply the speed function by a coefficient (1-d_mean / (R0·θ_total)), where θ_total is the total rotation angle in one revolution. Both methods essentially involve a comprehensive stretching or compression of the speed function.

[0086] After correction, recalculate the vibration displacement of all blades and verify again whether the statistical average has fallen within the static balance threshold range. If it still does not meet the requirement, the above correction steps can be repeated once or twice; usually, one correction is sufficient. It is important to note that this translation correction is only applicable to systematic, consistent low-frequency drift errors affecting all blades and should not affect the random vibration components of individual blades. Therefore, the relative vibration differences of each blade remain unchanged after correction; only the overall baseline is pulled back to zero.

[0087] The technical effects of this embodiment include: effectively eliminating DC bias and extremely low frequency drift in the rotation speed signal, improving the absolute accuracy of vibration displacement measurement; and avoiding false vibration amplitude caused by rotation speed function errors.

[0088] In a preferred embodiment, after outputting the frequency and amplitude of the blade vibration using the adjusted parameters, the method further includes the following steps: At the same rotational speed, the steps are repeated for all blades on the same stage rotor to obtain the vibration frequency and vibration amplitude of each blade; Calculate the average frequency and standard deviation of all blade vibration frequencies, and the average amplitude and standard deviation of all blade vibration amplitudes; For each blade, the absolute value of the frequency deviation between its vibration frequency and the average frequency is calculated. When the absolute value of the frequency deviation is greater than a preset first multiple multiplied by the frequency standard deviation, the blade is determined to be a frequency abnormal blade. For each blade, calculate the absolute value of the amplitude deviation between its vibration amplitude and the average amplitude. When the absolute value of the amplitude deviation is greater than a preset second multiple multiplied by the amplitude standard deviation, the blade is determined to be an amplitude abnormal blade. Output the number and type of all abnormal blades.

[0089] Specifically, this implementation method is used to perform batch analysis of all blades on the same stage of a rotor at the same rotational speed. It identifies blades with abnormal frequencies and amplitudes using statistical methods and outputs the blade numbers and abnormality types, providing direct evidence for condition monitoring and fault diagnosis of rotating machinery. In rotating machinery (such as aero engines, gas turbines, steam turbines, compressors, etc.), all blades on the same stage of a rotor theoretically have the same design parameters. Therefore, under the same rotational speed conditions, the vibration frequency and amplitude of each blade should be basically consistent. If the vibration frequency of a certain blade deviates too much from the group average, it may indicate that the blade has cracks, loosening, mass imbalance, or damping changes. If the vibration amplitude of a certain blade is significantly greater than that of other blades, it may indicate that the blade is subjected to additional aerodynamic excitation or has a resonance risk. Through statistical comparison in this implementation method, abnormal blades can be quickly located.

[0090] The specific implementation steps are as follows. First, under the same stable speed condition (e.g., cruise speed or rated speed), for each blade on the same stage rotor, all steps are repeated, including signal acquisition, time difference of arrival extraction, pulse width compensation, vibration displacement and vibration velocity calculation, resampling, multi-measurement vector sparse reconstruction, and adaptive block sparse Bayesian learning solution, ultimately obtaining the vibration frequency and vibration amplitude of each blade. For N blades, a total of N frequency values ​​f_i and N amplitude values ​​A_i are obtained, where i=1,2,...,N.

[0091] Next, calculate the average frequency (f_mean) and standard deviation (f_std) of all blade vibration frequencies. The average frequency is the arithmetic mean of the frequencies of each blade, and the standard deviation reflects the dispersion of frequencies within the group. Similarly, calculate the average amplitude (A_mean) and standard deviation (A_std) of all blade vibration amplitudes. It should be noted that vibration amplitude usually refers to the amplitude of the dominant frequency component (such as the first-order bending vibration of the blade). Multiple frequency components can be analyzed as needed, but the dominant frequency is generally used as the criterion.

[0092] Next, the abnormal blades are identified. For each blade i, the absolute value of the frequency deviation between its vibration frequency f_i and the average frequency f_mean is calculated, i.e., |f_i - f_mean|. If this absolute value of deviation is greater than a preset first multiple multiplied by the frequency standard deviation f_std, the blade is determined to be a frequency abnormal blade. The value of the first multiple is set based on engineering experience, usually between 3 and 5. Taking 3 times the standard deviation corresponds to the Raida criterion (3σ criterion) in statistics. Under the assumption of normal distribution, normal blades have a 99.7% probability of falling within 3 times the standard deviation; anything outside this range is considered abnormal. Taking 4 or 5 times the standard deviation is more stringent and suitable for situations with higher requirements for false alarm rate. Similarly, for each blade, the absolute value of the amplitude deviation between its vibration amplitude A_i and the average amplitude A_mean is calculated, i.e., |A_i - A_mean|. If this absolute value of deviation is greater than a preset second multiple multiplied by the amplitude standard deviation A_std, the blade is determined to be an amplitude abnormal blade. The second multiple can be the same as the first multiple, or they can be set separately according to actual needs. The amplitude fluctuations are usually greater than the frequency, so the second multiple can be appropriately relaxed, for example, by 4 or 5 times.

[0093] It is important to note that a blade may be simultaneously identified as having both frequency and amplitude anomalies; in such cases, both anomaly types should be recorded. After anomaly identification is complete, the numbers of all anomalous blades (i.e., their circumferential position numbers or physical numbers on the rotor) and their corresponding anomaly types (frequency anomaly, amplitude anomaly, or both) should be compiled and output in list or report form. The output can be displayed on the monitoring interface, stored in a log file, or sent to a host computer or maintenance system via a communication interface.

[0094] To further improve the reliability of the judgment, the above process can be repeated under multiple different speed conditions to observe whether the appearance of abnormal blades is consistent. For example, if a certain blade shows frequency deviation at multiple speeds, it is more certain that there is a structural problem. In addition, trend analysis can be performed by combining historical data, such as tracking the frequency and amplitude changes of each blade over time, and issuing an early warning when a gradually deteriorating trend appears.

[0095] The technical advantages of this implementation include: it can quickly and automatically screen out potentially faulty blades from a large number of blades, greatly reducing the workload of manual analysis; the output results are intuitive, making it easy for on-site maintenance personnel to directly locate the problematic blades, thereby improving the operational safety and maintenance efficiency of rotating machinery.

[0096] The above description is merely a preferred embodiment of this application and an explanation of the technical principles employed. Those skilled in the art should understand that the scope of the invention involved in this application is not limited to technical solutions formed by specific combinations of the above-described technical features, but should also cover other technical solutions formed by arbitrary combinations of the above-described technical features or their equivalents without departing from the inventive concept. For example, technical solutions formed by substituting the above features with (but not limited to) technical features with similar functions disclosed in this application.

Claims

1. A method for testing blade vibration parameters based on blade tip timing signals, characterized in that, Includes the following steps: During the variable speed operation of rotating machinery, a blade tip timing sensor installed on the casing is used to collect the pulse waveform signal when the blade sweeps across, and the measured arrival time and measured pulse width of each pulse are extracted from the pulse waveform signal. The speed sensor signal is acquired in real time, and a speed function describing the change of rotor speed with time is dynamically constructed based on the speed sensor signal; the theoretical arrival time of each pulse is calculated based on the speed function; the measured arrival time is compared with the theoretical arrival time to obtain the arrival time difference caused by blade vibration; the vibration displacement of the blade under the current speed fluctuation condition is calculated based on the arrival time difference. Obtain the pre-calibrated mapping relationship between pulse width and instantaneous speed; calculate the theoretical instantaneous speed of the blade when it passes the sensor based on the rotational speed function, and obtain the corresponding theoretical pulse width using the mapping relationship; The measured pulse width is compared with the theoretical pulse width to obtain the pulse width deviation caused by blade vibration. The vibration velocity of the blade under the current rotational speed fluctuation condition is calculated based on the pulse width deviation. The vibration displacement and the vibration velocity are integrated into a multi-measurement vector sparse reconstruction model, in which the vibration displacement and vibration velocity share the same sparse frequency support set in the frequency domain. The block sparse Bayesian learning method is used to jointly solve the multi-measurement vector sparse reconstruction model to recover the sparse frequency support set and obtain the frequency parameters of the blade vibration. Based on the current speed fluctuation amplitude and signal-to-noise ratio level, the regularization parameters in the block sparse Bayesian learning method are adaptively adjusted, and the frequency and amplitude of blade vibration are output using the adjusted parameters.

2. The method for testing blade vibration parameters based on tip timing signals according to claim 1, characterized in that, Before comparing the measured pulse width with the theoretical pulse width, a gap compensation step is also included: Obtain the pre-calibrated mapping relationship between the tip clearance and the pulse width compensation amount; Based on the current rotational speed and the pre-established rotor thermal deformation model, estimate the change in blade tip clearance under the current operating conditions; Using the mapping relationship, calculate the pulse width compensation amount corresponding to the change in blade tip clearance; Subtract the pulse width compensation amount from the measured pulse width to obtain the compensated measured pulse width.

3. The method for testing blade vibration parameters based on tip timing signals according to claim 1, characterized in that, Before integrating the vibration displacement and the vibration velocity into a multi-measurement vector sparse reconstruction model, the following steps are also included: Vibration displacement data and vibration velocity data are used as non-uniform sampling sequences to be processed; Based on the circumferential installation angle of each blade tip timing sensor and the rotational speed function, the actual sampling time corresponding to each pulse is calculated to form a non-uniform sampling time grid. A uniform sampling time grid is constructed, and the sampling interval of the uniform sampling time grid is adaptively determined according to the current speed fluctuation amplitude; The vibration displacement data and vibration velocity data on the non-uniform sampling sequence are mapped onto the uniform sampling time grid using the distance-weighted interpolation method to obtain the resampled uniform displacement sequence and uniform velocity sequence. The process of integrating the vibration displacement and the vibration velocity into a multi-measurement vector sparse reconstruction model includes the following steps: The resampled uniform displacement sequence and uniform velocity sequence are used as inputs to the multi-measurement vector sparse reconstruction model.

4. The method for testing blade vibration parameters based on tip timing signals according to claim 1, characterized in that, The method of jointly solving the multi-measurement vector sparse reconstruction model using block sparse Bayesian learning includes the following steps: Initialize the noise variance parameters for the displacement and velocity channels; During the expectation maximization iteration of block sparse Bayesian learning, the noise variance estimates of each channel are updated using displacement residuals and velocity residuals respectively. Based on the currently updated noise variance estimate, adaptive weights for the displacement and velocity channels are calculated. If the noise variance of the displacement channel is less than that of the velocity channel, the displacement channel receives a weight greater than that of the velocity channel; otherwise, the velocity channel receives a weight greater than that of the displacement channel. When updating the sparse frequency support set in each iteration, the likelihood terms of the displacement measurement vector and the velocity measurement vector are weighted using the adaptive weights so that the channels with small noise variance dominate in the joint solution. Iterate until convergence, and output the final sparse frequency support set.

5. The method for testing blade vibration parameters based on tip timing signals according to claim 2, characterized in that, The calculation of the blade's vibration displacement under the current rotational speed fluctuation condition based on the arrival time difference includes the following steps: The system acquires multiple measured pulse widths of the same blade within the same rotor revolution by using multiple blade tip timing sensors installed at different circumferential angles on the casing. Gap compensation is performed on each measured pulse width to obtain multiple compensated pulse widths; Obtain the pre-calibrated mapping relationship between pulse width and instantaneous velocity, and convert each compensated pulse width into the instantaneous linear velocity of the blade at the corresponding sensor. Based on the instantaneous linear velocity and current rotational speed function at each sensor, calculate the first tip radius at each sensor, where the first tip radius is equal to the instantaneous linear velocity divided by the tip angular velocity at the current rotational speed; The first blade tip radius obtained from multiple sensors is fused to obtain the dynamic blade tip radius at the current rotation speed; The blade tip linear velocity is calculated using the dynamic blade tip radius, and then the arrival time difference is converted into vibration displacement using the blade tip linear velocity.

6. The method for testing blade vibration parameters based on tip timing signals according to claim 1, characterized in that, The speed sensor is a sensor mounted on a toothed disc, which has multiple teeth evenly distributed around its circumference. Before dynamically constructing the speed function describing the rotor speed change over time based on the speed sensor signal, the method further includes the following steps: The arrival time of each pulse is extracted from the speed sensor signal, and the pulse interval value corresponding to each tooth is calculated based on the difference in arrival times of adjacent pulses. The calculated pulse interval values ​​are sequentially filled into the circular buffer in the circumferential order of the toothed disk. The length of the circular buffer is equal to the total number of teeth on the toothed disk. Each time a new pulse interval value is entered, the standard deviation of all pulse interval values ​​in the circular buffer is calculated, and the sum of the absolute values ​​of the differences between pulse interval values ​​at adjacent positions in the circular buffer is also calculated. The step of dynamically constructing a speed function describing the change of rotor speed over time based on the speed sensor signal includes the following steps: If the standard deviation is greater than or equal to the first threshold, or the sum of the absolute values ​​of the differences is less than or equal to the second threshold, then the original pulse interval value is directly used to construct the rotational speed function through interpolation.

7. The method for testing blade vibration parameters based on tip timing signals according to claim 6, characterized in that, After calculating the standard deviation of all pulse interval values ​​within the circular buffer, and simultaneously calculating the sum of the absolute values ​​of the differences between pulse interval values ​​at adjacent positions within the circular buffer, the method further includes the following steps: If the standard deviation is less than the first threshold and the sum of the absolute values ​​of the differences is greater than the second threshold, it is determined that there is a periodic measurement error in the speed sensor signal caused by the eccentricity of the gear plate. The pulse interval value in the cyclic buffer is corrected position by position, and the speed function is constructed by interpolation using the corrected pulse interval value.

8. The method for testing blade vibration parameters based on tip timing signals according to claim 7, characterized in that, The step-by-step correction of the pulse interval value within the circular buffer includes the following steps: During variable speed operation, whenever a new pulse interval value is filled into the cyclic buffer, the pre-built eccentricity error table is retrieved according to the tooth number of the corresponding tooth, and the eccentricity error compensation value corresponding to the tooth number is read. The pulse interval value is then subtracted from the read eccentricity error compensation value to obtain the corrected pulse interval value. The eccentricity error table is constructed through the following steps: During the startup phase or calibration condition after maintenance of rotating machinery, the machine is run at a constant speed at one-third lower than the first-order bending vibration natural frequency of the blade. The original pulse interval value of each tooth is collected during one complete rotation of the toothed disc. The original pulse interval value of each tooth is subtracted from the arithmetic mean of the pulse interval values ​​of all teeth in that rotation to obtain the eccentricity error compensation value of each tooth. The eccentricity error compensation value is then stored in an eccentricity error table according to the circumferential order of the teeth.

9. The method for testing blade vibration parameters based on tip timing signals according to claim 6, characterized in that, After dynamically constructing a speed function describing the rotor speed change over time based on the speed sensor signal, and before calculating the theoretical instantaneous speed of the blade when it passes the sensor based on the speed function, the following steps are also included: Obtain the vibration displacement data of all blades within the same rotor rotation circle, and calculate the statistical average value of the vibration displacement of all blades; If the statistical average value is not within the static balance threshold range, it is determined that the current rotational speed function has a systematic low-frequency drift error. Based on the sign and magnitude of the statistical average, the rotational speed function is subjected to overall translation correction so that the statistical average of all blade vibration displacements after correction is within the range of the static balance threshold.

10. The method for testing blade vibration parameters based on tip timing signals according to claim 1, characterized in that, After outputting the frequency and amplitude of blade vibration using the adjusted parameters, the method further includes the following steps: At the same rotational speed, the steps are repeated for all blades on the same stage rotor to obtain the vibration frequency and vibration amplitude of each blade; Calculate the average frequency and standard deviation of all blade vibration frequencies, and the average amplitude and standard deviation of all blade vibration amplitudes; For each blade, the absolute value of the frequency deviation between its vibration frequency and the average frequency is calculated. When the absolute value of the frequency deviation is greater than a preset first multiple multiplied by the frequency standard deviation, the blade is determined to be a frequency abnormal blade. For each blade, calculate the absolute value of the amplitude deviation between its vibration amplitude and the average amplitude. When the absolute value of the amplitude deviation is greater than a preset second multiple multiplied by the amplitude standard deviation, the blade is determined to be an amplitude abnormal blade. Output the number and type of all abnormal blades.

Images