Rotating blade disc dynamic response field reconstruction method based on blade tip timing test and Kalman filtering

Abstract

The invention belongs to the field of impeller mechanical vibration analysis and response reconstruction, and provides a rotating blade disc dynamic response field reconstruction method based on blade tip timing testing and Kalman filtering. Comprising the following steps: firstly, establishing a blade disc reduced order kinetic model, and identifying and updating detuning parameters on line by using blade tip timing data; thirdly, constructing an augmented state space model, reconstructing a global dynamic response field under unknown excitation by adopting augmented Kalman filtering, and estimating a noise covariance online in combination with a moving window method so as to improve the anti-noise capability; experiments verify that the actual strain error of the method is less than 10%, which indicates that the method can effectively suppress noise, is suitable for complex working conditions, and is suitable for blade disc vibration monitoring and safety evaluation. According to the method, the defect that only local vibration characteristics are reflected in a test can be effectively overcome, reconstruction of a dynamic displacement field and a dynamic strain field is achieved, the vibration characteristics of the blisk can be effectively evaluated from the angle of the blisk, and effective support is provided for fatigue reliability and fault early warning of the blisk.

Description

Technical Field

[0001] This invention belongs to the field of turbomachinery vibration analysis and response reconstruction, specifically a method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering. Background Technology

[0002] In the field of structural dynamics and condition monitoring, accurately acquiring the full-field dynamic response of complex structures is a key challenge. Response reconstruction technology, by fusing finite measurement point data with a high-fidelity finite element model, constructs a mapping from local to global data, achieving a high-fidelity virtual reproduction of the structure's dynamic behavior. With advancements in sensing and computing technologies, this method shows great promise for applications in health monitoring, vibration control, and fatigue assessment in aerospace and other fields.

[0003] Modal extension, a classic reconstruction method, relies on modal coordinate identification and mode shape extension. Its effectiveness is limited by model accuracy and measurement signal-to-noise ratio. However, in the monitoring of high-speed rotating machinery such as aero-engine bladed disks, this method faces a dual challenge: First, the "detuning" of blades caused by manufacturing and wear alters the structural dynamics, leading to vibration localization and frequency shifts. This results in inherent errors in finite element analysis based on ideal models, which, when directly used for reconstruction, will produce distortion. Second, under actual operating conditions, timing measurements at the blade tip are susceptible to speed fluctuations and signal noise interference, affecting the accuracy of modal coordinate estimation. The errors, amplified by the mode shape matrix, further reduce the accuracy of reconstruction and may even mask detuning characteristics, thus severely impacting the reliability of fault diagnosis and life prediction.

[0004] Therefore, traditional modal extension methods struggle to meet engineering accuracy requirements under the combined effects of model detuning and measurement noise. Developing a response reconstruction method capable of simultaneously correcting model errors, suppressing noise interference, and adapting to complex operating environments has become an urgent requirement for improving the condition monitoring capabilities of rotating machinery. Summary of the Invention

[0005] To address the problem that existing reconstruction methods lack sufficient accuracy in reconstructing the dynamic response field under complex operating conditions such as variable speed, detuning, and strong noise, this invention proposes a method for reconstructing the dynamic response field of a rotating bladed disk based on inter-blade timing testing and Kalman filtering.

[0006] This invention includes the following technical solutions:

[0007] A method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering includes the following steps:

[0008] Step 1: Variable Rotational Mistuning Order Reduction Modeling Method

[0009] To address the vibration analysis requirements under the combined effects of rotational and detuning effects, a reduced-order dynamic model of a rotating detuned integral bladed disk applicable to arbitrary speeds is established by combining parametric methods and substructure modal detuning methods. The specific process is as follows: First, based on the stiffness matrices of the harmonic bladed disk at three characteristic speeds, a quadratic polynomial is fitted to express the change in stiffness matrix with speed, enabling rapid calculation of the stiffness matrix K0(Ω) at arbitrary speeds. Second, the system degrees of freedom are reduced to modal coordinates through modal transformation, and Λ0(Ω) at arbitrary speeds is obtained by interpolation using the eigenvalue matrices Λ0 of the harmonic bladed disk at the three characteristic speeds to characterize the rotational effect. Finally, neglecting the Coriolis matrix, a steady-state vibration equation in the reduced-order modal space is established, incorporating rotational stiffness changes, detuning effects, and Rayleigh damping, forming a reduced-order dynamic model of the rotating detuned integral bladed disk applicable to arbitrary speeds.

[0010] Furthermore, in step 1, the dynamic equation of the integral bladed disk rotating at any speed can be expressed as:

[0011]

[0012] In the formula, D = C + G, where C and G are the damping matrix and Coriolis matrix of the bladed disk, respectively, and x(t) is the vibration displacement response of the bladed disk. , Let x(t) be the second and first derivatives of x(t), respectively, F(t) be the excitation force acting on the bladed disk, and M and K be the mass matrix and stiffness matrix of the bladed disk, respectively.

[0013] In the above equation, the mass matrix M, damping matrix C, and stiffness matrix K of the rotating bladed disk can be expressed as:

[0014]

[0015] In the formula, M0, C0, and K0 are the mass matrix, damping matrix, and stiffness matrix of the tuned bladed disk, respectively, and ΔM, ΔC, and ΔK represent the changes in the mass matrix, damping matrix, and stiffness matrix caused by detuning, respectively.

[0016] For an integral bladed disk rotating at any speed, the stiffness matrix K0 is expressed as:

[0017]

[0018] In the formula, , and The coefficient matrix is ​​a polynomial, where Ω1, Ω2, and Ω3 represent three different rotational speeds of the bladed disk.

[0019] Parameterize the rotational stiffness matrix to obtain

[0020]

[0021] In the formula, K is the stiffness matrix of the bladed disk. This refers to the maximum rotational speed of the impeller; for speeds between 0 and Ω... max Between the rotating bladed disks, select three different rotational speeds: Ω0, Ω1, and Ω2, and let... .

[0022] The parameterization method described above only requires calculating the polynomial coefficient matrix using the stiffness matrix of the harmonic bladed disk at three different rotational speeds. , and Thus, the stiffness matrix K0 of the harmonic bladed disk at any rotational speed is obtained.

[0023] For small-scale detuning, the component mode mistuning (CMM) method is used to introduce detuning. The mode shapes of the detuned bladed disk are approximated by a linear superposition of the frequency-dense mode shapes of the harmonic bladed disk. Therefore, the first r-order mode shapes of the harmonic bladed disk are selected. Perform the following transformation

[0024]

[0025] In the formula, x is the physical coordinate and u is the modal coordinate.

[0026] Substitute the above equation into the dynamic equation and multiply both sides of the equation by the transpose of Φ0. The basic equations of the reduced-order model are obtained as follows:

[0027]

[0028] In the formula, Λ0(Ω) is the eigenvalue matrix of the harmonic bladed disk, and I is the identity matrix. , , and This represents the modal mass matrix, modal damping matrix, modal Coriolis matrix, and modal stiffness matrix in modal coordinates. The excitation force acting on the bladed disk in modal coordinates. , They are respectively The first and second derivatives.

[0029] Assuming the mode shape is independent of rotational speed, Λ0(Ω) at any rotational speed is obtained by interpolation of the eigenvalue matrix Λ0 of the harmonic bladed disk at three rotational speed points:

[0030]

[0031] In the formula, Λ0(0) and Λ0(Ω) max / 2), Λ0(Ω) max) represent rotational speeds of 0 and Ω respectively. max / 2 and Ω max The eigenvalue matrix of the harmonic bladed disk at that time.

[0032] Assuming the Coriolis matrix G(Ω) caused by rotation is negligible, we can further obtain the final simplified reduced-order model of the rotating detuned global bladed disk dynamics applicable to any rotational speed:

[0033]

[0034] In the formula, i is the imaginary unit. ω represents the excitation force frequency, which is also the vibration frequency of the bladed disk; I is the identity matrix; α and β are the Rayleigh damping coefficients. and Let be the stiffness mistuning parameter and the mass mistuning parameter of the nth blade, respectively. and These are the modal participation factors corresponding to mass mistuning and stiffness mistuning, respectively. and This represents its transpose. and These are the modal stiffness contribution matrix and modal mass contribution matrix of the nth blade in modal coordinates, respectively.

[0035] Step 2: Timed leaf tip testing and test parameter identification

[0036] Tip timing tests were conducted: First, multiple tip timing probes (BTT probes) were arranged circumferentially around the aero-engine casing, and a key phase signal sensor was installed on the shaft as a time reference. After the tip timing system was running, the tip vibration displacement was calculated by collecting the time difference sequence between the actual arrival time and the theoretical arrival time of each blade at each probe. Subsequently, the multi-probe data was processed using the circumferential Fourier fitting method to identify key vibration parameters such as vibration amplitude, resonant frequency, initial phase, and constant deviation of each blade. Combined with N sets of frequency response data obtained experimentally (N being the number of blades), accurate input was provided for subsequent mistuning identification and model updates.

[0037] Step 3: Mistuning Identification and Model Update

[0038] First, based on the reduced-order dynamics model of the rotating detuned integral bladed disk applicable to any rotational speed established in Step 1, by selecting two sets of excitation frequencies near the resonance frequency and constructing difference equations to eliminate the influence of unknown excitation forces, a linear identification equation for the mass and stiffness detuning parameters, i.e., the detuning identification equation, is established. Second, the detuning identification equation is rewritten as follows: By establishing the linear form of the detuning parameters and the corresponding deviations, and by substituting the N sets of frequency response data obtained from the experiment in step 2, the mass detuning parameters and stiffness detuning parameters of all blades can be identified. This enables the model update of the reduced-order dynamics model of the rotating detuned integral bladed disk applicable to any rotational speed, as well as the model update of the dynamics equations of the integral bladed disk rotating at any rotational speed.

[0039] Furthermore, in step 3: the reduced-order model of the rotating detuned integral bladed disk dynamics proposed in step 1, applicable to arbitrary rotational speeds, is rewritten as follows:

[0040]

[0041] In the formula, the unknown quantity and These are the detuning parameters that need to be identified.

[0042] Two sets of excitation frequencies, ω, are selected near the resonance frequency. j and ω k , and Substitute the corresponding modal coordinates into the above equation, and subtract the two equations to obtain the harmonic identification equation.

[0043]

[0044] To simplify the expression, let

[0045]

[0046]

[0047] The detuning identification equation is rewritten in the following form.

[0048]

[0049] In the formula

[0050]

[0051]

[0052]

[0053] In the formula, L is the system matrix and T is the response vector; the parameters in the above formula are the overall bladed disk mistuning parameters to be determined, which are calculated by the following formula, i.e.

[0054]

[0055] In the formula, the superscript "†" indicates the Moore-Penrose generalized inverse of the matrix, and the number of unknowns is N (number of blades). Therefore, only the N sets of frequency response data measured in step 2 (N is the number of blades) are needed to solve for the detuning parameters of the overall bladed disk, thereby realizing the updated dynamics model of the rotating detuned overall bladed disk applicable to any speed and the updated dynamics equations of the overall bladed disk rotating at any speed.

[0056] Step 4: Establishing the augmented state-space model

[0057] First, the updated dynamic equations of the integral bladed disk rotating at arbitrary speeds are reduced to decrease the system's degrees of freedom through modal transformation, and the system is simplified using modal orthogonality. The second-order equations are further expressed in state-space form by defining state vectors and then discretized. Second, to handle the problem of unknown loads in practice, the system state vector and the unknown load vector are combined into an augmented state vector, based on which an augmented state-space model is established. Finally, based on the augmented Kalman filter algorithm, the time prediction step and measurement update step are calculated, providing the core state input for subsequent dynamic response field reconstruction.

[0058] Furthermore, in step 4: the updated dynamic equations of the integral bladed disk rotating at arbitrary speeds obtained in step 3 are written as follows:

[0059]

[0060] In the formula, M, D, K∈ These represent the system mass, damping, and stiffness matrices, respectively, n s Let x(t) be the number of degrees of freedom of the system. Let be the system vibration displacement vector. , Let F(t) be the first and second derivatives of x(t), respectively. Let q(t) be the load vector corresponding to all degrees of freedom of the system. Let n be the load vector corresponding to the degrees of freedom of the excitation position. p S represents the number of degrees of freedom of the system under excitation. p ∈ The excitation location selection matrix consists of 0s and 1s.

[0061] Using mode transformation to reduce the number of system degrees of freedom, i.e. Select the top n dominant system responses m First mode Taking advantage of the orthogonality of the system mode shapes with respect to the mass matrix M, the updated dynamic equations of the integral bladed disk rotating at arbitrary speeds are expressed as follows:

[0062]

[0063] In the formula, p(t)∈ Let be the vibration displacement vector of the system in modal coordinates. Here is the modal damping matrix. Here is the modal stiffness matrix. This is the modal force vector.

[0064] Define state vector The above equation can then be further expressed in state-space form. After discretization, a discrete model for filtering is obtained, i.e.

[0065]

[0066] In the formula, and These represent the measured quantity of the blade disk vibration displacement and the system state quantity of the rotating blade disk tip timing measurement system, respectively. , w k Let be the model (process) noise vector, and its covariance matrix be... v k To measure the noise vector, its covariance matrix is: The model noise and measurement noise are uncorrelated Gaussian white noise processes, i.e. ; For continuous-time system matrices; t is the time step; For continuous-time input matrices; To measure the noise covariance matrix, Let be the process noise covariance matrix.

[0067] The system's state vector z in modal coordinates k and load vector f k Combined into an augmented state vector

[0068] ,Right now

[0069]

[0070] The unknown load (input) acting on the system is represented in incremental form, i.e.

[0071]

[0072] In the formula, Let be a zero-mean random process, and its covariance matrix be... , For the random fluctuation of the l-th unknown parameter or detuned parameter.

[0073] The augmented state-space model is then...

[0074]

[0075] In the formula,

[0076]

[0077] Augmented noise vector The covariance matrix is

[0078]

[0079] In the formula, Let be the covariance parameter corresponding to the l-th noise component.

[0080] According to the augmented Kalman filter method, state estimation can be completed in two steps.

[0081] The first step is the time prediction step, namely:

[0082]

[0083]

[0084] In the formula, This is the predicted value of the state at time k based on the information at time k-1; This is the prior estimate of the error covariance matrix.

[0085] Then comes the measurement update step, i.e.

[0086]

[0087]

[0088]

[0089] In the formula, The Kalman gain matrix; and Let be the state prediction error covariance matrix and the state update error covariance matrix.

[0090] Step 5: Dynamic Response Field Reconstruction

[0091] Based on the augmented state-space model and augmented Kalman filter method established in step 4, the augmented state vector in modal coordinates is first transformed to physical coordinates, and the displacement response state vector is extracted to reconstruct the dynamic displacement field of the bladed disk. Then, based on the finite element strain-displacement relationship, the element strain is calculated by shape function and Jacobian matrix and integrated into the overall strain field. Finally, the moving window method is used to estimate the measurement noise characteristics online, and the fixed noise covariance is replaced in the Kalman filter update step to achieve adaptive noise suppression, thereby improving the response reconstruction accuracy.

[0092] Furthermore, in step 5: based on the solution result in step 4, the augmented state vector in modal coordinates is... Transforming to physical coordinates yields the displacement response state vector x. k Take only x k This yields the dynamic displacement field of the rotating bladed disk.

[0093] According to finite element theory, the element strain ε e With nodal displacement δ e The relationship is

[0094]

[0095] In the formula, B is the strain-displacement matrix.

[0096] Based on the shape function and Jacobian matrix, B can be obtained and the strain of each element can be calculated. Finally, the strain field of the entire structure is obtained by integration.

[0097] The measurement noise characteristics are estimated online using the moving window method. The measurement noise covariance matrix at step k is... It can be calculated by the following formula

[0098]

[0099] In the formula, j represents the time step index within the moving window; To measure the noise substitution value, For a smoothed estimate of the true value, y k Let λ be the measurement value at step k, λ be the weighting factor, and Nr be the moving window length. and These are the measured noise surrogate values. The moving mean and moving covariance matrix.

[0100] In the Kalman filter update step, using It replaces the fixed measurement noise covariance R to achieve adaptive noise suppression and improve reconstruction accuracy.

[0101] Compared with the prior art, the present invention has the following advantages:

[0102] 1. Efficient Modeling for Variable Speeds: A matrix parameterization method is proposed, which can quickly fit the stiffness matrix at any speed with only three characteristic speed data, enabling efficient analysis over a wide speed range.

[0103] 2. Online mistuning identification: Combining tip timing measurement and reduced-order model, the blade mistuning parameters can be identified in real time while rotating, without the need for shutdown or additional excitation devices.

[0104] 3. Noise-robust reconstruction: The augmented Kalman filter combined with moving window noise estimation significantly improves the accuracy and stability of dynamic response field reconstruction in environments with strong measurement noise.

[0105] 4. Localized vibration characterization: The updated model can accurately predict the vibration localization phenomenon caused by detuning, providing a more realistic dynamic basis for the fatigue life assessment of bladed disks.

[0106] 5. Accuracy-efficiency balance: While maintaining high accuracy (error <10%), computational efficiency is greatly improved through order reduction and parameterization strategies, making it suitable for real-time monitoring and offline analysis in engineering. Attached Figure Description

[0107] To more clearly illustrate the specific embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the specific embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of the present invention. For those skilled in the art, other drawings can be obtained from these drawings without creative effort.

[0108] Figure 1 : Flowchart of the present invention.

[0109] Figure 2 The vibration table foundation excitation test includes: ① strain acquisition system, ② laser displacement sensor, ③ sensor support, ④ electromagnetic vibration table, ⑤ strain gauge, and ⑥ bladed disk.

[0110] Figure 3 : Time-domain response of blade tip vibration displacement under sweep frequency excitation.

[0111] Figure 4 Frequency response of blade tip vibration displacement for different blades.

[0112] Figure 5 Results of basic excitation mismatch identification.

[0113] Figure 6 Rotary integral bladed disk test bench, including ① servo motor, ② gas nozzle, ③ integral bladed disk, ④ air compressor, ⑤ BTT probe, and ⑥ data acquisition system.

[0114] Figure 7Vibration displacement of the blade tip of blade #1 as measured by different BTT probes.

[0115] Figure 8 Vibration displacement of the blade tip of blade #4 measured by different BTT probes.

[0116] Figure 9 : The frequency response of blade vibration displacement measured at the blade tip, where (a) before Savitzky-Golay filtering, (b) after Savitzky-Golay filtering.

[0117] Figure 10 : Detuning parameter identification results based on tip timing.

[0118] Figure 11 Reconstructed results of dynamic displacement and dynamic strain fields when the vibration of blade #1 reaches its peak (2630 rpm), where (a) dynamic displacement field and (b) dynamic strain field.

[0119] Figure 12 Reconstructed results of dynamic displacement and dynamic strain fields when the vibration of blade #4 reaches its peak (2810 rpm), where (a) dynamic displacement field and (b) dynamic strain field. Detailed Implementation

[0120] like Figure 1 As shown, this invention provides a method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering, comprising the following steps:

[0121] Step 1: Variable Rotational Mistuning Order Reduction Modeling Method

[0122] First, based on the stiffness matrix of the harmonic bladed disk at three characteristic rotational speeds, a quadratic polynomial is fitted to describe the change of the stiffness matrix with rotational speed, enabling rapid calculation of the stiffness matrix K0(Ω) at any rotational speed. Second, the system degrees of freedom are reduced to modal coordinates through modal transformation, and Λ0(Ω) at any rotational speed is obtained by interpolation using the eigenvalue matrix Λ0 of the harmonic bladed disk at the three characteristic rotational speeds to characterize the rotational effect. Finally, neglecting the Coriolis matrix, a steady-state vibration equation in the reduced-order modal space is established, incorporating rotational stiffness variation, detuning effect, and Rayleigh damping, forming a reduced-order dynamic model of the rotating detuned overall bladed disk applicable to any rotational speed.

[0123] Step 2: Timed leaf tip testing and test parameter identification

[0124] To address the issue of overall bladed disk mistuning identification, static shaking table experiments and tip-timing (BTT) rotating table experiments were conducted. Figure 2The vibration table excitation system shown applies a sweep frequency excitation to the integral bladed disk to identify its detuning parameters. This vibration table has a maximum sinusoidal excitation force of 6 kN, a maximum excitation frequency of 3500 Hz, and a maximum acceleration of 100 g, which meets the experimental requirements for excitation intensity and frequency range.

[0125] Static shaking table experiment: A 9-bladed integral bladed disk (e.g.) was selected. Figure 2 As shown in the upper left corner (⑥), a sweep frequency excitation (sweep speed 0.5 Hz / s, acceleration 1 g) is applied within the range of [400Hz, 500Hz] by rotating the integral bladed disk test bench. A high-precision laser displacement sensor (such as...) is used. Figure 2 (As shown in the upper left corner ②) The blade tip displacement of each blade was measured to obtain the time-domain diagram of the vibration displacement response of each blade to the sweep frequency excitation (as shown in the upper left corner ②). Figure 3 (As shown).

[0126] The inputs to the detuning identification equation include the natural frequency of the tuned bladed disk, the Rayleigh damping coefficient of the actual bladed disk, and the displacement response of the bladed disk at different frequencies. The natural frequency of the tuned bladed disk can be obtained through finite element simulation. The Rayleigh damping coefficient of the actual bladed disk can be obtained using the half-power bandwidth method. This method also obtains the frequency response of the bladed disk through frequency sweep excitation, and then obtains the resonant frequency and two half-power frequency points from the frequency response to calculate the damping of the structure. Therefore, the displacement frequency response of different blade tips was calculated, such as... Figure 4 As shown. The damping coefficient of the actual bladed disk, calculated using the half-power bandwidth method, is: , .

[0127] Since the bladed disk used in the experiment contains nine blades, the tip displacement response at nine frequency points is required according to the mistuning identification equation. Figure 4 Nine frequency points were selected in the resonance region of the disk, and their displacement responses were used as input information to identify the detuning of the actual bladed disk, resulting in the following: Figure 5 The results of the detuning parameter identification are shown. From... Figure 5 It can be seen that each blade exhibits varying degrees of mistuning.

[0128] Rotating leaf tip timed experiment Figure 6 The rotating test bench shown uses the same test blade as the static vibration table experiment. The test bench is driven by a servo motor with a rotational speed range of 1200-3000 rpm, and uses a 10-nozzle pneumatic system for excitation. Four BTT probes (angles of 0°, 36°, 70°, and 108°) and a key phase sensor are mounted circumferentially on the casing to collect blade tip vibration signals.

[0129] Figure 7 and Figure 8The blade tip vibration displacement signals of blades #1 and #4, collected by the BTT probes, show different vibration constant deflection values ​​due to the different installation angles of each BTT probe. When the rotor speed approaches the resonant speed, a significant increase in blade amplitude is observed. Since the casing has 10 gas nozzles, each blade experiences 10 periodic airflow excitations per revolution of the bladed disk. Figure 7 and Figure 8 The dominant excitation order is observed to be EO = 10. Multiplying the excitation order by the corresponding rotational frequency allows estimation of the resonant frequency range of the blade. For blades #1 and #4, the resonant frequency range is calculated using... Figure 7 and Figure 8 It can be preliminarily determined that their resonant frequencies are located between 420 Hz ~ 450 Hz and 450 Hz ~ 480 Hz, respectively.

[0130] The vibration amplitude of each blade at different rotational speeds was extracted from the BTT signal using the circumferential Fourier fitting method, and the data was smoothed and denoised using Savitzky-Golay filtering. Figure 9 (a) and Figure 9 (b) in the diagram represents the blade tip amplitude signals before and after filtering, respectively. From... Figure 9 Coupled vibrations between adjacent blades can be observed, which are shown in the figure as multiple resonance peaks appearing on a single blade.

[0131] Using the filtered signal as input, the detuning parameters of the rotating bladed disk are identified based on the reduced-order model and the tip displacement response. The identification results are as follows: Figure 10 As shown in the figure, due to factors such as manufacturing errors and operational wear, each blade exhibits varying degrees of mistuning.

[0132] Step 3: Mistuning Identification and Model Update

[0133] The key vibration parameters obtained from the two types of experiments in step 2 are combined with the reduced-order dynamics model of the rotating detuned integral bladed disk applicable to any rotational speed established in step 1 to construct the detuning identification equation. This detuning equation is then rewritten as follows: The linear form of the model is used to establish the relationship between the detuning parameters and the corresponding deviations. The N sets of frequency response data obtained from the experiment in step 2 are then substituted into the model to finally identify the detuning parameters of each blade. This enables the model update of the reduced-order dynamics model of the rotating detuned integral bladed disk applicable to any rotational speed, as well as the model update of the dynamics equations of the integral bladed disk rotating at any rotational speed.

[0134] Step 4: Establishing the augmented state-space model

[0135] First, the updated dynamic equations of the integral bladed disk rotating at arbitrary speeds are reduced to decrease the system's degrees of freedom through modal transformation, and the system is simplified using modal orthogonality. The second-order equations are further expressed in state-space form by defining state vectors and then discretized. Second, to handle the problem of unknown loads in practice, the system state vector and the unknown load vector are combined into an augmented state vector, based on which an augmented state-space model is established. Finally, based on the augmented Kalman filter algorithm, the time prediction step and measurement update step are calculated, providing the core state input for subsequent dynamic response field reconstruction.

[0136] Step 5: Dynamic Response Field Reconstruction Test

[0137] This step uses Figure 6 As shown in the rotating test bench, in step 2, the resonant frequency and detuning parameters of the rotating bladed disk have been obtained based on the blade tip timing data of the rotating test bench, and this step directly uses the results. Figure 11 The figure shows the dynamic displacement field and dynamic strain field reconstructed based on the Kalman filter method when the vibration of blade #1 reaches its peak at a rotating disk speed of 2630 rpm. Figure 12 The figure shows the dynamic displacement field and dynamic strain field reconstructed by the Kalman filter method when the vibration of blade #4 reaches its peak at a speed of 2810 rpm. The reconstruction results reflect the localized vibration characteristics of the detuned bladed disk.

[0138] Step 6: Experimental verification based on leaf tip timing data and Kalman filtering method

[0139] Tables 1 and 2 show the dynamic strain response amplitudes at the measuring points of blades #1, #2, #3, and #4, reconstructed using the Kalman filter method at speeds of 2630 rpm and 2810 rpm, respectively. These values ​​were compared with the measured values ​​from strain gauges, and the relative reconstruction error was calculated. It can be seen that the dynamic strain errors obtained using the Kalman filter method proposed in this invention are all below 10%. This is mainly because the Kalman filter method can consider the influence of process noise and measurement noise. In test environments with severe noise interference, adjusting the noise variance in the Kalman filter algorithm can effectively reduce the noise impact and improve the accuracy of the dynamic response field reconstruction.

[0140] Table 1. Comparison of dynamic strain response amplitudes of each blade at 2630 rpm.

[0141]

[0142] Table 2. Comparison of dynamic strain response amplitudes of each blade at 22810 rpm

[0143]

[0144] Verified by tip-timed tests, the dynamic response field reconstruction method based on tip-timed testing and Kalman filtering proposed in this invention exhibits significant noise resistance and engineering adaptability. In actual tests, the reconstructed dynamic response fields under different speed conditions accurately capture the localized vibration characteristics of the detuned bladed disk, with the error between the dynamic strain response amplitude and the measured value being less than 10%. The results show that this method can effectively suppress measurement noise interference, adapt to unknown excitation and complex time-varying conditions, and possesses good robustness and reconstruction accuracy, making it suitable for dynamic response monitoring and safety assessment of bladed disks in noisy environments in practical engineering.

Claims

1. A method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering, characterized in that, Includes the following steps: Step 1: Variable Rotational Mistuning Order Reduction Modeling Method First, based on the stiffness matrix of the harmonic bladed disk at three characteristic speeds, a quadratic polynomial is fitted to the stiffness matrix as a function of speed, thus realizing the stiffness matrix at any speed. K The system achieves rapid calculation of Ω; secondly, it reduces the system's degrees of freedom to modal coordinates through modal transformation and utilizes the eigenvalue matrix of the harmonic bladed disk at three characteristic rotational speeds. Λ 0 interpolation yields values ​​at any speed Λ 0(Ω) is used to characterize the rotational effect. Finally, after ignoring the Coriolis matrix, a steady-state vibration equation including the rotational stiffness change, detuning effect and Rayleigh damping is established in the reduced-order modal space, forming a reduced-order model of the rotational detuning overall bladed disk dynamics applicable to any rotational speed. Step 2: Timed leaf tip testing and test parameter identification Blade tip timing tests were conducted: First, multiple blade tip timing probes (BTT probes) were arranged circumferentially around the aero-engine casing, and a key phase signal sensor was installed on the shaft as a time reference. After the blade tip timing system was running, the blade tip vibration displacement was calculated by collecting the time difference sequence between the actual time and the theoretical arrival time of each blade at each probe. Subsequently, the multi-probe data was processed using the circumferential Fourier fitting method to identify key vibration parameters such as vibration amplitude, resonant frequency, initial phase, and constant deviation of each blade. These parameters were then combined with experimentally measured values. N The set of frequency response data provides accurate input for subsequent mistuning identification and model updates; Step 3: Mistuning Identification and Model Update First, based on the reduced-order dynamic model of the rotating detuned integral bladed disk applicable to any rotational speed established in step 1, by selecting two sets of excitation frequencies near the resonance frequency and constructing difference equations, the influence of unknown excitation force is eliminated, and a linear identification equation for the mass and stiffness detuning parameters is established, namely the detuning identification equation. Secondly, the mistuning identification equation is rewritten as Using the linear form of the detuning parameter, establish the relationship between the detuning parameter and the corresponding deviation, and use the experimentally measured data from step 2... N By inputting the set of frequency response data, the mass mistuning parameters and stiffness mistuning parameters of all blades can be identified, thereby enabling the model update of the reduced-order dynamics model of the rotating mistuned integral bladed disk applicable to any rotational speed, as well as the model update of the dynamics equations of the integral bladed disk rotating at any rotational speed. Step 4: Establishing the augmented state-space model First, the updated dynamic equations of the integral bladed disk rotating at arbitrary speeds are reduced by modal transformation to decrease the system's degrees of freedom, and the system is simplified using modal orthogonality. Second, by defining state vectors, the second-order equations are further expressed in state-space form and discretized. Third, to handle the problem of unknown loads in practice, the system state vector and the unknown load vector are combined into an augmented state vector, based on which an augmented state-space model is established. Finally, based on the augmented Kalman filter algorithm, the time prediction step and measurement update step are calculated, providing the core state input for subsequent dynamic response field reconstruction. Step 5: Dynamic Response Field Reconstruction Based on the augmented state-space model and augmented Kalman filter method established in step 4, the augmented state vector in modal coordinates is first transformed to physical coordinates, and the displacement response state vector is extracted to reconstruct the dynamic displacement field of the bladed disk. Then, based on the finite element strain-displacement relationship, the element strain is calculated by shape function and Jacobian matrix and integrated into the overall strain field. Finally, the moving window method is used to estimate the measurement noise characteristics online, and the fixed noise covariance is replaced in the Kalman filter update step to achieve adaptive noise suppression, thereby improving the response reconstruction accuracy.

2. The method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering according to claim 1, characterized in that, In step 1, the dynamic equation of the integral bladed disk rotating at any speed can be expressed as: ； In the formula, D = C + G, where C and G are the damping matrix and Coriolis matrix of the bladed disk, respectively, and x(t) is the vibration displacement response of the bladed disk. , Let x(t) be the second and first derivatives of x(t), respectively, and F(t) be the excitation force acting on the bladed disk. M and K These are the mass matrix and stiffness matrix of the impeller, respectively; The mass matrix of the rotating bladed disk in the above formula M Damping matrix C and stiffness matrix K It can be represented as ； In the formula, M 0、 C 0、 K 0 represents the mass matrix, damping matrix, and stiffness matrix of the harmonic bladed disk, respectively, and Δ represents the stiffness matrix. M Δ C and Δ K These represent the changes in the mass matrix, damping matrix, and stiffness matrix caused by detuning, respectively. For an integral bladed disk rotating at any speed, the stiffness matrix K 0 Represented as ； In the formula, , and The coefficient matrix is ​​a polynomial, where Ω1, Ω2, and Ω3 represent three different rotational speeds of the bladed disk. Parameterize the rotational stiffness matrix to obtain ； In the formula, K Let be the stiffness matrix of the bladed disk. This refers to the maximum rotational speed of the impeller; for speeds between 0 and Ω... max Between the rotating bladed disks, select three different rotational speeds: Ω0, Ω1, and Ω2, and let... ; The parameterization method described above only requires calculating the polynomial coefficient matrix using the stiffness matrix of the harmonic bladed disk at three different rotational speeds. , and This allows us to obtain the stiffness matrix of the harmonic bladed disk at any rotational speed. K 0; For small-scale detuning, the component mode mistuning (CMM) method is used to introduce detuning; the mode shapes of the detuned bladed disk are approximated by the linear superposition of the mode shapes of the frequency-dense harmonic bladed disk; therefore, the harmonic bladed disk front is selected. r First-order mode shape Perform the following transformation ； In the formula, x For physical coordinates, u Modal coordinates; Substitute the above equation into the dynamic equation and multiply both sides of the equation by the left side. Φ transpose of 0 The basic equations of the reduced-order model are obtained as follows: ； In the formula, Λ 0(Ω) is the eigenvalue matrix of the harmonic bladed disk. I It is the identity matrix. , , and This represents the modal mass matrix, modal damping matrix, modal Coriolis matrix, and modal stiffness matrix in modal coordinates. The excitation force acting on the bladed disk in modal coordinates. , They are respectively The first and second derivatives; Assuming the mode shape is independent of rotational speed, the eigenvalue matrix of the harmonic bladed disk at three rotational speed points is obtained. Λ 0 interpolation yields values ​​at any speed Λ 0(Ω): ； In the formula, Λ 0(0), Λ 0(Ω max / 2) Λ 0(Ω max ) represent rotational speeds of 0 and Ω respectively. max / 2 and Ω max The eigenvalue matrix of the harmonic bladed disk at that time; Let the Coriolis matrix be induced by rotation. G Neglecting (Ω), we further obtain the final simplified reduced-order model of the rotating detuned integral bladed disk dynamics applicable to any rotational speed: ； In the formula, i is the imaginary unit. ω represents the excitation force frequency, which is also the vibration frequency of the bladed disk. I It is the identity matrix; α And β is the Rayleigh damping coefficient, and Let be the stiffness mistuning parameter and the mass mistuning parameter of the nth blade, respectively. and These are the modal participation factors corresponding to mass mistuning and stiffness mistuning, respectively. and This represents its transpose. and These are the modal stiffness contribution matrix and modal mass contribution matrix of the nth blade in modal coordinates, respectively.

3. The method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering according to claim 2, characterized in that, In step 3: the reduced-order dynamics model of the rotating detuned integral bladed disk applicable to any rotational speed proposed in step 1 is rewritten as follows: ； In the formula, the unknown quantity and These are the detuning parameters that need to be identified; Two sets of excitation frequencies were selected near the resonance frequency. ω j and ω k , and Substitute the corresponding modal coordinates into the above equation, and subtract the two equations to obtain the harmonic identification equation. ； To simplify the expression, let ； ； The detuning identification equation is rewritten in the following form. In the formula ； ； ； In the formula, L is the system matrix and T is the response vector; the parameters in the above formula are the overall bladed disk mistuning parameters to be determined, which are calculated by the following formula, i.e. ； In the formula, the superscript "†" denotes the Moore-Penrose generalized inverse of the matrix, and the number of unknowns is . N (Number of blades) Therefore, only the N sets of frequency response data measured in step 2 are needed to solve for the detuning parameters of the overall bladed disk, thereby realizing the updated dynamics model of the rotating detuned overall bladed disk applicable to any speed and the updated dynamics equations of the overall bladed disk rotating at any speed.

4. The method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering according to claim 3, characterized in that, In step 4: the updated dynamic equations of the integral bladed disk rotating at arbitrary speeds obtained in step 3 are written as follows: ； In the formula, M , D , K ∈ These represent the system mass, damping, and stiffness matrices, respectively. n s Let be the system's degrees of freedom. x ( t ) ∈ Let be the system vibration displacement vector. , They are respectively x ( t The first and second derivatives of ) F ( t ) ∈ The load vector corresponds to all degrees of freedom of the system. q ( t ) ∈ Let the load vector correspond to the degrees of freedom of the excitation position. n p The number of degrees of freedom of the system under stimulus. S p ∈ The excitation location selection matrix consists of 0s and 1s. Using mode transformation to reduce the number of system degrees of freedom, i.e. Select the dominant system response n m First mode Using the system modal vibration modes with respect to the mass matrix M Due to the orthogonality of the equations, the updated dynamic equations of the integral bladed disk rotating at arbitrary speeds are expressed as follows: ； In the formula, p ( t )∈ Let be the vibration displacement vector of the system in modal coordinates. Here is the modal damping matrix. Here is the modal stiffness matrix. This is the modal force vector; Define state vector The above equation can then be further expressed in state-space form. After discretization, a discrete model for filtering is obtained, i.e. ； In the formula, and These represent the measured quantity of the blade disk vibration displacement and the system state quantity of the rotating blade disk tip timing measurement system, respectively. , , w k Let be the model noise vector, and its covariance matrix be... , v k To measure the noise vector, its covariance matrix is: ; The model noise and measurement noise are uncorrelated Gaussian white noise processes, i.e. ; For continuous-time system matrices; t For time step; Input matrix in continuous time; To measure the noise covariance matrix, The process noise covariance matrix; The state vector of the system in modal coordinates z k and load vector f k Combined into an augmented state vector ,Right now ； The unknown load on the system is represented in incremental form, i.e. ； In the formula, Let be a zero-mean random process, and its covariance matrix be... , For the first l Random fluctuations of an unknown or detuned parameter; The augmented state-space model is then... ； In the formula, ； Augmented noise vector The covariance matrix is ； In the formula, For the first l Covariance parameters corresponding to each noise component; According to the augmented Kalman filter method, state estimation can be completed in two steps; The first step is the time prediction step, namely: ； ； In the formula, This is the predicted value of the state at time k based on the information at time k-1; The prior estimate is the error covariance matrix; Then comes the measurement update step, i.e. ； ； ； In the formula, The Kalman gain matrix; and Let be the state prediction error covariance matrix and the state update error covariance matrix.

5. The method for reconstructing the dynamic response field of a rotating bladed disk based on tip timing testing and Kalman filtering according to claim 1, characterized in that, In step 5: Based on the solution results in step 4, the augmented state vector in modal coordinates is... Transforming to physical coordinates yields the displacement response state vector. x k Only take x k That is, the dynamic displacement field of the rotating bladed disk is obtained; According to finite element theory, element strain ε e With nodal displacement δ e The relationship is ； In the formula, B This is the strain-displacement matrix; Based on the shape function and the Jacobian matrix, it can be obtained that... B The strain of each element is calculated, and the strain field of the entire structure is finally obtained by integration. The moving window method is used to estimate the measurement noise characteristics online. k Step measurement noise covariance matrix It can be calculated by the following formula ； In the formula, j represents the time step index within the moving window; To measure the noise substitution value, For a smoothed estimate of the true value, y k Let λ be the measurement value at step k, λ be the weighting factor, and Nr be the moving window length. and These are the measured noise surrogate values. The moving mean and moving covariance matrix; In the Kalman filter update step, using Alternative fixed measurement noise covariance R This enables adaptive noise suppression and improves reconstruction accuracy.

Images