Research Objectives

The objective of this research is to, first, enhance the Covariance-driven Stochastic Subspace Identification method (SSI-COV) to the named “Singular Spectrum Analysis–

Covariance driven Stochastic Subspace identification method” (SSA-SSI-COV), validated both by numerical simulation and the application in system identification of Canton Tower, a benchmark problem for structural health monitoring of high-rise slender structures.

Second, develop the recursive Singular Spectrum Analysis method (rSSA), and in conjunction with the recursive Covariance-driven Stochastic Subspace Identification method to construct the named “recursive Singular Spectrum Analysis – Covariance driven Stochastic Subspace Identification method” (rSSA-SSI-COV), through a moving window approach. The method will be validated firstly by numerical simulation and later by application in the damage detection and health monitoring of laboratory experiments.

The organization of this thesis is briefly described as follows:

Chapter 2: The basic methodology of subspace identification algorithm is recalled through, firstly, the introduction of the dynamic model of a linear system, followed by the formulation of SSI-COV and SSI-DATA method, and finally, the Singular Spectrum Analysis (SSA) procedure will be described.

In system identification algorithms, it is important to distinguish the structure modes from the spurious modes because the order of the real system is always unknown. The alternatives to build the stabilization diagram will be introduced and compared one to

another. A comparison benefit-drawback and implementation issues will be discussed through a numerical simulation example and the identification of a laboratory test.

Chapter 3: A comprehensive numerical study and comparison between different SSI algorithms is carried out. Measurement noise effect and the addition of a noise which violates the SSI assumption is discussed. Identification of the simulated nonlinear signals, signals with time-varying frequency, signals with closely-spaced frequencies mixed with white noise is done to understand the performance of SSI algorithms under different scenarios of assumption violation and the mechanism to overcome this difficulties. The SSA-SSI-COV algorithm is introduced in this chapter to solve the identification problem of closely-spaced frequencies with added white noise.

Chapter 4: Application of SSI algorithms in system identification of the Canton Tower is discussed. The order determination procedure through the SSA-SSI-COV algorithm will be described. Comparison between different SSI approaches is made in this chapter.

The procedure called decimation is although studied and applied to increase the convergence speed of the stabilization diagram.

Chapter 5: the derivation of Covariance-driven Recursive Stochastic Subspace Identification algorithm (RSSI-COV) can be found in this chapter. The Projection Approximation Subspace Tracking algorithm (PAST) and its Instrumental Variable extensions (EIV-PAST) is also described and implemented to RSSI-COV. To consider the noise contaminated data, a recursive pre-processing technique called recursive singular spectrum analysis technique (rSSA) is derived to enhance the accuracy and stability in the online tracking capability.

8

Chapter 6: the RSSI-COV method and the proposed rSSA-SSI-COV algorithm through a moving window approach are validated in this chapter by means of numerical simulation of a 6 DOF system, cases with sudden reduction and slow decreasing in system stiffness are studied. The effects of the selected RSSI model parameters in the online modal analysis, and the influence of time-varying frequencies in the selection of system order are also discussed.

Chapter 7: the RSSI-COV method and the proposed rSSA-SSI-COV algorithm through a moving window approach are applied to the monitoring and damage detection of, first, shaking table test of a 3-story steel structure with instantaneous stiffness reduction.

Second, the shaking table test of a 1-story 2-bay reinforced concrete frame subjected to earthquake excitations with increasing intensity. Finally, application to the monitoring of a three pier and four span steel bridge under continuous scour is carried out.

Chapter 8: Summaries and suggestions for the use of the proposed algorithms will be given here. The potential research topics are indicated at the end.

Chapter 2

Stochastic Subspace Identification Methods

2.1 Introduction

In output-only characterization, the ambient response of a structure is recorded during ambient influence (i.e. without artificial excitation) by means of highly-sensitive velocity or acceleration sensing transducers. The Stochastic Subspace Identification (SSI) technique is a well known multivariate identification technique for output-only measurements. It was proved by several researchers to be numerically stable, robust to noise perturbation and suitable for conducting non-stationarity of the ambient excitations although its stationary assumption is violated [5, 37, 53].

The SSI-DATA algorithm was fully enhaced by Van Overschee and De Moor [47], while SSI-COV algorithm has as its antecedent the Eigensystem Realization Algorithm [25] for the free response of a structure, which are applied along with the Natural Excitation Technique (NExT) or Random Decrement (RD) functions. This chapter will begin with the introduction of the dynamic model of structures, followed by the stochastic properties and the system realization methods of each subspace algorithm.

2.2 Models of vibrating structures

2.2.1 Continuous-time state-space model

The Finite Element model of a linear time-invariant dynamic system can be expressed as:

( )

^{t q}

( )

^{t q}

( )

^{t F(t) u(t)}

q C K L

M&& + ₂& + = = _(2.1)

10

q&& is the acceleration vector with the same dimension as the displacement vector.

(t)∈ℜn

L is the input location matrix.

(t)∈ℜm

u is the vector describing m inputs as a function of time t.

n is the number of DOFs and m is the number of inputs. input matrix in the state equation. A_c and B_c are arranged as follows:

n

solution [57]:

( )

^{t =e (t−t)}

( )

^{t +}

∫

_t^{te ( )t−}

( )

^d

0 0

0 ^{τ u}τ τ

x

x ^{Ac Ac} B_c (2.4)

where the 1^st term is the free vibration solution given an initial condition x(t₀), and the 2^nd term is a typical convolution integral. Through an eigen-analysis of the system matrix Ac, the state equation can be decoupled through a coordinate transformation using the obtained complex eigenvectors.

A_c=ΨΛ_cΨ⁻¹ , ^x

( )

^t =^Ψη

( )

^t (2.5)

where^η

( )

^t is the generalized coordinate. Λc∈ℜ^{2n 2× n}is a diagonal matrix containing complex eigenvalues λ_i in the diagonal which appear in conjugate pairs, Ψ∈ℜ^{2n 2× n} are the complex eigenvectors. From the eigen-analysis A_cΨ=ΨΛ_c, one may find that they have the following structure:





 



= _* Λ

Λ_c Λ = , 









= _*^*_* Λ Θ ΘΛ

Θ

Ψ Θ _(2.6)

In fact, it can be easily verified that Λ are the same eigenvalues and Θ the same eigenvectors, i.e., mode shapes, than those obtained by conducting eigen-analysis directly in the unforced equation of motion (2.1), but they cannot be used to decouple the equation of motion unless it is a proportionally damped system.

Then, the decoupled state equation can be written as follows:

^η&

( )

^t =Λc^η

( )

^t +Ψ⁻¹Bc^u

( )

^t (2.7)

Furthermore, to relate the obtained complex eigenvalues to a physical interpretation, a Taylor Expansion is required to decouple the free vibration term e^Act in (2.4), which is a matrix exponential:

12 term in (2.4) and having in mind that the complex eigenvalue has its real and imaginary part: λ_i =α_i+ jβ_i, solution to the i-th mode free vibration is: free vibration solution. Comparing (2.9) with the well-known free vibration solution of a SDOF system, the so-called i-th effective modal frequency ω_i^′ and effective damping ratio ζ_i^′ can be realized:

ζi^′ is related to the phase. Hence, when a structural system is changed due to damage, the migration of system poles will be directly reflected by the computed effective modal frequency and damping ratio, which the term “effective” will be omitted hereafter.

The observation equation:

If only subsets of the n DOF can be measured, and considering that measurements are taken at l locations and the sensors can be either accelerometers, velocity or displacement transducers, the observation equation can be defined as:

( )

^{t q}

( )

^{t q}

( )

^{t q}

( )

^t

y =C_a&& +C_v& +C_d (2.11)

where ^y

( )

^{t ∈ℜl} represents the l outputs. C_a, C_v and C_d∈ℜ^l×n are the output location matrices corresponding to acceleration, velocity and displacement respectively.

To relate the output y(t) to the system state x(t), the equation of motion (2.1) can be used to eliminate ^q&&

( )

^t , and by arranging and grouping location matrices, the observation equation become:

( )

^t =Cc^x

( )

^t +Dc^u

( )

^t

y (2.12)

where ^Cc=

(

^{Cd−CaM−1K Cv−CaM−1C2}

)

^{∈ℜl 2× n} is the output matrix, and

l×m

− ∈ℜ

=C M L

D_{c a} ¹ is the direct transmission matrix.

Although the eigenvectors of system matrix A contains mode shapes information as that shown in (2.6), however, there is no knowledge about the location of each DOF when the matrix A is identified, moreover, usually the number of modes, i.e., order of the system extracted from measurement data is different than the number of sensors, thus, the system eigenvectors should be mapped to the sensor locations through the output location matrix C_c:

Ψ C

V_c = _c (2.13)

where V_c are the observed mode shapes.

14

2.2.2 Discrete-time state-space model

Since all data is sampled in discrete time, the above continuous time state-space model can be converted into a discrete time state-space model. By gathering together the state and observation equation:

k measurement. A is the system matrix, B is the input matrix, C is the observation matrix and D, the direct transmission matrix, all in discrete-time. The relationships between these matrices in discrete-time and continuous time are the following [24]:

e ^∆t

The eigenvalues µi of the discrete-time system matrix A can be, therefore, related to the continuous-time eigenvalues by

( )

Then, frequencies and damping ratios can be computed as mentioned before. Both the observation matrix and complex eigenvectors are not affected by the discretization in time, the above-mentioned equations can be used without any change.

This model is called the deterministic state-spaced model since both input and output are known.

2.2.3 Stochastic state-space model

Considering that there is always noise and perturbations both in the system (due to modelling inaccuracies) as in the measurementes, therefore (2.14) can be modified to its combined deterministic-stochastic state-space model:

k unknows. Both the input terms and the noise terms are assumed to be a spatially white noise and they can be combined together as the process noise . Therefore, the discrete-time sthochastic state-space model can be simply stated as:

k and with the following covariance matrices:

( )

T ^pq otherwise δ_pq=0). p and q are arbitrary time instants.

One should note that if acceleration measurement is used, the direct transmission

16

term Du_k is also considered by the stochastic model as a process noise. Reviewing expressions (2.11) and (2.12), if only velocity or displacement transducers are used, Cc

= ( 0 Cv ) or Cc = ( Cd 0 ), and Ca = 0, the direct transmission matrix vanishes, therefore, theoretically the external excitation will not be measured.

In the case of structures subjected to ground motions such as earthquakes, L is an identity matrix, and u

( )

^t =−Mq&&g

( )

^t , hence, quantities in the state vector shown in (2.2) will be relatives. If accelerometers are used and since it measures absolute accelerations, the ground acceleration should be added to the observation equation (2.12):

( )

^t x

( )

^t u

( )

^t x

( )

^t

[

^Mq

( )

^t

]

qg x

( )

^t

y =C_c +D_c =C_c +C_aM⁻¹− &&_g +C_a&& =C_c

The Dcu(t) term will be cancel out with the ground acceleration, i.e., in the case of base excited structures, the acceleration measurements are free from the external noise contributed by D_c term. But this is not the case for structures excited by wind or other sources acting directly in the body of the structure, this externally imposed acceleration will be transmitted in the measurements as a measurement noise.

Properties of the stochastic state-space model is summarized in the following chart [37]:

The system state is a stochastic process and assumed to be stationary with zero mean:

=Σ ] [ _{k k}^T

E x x , E

[ ]

^xk ⁼⁰ (2.20) The noise terms are zero mean white noise and uncorrelated with the current system state:

0 w x ]= [ _k ^T_k

E , E[x_kv^T_k]=0 (2.21)

The output covariance matrices R_i∈ℜ^l×l of arbitrary time lag i are defined as: From the stationarity, spatially white noise assumption and the previous definitions, following properties can be deduced:

Q From the stochastic state-space model and applying stochastic properties shown in (2.24), the most important property can be deduced:

G CA

Ri = ⁱ⁻¹ (2.25) This last will be the key property to derive the Covariance-driven SSI algorithm.

2.3 Covariance-driven Stochastic Subspace Identification (SSI-COV)

The SSI-COV stems from the need to solve the problem through identifying a stochastic state-space model (matrices A and C) from output-only data. The first step is to gather the measurement vectors in a Hankel Data matrix:



18

where Y_p denotes the past measurements and Y_f denotes for the future measurements. It can be easily find that the block Toeplitz matrix can be obtained by a multiplication between future and transpose of past measurements:

( )

p ^T

where Ri is the block output covariance with time lag i defined in (2.22).

Through the stochastic property in (2.25), the Toeplitz matrix can be factorized into the extended observability matrix O_i∈ℜ^li×2n and the reversed extended stochastic controllability matrix Γ_i∈ℜ^2n×li , as shown below: constitutes. Singular Value Decomposition (SVD) is the tool used to perform the above mentioned factorization:

( )

_T^{T T} and C) can be computed by splitting the SVD in two parts:

2 MATLAB notation, the C matrix is just the first block of Oi:

( )

^:l,:

i1 O

C= (2.31)

System matrix A can be computed by exploiting the shift structure of the extended observability matrix Oi:

A

where ( )⋅^Ddenotes pseudo-inverse. Then, by conducting eigenvalue decomposition on the system matrix A, after the eigenvalues are converted to a continuous-time poles with (2.16), modal frequencies and damping ratios can be computed with (2.10). The observed mode shapes can be obtained by applying (2.13).

2.4 Data-driven Stochastic Subspace Identification (SSI-DATA)

As opposite to SSI-COV, the Data-driven algorithm (SSI-DATA) avoids the calculation of covariance. Instead, the data reduction step is accomplished by projecting the row space of the future outputs into the row space of past outputs. Covariance and orthogonal projection are closely related, in that they are both intended to eliminate uncorrelated noise contributions. From the data structure shown in (2.26) the orthogonal

20

projection can be defined as follows:

(

^{p Tp}

)

^{p i} orthogonal projection matrix. The main theorem of stochastic subspace identification [47] implies that the extended observability matrix Oi can be found from the result of the estimates from the forward non-steady state Kalman filter [47].

(

^{1 ... 2 1}

)

ˆi = xi xi+ xi+j− xi+j−

X (2.36)

Instead of (2.34), the orthogonal projection can be performed by a numerically robust and stable tool called LQ decomposition (this is, actually, the transpose version of the well-known QR decomposition), which is applied directly on the Hankel data

where L is a lower triangular matrix, and Q is an orthogonal matrix. Lij are partitions of the lower triangular matrix and Qij are partitions of the Q matrix.

j

space of Oi can be obtained directly from the column space _

why these algorithms are called “subspace” identification algorithms. They retrieve system matrix as the subspace of the projection matrix. If only system matrices (A,C)

singular values) from the noise subspace, which corresponds to vanishing singular values, and from now on, everything can be computed as that outlined in SSI-COV.

Estimating the noise covariances: Q, R and S

The Kalman filter state sequence Xˆ_i can be determined from the projection matrix by:

i i

i O P

Xˆ = ^D _(2.39)

where Oi is obtained by applying (2.29) and (2.30), ( )⋅^Ddenotes for pseudo-inverse.

The Hankel data matrix can be split in a different way:



to the past outputs as its last block row. Similar to the main theorem of the orthogonal projection [47], it can be realized that:

( )

^{1 1}

Similarly to (2.38), Oi-1 computed in a numerically stable way is in terms of the LQ

22

Since O_i-1 can be computed just deleting the last l rows, P_i-1 is calculated by the LQ decomposition, then, the shift Kalman state sequence can be determined by doing pseudo-inverse to Oi-1:

1 (A,C) are determined, the residual sequence, i.e., the process noise sequence ρw and ρv

can be easily calculated:

i algorithm guarantees the positive realness of the identified error covariances [47]:

( )

^

2.5 Singular Spectrum Analysis (SSA)

SSA [65] is a novel non-parametric technique used in the analysis of time series based on multivariate statistics. This method was firstly applied to extract tendencies and harmonic components in meteorological and geophysical time series [3]. Except the extraction of tendency, SSA can be applied to eliminate noise effect, or to detect the singularities, e.g., to extract structural residual deformation [32]. Basically, this method is capable of decomposing the original series into a summation of principal components, so that each component in this sum can be identified as a tendency, periodic components (stationary), nontationary signal or noise.

The SSA procedure starts from: (1) Embedding: generate a Hankel matrix from the time series itself by sliding a window that is shorter in length than the original series.

Firstly, let Y =_(y1,y2,K_,yN) be the time series of length N. And let L be the window length or number of block rows of the Hankel matrix, which is an integer in 1< L<N.

Each sliding window vector X_j with length of L would then be derived: X_j ={ y_j, y_j+1 , …,

(2) Singular Value Decomposition: the Hankel matrix can be represented in the form: X

= E1+ E2+…+ Ed, where d is the number of non-zero eigenvalues of theL×Lsample covariance matrix SCOV = X．X^T. The i-th elementary matrix, or called i-th eigentriple, are given by E_k =

λk u_k v_k^T = s_k u_k v_k^T where_λ1,_λ2,..._λd are the non-zero eigenvalues of

24

(3) Grouping: the plot of the singular values in descending order is called the singular spectrum and is essential in deciding the index from where to truncate the summation.

Finally, decide a parameter r to reconstruct an approximate matrix of X, i.e.

is reconstructed, these may become slightly different, therefore, these entries are averaged to reconstruct signal:

 approximate matrix X~

, makes minimum the Frobenius norm of the error between X~ and the Hankel matrix assembled by the reconstructed signal, i.e., the averaging leads to the optimal signal reconstruction in terms of the principal components.

2.6 Pole discrimination: the stabilization diagram

2.6.1 Alternatives to build the stabilization diagram

In real world applications, noise and perturbations are always present at any measurement, and there is no prior information about the number of modes can be extracted from the data, i.e., there are always uncertainties in the determination of system order. Therefore, a stabilization diagram is used to discriminate between noise or spurious poles and true system poles. Based on the procedure of SSI-COV, there are several ways to build the stabilization diagram:

1^st version: Decide first the maximum dimension of the Toeplitz matrix shown in (2.27), perform SVD, and let the order of system matrix A increases from a lower value till reaching the maximum dimension of the Toeplitz matrix defined by the user.

The advantage of this version is that only once has to be done the SVD; less time is consumed in the construction of the stabilization diagram. The drawback is that, there is no clear criterion to ensure that the chosen maximum dimension is sufficient or not to reveal true system information. While the order of system matrix A is increasing, more noise information will be included in the system matrix A, consequently, more noise or spurious poles will appear on the diagram. For this purpose, modal transfer norm [41]

was introduced in addition, to clear out the large number of spurious poles at higher orders thus clarifying the stabilization diagram. But again, a threshold level has to be defined for the modal transfer norm. The concept behind this version to construct the stabilization diagram is that, even including more spurious poles in the system matrix A, the true modes (frequency, damping, mode shape) extracted by eigen-decomposition will remain stable.

26

2^nd version: Determine the order of system matrix A by observing the variation of the singular values, and then, increase the size of the Toeplitz matrix both rows and columns holding the order of system matrix unchanged. It is important to note that if full sensors are used to compute the covariance as usually does, the formed covariance block is a square matrix, for convenience, this version will be called the “square Toepliz matrix” because the shape of Toeplitz matrix remains squared. The main drawback of this alternative is, first, more time consuming, and second, the system matrix order must be defined previously. For field measurement data, generally there is no clear gap on the distribution of singular values as that appearing in numerical simulation. The advantage of this version is that, one do not have to try at the beginning the maximum Toeplitz matrix size, since the required size to achieve good results may vary from case to case.

Increase of the Toeplitz matrix dimension means a larger subspace dimension and also more data to extract the orthonormal base which spans the system-related

Increase of the Toeplitz matrix dimension means a larger subspace dimension and also more data to extract the orthonormal base which spans the system-related

在文檔中 協方差型隨機子空間識別法之應用 (頁 27-0)