# A subspace method for modal identification of bladed assemblies using blade tip-timing data

Bendali Salhi, Marc Berthillier, Joseph Lardies, Philippe Voinis, Charles Bodel

► To cite this version:

Bendali Salhi, Marc Berthillier, Joseph Lardies, Philippe Voinis, Charles Bodel. A subspace method for modal identification of bladed assemblies using blade tip-timing data. ASME Turbo Conference 2007 on Controls, Diagnostics and Instrumentation, May 2007, Montréal, Canada. pp.839-848, 10.1115/GT2007-28151. hal-00330035

HAL Id: hal-00330035

https://hal.science/hal-00330035v1

Submitted on 3 May 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# A SUBSPACE METHOD FOR MODAL IDENTIFICATION OF BLADED ASSEMBLIES USING BLADE TIP-TIMING DATA

**Bendali SALHI**

PhD Student
bendali.salhi@univ-fcomte.fr

**Marc BERTHILLIER**

Professor
marc.berthillier@univ-fcomte.fr

**Joseph LARDIES**

Professor
joseph.lardies@univ-fcomte.fr

Department of Applied Mechanics Raymond Chaléat
FEMTO-ST Institute
24, chemin de l'Epitaphe, 25000 Besançon, France

**Philippe VOINIS**

Chief Engineer
philippe.voinis@edf.fr

**Charles BODEL**

Projects Engineer
charles.bodel@edf.fr

Department of research and development (AMA)
Electricité de France (EDF firme)
1, avenue du Général de Gaulle, 92141 Clamart, France

## ABSTRACT

Unknown excitation forces are applied to bladed disk assemblies, such as turbines blades, leading to forced vibration responses. Non contact measurement of such vibrations using blade tip-timing data has become an industrial standard procedure and current research focuses on analysis methods for interpretation of measured vibrations. Our purpose is to develop a method for identification of the blade's natural frequencies and damping ratios using blade tip-timing data. The method is based on a subspace analysis. Its performance is compared to the traditional Fast Fourier Transform analysis. A detailed description of these methods and results are presented.

**Keywords:** blade vibration analysis, blade tip-timing method, subspace method.

## INTRODUCTION

Bladed assembly vibration measurements are crucial to the prediction of blade life. Generally, blade vibrations in compressor or turbine stages during operation are difficult to

measure and require many sensors in order to be fully detected [3]. Two standard measuring methods are used to measure blade vibrations: strain gages and frequency-modulated Grid Systems. The first method consists of installing strain gages on some of the blades. This is a very expensive process since it requires slip rings or high-quality telemetry and long installation times. For the second standard method, the frequency modulated Grid Systems (FM grid), permanent magnets are fitted on the tips of some blades and a specially formed wire is installed in the compressor casing above the magnets. As a blade passes a wire, the magnet produces an alternating current (AC) in the wire. Blade vibrations are interpreted from the voltage modulation of this AC signal. No signal transmitter between the rotor and the stator is needed, but the wire requires a special casing and the number of blades monitored at a time is limited to three per rotor stage. An alternative method to the two standard measurement methods is the blade tip-timing data (BTT) or non intrusive stress measurement which is based on optical, capacitive or magnetic probe technology.
This paper describes the use of subspace modeling for identification of blade vibration frequencies and their damping ratios using blade tip-timing data. This method is applied to both: a lumped parameter model which produces asynchronous responses and to real measurements of a turbine in operating conditions.

## BLADE TIP-TIMING METHOD

The basic concepts of this technique have already been known for many years [1-3]. Several probes are installed in the engine casing above the rotor, and the blade transit times between the probes are measured. Figure 1 summarizes the basic principle of this technique. Two sensors are placed around the blades. They measure the passage time of each blade in front of the sensors. The angle between these two sensors is of 180°. A third sensor is placed in the proximity of the engine shaft. It measures the once-per-revolution signal.

<figure><img src="image_2.png" /><figcaption>Fig.1. A SIMPLE BLADE TIP-TIMING SYSTEM</figcaption></figure>

When there are no blade vibrations, the blade transit times are a function of rotor speed, rotor radius, and circumferential probe position. In the case of blade vibrations, the blade transit times deviate from those in the undisturbed condition, with the blades passing the probes earlier or later than normal, depending on their instantaneous deflection. Compared with the two standard methods cited previously, this technique provides notable advantages: 1) it is non-contacting; 2) it senses all blades; and 3) it reduces costs because it eliminates the need for rotor instrumentation and telemetry for signal transmission. This technique however, has some drawbacks. One such drawback is that the measurement sampling frequency is completely dependent of the rotor speed and the number of measurement probes installed.

The signals presented in figure 2 are the once-per-revolution signal and the passage time of seventy-seven blades in front of the two sensors. The difference between the actual time of arrival TOA (the blade arrival time at the front of the sensor) of the vibrating blade and its computed TOA (when it

has not been in vibration) provides the raw data from which the instantaneous blade displacement is found. These displacements are found in the following way: in each revolution the values of the once-per-revolution data and the average value of the rotation speed are removed from raw data to find blade vibration times. To obtain blade displacements, the vibration time of each blade is multiplied by the value of rotation speed.

<figure><img src="image_3.png" /><figcaption>Fig.2. RAW BLADE TIP-TIMING DATA</figcaption></figure>

Methods for analysing blade tip-timing data have been developed since 1970. However, given a limited number of probes, there are still no standard approaches that can identify synchronous response resonance frequencies with adequate precisions. Let us recall that synchronous response occurs when a blade's vibration frequency is an integer multiple of the engine rotation speed and asynchronous response occurs when the blades frequency response is a non-integer multiple of the assembly rotation speed. Synchronous vibration can be caused by mechanical effects such as residual unbalance of the rotors and non-concentric casing, as well as aerodynamic effects such as irregular pressure distributions within the airflow due to the engine intake geometry and wakes produced by upstream stators. Asynchronous vibrations are typically caused by flutter, rotating stall, or by acoustic resonance.

The analysis of blade tip-timing data can be performed using on-line data analysis which provides near real-time information updated at least once per second or by off-line data analysis. In this article we consider off-line analysis to obtain more detailed information on individual blades.

## THE SUBSPACE METHOD

The subspace method is based on the use of a state space model which is obtained by transforming equation of motion for the bladed assemblies: $M\ddot{q} + C\dot{q} + Kq = f(t)$ into a
state equation and an observation equation. The discrete time state space model is given by [4-6]:

$$z_{k+1} = Az_k + Be_k \quad \text{state equation} \quad (1)$$

$$y_k = Cz_k + e_k \quad \text{observation equation} \quad (2)$$

where $y_k$ is the (m×1) vector of observation, $e_k$ the vector of stochastic innovations, $z_k$ the (n×1) vector of states, A the (n×n) transition matrix, B the (n×m) Kalman gain matrix, and C the (m×n) observation matrix. The matrices A, B, C and the order n of the multivariate stochastic process are parameters to be estimated, from output data only.

In order to estimate these parameters, the shifted observability matrix method presented here is used. Let $y_k^+$ be the (mf×1) future data vectors and $\bar{y}_k$ the (mp×1) past data vectors such that:

$$y_k^+ = (y'_k, y'_{k+1}, \dots, y'_{k+f-l})' \quad (3)$$

$$y_k^- = (y'_k, y'_{k-1}, \dots, y'_{k-p+1})' \quad (4)$$

The covariance matrix between the future and the past observations is the (mf×mp) Block Hankel matrix given by:

$$H = E[y_k^+ y_{k-1}^-] = \begin{bmatrix} R_1 & R_2 & \cdots & R_p \\ R_2 & R_3 & \cdots & R_{p+1} \\ \vdots & \vdots & \ddots & \vdots \\ R_f & R_{f+1} & \cdots & R_{f+p-1} \end{bmatrix} \quad (5)$$

where E denotes the expectation operator, H is the Block Hankel matrix (a block band counter diagonal matrix) formed with the (m×m) individual covariance matrices:

$$R_i = E[y_{k+i} y_k^+] = CA^{i-1}G \quad (6)$$

with $G = E[z_{k+1} y_k^+]$. The sample covariance matrices are computed from the data $y_k, k=1, 2...T$

$$\hat{R}_i = T^{-1} \sum_{k=1}^{T-1} y_{k+i} y_k^+; \quad i=0, 1...p+f \quad (7)$$

and are used to estimate a sample Block Hankel matrix $\hat{H}$. In order to identify the model coefficients, two factorizations of the Block Hankel matrix are considered. The first factorization uses the orthogonal-triangular (QR) decomposition of $\hat{H}$ as: $\hat{H} = QR$ with R an upper triangular matrix of the same dimension as $\hat{H}$ and Q an identity matrix. The second

factorization of the Block Hankel matrix considers its (mf×n) observability and (n×mp) controllability matrices, O and K as [4, 5]:

$$H = \begin{bmatrix} C \\ CA \\ \vdots \\ CA^{f-l} \end{bmatrix} \begin{bmatrix} G & A & A^{p-l} & G \end{bmatrix} = OK \quad (8)$$

The two factorizations of the Block Hankel matrix are for a given order n:

$$\hat{H}_n = \hat{Q}_n \hat{R}_n = \hat{O}_n \hat{K}_n \quad (9)$$

To estimate A, it is necessary to introduce the following shifting matrices

$$O^{\downarrow} = \begin{bmatrix} C \\ CA \\ \vdots \\ CA^{f-2} \end{bmatrix} \quad \text{and} \quad O^{\uparrow} = \begin{bmatrix} CA \\ CA^2 \\ \vdots \\ CA^{f-1} \end{bmatrix} \quad (10)$$

where, $O^{\downarrow}$ is the (m(f-1)×n) matrix obtained by deleting the last block line matrix of the observability matrix O and $O^{\uparrow}$ is the (m(f-1)×n) matrix obtained by deleting the first block line matrix of the observability matrix O. It is easy to show that the transition matrix A is related to these two matrices by the following expression:

$$O^{\uparrow} = O^{\downarrow} A \quad (11)$$

Using generalized inverses we deduce the transition matrix as:

$$A = [O^{\downarrow}]^T O^{\uparrow} \quad (12)$$

It is easy to show that the observation matrix C is obtained directly from the observability matrix O by considering the first block line of O:

$$C = O_{IL} = \{\text{first block line of O}\} \quad (13)$$

Once the transition matrix A has been estimated its eigenvalues $\lambda_i$ and eigenvectors $\phi_i$ are used to identify the modal characteristics of the system: eigenfrequencies $F_i$, damping ratios $\xi_i$ and mode shapes $\Phi_i$:
$$F_i = \frac{1}{2\pi\Delta t} \sqrt{\frac{\left[\ln(\lambda_i \lambda_i^*)\right]^2}{4} + \left[\cos^{-1}\left(\frac{\lambda_i + \lambda_i^*}{2\sqrt{\lambda_i \lambda_i^*}}\right)\right]^2} \quad (14)$$

$$\xi_i = \sqrt{\frac{\left[\ln(\lambda_i \lambda_i^*)\right]^2}{\left[\ln(\lambda_i \lambda_i^*)\right]^2 + 4 \left[\cos^{-1}\left(\frac{\lambda_i + \lambda_i^*}{2\sqrt{\lambda_i \lambda_i^*}}\right)\right]^2} \quad (15)$$

$$\Phi = C\phi \quad (16)$$

## VALIDATION METHODOLOGY

This section describes the methodology used to develop and verify efficiency of the subspace technique. Figure 3 shows the two degrees of freedom per sector model. This model simulates the dynamic behavior of specific bladed assemblies which has seventy-seven blades, and provides a temporal forced response of each blade. The mechanical and aerodynamic parameters used in this model were obtained by a theoretical mode shape analysis. Afterwards, the subspace method [4-6] was employed to identify the modal parameters of the structure: eigenfrequencies and damping ratios by using the temporal forced response of each blade. To verify the efficiency of the identification technique, the results obtained will be compared with the model's exact data found by a simple modal computation. Figure 4 presents the flowchart of the proposed methodology.

<figure><img src="image_4.png" /><figcaption>Fig.3. TWO DEGREES OF FREEDOM PER SECTOR MODEL</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig.4. THE METHODOLOGY USED IN THIS PAPER</figcaption></figure>

## MECHANICAL MODEL

The mechanical model used in this study and shown in figure 1, has blade (d o f) and disk (d o f) for each sector. Each disk (d o f) is coupled with neighbouring sectors by springs of stiffness $K_c$. This model simulates the dynamic behavior of a bladed assembly which has seventy-seven blades. The mechanical and aerodynamic parameters used in this model were obtained from a theoretical mode shape analysis. The external efforts are simulated by a random excitation and structural damping is assumed to be present.

Mistuning is modeled by adding stiffness variations to nominal blade stiffness. The i^th^ mistuned stiffness is expressed as:

$$k_{b,i} = k_b(1 + \delta_i) \quad (17)$$

where, $\delta_i$ is the i^th^ mistuning value. In this study $\delta_i$ is given using experimental data.

$$\delta_i = \left( \frac{f_i}{f_{mean}} \right)^2 - 1 \quad (18)$$

where $f_i$ is the frequency of the i^th^ blade and $f_{mean}$ is the mean frequency of all blades.

Aeroelastic coupling has been introduced, in the way of circular matrix coefficients, considering only five different coefficients. The equations of motion for the model are:

$$M\ddot{q} + C\dot{q} + Kq = f(t) \quad (19)$$
where, q is the vector of displacements, M is the mass matrix, C is the damping matrix and K is the stiffness matrix. Matrix details are provided in the appendix A.

The eigenvalues of the tuned and mistuned model are presented in figure 5. These computed eigensolutions will be called exact solutions in the rest of the paper. The tuned system is unstable, due to the presence of negative damping coefficients, while the mistuned system is stable. Mistuning couples all the tuned modes and stabilizes the modes unstable.

<figure><img src="image_6.png" /><figcaption>Fig.5. STABILITY LOOP FOR THE TUNED AND MISTUNED MODEL</figcaption></figure>

## MATHEMATICAL MODEL

The system described by the equation (19) was presented in the state space of 2N dimension (N is the number of (d o f)).

$$ \begin{bmatrix} C & M \\ M & 0 \end{bmatrix} \begin{bmatrix} \dot{q} \\ \ddot{q} \end{bmatrix} = \begin{bmatrix} -K & 0 \\ 0 & M \end{bmatrix} \begin{bmatrix} q \\ \dot{q} \end{bmatrix} + \begin{bmatrix} f(t) \\ 0 \end{bmatrix} \quad (20) $$

The temporal forced response is then given by:

$$ x(t) = Xe^{S(t-t_0)} \bar{Y}^T U x_0 + \int_{t_0}^{t} X e^{S(t-\tau)} \bar{Y}^T g(\tau) d\tau \quad (21) $$

where Y, X and S are respectively the left eigenvector matrix, the right eigenvector matrix, and the diagonal matrix of eigenvalues (using the usual mathematical notations).

This method of computation gives the vibration signal of each blade in a negligible computing time compared to methods of temporal integration (Newmark and Runge Kuta tested on the same system). The second advantage is that this

method gives exact results because the differential equation (20) is solved without approximation. Figure 6 presents a simulated temporal forced response of one blade. The simulation was made over 24 minutes with a sampling step 0.007 seconds.

<figure><img src="image_7.png" /><figcaption>Fig.6. TEMPORAL FORCED RESPONSE OF ONE BLADE</figcaption></figure>

Our purpose is to identify the eigenfrequencies and damping ratios of the blades using only such temporal responses.

## SIGNAL RECONSTRUCTION

The subspace method requires signals without aliasing to identify modal parameters, whereas tip timing data gives signals with aliasing. To solve this problem and find a link between the subspace method and tip-timing data, we propose in this paragraph an interpolation technique based on the Shannon sampling theorem [8]. This theorem is given by the following expression:

$$ s(t) = \sum_{k=-\infty}^{+\infty} s_{\text{analytic\_complex}} \left( \frac{k}{v_s} \right) \sin c(v_s t - k) \quad (22) $$

where k is the number of samples, sinc is the sinus cardinal function, and $s_{\text{analytic\_complex}}$ is the complex signal of the real signal $s_{\text{aliasing}}$ (signal with aliasing) given by tip-timing data. The Hilbert Transform is used to compute $s_{\text{analytic\_complex}}$ from $s_{\text{aliasing}}$. $s_{\text{analytic\_real}}$ is the interpolated signal of bandwidth B and center frequency $f_0$. Its expression is given by the following relation:

$$ s_{\text{analytic\_real}} = \text{real}(s(t) \times e^{i2\pi\phi_0 t}) \quad (23) $$

An analytic band-pass signal with bandwidth B and center frequency $f_0$ can be perfectly reconstructed from samples of the corresponding complex envelope at a sub-sampling rate $v_s > B$. Hence the sampling rate $v_s$ can be much smaller than twice the
highest frequency ($f_0 + B/2$). The following formula describes the reconstruction of this signal:

$$s_{\text{analytic\_real}} = \text{real} \left[ \sum_{k=-\infty}^{+\infty} s_{\text{analytic\_complex}} \left( \frac{k}{v_s} \right) \times e^{i2\pi(f_0/v_s)(v_s t - k)} \times \sin c(v_s t - k) \right] \quad (24)$$

Using this formula, the original high-frequency signals $s_{\text{analytic\_real}}$ can be reconstructed at an arbitrary time samples, which is important for many applications, such as tip-timing data.

However, formula of this type are known to be impractical, since the sinc function decays very slowly, so that many terms of the series have to be evaluated for accurate reconstruction. A well known approach to avoid this problem in the complex envelope or low-pass case is to use reconstruction kernels K that are time-limited so that only a finite number of samples have to be considered. The general reconstruction formula, for $s_{\text{analytic\_real}}$ sampled at rate $v_s$, is given by:

$$K_{v_s} s_{\text{analytic\_real}}(t) := \text{real} \left[ \sum_{k=-\infty}^{+\infty} s_{\text{analytic\_complex}} \left( \frac{k}{v_s} \right) \times e^{i2\pi(f_0/v_s)(v_s t - k)} \times K(v_s t - k) \right] \quad (25)$$

where, K is piecewise continuous function.

The optimum kernel for an arbitrary order m is given in the time domain by:

$$K_{\text{opt}}(m)(t) = \sum_{l=0}^{[(m-l)/2]} (-1)^l b_{lm} D^{2l} B_m(t) \quad (26)$$

where, $B_n$ is the B-spline function of order n. It is a piecewise polynomial of degree n-1 with finite support in [-n/2, n/2], defined recursively by:

$$B_1(t) = \frac{1}{(n-1)!} \sum_{k=0}^{[n/2-t]} (-1)^k \binom{n}{k} \left( \frac{n}{2} - |t| - k \right)^{(n-1)} \quad (27)$$
$$\times t \in [-n/2, n/2]$$

and

$$B_n(t) = (B_{n-1} * B_1)(t) \quad (28)$$

is the n-fold convolution of the rect-function, which is equal to $B_1$, with itself. The number $b_{lm}$ is defined by:

$$\left( \frac{t/2}{\sin(t/2)} \right)^m = \sum_{l=0}^{\infty} b_{lm} t^{2l} \quad (29)$$

The appendix B presents the optimum kernels used for the reconstruction of the sub-sampling signals.

Figure 7 presents the results obtained using the reconstruction formula (25).

<figure><img src="image_8.png" /><figcaption>Fig.7. RECONSTRUCTION SIGNAL</figcaption></figure>

In this figure the reconstructed signal (obtained from the signal with aliasing) is compared to the original signal (the same signal without aliasing). The reconstruction error is given by the 2-norm of difference value between reconstructed and original signal. The following relation describes how this error is computed:

$$\text{Error} = \| S_{\text{reconstructed}} - S_{\text{original}} \|_2 \quad (30)$$

for $m=1$, the error of reconstruction is equal to $6 \cdot 10^{-3}$.

## IDENTIFICATION OF MODAL PARAMETERS

Two approaches are proposed in the identification process using the subspace method described previously: global and local approach.

The global approach uses all temporal blade responses. This approach allows us to identify the mode shapes of the bladed assemblies. The local approach uses the temporal response of each blade alone. The results obtained, using the local approach, are presented here.

Figure 8 gives a comparison between the identified and the exact eigenfrequencies. The signal used in the identification process respects the Shannon theorem (no aliasing).
<figure><img src="image_9.png" /><figcaption>Fig.8. COMPARISON BETWEEN IDENTIFIED AND EXACT EIGENFREQUENCIES</figcaption></figure>

Figure 9 presents the relative errors between the identified and the exact eigenfrequencies.

<figure><img src="image_10.png" /><figcaption>Fig.9. ERRORS IN EIGENFREQUENCIES</figcaption></figure>

The average quadratic error in the identified eigenfrequencies is 3.91. It is given by the following expression:

$$J = \sum_{j} |\omega_{1j} - \omega_{2j}|^2 \quad (31)$$

To check the efficiency of the subspace method, the results found are compared to those obtained through a standard discrete Fourier Transform.

Figures 10 and 11 present the results obtained by the standard discrete Fourier transform (FFT).

<figure><img src="image_11.png" /><figcaption>Fig.10. BLADES FREQUENCIES SPECTRUM</figcaption></figure>

Figure 11 gives comparison between the identified eigenfrequencies, found by using FFT method, and the exact eigenfrequencies.

<figure><img src="image_12.png" /><figcaption>Fig.11. COMPARISON BETWEEN THE IDENTIFIED FREQUENCIES, OBTAINED BY THE FFT METHOD, AND THE EXACT FREQUENCIES</figcaption></figure>

Figure 12 presents the relative errors between the identified and the exact eigenfrequencies (the exact eigenfrequency for the i-th blade is the eigenfrequency of the mode shape for which maximum displacement is obtained for this i-th blade).
<figure><img src="image_13.png" /><figcaption>Fig.12. ERRORS IN EIGENFREQUENCIES</figcaption></figure>

The average quadratic error in the identified eigenfrequencies obtained by a standard Fourier transform is 44.78.

These results indicate that the subspace method is better suited for identifying closely spaced blade vibration frequencies than the FFT method.

The second advantage of the subspace technique is identification of the damping ratios. The random decrement technique [7] is used to extract the damped vibration from random responses of each blade. Afterwards the subspace method is used to identify the damping ratios of each blade. Figure 13 gives the damped vibration of blade n° 77.

<figure><img src="image_14.png" /><figcaption>Fig.13. DAMPED VIBRATION OF BLADE N°77</figcaption></figure>

Figure 14 gives the relative errors between the identified and the exact damping ratios (the exact damping ratio for the i^th^ blade is the damping ratio of the mode shape for which maximum displacement is obtained for this i^th^ blade).

<figure><img src="image_15.png" /><figcaption>Fig.14. ERRORS IN DAMPING RATIOS</figcaption></figure>

We observe that the level of relative error is very high for 12 blades (more than fifty per cent). We are currently working to improve the identification of damping coefficients using a weighted subspace method.

CONCLUSION

A new analysis method of blade tip-timing data is proposed in this article. This technique is based on the subspace method using QR decomposition. The identified eigenfrequencies have been compared to those obtained by an FFT standard method. The comparison between the average quadratic errors of the two methods shows that the subspace method is more efficient than the standard FFT method: the average quadratic error in identified eigenfrequencies using the subspace method is 3.91 while it is of 44.78 using the standard FFT method.

It can be noted that the subspace method gives easily the damping ratios; these coefficients can not be identified by FFT method. However, the level of relative errors in damping ratios is remarkably high for certain blades. We are currently working to improve identification of the damping coefficients using a weighted subspace method.

Another advantage of the subspace method is the identification of the structure mode shapes. This identification method gives good results only if all temporal forced responses of the bladed assembly are simultaneously treated. However, identification using the subspace method requires long signals in the case of eigenfrequencies closely spaced. Consequently, computing time will be bulky and identification of mode shapes becomes impossible for bladed assemblies which have more than thirty blades. We are currently working to improve identification of mode shapes by a new technique based on a weighted subspace method. This technique allows identification of mode shapes in a reasonable computing time.

Despite the good results of the subspace method, in particular the identification of closely spaced eigenfrequencies, this method suffers from some disadvantages. One such disadvantage is that the subspace method, like others, is robust only if the signals are not aliased. To surpass this problem a
reconstruction technique of the sub-sampled signal has been proposed. This technique provides good results and allows identification using the subspace method.

## REFERENCES

[1] Carrington I.B., Wright J.R., Cooper J.E., and Dimitriadis G., A comparison of blade tip timing data analysis methods, *Proc Instn Mech Engrs* Vol 215 Part G, 2001.

[2] Craig P., Lawson P., and Ivey C., Turbo-machinery blade vibration amplitude measurement through tip timing with capacitance tip clearance probes, *Sensors and actuators A* 118, pp 14-24, 2005.

[3] Dimitriadis G., Carrington I.B., Wright J.R., and Cooper J.E., Blade tip timing measurement of synchronous vibrations of rotating bladed assemblies, *Mechanical Systems and Signal Processing*, 16 (4), pp. 599-622, 2002.

[4] Lardies J., Estimation of parameters and model order in state space innovation forms, *Inverse Problem in Engineering*, Vol. 8, pp. 75-92, 2000.

[5] Lardies J., A stochastic realisation algorithm with application to modal parameter estimation, *Mechanical Systems and Signal Processing*, Vol. 15(2), pp. 275-285, 2000.

[6] Zhang Y., Zhang Z., Xu X., and Hua H., Modal parameter identification using response data only, *Journal of sound and vibration*, Vol. 282, pp. 367-380, 2005.

[7] Ibrahim S.R., Random Decrement Technique for Modal Identification of Structures, *Journal of Spacecraft and Rockets*, Vol. 14 (11), pp. 696-700, 1977.

[8] Ries S., Digital time-delay beam forming with interpolated signals, *Signal Processing* 84, pp. 2403-2423, 2004.

## APPENDIX A – C, M, & K MATRICES

The mass matrix M is given by the following expression:

$$ M = \begin{bmatrix} M_a & 0 & \cdots & 0 \\ 0 & \ddots & \ddots & \vdots \\ \vdots & \ddots & \ddots & 0 \\ 0 & \cdots & 0 & M_a \end{bmatrix} \quad (32) $$

where $M_a = \begin{bmatrix} m1 & 0 \\ 0 & m2 \end{bmatrix}$, m1 is equal to 4.78 kg, and

m2 is equal to 30.79 kg.

The structural stiffness matrix K is given by the following expression:

$$ K_{struct} = \begin{bmatrix} \hat{K} & K'' & 0 & K'' \\ K'' & \ddots & \ddots & 0 \\ 0 & \ddots & \ddots & K'' \\ K'' & 0 & K'' & \hat{K} \end{bmatrix} \quad (33) $$

where $\hat{K} = \begin{bmatrix} k1 & -k1 \\ -k1 & k1+k2+2k3 \end{bmatrix}$ and $K'' = \begin{bmatrix} 0 & 0 \\ 0 & -k3 \end{bmatrix}$

$k1$ is equal to $6.7206 \times 10^5 \text{ N/m}$, $k2$ is equal to $3.378 \times 10^7 \text{ N/m}$, and $k3$ is equal to $2.158 \times 10^9 \text{ N/m}$.

The structural damping matrix is given by

$$ C_{struct} = P^{-T} \text{diag}(2\xi\omega_n) P^{-1} \quad (34) $$

where, $PP^{-1} = I$ ($P$ is the matrix of conservative mode shapes), and the damping coefficient $\xi$ is equal to 0.02 %.

The damping matrix with aeroelastic effects is given by the following relation:

$$ C = C_{struct} - F_{aero}^{R} \quad (35) $$

where $F_{aero}^R$ is:

$$ F_{aero}^R = \begin{bmatrix} C_1^a & C_2^a & C_3^a & \cdots & C_{n-1}^a & C_{n}^a \\ C_1^a & C_2^a & C_3^a & \cdots & C_{n-1}^a & C_{n}^a \\ C_1^v & C_2^v & C_3^v & \cdots & C_{n-1}^v & C_{n}^v \\ C_1^v & C_2^v & C_3^v & \cdots & C_{n-1}^v & C_{n}^v \\ 0 & 0 & 0 & \cdots & 0 & 0 \\ C_1^k & C_2^k & C_3^k & \cdots & C_{n-1}^k & C_{n}^k \\ C_1^k & C_2^k & C_3^k & \cdots & C_{n-1}^k & C_{n}^k \end{bmatrix} \quad (36) $$

The stiffness matrix with aeroelastic effects is given by the following expression:

$$ K = K_{struct} + \omega F_{aero}^I \quad (37) $$