Recursive Singular Spectrum Analysis (RSSA)

Chapter 5 Recursive Stochastic Subspace Identification algorithms

5.2 Recursive Singular Spectrum Analysis (RSSA)

+− +

+₁ = H^* ₁ − y ₁ y ₁

H µ ^(5.28i)

Since the updating task is done for the time-varying column subspace U1(t) from the formulas shown above, the system information can be then extracted from U1(t) as previously discussed in section 2.3, for each time instant.

5.2 Recursive Singular Spectrum Analysis (rSSA)

To be able to apply SSA in online filtering of vibration measurements, an on-line version of the algorithm that describe the current signal structure at each time instant is required. As mentioned before in section 2.7, the first step is to assemble the measurement data in a Hankel data matrix of N data points X(N) as shown in (2.47).

Although a forgetting factor 0<λ<1 can be applied to gives different weights to the data in terms of its age, through several simulation studies have been carried out, one conclude again that for subspace-based algorithms, the moving window approach should be adopted. The number of block rows is i’ and is kept constant meanwhile a new data point is added as a new column appended to the moving window Hankel matrix:

72 i.e., number of block rows of the Hankel data matrix; L’ is the length of moving window, K’=L’–i’+1 is the number of columns, with K’ > i’. For convenience, the subscript notation for the sliding window vector is different than that used in off-line SSA, which becomes now X_{N-L’+i’+j} , j =0,1,2, …, K’.

Similar to that shown in (5.24), by applying SVD to the Hankel data matrix

(

N+¹

)

X , the left singular vectors U can be computed via Eigen-Decomposition (ED) of the covariance of sliding window vectors:

where C_SSA(N+1) is the covariance matrix of the sliding window vector X_{N-L’+i’+j}.

To make rSSA algorithm possible, only the left singular vectors corresponding to the non-zero singular values of X(N+1) will be used and updated as on-line filters: since they correspond to the column subspace which span the range of X(N+1), i.e., any sliding window vector X_{N-+L’+i’+j} can be expressed as a linear combination of the computed column subspace from X(N+1).

Since the eigenvectors of the covariance matrix C_SSA(N+1) correspond to the desired column subspace, this is actually a typical rank-two modification of the symmetric eigen-problem, i.e., meanwhile a new data column is appended to the Hankel matrix (i.e., rank-one modification), an old data column is subtracted (i.e., rank-two

modification). The above-mentioned PAST algorithm is suitable to be implement to rSSA, because it is able to update in a recursive fashion the dominant eigenvectors of the signal covariance C_SSA(N+1), i.e., the subspace of X(N+1).

Adaptation of PAST to rSSA

From (5.8), the random vector z(k) can be substituted by the sliding window vector X_k:

( )

By introducing the same approximation:

( ) ( )

H k

k X

h′ = W −1 (5.32)

the original cost function is converted into a quadratic criterion:

( )

When the matrix inversion lemma is applied to (5.34), the well-known RLS algorithm can be easily derived, which is shown in the chart below.

1. From an initial SVD the recursive algorithm can be initialized with U′1

( )

^N , later

h can be updated using the following algorithm:

Updating:

The recursive SSA is different than the off-line SSA in the reconstruction step. The recursive SSA is applied here as a preprocessing tool to filter out the undesired measurement noise keeping the system related information, thus, the same “orthogonal projection” concept used to derive PAST can be directly applied for rSSA. After the

Hence, for each new incoming data, a new sliding vector column XN+1 is appended, the column subspace is updated to U1′

(

N+¹

)

, and the reconstructed data vector X^~ _N₊₁

can be obtained by procedure shown above. Finally the reconstructed data vector

X N is placed in the corresponding location of the reconstructed Hankel data matrix

(

)

~ N+

X , then, elements of the same time instant (in the anti-diagonal direction) can be averaged to reconstruct the signal:

( )

somehow different than the algorithm used in SSA, more precisely, this is a combined

approach of both SSA and data compression based on the Karhunen-Loéve (KL) transformation [56], where a sequence of data vectors is coded by their principal components.

Chapter 6 Simulation study of RSSI-COV and rSSA-SSI-COV

To validate the adaptation of EIV-PAST algorithm to RSSI-COV, and the proposed rSSA solved by PAST, together with its applicability and stability in the on-line tracking of system modal parameters, through a moving window and recursive approach, simulation study was carried out firstly.

Consider a simulated linear 6-DOF system consists in a lumped mass model and a shear building type stiffness matrix, which has the following modal frequencies and damping ratios in its original state:

The system natural frequencies are:

f = [0.9972 , 2.9254 , 4.6600 , 6.0905 , 7.1353 , 7.7556 ] Hz

Rayleigh damping was assumed for the derivation of damping matrix, the assumed damping ratio for each mode are:

ξ = [0.03 , 0.03 , 0.01 , 0.01 , 0.02, 0.02 ],

Response is generated using discrete time deterministic state-space model having a spatially white noise as the input, outputs are measured at each DOF. Measurement noise can be added after the system response is obtained. The sampling rate is 200 Hz and the total generated data length is 20000 points, which equals to 100 sec.

6.1 Implementation of the RSSI-COV and rSSA-SSI-COV algorithm

In conducting RSSI method, considering the recursive SVD-updating task carried out by EIV-PAST algorithm updates only the dominant singular vectors, a procedure is proposed to avoid the accumulation of error during the subspace computation:

1. Select a length for the statistic moving window WLstatistic which samples the identified modal frequencies along the time axis.

2. Calculate the mean frequency over data sampled by the moving window for each modal frequency.

3. Calculate the Euclidean norm of the standard deviation of each modal frequency, i.e., take the square root over the sum of squares for the computed standard deviations over the moving window of each modal frequency.

4. Repeat 2 and 3 by moving forward the statistic window over a specific frequency data points (d_statistic).

5. Compare the Euclidean norm of displacement “k” and “k+1”, and calculate their percentage of difference. If the percentage between the two segments is greater than the specific criteria (e.g. 50%), the traditional SVD will be computed at that given time instant and serves as a restart for EIV-PAST.

6. Continues the recursive computation of modal frequencies, the statistic window moves a step forward until the displacement length of the statistic moving window dstatistic has been reached, and repeat the steps 2, 3, 4, 5.

The above procedure is summarized in the flow chart shown in Figure 6-1.

6.2 Simulation study 1: time invariant 6-DOF system

The simulation study starts with the case where the system is time-invariant. Thus, the most basic requirement to be meet is the sufficient stability to track the modal parameters in a recursive manner, when these latter have not any change.

Noise free

Figure 6-2 shows the recursive tracking results for the response generated by the 6-DOF system subjected to a spatially white noise input excitation, without adding any kind of noise. From now on, the same number of block rows i will be used for block columns, i.e., a square Hankel Covariance matrix. Since this is a simulation example with known DOF, naturally the system order is defined to be 12. Comparing a) and b), one can note that by increasing the moving window length L used to form the Hankel Covariance matrix, the stability is better. This matches the assumption for Covariance for which the data length must be theoretically tending to the infinity.

Figure 6-3 shows the result for the tracking of damping ratio. From a), the use of a short moving window (1500 points equals to 7.5 sec) leads to a very unstable damping ratio and wrong trace for the 1^st mode. By increasing the window length to 3000 points, although the deviation from correct answers is still large as it always occurs in damping ratio estimates due to its high sensitivity to any perturbation, at least, stable results were achieved. The two largest damping ratios correspond to the 1^st and 2^nd mode.

Adding noise correlated with output (SSI assumption violated)

Consider the addition of the input acceleration to the measurement as a noise correlated with output, which occurs in acceleration measurements as that discussed in section

2.2.3, the purpose is to test the robustnss of RSSI-COV when its assumptions were violated. In this case, the moving window length is fixed to 2000 points, system order is 12. Frequency tracking results is shown in Figure 6-4. The effect of increase number of block rows i is shown by comparing a) and b). The 1^st mode shown in a) was perturbed by noise; the higher the number of block rows i, the better the noise separation, consequently a more stable tracking result was obtained in b).

The recursive tracking of damping ratio is shown in Figure 6-5 for the same case when the added noise violate the assumptions. If the number of block rows i is insufficient to separate noise from the system information, negative damping appears due to noise perturbation as shown in Figure 6-5 a). Although negative damping was corrected by increasing the number of block rows in Figure 6-5 b), the trace of damping is still unstable and several significant changes occurs in damping tracking meanwhile there is not any variation in the correct damping ratio answer. Given this situation, damping ratio cannot be used as an reliable index of system change or damage.

6.3 Simulation study 2: time varying 6-DOF system with sudden stiffness reduction.

To simulated drastic system change scenarios, sudden stiffness reduction and damping ratio changes were introduced in the 6-DOF system. The loss of stiffness was introduced in the 1^st DOF of the simulated shear-building type system, the resulting changes in modal frequency are presented in Table 6-1. Comparing with the precision level of SSI shown in section 3.1, the resulting change in the 1^st modal frequency become significant (more than 0.1 Hz for the 1^st mode) until the stiffness reduction in 1^st DOF is more than 50%.

Noise free

Firstly, the RSSI-COV algorithm is applied to trace these sudden drops in modal frequencies without adding noise to the simulated measurement. Result for modal frequencies is shown in Figure 6-6. The number of blocks i is 70, system order is 12;

the statistic window length WL_statistic is 400 frequency points, displacing by every 100 points (as explained in Figure 6-1). The restart criterion is 50% of difference between this and the previous standard deviation Euclidean norm. Similar to the results obtained in 6.2, the tracking stability was enhanced by increasing the moving window length.

However, as that shown in Figure 6-7 b) and Figure 6-8 b), variations in damping ratio cannot be traced correctly, no clear tendency can be observed even increasing the moving window length or system order. Since RSSI algorithm only shows good frequency tracking capability, hereafter the results for damping ratio will be omitted.

Between 7000 and 9000 data points, the 1^st mode damping ratio shown in Figure 6-7 b) has increased to an abnormal level, which indicates existence of problems in the identification. Until the system order is increased to 16, 2 pairs of SV more than that required theoretically, as that shown in Figure 6-8 a) for frequency and b) for damping ratio, a more reasonable range for the 1^st mode damping ratio is reached (the negative damping ratio corresponds to spurious poles) although it is not accurate at all.

Based on the previous experience adquired in off-line system identification, since the system in no longer time-invariant, more orthogonal components is required to span the system subspace, therefore, if insufficient system order is defined, information loss in the recovery of system matrix A will lead to an unstable tracking on the 1^st modal frequency or damping ratio. The use of more system order, however, introduces spurious poles in the time-frequency plot.

Figure 6-9 shows the six mode shape identification results for several selected points, which are obtained with a system order of 12. Except mode 1, all the identified mode shapes are almost identical to the correct answer. This phenomena is similar to that we find in modal frequency and damping ratio tracking, the 1^st mode seems to be the most affected by any violating any assumption of SSI, because we are trying to use a recursive linear identification algorithm, with a certain moving window length, to identify time-varying systems.

Adding noise correlated with output (SSI assumption violated)

The combine effect of system with sudden changes in frequency and the addition of noise correlated with output is presented in Figure 6-10, both violating SSI assumptions. The moving window length is 5000 points for both a) and b). With a system order of 12 and 150 block rows in a), bad tracking results were obtained and the 1^st mode almost disappeared. The combination of noise and time-varying signal increased the required column subspace to describe the system information, although the system order is 12 theoretically as that shown in Figure 6-10 b): by increasing the system order to 18, the 1^st modal frequency can be traced correctly, but again, spurious poles appears in the time-frequency plot.

6.4 Simulation study 3: time-varying 6-DOF system with gradual stiffness reduction

A gradual reduction in the 1^st DOF stiffness from point 4001 to point 16000 is introduced in the system matrix to simulate slow varying frequencies. There are totally 20000 points in the simulation example. The system natural frequencies and the corresponding damping ratios from 1 to 4000 points are:

f initial = [0.9972 , 2.9254 , 4.6600 , 6.0905 , 7.1353 , 7.7556 ] Hz

ξ initial = [0.03 , 0.03 , 0.01 , 0.01 , 0.02, 0.02 ]

After reaching 16000 points the frequencies and damping ratios become finally:

f final = [0.5380 , 2.2141 , 4.0478 , 5.6599 , 6.9051 , 7.6898 ] Hz

ξ final = [0.10 , 0.07 , 0.051 , 0.02 , 0.03, 0.05 ]

Between 4000 and 16000 data point, the natural frequencies and damping ratios are interpolated linearly based on the initial and final state of the modal frequency and damping ratio.

Noise free

Firstly the ability of RSSI-COV to trace slow time-varying frequencies is shown in Figure 6-11, a) for the case of using 2500 points for moving window and b) for 4000 points of moving window. For both cases the number of block rows i is 70, system order is defined as 12. The same as that defined in the previous section, the statistic window length WL_statistic is 400 frequency points, displacing by every 100 points; the restart criterion is 50%. Again, the increase of the moving window length has enhanced the tracking stability. However, similar to that found previously, the 1^st mode seems to be very sensitive and lost the stability in the final points. If the moving window length is increased to 5000 points as shown in Figure 6-11 c) keeping the system order as 12, certain stability problems still remains; however, similar to what is done in section 6.3, if the system order is increased to 18 keeping intact the remaining parameters, stable results can be achieved for the final points.

In fact, this outcome has a similar sense than that obtained in the simulation study

of section 3.3, where off-line identification task is carried out for two closely-spaced frequencies blended with time-varying frequencies: larger the number of block rows i, the time-varying frequencies will be spanned by more orthogonal vectors after being decomposed by SVD, thus, higher system order is required to prevent the loss of certain mode information due to the order consumption in covering the most time-varying component. The fact is ilustrated in Figure 6-12, where the number of block rows has been incremented from 100 to 130, moving window length fixed to 5000 points, and different system orders is considered.

Adding noise correlated with output (SSI assumption violated)

Consider the addition of a noise correlated with output, to test the robustnss of RSSI-COV in the case of slow-varying systems when its assumption were violated.

Firstly, two different number of block rows i is compared, and system order is kept as 12. Frequency tracking results is shown in Figure 6-13. The effect of the assumption violation impacts mostly in the recursive identification of the 1^st mode; again, due to the time-varying characteristic of the system, higher system order is needed for larger number of block rows i, otherwise, the 1^st mode information is not covered. This explains why the outcome in Figure 6-13 a) using fewer block rows is better than b).

Consider the use of the Recursive SSA-SSI-COV algorithm, there are two sets of parameters to be determined, one for recursive SSA (rSSA): moving window length L’

and number of block rows i’; and another for RSSI-COV: L and i. Regarding to the selection of number of block rows i and i’, 100 block rows is a good choice from previous simulation experiences, this is actually a tradeoff between computation effort and accuracy. 5000 points is selected for the moving window length of RSSI-COV, because the larger the window length, more stable the result. To investigate the

influence of the moving window length L’ and the number of principal components to be selected in rSSA step, nine different combinations for model parameters are considered. For these cases, firstly the singular spectrum in rSSA step is shown in Figure 6-14. To determine the best choice of number of SV in rSSA step, implementation algorithm shown in section 4.1.3.1 can be applied. The singular spectrum in RSSI-COV step is shown in Figure 6-15; they are computed from the first 5000 data points which corresponds to the moving window length of RSSI-COV.

From the 9 combinations shown in Figure 6-15, the best choice of model parameters for rSSA seems to be a window length L’ of 1000 points and 20 SV, for the reason that this choice leads to the best separation of the 6 pairs of singular values (corresponding to the system order). This can be verified by selecting various cases for comparison as shown in Table 6-2, the frequency tracking results is shown in Figure 6-16.

Although the 1^st modal frequency were miss for the last points, in Figure 6-16 d), if the system order is increased to 16, the final segment stabilizes as shown in Figure 6-17 a). This is due to the time-varying property of the signal and the reason explained above. However, if only RSSI-COV is used, even increasing the system order to 30, this final segment is still a little bit unstable. This is shown in Figure 6-17 b).

After the slow time-varying simulation example is studied, several conclusions about the effect of time-varying system in the model parameters are obtained:

The required number of orthogonal components (i.e., system order to be chosen ) to span the system subspace is more than the actual system order, when it comes to the signal processing of time-varying systems.

Although a better noise elimination can be achieved using larger number of block columns i, however, the required number of orthogonal components (i.e., system order to be chosen) to span the system subspace is also greater. If the selected system order is lower than the required, spurious modes appear instead of the correct frequency.

The larger the moving window length L in RSSI-COV step, more stable the tracking results. On the contrary, the shorter is the moving window length L’ in rSSA, better the filtering result and, thus, the identification quality.

Chapter 7 Application of recursive SSI algorithms in damage detection and early warning

The RSSI-COV and rSSA-SSI-COV algorithms has been validated through numerical examples, sensitivity study of subspace model parameters and their selection criteria were also conducted in chapter 6. In the following, these two algorithms will be applied to trace modal parameters for several cases: firstly the shaking table test of a 3-story steel structure with instantaneous stiffness reduction, followed by the shaking table experiment of a 2-bay reinforced-concrete frame, both subjected to earthquake ground motions, and finally to the bridge scouring experiment.

7.1 Application: shaking table test of a 3-story steel structure with instantaneous stiffness reduction

A 3-story full-scale steel frame is designed and constructed at the National Center for Research on Earthquake Engineering (NCREE) in Taipei, Taiwan, in 2007/01/17.

As shown in Figure 7-1, the structure consists of a single bay with a 3m by 2m floor area and 3m tall stories. The structure is constructed using H150x150x7x10 steel I-beam elements with each beam-column joint designed as a bolted connection. To apply additional dead load upon each floor, concrete blocks are fastened to the floor

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