Abstract
This study proposes a novel and efficient approach for locating the floors of a building whose stiffnesses change after being subject to a strong earthquake. The floors that may be damaged are determined by comparing the unitary stiffness matrix in different stages in the life cycle of a building. To evaluate the coefficient matrices of a statespace model, the proposed procedure applies a subspace approach in conjunction with the short time Fourier analysis. The dynamic characteristics of a structure are determined from the coefficient matrices. Next, the unitary stiffness matrix is constructed by identifying natural frequencies and mode shapes. The effectiveness of the proposed procedure is verified using the measured earthquake acceleration responses of a threestorey structure that sustains damage on one or two floors and an eightstorey steel frame under a 200 Gal and a 1200 Gal earthquake, such as the Chichi earthquake that shook Taiwan on September 1999. The proposed scheme is compared to mode shape based approaches in identifying damaged floors, and is demonstrated to be superior to both MAC (Modal Assurance Criterion) and COMAC (Coordinate Modal Assurance Criterion).
1. Introduction
Using the identified modal parameters of an existing structure or bridge from measured responses to locate damaged locations is essential in the field of structural health monitoring. In literature, there are vast algorithms for identifying modal parameters of a structure from its dynamic responses. The approaches of system identification are divided into broad categories, parametric and nonparametric identification approaches.
– The parametric approach methods uses a chosen time series model and estimates the model parameters from a regression process of the measured data. Autoregressive (AR), autoregressive exogenous (ARX), autoregressive moving average (ARMA), autoregressive moving average exogenous (ARMAX) methods [14], stochastic subspace approaches [5, 6] and Ibrahim time domain schemes [7] have been frequently used to identify modal parameters of civil structures.
– The nonparametric approach is performed by estimating the dynamic characteristics of a structure without necessarily using a given parameterized model. Correlation and spectral analyses are often used to estimate impulse response function and frequency response functions, respectively, in the nonparametric approach [8]. However, they do not work well for structural systems with heavy damping and strong modal interference.
The method of damage identification based on examining change in measured vibration responses is the vibrationbased damage detection method. Most vibrationbased damage detection methods can be classified by identifying their dynamic properties [9, 10]:
– Natural frequency based methods: applying identified natural frquency from measured structural responses for detcting structural damage locations [1113].
– Mode shape based methods: developing a damage index from identified modal shapes to locate floors with structural damage [1416].
– Stiffness matrix based methods: comparing the difference of the identified stiffness matrix between the damaged and undamaged structure to assess structural damage [1719].
– Flexibility matrix based methods: using identified modal parameters and given mass properties to construct a flexibility matix to identify damaged locations [20, 21].
– Frequency response function based methods: frequency response function is a transfer function in the frequency domain, with dissimilarities between frequency responses functions for damaged and undamaged states used to identiy the damage [22, 23].
This work presents a simple and efficient stiffnessbased approach to determine which floors of a structure exhibit a change in stiffness. First, the statespace model of structure is established and a subspace approach is used to estimate modal parameters. The measured responses are in terms of acceleration or velocity, and the acceleration responses are used in this approach. Next, the unitary stiffness matrix can be constructed via the identified modal parameters. The unitary stiffness matrix in the reference stage (without damage) are then compared with those in the current state (possibly damaged) and, by using the singular value decomposition (SVD), the location of any damaged floors can be easily determined. The statespace models are developed via the short time Fourier transform.
The proposed procedure is validated using the measured earthquake acceleration responses of a threestorey steel frame and an eightstorey steel frame. A series of cases, involving singlesite and dualsite damage, is examined for the threestorey frame. In another, the eightstorey frame was subjected to 200 Gal and 1200 Gal earthquake forces. For comparison, the MAC and COMAC is also adopted for damage detection.
2. Methodology
2.1. Subspace approach in timefrequency domain
The dynamic responses of a linear structure satisfy the equation of motion:
where $\mathbf{M}$, $\mathbf{C}$ and $\mathbf{K}$ are mass, damping and stiffness matrices, respectively; $\ddot{\mathbf{x}}$, $\dot{\mathbf{x}}$ and $\mathbf{x}$ are the acceleration, velocity, and displacement responses vectors of the system, and $\mathbf{f}$ is the input force vectors. Usually, not all degrees of freedoms of system are measured in a field experiment, for reasons of economy. Only some parts of $\ddot{\mathbf{x}}$ or $\dot{\mathbf{x}}$ are measured. Consequently, the measured response vector $\mathbf{y}$, which can be velocity or acceleration responses, can be described by the statespace model. Generally, the statespace model considered in the following will be presented as [5]:
where $\mathbf{A}$, $\mathbf{B}$, $\mathbf{E}$ and $\mathbf{D}$ are system matrices that relate to $\mathbf{M}$, $\mathbf{C}$ and $\mathbf{K}$; $\mathbf{z}$ is statevariable $\mathbf{z}={\left[\begin{array}{ll}{\mathbf{x}}^{T}& {\dot{\mathbf{x}}}^{T}\end{array}\right]}^{T}$.
From Eq. (2), one can construct:
In the frequency domain, one can arbitrarily remove uninteresting frequencybands to reduce the effects of noise in accurately estimating the modal parameters of a linear time invariant system. Consequently, Fourier transform is further introduced into Eq. (3). Treating the columns of ${\mathbf{z}}_{k}$, ${\mathbf{y}}_{k}$ and ${\mathbf{f}}_{k}$ as vector functions and applying the Fourier transform to Eq. (3) yields:
Via Eq. (4), one can further construct:
where:
where the two sets $\mathcal{F}{\left\{\mathbf{y}\right\}}_{k,s}$ and $\mathcal{F}{\left\{\mathbf{f}\right\}}_{k,s}$ are the shorttime Fourier transform coefficients of ${\mathbf{y}}_{k}$ and ${\mathbf{f}}_{k}$, respectively. From Eq. (5), the following relation can be established:
where:
For further details of subspace system identification approaches, see Huang and Lin [5], and the references therein. As a brief summary, the subspace approach procedure is as follows:
1) Define an orthogonal projection matrix ${\prod}_{f}^{\perp}=\mathbf{I}{\widehat{\mathbf{F}}}_{k,N}^{\mathbf{T}}{\left({\widehat{\mathbf{F}}}_{k,N}{\widehat{\mathbf{F}}}_{k,N}^{\mathbf{T}}\right)}^{1}{\widehat{\mathbf{F}}}_{k,N}$, onto the nullspace of ${\widehat{\mathbf{F}}}_{k,N}$.
2) Introduce instrumental variables $\mathbf{P}={\left[\begin{array}{ll}{\widehat{\mathbf{F}}}_{p,N}^{T}& {\widehat{\mathbf{Y}}}_{p,N}^{T}\end{array}\right]}^{T}$.
3) Calculate weighting matrices ${\mathbf{W}}_{r}=\mathbf{I}$ and ${\mathbf{W}}_{c}={\left(\mathbf{P}{\mathrm{\Pi}}_{f}^{\perp}{\mathbf{P}}^{T}/N\right)}^{1/2}$.
4) Define a matrix $\stackrel{}{\mathbf{H}}={\mathbf{W}}_{r}{\widehat{\mathbf{Y}}}_{k,N}{\prod}_{f}^{\perp}\mathbf{P}{\mathbf{W}}_{c}/N$.
5) Apply singular value decomposition on $\stackrel{}{\mathbf{H}}\cong {\mathbf{Q}}_{\stackrel{}{n}}{\mathbf{D}}_{\stackrel{}{n}}{\mathbf{V}}_{\stackrel{}{n}}^{T}$.
6) One has ${\mathbf{\Gamma}}_{\alpha}={\stackrel{}{\mathbf{\Gamma}}}_{\alpha}{\mathbf{T}}_{\stackrel{}{n}}$, where ${\stackrel{}{\mathbf{\Gamma}}}_{\alpha}={\mathbf{W}}_{r}^{1}{\mathbf{Q}}_{\stackrel{}{n}}{\mathbf{D}}_{\stackrel{}{n}}^{1/2}$; ${\mathbf{T}}_{\stackrel{}{n}}={\mathbf{D}}_{\stackrel{}{n}}^{1/2}{\mathbf{V}}_{\stackrel{}{n}}{\left({\widehat{\mathbf{Z}}}_{k,N}{\prod}_{f}^{\perp}\mathbf{P}{\mathbf{W}}_{c}/N\right)}^{1/2}$.
7) Rewrite Eq. (5) as $\mathcal{F}{\left\{\mathbf{y}\right\}}_{k,s}={\stackrel{}{\mathbf{\Gamma}}}_{\alpha}{\left\{\stackrel{}{\mathbf{z}}\right\}}_{k}+{\mathbf{\Phi}}_{\alpha}\mathcal{F}{\left\{\mathbf{f}\right\}}_{k,s}$, where ${\left\{\stackrel{}{\mathbf{z}}\right\}}_{k}={\mathbf{T}}_{\stackrel{}{n}}\mathcal{F}{\left\{\mathbf{z}\right\}}_{k}$.
8) ${\stackrel{}{\mathbf{\Gamma}}}_{\alpha}$ can be established by ${\stackrel{}{\mathbf{\Gamma}}}_{\alpha}={\left[\begin{array}{llll}{\stackrel{}{\mathbf{E}}}^{T}& {\left(\stackrel{}{\mathbf{E}}\stackrel{}{\mathbf{A}}\right)}^{T}& \cdots & {\left(\stackrel{}{\mathbf{E}}{\stackrel{}{\mathbf{A}}}^{s1}\right)}^{T}\end{array}\right]}^{T}$.
9) Define:
10) Calculate $\stackrel{}{\mathbf{A}}$ by solving linear equation ${\stackrel{}{\mathbf{\Gamma}}}_{\alpha 2}={\stackrel{}{\mathbf{\Gamma}}}_{\alpha 1}\stackrel{}{\mathbf{A}}$.
2.2. Determination of modal parameters
When the equation of motion is expressed in terms of the spacestate variable, it is well known that the dynamic characteristics of the structural system are determined by the eigenvalues and eigenvectors of $\mathbf{A}$ in Eq. (2a). However, $\mathbf{A}$ cannot be determined from the preceding derivation. Thus, to find the modal parameters, $\stackrel{}{\mathbf{A}}$ is employed:
where ${\lambda}_{j}$ and ${\stackrel{}{\lambda}}_{j}$ are the $j$th eigenvalue of $\mathbf{A}$ and $\stackrel{}{\mathbf{A}}$, respectively, while ${\mathbf{\varphi}}_{j}$ and ${\stackrel{}{\mathbf{\varphi}}}_{j}$ are the corresponding eigenvectors. Since ${\mathbf{T}}_{\stackrel{}{n}}$ remains unknown, the eigenvectors of $\mathbf{A}$ from Eq. (9) cannot be evaluated. Nevertheless, usually, one only needs to determine the modal shape corresponding to the observed degrees of freedom, ${\mathbf{\varphi}}_{j,y}$. From the procedure in Section 2.1, the following relationship can be determined:
The eigenvalues are complex numbers. Let:
From Eq. (9), one has:
Then, the pseudoundamped circular natural frequency and the modal damping ratio for the system are:
In a proportionally damped system, this is equivalent to its undamped circular natural frequency.
2.3. Damage assessment
Consider the response of a structure described by equation of motions such as Eq. (1). Assuming proportional damping, the orthogonality property of the mode shapes with respect to the mass and stiffness matrix leads to [24]:
in which ${\mathbf{M}}_{D}$ and ${\mathbf{K}}_{D}$ are the diagonal modal mass matrix and diagonal modal stiffness matrix, respectively. The equation used to express the square of modal frequencies as a matrix is:
Here, the unitary stiffness matrix were defined as:
The unitary stiffness matrix at the sensor locations constructed from measured data before and after damage are denoted as ${\widehat{\mathbf{K}}}_{u}$ and ${\widehat{\mathbf{K}}}_{d}$, respectively. Applying singular value decomposition (SVD) to ${\widehat{\mathbf{K}}}_{u}{\widehat{\mathbf{K}}}_{d}$ yields:
The singular values $\mathbf{D}$ are the difference in potential energy between ${\widehat{\mathbf{K}}}_{u}$ and ${\widehat{\mathbf{K}}}_{d}$ under each state $\mathbf{V}$ and $\mathbf{W}$. The location of maximum relative displacement in each $\mathbf{V}$ and $\mathbf{W}$ that lead to the difference of potential energy implies the maximum difference between the damaged and undamaged systems. Therefore, the possible damaged floors can be identified by nonzero singular values and corresponding states $\mathbf{V}$ and $\mathbf{W}$. The damage index can be obtained from the following steps:
1) Find the location (${l}_{i}$) with maximum relative displacement in ${\mathbf{V}}_{i}$ or ${\mathbf{W}}_{i}$_{.}
2) Define damage index ($DI$):
3. Applications
To demonstrate the feasibility of the proposed procedure in processing data in real applications, the procedures were applied to identify the modal parameters of structures from their measured response data obtained from shaking table tests.
3.1. Threestorey steel frame
Shaking table tests are often employed in a laboratory to investigate the behaviors of structures under earthquake conditions. To generate a set of earthquake response data for a benchmark model of a 1/8scaled threestorey steel structure, CVNCTU (Department of Civil Engineering, National Chiao Tung University) performed a series of shaking table tests on three frames (Fig. 1). One steel frame was represented as “std”, where each floor weighed about 3.8 kg and each column had a crosssectional area of 80 mm^{2} and was 440 mm in height. The second frame and third frame, which were structurally irregular with respect to stiffness and denoted “R_k1”, were identical to “std” except that a stiffening central column was removed from the first story. The third frame, denoted as “R_k12”, was identical to “R_k1” except that the stiffening central column was removed from the second story. These frames were subjected to base excitations specified by records of real earthquakes such as the 1995 Kobe earthquake. The responses were sampled at 100 Hz. Fig. 1 depicts the acceleration responses of each floors in the weak axis direction for frame “std”, subject to 20 % of the strength of the Kobe earthquake.
Fig. 1A photo of frame “std” on shaking table and its responses
Table 1Comparison of identified modal parameters for different frames
Frame  Mode  ${f}_{n}$ (Hz)  $\xi $(%)  Mode shapes  
std  1  1.65  1.81  1.000  0.803  0.454 
2  4.83  1.12  –0.890  0.395  1.000  
3  7.22  0.94  –0.513  1.000  –0.789  
R_k1  1  1.44  1.16  1.000  0.843  0.563 
2  4.53  0.84  –0.912  0.264  1.000  
3  7.16  0.57  –0.525  1.000  –0.721  
R_k12  1  1.36  0.83  1.000  0.865  0.484 
2  4.38  0.60  –0.654  0.115  1.000  
3  6.62  0.91  –0.704  1.000  –0.530 
Table 1 summarizes the identified dynamic characteristics of the three frames obtained using the acceleration responses of all floors and the input excitation in the weak axis direction. The results were obtained from the responses of all floors at $t=$ 525 seconds. The natural frequencies of “R_k1” and “R_k12” differ significantly from those of “std”. As expected, frame “R_k1” has lower natural frequencies than frame “std”, whereas frame “R_k12” has the lowest natural frequencies in all frames. To evaluate the correlation of mode shapes obtained from different methods and frames, the index of modal assurance criterion (MAC) [15] and coordinate modal assurance criterion (COMAC) [16] was computed to indicate the correlation between any two mode shapes of interest, and which are defined as:
To compare the identified mode shapes of undamaged frames in “std” with damaged frames in “R_k1” and “R_k12”, ${\phi}_{ij,s}$ and ${\phi}_{ij,R}$ respesent the deformation of the $i$th degree of freedom in the $j\mathrm{t}\mathrm{h}$ identified mode shapes, respectively.
Fig. 2 compares the results obtained using the various methods. The frame “std” is treated as a reference structure. Since the MAC and COMAC values are always between zero and unity, the maximum values of the results from proposed approach were normalized to unity for comparison. The results denoted by “present” in Fig. 2 normalize the difference of potential energy between damaged and undamaged systems. Fig. 2 clearly demonstrates that the proposed approach is superior to MAC and COMAC in identifying the floor of the structure whose properties differ from those of the corresponding floor in the reference structure.
Fig. 2Normalized index
a) R_k1
b) R_k12
3.2. Eightstorey steel frame
To generate a set of earthquake response data for a benchmark model of an eightstorey steel frame, the National Center for Research on Earthquake Engineering (NCREE) in Taiwan performed a series of shaking table tests on this model (Fig. 3). The eightstorey steel frame analyzed in this test was 1.8 m long, 1.2 m wide and 8.5 m high. Lead plates were piled on each floor such that the total mass of the steel frame was approximately 4.519 tons. In the shaking table tests, two accelerometers and two linear displacement transducers were installed in the longspan direction on the two edges of each floor to measure acceleration and displacement responses of the floor, respectively. The frame was subjected to base excitations equivalent to the Chichi earthquake at different reduced levels, including 200 Gal and 1200 Gal of base excitation levels. Data were recorded at a sampling rate of 200 Hz. The acceleration responses of the base and all floors at $t=$ 535 seconds were used in calculating modal parameters of the frame. Fig. 4 shows the base excitations and the responses of all floors in the longspan direction of the frame, subjected to 1200 Gal, equivalent to the Chihi earthquake. The strain records shown in Fig. 5 for different levels of earthquake input confirm that the columns of the first to third floors yielded when the frame was subjected to 1200 Gal base excitation.
Fig. 3A photo of frame on shaking table
Fig. 6 shows the results of the identified modal parameters of eight modes for the steel frame under 200 Gal and 1200 Gal of the Chichi earthquake force. Since the frame is symmetric with eight storeys, it is expected to have eight modes in the longspan direction of the frame. Damage that is caused by a strong earthquake will change the natural frequencies of the structure.
Fig. 7 compares the results obtained using the various methods. When the frame was subjected to 200 Gal of the Chichi earthquake, no nonlinear behaviors were observed, and it was treated as a reference structure. Again, the maximum values of the results from proposed approach were normalized to unity for comparison. The results denoted by “present” in Fig. 7 normalize the difference of potential energy between damaged and undamaged systems. Fig. 7 clearly indicates that the stiffness of the first to third floors significantly changed when the frame was under 1200 Gal of the Chichi earthquake, and these floors were somewhat damaged. Notably, Fig. 5 shows the measurements obtained by strain gauges installed on the columns of the first to third stories, which also indicated that the columns yielded during such base excitation. The success of the proposed procedure, as revealed in capturing the differences in properties of experimental frames, demonstrates the practical applicability of this procedure on an actual building with symmetry.
Fig. 4Responses of frame under 1,200 Gal of the Chichi earthquake
Fig. 5Strain responses at first, second and third floor columns
4. Conclusions
The work developed a simple and efficient approach for identifying the floors of a structure that sustain earthquake damage, based on the fact that such damage to a particular floor changes the stiffness matrix of the structure. This procedure was established through a statespace model to describe the measured responses. The coefficient matrices related to the dynamic characteristics were evaluated by a subspace method in conjunction with the concept of an instrumental variable. Next, the modal parameters were evaluated from the eigenvalues and eigenvectors of a coefficient matrix. The unitary stiffness matrix was constructed from identified modal parameters. The coefficient matrix was established from the acceleration or velocity responses of structure using the short time Fourier transform. Comparing the unitary stiffness matrix in the current state with those in the undamaged state enables damaged floors to be accurately located.
Fig. 6Identified natural frequency and model shape of eightstorey steel frame
Fig. 7Normalized index
To demonstrate the feasibility of proposed approach for actual applications, this procedure was applied to process in situ earthquake response measurements of a threestorey and an eightstorey steel frame. The threestorey steel frame with singlesite or dualsite damage was considered, and the eightstorey was subjected to 200 Gal and 1200 Gal earthquake forces. The proposed approach was validated by successfully identifying damaged floors through processing measured responses. Comparing the results obtained by the proposed approach with MAC and COMAC indexes revealed that the present approach is substantially superior to both in identifying damaged floors. The success of the proposed approach when applied to the experimental responses demonstrates its practical applicability to an actual symmetrical building.
References

Wang Z. N., Fang T. A time domain method for identifying modal parameters. Transaction of ASME, Journal of Applied Mechanics, Vol. 53, Issue 1, 1986, p. 2832.

Loh C. H., Wu T. S. Identification of FeiTsui arch dam from both ambient and seismic response data. Soil Dynamics and Earthquake Engineering, Vol. 15, Issue 7, 1996, p. 465483.

Huang C. S. Structural identification from ambient vibration measurement using the multivariate AR model. Journal of Sound and Vibration, Vol. 241, Issue 1, 2001, p. 337359.

Lardies J. Modal parameter identification based on ARMAV and statespace approaches. Archive of Applied Mechanics, Vol. 80, Issue 1, 2010, p. 335352.

Huang C. S., Lin H. L. Modal identification of structures from ambient vibration, free vibration, and seismic response data via a subspace approach. Earthquake Engineering and Structural Dynamics, Vol. 30, Issue 12, 2001, p. 18571878.

Ali M. R., Okabayashi T. System identification of highway bridges from ambient vibration using subspace stochastic realization theories. Earthquake and Structures, Vol. 2, Issue 2, 2011, p. 189206.

Ibrahim S. R., Mikulcik E. C. The experimental determination of vibration parameters from time responses. Shock and Vibration Bulletin, Vol. 46, Issue 5, 1976, p. 187196.

Bendat J. S., Piersol A. G. Engineering Applications of Correlation and Spectral Analysis, 2nd Ed. John Wiley and Sons, New York, U.S., 1993.

Carden E. P., Fanning P. Vibration based conditioning monitoring: a review. Structural Health Monitoring, Vol. 3, Issue 4, 2004, p. 355377.

Montalvao D., Maia N., Ribeiro A. A review of vibrationbased structural health monitoring with special emphasis on composite materials. Shock and Vibration Digest, Vol. 38, Issue 4, 2006, p. 295324.

Adams R. D., Cawley P., Pye C. J., Stone B. J. A vibration technique for nondestructively assessing the integrity of structures. Journal of Mechanical Engineering Science, Vol. 20, Issue 2, 1978, p. 93100.

Cawley P., Adams R. D. The location of defects in structures from measurements of natural frequencies. Journal of Strain Analysis, Vol. 14, Issue 2, 1979, p. 4957.

Salawu O. S. Detection of structural damage through changes in frequency: a review. Engineering Structures, Vol. 19, Issue 9, 1997, p. 718723.

Alampalli S., Fu G., Dillon E. W. Signal versus noise in damage detection by experimental modal analysis. Journal of Structural Engineering, Vol. 123, Issue 2, 1997, p. 237245.

Allemang R. L., Brown D. L. A correlation coefficient for modal vector analysis. Proceedings of the 1st International Modal Analysis Conference, 1983, p. 110116.

Lieven N. A. J., Ewins D. J. Spatial correlation of modespaces: the coordinate modal assurance criterion (COMAC). Proceedings of the 6th International Modal Analysis Conference, Kissimmee, Florida, USA, 1988, p. 10631070.

He J., Ewins D. J. Analytical stiffness matrix correction using measured vibration modes. Modal Analysis: The International Journal of Analytical and Experimental Modal Analysis, Vol. 1, Issue 1, 1986, p. 914.

Salawu O. S., Williams C. Structural damage detection using experimental modal analysis – a comparison of some methods. Proceedings of 11th International Modal Analysis Conference, 1993, p. 254260.

Peterson L. D., Alvin K. F., Doebling S. W., Park K. C. Damage detection using experimentally measured mass and stiffness matrices. Proceedings of 34th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, 1993, p. 15181528.

Pandey A. K., Biswas M. Damage detection in structures using changes in flexibility. Journal Sound and Vibration, Vol. 169, Issue 1, 1994, p. 317.

Bernal D. Load vectors for damage localization. Journal of Engineering Mechanic, Vol. 128, Issue 1, 2002, p. 714.

Sampaio R. P. C., Maia N. M. M., Silva J. M. M. Damage detection using the frequency response function curvature method. Journal of Sound and Vibration, Vol. 226, Issue 5, 1999, p. 10291042.

Lee U., Shin J. A frequency response functionbased structural damage identification method. Computers and Structures, Vol. 80, Issue 2, 2002, p. 117132.

Tedesco J. W., McDougal W. G., Ross C. A. Structural Dynamics: Theory and Applications. AddisonWesley, Menlo Park, Calif, 1999.
About this article
The authors would like to thank the National Science Council of the Republic of China, Taiwan, for financially supporting this research under Contract No. MOST 1042622M492001CC2. The appreciation is also extended to the National Center for Research on Earthquake Engineering for providing shaking table test data.