# Complete Inverse Method Using Ant Colony Optimization Algorithm for Structural Parameters and Excitation Identification from Output Only Measurements

Authors: Jun Chen; Xin Chen; Wei Liu
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

In vibration-based structural health monitoring of existing large civil structures, it is difficult, sometimes even impossible, to measure the actual excitation applied to structures. Therefore, an identification method using output-only measurements is crucial for the practical application of structural health monitoring. This paper integrates the ant colony optimization (ACO) algorithm into the framework of the complete inverse method to simultaneously identify unknown structural parameters and input time history using output-only measurements. The complete inverse method, which was previously suggested by the authors, converts physical or spatial information of the unknown input into the objective function of an optimization problem that can be solved by the ACO algorithm. ACO is a newly developed swarm computation method that has a very good performance in solving complex global continuous optimization problems. The principles and implementation procedure of the ACO algorithm are first introduced followed by an introduction of the framework of the complete inverse method. Construction of the objective function is then described in detail with an emphasis on the common situation wherein a limited number of actuators are installed on some key locations of the structure. Applicability and feasibility of the proposed method were validated by numerical examples and experimental results from a three-story building model.


---

## Body

## 1. Introduction

Structural health monitoring (SHM) has remained an active research topic in structural engineering since the 1970s. SHM identifies the occurrence, location, and severity of structural damage via significant adverse changes of structural parameters or properties. Thus, the core part of an SHM system is the algorithm used to accurately identify the structural parameters. A number of methods are available nowadays for us to accomplish this identification task. Housner [1] presented an extensive summary on the state-of-the-art methods in the vibration control and health monitoring of civil engineering structures. Zou et al. [2] summarized the methods of vibration-based damage detection and health monitoring for composite structures. Very recent reviews on identification methods used in SHM can be found in Ou and Li [3] and Fan and Qiao [4] among several others. Traditional identification algorithms are generally based on the assumption that the system’s input (excitations) and output (responses) are completely known (measured). However, in vibration-based structural health monitoring of existing large civil structures, it is difficult, sometimes even impossible, to measure the actual excitation applied to the structures. Therefore, output-only structural parameter identification methods are of great significance for practical application of SHM.In contrast to a number of publications on structural parameter identification with complete output and input information, there is a paucity of publications addressing the identification method using responses only. Since the input/excitation is unknown, assumptions on the input are necessary in order to make the traditional identification method applicable. The most commonly adopted assumption is treating the input as a white-noise process, whose power spectrum is theoretically known. The classical Ibrahim method plus random decrement technique, for instance, is based on this assumption. This assumption, however, is not tenable for excitation like earthquakes or strong winds that are in fact a nonstationary random process. Moreover, the time history of the real input cannot be directly identified. The second common way is assuming that the input has some special features. For instance, Toki et al. [5] assumed that the coda of the measured structural responses during an earthquake could be treated as free vibration responses from which the structural parameters could be easily identified. Wang and Haldar [6] identified the unknown earthquake input and structural parameters using output-only measurements through a recursive identification procedure consisting mainly of three steps. Step 1: The unknown structural parameters were identified by assuming that the input excitations were zero at the beginning building up to four time instants, say t0 to t4. Step 2: The input excitation was conjectured by the measured structural responses and the parameters estimated in the first step. Step 3: Replace excitation at t0 to t4 in step 1 with new values obtained at step 2, and then repeat step 1 and 2 until the identified input excitation at t0 to t4 converged within a preset tolerance limit. For unknown wind load, Law et al. [7] assumed the wind load model was known even though the time histories of the wind load had not been recorded. Based on this assumption, the time history of wind load and the structural parameters could be identified using the response measurements. More recently, Yang et al. [8] suggested a recursive least squares estimation with unknown inputs to identify the stiffness, damping, and other nonlinear parameters at element level. The locations of the unknown excitations are assumed known in their approach. The ASCE structural damage benchmark structure was used to show the feasibility of the method. Yang and Huang [9] further extended this method to situations where external excitations and some acceleration responses are not measured.To deal with the unknown input identification problem, we have proposed the concept of complete inverse problem, which means identification of structural parameters and the input’s time histories simultaneously from the output-only measurements as per Chen and Li [10, 11], Li and Chen [12, 13], and Zhao et al. [14, 15]. Within the framework of complete inverse problem, we have suggested a series of identification methods named as complete inverse methods (CIM) addressing different types of unknown excitations. The core idea of CIM is to convert any additional information of the excitations (whose time histories are unknown) into mathematical constraint conditions that can be further integrated into an iteration identification procedure based on least square technique. For instance, when the locations of the excitations are known, the spatial information of the external excitation can be used in the parameter identification method [10, 11]. For ground motion excitation like that resulting from an earthquake, its mechanical features indicate that the inertial force proportional to the mass can be introduced into the interaction procedure as a mathematical constraint to ensure a stable and unique solution [12, 13]. For a proportional-type excitation, like wind loads, the ratio of forces at different structural heights can be used as a mathematical condition to identify the structural parameters and inputs. For shear-type building under earthquake excitation and with limited response measurements, we proposed a hybrid identification method where the unknown structural parameters for the first floor are identified using measured modal shapes, and parameters of all the other stories are identified using measured acceleration responses [14, 15]. We have also provided strict mathematical proofs from Chen and Li [11] and Li and Chen [13] to demonstrate the unconditional convergence features of the proposed complete identification methods.Despite CIM’s success, it does have several limitations. This paper thus tries to improve the CIM in two directions. The first arises from replacing the least square method with another efficient and more robust optimization method. The recently emerged ant colony optimization algorithm (ACO) has been adopted in this paper to identify the structural parameter. The second improvement aims to validate the proposed method by experimental examples. The feasibility and effectiveness of all the aforementioned identification methods have already been demonstrated by different kinds of numerical examples. However, few experimental investigations have been conducted to assess the practical application of these methods. The effect of measurement noise and modeling error cannot be fully investigated in a numerical simulation. However, time-domain identification methods are known to be sensitive to measurement noise.To this end, the paper first presents the principle of ACO for solving both discrete optimization problems and continuous optimization problems as well. Then, the ACO is integrated into the CIM methods by constructing the objective function according to the type of excitation. After that, the proposed method is applied to numerical model and an experimental model. The results show that the CIM+ACO algorithm performs very well for a noise-free and noise-polluted case and has good identification accuracy in parameters and inputs.


## 2. Ant Colony Optimization

Inspired by the ants’ foraging behavior, Dorigo [16] proposed the ant colony optimization algorithm (ACO). ACO emulates the behavior of a group of ants in searching food from their nest to the food sources. Every ant leaves an amount of pheromone on the path that it passes and chooses the path with more pheromone left on it from the previous ants. Then, after more and more ants pass, the path having the maximum pheromone will be the best (shortest) way from the nest to the food source. ACO algorithms have already been successfully applied for solving combinatorial optimization problems, including the traveling salesman problem (TSP) [17], the routing problem in a computer network [18], the quadratic assignment problem (QAP) [19], and structural health monitoring problems [20, 21]. ACO algorithms for continuous optimization have been proposed in the literature [22–24]. All the above work has proven ACO to be an efficient and versatile tool for solving various continuous optimization problems. The ACO algorithm has already been well established so far. Principal and application procedures of ACO are briefly summarized in Section 2.1.


### 2.1. Ant System Model [22]

As mentioned earlier, ACO emulates the behavior of a group of ants in searching for food. Modeling of the ant system and the pheromone left by ants on the path are of great importance for application of ACO. When ACO is used for discrete optimization problems, ants construct solutions incrementally. That is, each ant starts with an empty solutionS0 and a component of the solution is added at each construction step. If Ci denotes all the available solutions at step i, we need to choose one best solution from Ci and add the solution to previous solution Si-1 leading to the new solution as Si. The definition of the available solution componentdepends on the problem tackled. For instance, in the popular travelling salesman problem (TSP), a component of the solution is a city that is added to a tour. The solution components may be defined differently for other problems.A probability-based strategy is adopted to choose the best solution components fromCi for constructing the current partial solution Si. This decision is usually influenced by amount of pheromone τ associated with available choices and by heuristic information about the problem. To avoid the loss of generality, in the following problem, no heuristic information is used. Assuming that the partial solution constructed is Si so far, the probability p is usedto choose a solution component c∈Ci at step i in iteration t as calculated by the following;
(1)pSic(t)=τSic(t)∑j∈CiτSij(t),
where τ is the pheromone, the subscript S denotes the solution domain, and i means the ith iteration step. Hence, in case of discrete optimization problems, the ants make a probabilistic decision according to some discrete probability distribution at each construction step. For continuous optimization problems, the domain changes from discrete to continuous. The logical adaptation also would be moved from using the discrete probability distribution to a continuous one—the probability density function (PDF). Instead of choosinga component c∈Ci at step i, the ants would generate a random number according to a certain PDF P(x).


### 2.2. ACO for Continuous Domain:AC⁡OR [22, 23]

The original ACO algorithm applies only to discrete domain and cannot be directly introduced into continuous domain optimization problems. Structural identification in nature is an optimization problem that aims to find the best parameters for a given objective function defined in a continuous domain. The ACO for continuous domain is, therefore, necessary. We adopted the method suggested by Socha [22] and Socha and Dorigo [23], which is an ACO for continuous domain (ACOR) based on Gaussian probability density functions. Application of ACOR is briefly summarized here for the purpose of completion. More technical details of ACOR can be found in Socha [22] and Socha and Dorigo [23].A Gaussian kernel PDF is used inACOR to account for the multiple-peaks (multiple optimization solutions) domain. Suppose the optimization problem on hand is i=1,2,…,n dimensions; for the ith dimension, the Gaussian kernel PDF Gi(x) is defined as a weighted sum of several one-dimensional Gaussian functions as follows:
(2)Gi(x)=∑l=1kωlgli(x)=∑l=1kωl1σli2πexp⁡(-(x-μli)22(σli)2),
where, i=1,2,…,N is the number of dimensions of the problem; l=1,2,…,k is the number of single Gaussian functions constituting the Gaussian kernel PDF; ωl is the weight associated with the lth individual Gaussian function; μli and σli are the mean and standard deviation for the lth function in the ith dimension. They can also be expressed in vector form as {μi} and {σi} whose cardinality is equal to k.For anN-dimension problem, an ant constructs a solution in N steps. At each step, an ant gets a value for the unknown variable xi. For each solution sl, ACOR will store the current results of N as unknown variables and the value of the objective function f(sl) in an archive T. All the solutions are ordered in T according to their qualities, that is, f(s1)≤f(s2)≤⋯≤f(sl)≤⋯≤f(sk). Each solution has an associated weight ω that is proportional to the solution quality, ω1≥ω2≥⋯≥ωl≥⋯≥ωk. The PDF Gi is constructed using only the ith coordinates of all k solutions from the archive.For each dimensioni=1,2,…,N of the problem, there is a different Gaussian kernel PDF Gi(x) defined. For each Gi(x), the values of the ith variable of all the solutions in the archive T become the elements of the vector μi:
(3)μi={μ1i,…,μki}={s1i,…,ski}.Before each solution is added to the archiveT, it must be evaluated and ranked. The solutions in the archive are sorted according to their rank—that is, solution sl has rank l. Better solutions will have a higher weight. The weight ωl of the solution sl is calculated by the following:
(4)ωl=1q2πexp⁡(-(l-1)22q2).The weight is for the Gaussian function with argumentl, mean 1.0, and standard deviation q. In the ACOR, q is in fact a parameter to balance between diversification and intensification. When q is small, the best-ranked solutions are strongly preferred, and when it is large, the probability becomes uniform.In practice, generating the Gaussian kernel PDF is accomplished as follows. First, choose one of the individual Gaussian functions that compose the Gaussian kernel with probabilitypl given by (5). Then generate the chosen Gaussian function
(5)pl=ωl∑r=1kωr.We consider the chosen Gaussian function given in (5) with μli a standard solution for other ant solutions to explore. To establish the value of the standard deviation σli at step i, we calculate the average distance from μli to all solutions in the archive T, and multiply it by the parameter ξ:
(6)σli=ξ∑e=1k|sei-sli|k-1,
where the parameter ξ>0 has a similar effect to the pheromone evaporation rate in traditional ACO. The higher the value of ξ, the lower the convergence speed of the algorithm.At the start of the algorithm, the solution archiveT is initialized generating k solutions by uniform random sampling. Pheromone update is accomplished by adding the set of newly generated solutions to the solution archive T and then removing the same number of the worst solutions so that the total size of the archive does not change. This process ensures that only the best solutions are kept in the archive, so that they effectively guide the ants in the search process. When the algorithm is completed, the solutions are ordered in the archive according to their quality, and the best solution is s1={s11,s12,…,s1l,…,s1N}.


## 2.1. Ant System Model [22]

As mentioned earlier, ACO emulates the behavior of a group of ants in searching for food. Modeling of the ant system and the pheromone left by ants on the path are of great importance for application of ACO. When ACO is used for discrete optimization problems, ants construct solutions incrementally. That is, each ant starts with an empty solutionS0 and a component of the solution is added at each construction step. If Ci denotes all the available solutions at step i, we need to choose one best solution from Ci and add the solution to previous solution Si-1 leading to the new solution as Si. The definition of the available solution componentdepends on the problem tackled. For instance, in the popular travelling salesman problem (TSP), a component of the solution is a city that is added to a tour. The solution components may be defined differently for other problems.A probability-based strategy is adopted to choose the best solution components fromCi for constructing the current partial solution Si. This decision is usually influenced by amount of pheromone τ associated with available choices and by heuristic information about the problem. To avoid the loss of generality, in the following problem, no heuristic information is used. Assuming that the partial solution constructed is Si so far, the probability p is usedto choose a solution component c∈Ci at step i in iteration t as calculated by the following;
(1)pSic(t)=τSic(t)∑j∈CiτSij(t),
where τ is the pheromone, the subscript S denotes the solution domain, and i means the ith iteration step. Hence, in case of discrete optimization problems, the ants make a probabilistic decision according to some discrete probability distribution at each construction step. For continuous optimization problems, the domain changes from discrete to continuous. The logical adaptation also would be moved from using the discrete probability distribution to a continuous one—the probability density function (PDF). Instead of choosinga component c∈Ci at step i, the ants would generate a random number according to a certain PDF P(x).


## 2.2. ACO for Continuous Domain:AC⁡OR [22, 23]

The original ACO algorithm applies only to discrete domain and cannot be directly introduced into continuous domain optimization problems. Structural identification in nature is an optimization problem that aims to find the best parameters for a given objective function defined in a continuous domain. The ACO for continuous domain is, therefore, necessary. We adopted the method suggested by Socha [22] and Socha and Dorigo [23], which is an ACO for continuous domain (ACOR) based on Gaussian probability density functions. Application of ACOR is briefly summarized here for the purpose of completion. More technical details of ACOR can be found in Socha [22] and Socha and Dorigo [23].A Gaussian kernel PDF is used inACOR to account for the multiple-peaks (multiple optimization solutions) domain. Suppose the optimization problem on hand is i=1,2,…,n dimensions; for the ith dimension, the Gaussian kernel PDF Gi(x) is defined as a weighted sum of several one-dimensional Gaussian functions as follows:
(2)Gi(x)=∑l=1kωlgli(x)=∑l=1kωl1σli2πexp⁡(-(x-μli)22(σli)2),
where, i=1,2,…,N is the number of dimensions of the problem; l=1,2,…,k is the number of single Gaussian functions constituting the Gaussian kernel PDF; ωl is the weight associated with the lth individual Gaussian function; μli and σli are the mean and standard deviation for the lth function in the ith dimension. They can also be expressed in vector form as {μi} and {σi} whose cardinality is equal to k.For anN-dimension problem, an ant constructs a solution in N steps. At each step, an ant gets a value for the unknown variable xi. For each solution sl, ACOR will store the current results of N as unknown variables and the value of the objective function f(sl) in an archive T. All the solutions are ordered in T according to their qualities, that is, f(s1)≤f(s2)≤⋯≤f(sl)≤⋯≤f(sk). Each solution has an associated weight ω that is proportional to the solution quality, ω1≥ω2≥⋯≥ωl≥⋯≥ωk. The PDF Gi is constructed using only the ith coordinates of all k solutions from the archive.For each dimensioni=1,2,…,N of the problem, there is a different Gaussian kernel PDF Gi(x) defined. For each Gi(x), the values of the ith variable of all the solutions in the archive T become the elements of the vector μi:
(3)μi={μ1i,…,μki}={s1i,…,ski}.Before each solution is added to the archiveT, it must be evaluated and ranked. The solutions in the archive are sorted according to their rank—that is, solution sl has rank l. Better solutions will have a higher weight. The weight ωl of the solution sl is calculated by the following:
(4)ωl=1q2πexp⁡(-(l-1)22q2).The weight is for the Gaussian function with argumentl, mean 1.0, and standard deviation q. In the ACOR, q is in fact a parameter to balance between diversification and intensification. When q is small, the best-ranked solutions are strongly preferred, and when it is large, the probability becomes uniform.In practice, generating the Gaussian kernel PDF is accomplished as follows. First, choose one of the individual Gaussian functions that compose the Gaussian kernel with probabilitypl given by (5). Then generate the chosen Gaussian function
(5)pl=ωl∑r=1kωr.We consider the chosen Gaussian function given in (5) with μli a standard solution for other ant solutions to explore. To establish the value of the standard deviation σli at step i, we calculate the average distance from μli to all solutions in the archive T, and multiply it by the parameter ξ:
(6)σli=ξ∑e=1k|sei-sli|k-1,
where the parameter ξ>0 has a similar effect to the pheromone evaporation rate in traditional ACO. The higher the value of ξ, the lower the convergence speed of the algorithm.At the start of the algorithm, the solution archiveT is initialized generating k solutions by uniform random sampling. Pheromone update is accomplished by adding the set of newly generated solutions to the solution archive T and then removing the same number of the worst solutions so that the total size of the archive does not change. This process ensures that only the best solutions are kept in the archive, so that they effectively guide the ants in the search process. When the algorithm is completed, the solutions are ordered in the archive according to their quality, and the best solution is s1={s11,s12,…,s1l,…,s1N}.


## 3. Integration ofAC⁡OR with CIM

### 3.1. Introduction of CIM

The equation of motion of an-DOFs system can be expressed as follows:
(7)MX¨(t)+CX˙(t)+KX(t)=F(t),
where M, C, and K represent, respectively, mass, damping, and stiffness matrices of the structures, X(t) and X˙(t), while X¨(t) represents, respectively, the displacements, velocities, and accelerometers response vector of the structure; F(t) is the external excitation on the structures. Equation (7) can be rewritten as (8), which can be further rearranged as (9) at time instant ti [10]:
(8)CX˙(t)+KX(t)=F(t)-MX¨(t),(9)H(ti)θ=Y(ti),
where vector θ contains all the unknown parameters to be identified; matrix H consists of the measured structural responses; and vector Y is the system input. To assemble (9) at all the sampling time instants ti, i=1,2,…,Q together, will give
(10)HQθ=YQ.The components of matrixH and Y depend on the type of structure. For a shear-type structure, the expressions for H, Y, and θ are(11)θ=[c1,k1,c2,k2,…,cn,kn]T,(12)HQ=[H(t1),H(t2),…,H(tQ)]T,(13)H(ti)=[x˙1(ti)x1(ti)x˙1(ti)-x˙2(ti)x1(ti)-x2(ti)⋯000x˙2(ti)-x˙1(ti)x2(ti)-x1(ti)x˙2(ti)-x˙3(ti)x2(ti)-x3(ti)⋯⋯⋯⋯⋯⋯0⋯⋯⋯x˙n(ti)-x˙n-1(ti)xn(ti)-xn-1(ti)],(14)YQ=[Y(t1),Y(t2),…,Y(tQ)]T,(15)Y(ti)=[f1(ti)-m1x¨1(ti)f2(ti)-m2x¨2(ti)⋯fn(ti)-mnx¨n(ti)]T,wherek1,k2,…,kn, and c1,c2,…,cn represent, respectively, the stiffness and damping coefficients of the structure for each story, and xj(ti) and fj(ti) are the displacement response and external excitation force of the jth DOF (j=1,2,…,n) at the time instant ti.In traditional calculations, the parameters of the structure can be identified from (10) by the least-squares technique as
(16)θ=[HQTHQ]-1HQTYQ.However, (16) cannot be easily solved to determine the stiffness and damping parameters since there are unknown quantities involved in calculating HQ and YQ.As mentioned earlier, we have suggested the complete inverse method to tackle the unknown input situation in structural identification. Following the framework of CIM, Chen and Li [11] suggested the total compensation method. The total compensation method rests on the assumption that the locations of the external forces are known even though their time histories are unknown, and the number of DOF with applied (unknown) forces is less than the number of DOF whose responses are measured. This assumption reflects the situation of a forced vibration survey of structure where a limited number of one or several actuator(s) are installed on key locations of the structure to excite it. In this case, the input excitations in (15) can be further expressed into two parts:
(17)yu(ti)=fu(ti)-mu(ti)x¨u(ti),i=1,2,…,Q,yl(ti)=fl(ti)-ml(ti)x¨l(ti)=-ml(ti)x¨l(ti),yl(ti)=fl(ti)-ml(ti)x¨l(ti)wiii=1,2,…,Q,
where yu(ti) denotes those DOFs with unknown external excitation and yl(ti) stands for those DOFs without applied force, that is, fl(ti)≡0.


### 3.2. Objective Function

In order to identify structural parameters and the input time history from output-only measurements, an objective function is defined as minimizing discrepancies between the forceyl(t) and the calculated force yl⋄(t). The minimization of the objective function is expressed as a bound-constrained nonlinear least squares problem:
(18)min⁡⁡∑ti=t1tQz(ti)=∑ti=t1tQ∑l∈L∥yl(ti)-yl⋄(ti)∥2,i=1,2,…,Q,
where L is the set of l and the identified force yl⋄(ti) is given by
(19)yl⋄(ti)=Hl(ti)θ⋄,i=1,2,…,Q,
where Hl(ti) is lth row of the matrix H(ti), and θ⋄ is the result of identified parameters based on ACOR with the details as
(20)θ⋄=[c1⋄,k1⋄,c2⋄,k2⋄,…,cn⋄,kn⋄]T=[s11,s12,…,s12n-1,s12n]T,
where s1r, r=1,2,…,2n is the value of the first row of the archive T, which is the best solution for (18).Once the iterative procedure converges, the updated parameter vectors in (20) will give the final identification result of all the structural parameters, whilst the time history of the input F(t) can be easily determined by (7).


## 3.1. Introduction of CIM

The equation of motion of an-DOFs system can be expressed as follows:
(7)MX¨(t)+CX˙(t)+KX(t)=F(t),
where M, C, and K represent, respectively, mass, damping, and stiffness matrices of the structures, X(t) and X˙(t), while X¨(t) represents, respectively, the displacements, velocities, and accelerometers response vector of the structure; F(t) is the external excitation on the structures. Equation (7) can be rewritten as (8), which can be further rearranged as (9) at time instant ti [10]:
(8)CX˙(t)+KX(t)=F(t)-MX¨(t),(9)H(ti)θ=Y(ti),
where vector θ contains all the unknown parameters to be identified; matrix H consists of the measured structural responses; and vector Y is the system input. To assemble (9) at all the sampling time instants ti, i=1,2,…,Q together, will give
(10)HQθ=YQ.The components of matrixH and Y depend on the type of structure. For a shear-type structure, the expressions for H, Y, and θ are(11)θ=[c1,k1,c2,k2,…,cn,kn]T,(12)HQ=[H(t1),H(t2),…,H(tQ)]T,(13)H(ti)=[x˙1(ti)x1(ti)x˙1(ti)-x˙2(ti)x1(ti)-x2(ti)⋯000x˙2(ti)-x˙1(ti)x2(ti)-x1(ti)x˙2(ti)-x˙3(ti)x2(ti)-x3(ti)⋯⋯⋯⋯⋯⋯0⋯⋯⋯x˙n(ti)-x˙n-1(ti)xn(ti)-xn-1(ti)],(14)YQ=[Y(t1),Y(t2),…,Y(tQ)]T,(15)Y(ti)=[f1(ti)-m1x¨1(ti)f2(ti)-m2x¨2(ti)⋯fn(ti)-mnx¨n(ti)]T,wherek1,k2,…,kn, and c1,c2,…,cn represent, respectively, the stiffness and damping coefficients of the structure for each story, and xj(ti) and fj(ti) are the displacement response and external excitation force of the jth DOF (j=1,2,…,n) at the time instant ti.In traditional calculations, the parameters of the structure can be identified from (10) by the least-squares technique as
(16)θ=[HQTHQ]-1HQTYQ.However, (16) cannot be easily solved to determine the stiffness and damping parameters since there are unknown quantities involved in calculating HQ and YQ.As mentioned earlier, we have suggested the complete inverse method to tackle the unknown input situation in structural identification. Following the framework of CIM, Chen and Li [11] suggested the total compensation method. The total compensation method rests on the assumption that the locations of the external forces are known even though their time histories are unknown, and the number of DOF with applied (unknown) forces is less than the number of DOF whose responses are measured. This assumption reflects the situation of a forced vibration survey of structure where a limited number of one or several actuator(s) are installed on key locations of the structure to excite it. In this case, the input excitations in (15) can be further expressed into two parts:
(17)yu(ti)=fu(ti)-mu(ti)x¨u(ti),i=1,2,…,Q,yl(ti)=fl(ti)-ml(ti)x¨l(ti)=-ml(ti)x¨l(ti),yl(ti)=fl(ti)-ml(ti)x¨l(ti)wiii=1,2,…,Q,
where yu(ti) denotes those DOFs with unknown external excitation and yl(ti) stands for those DOFs without applied force, that is, fl(ti)≡0.


## 3.2. Objective Function

In order to identify structural parameters and the input time history from output-only measurements, an objective function is defined as minimizing discrepancies between the forceyl(t) and the calculated force yl⋄(t). The minimization of the objective function is expressed as a bound-constrained nonlinear least squares problem:
(18)min⁡⁡∑ti=t1tQz(ti)=∑ti=t1tQ∑l∈L∥yl(ti)-yl⋄(ti)∥2,i=1,2,…,Q,
where L is the set of l and the identified force yl⋄(ti) is given by
(19)yl⋄(ti)=Hl(ti)θ⋄,i=1,2,…,Q,
where Hl(ti) is lth row of the matrix H(ti), and θ⋄ is the result of identified parameters based on ACOR with the details as
(20)θ⋄=[c1⋄,k1⋄,c2⋄,k2⋄,…,cn⋄,kn⋄]T=[s11,s12,…,s12n-1,s12n]T,
where s1r, r=1,2,…,2n is the value of the first row of the archive T, which is the best solution for (18).Once the iterative procedure converges, the updated parameter vectors in (20) will give the final identification result of all the structural parameters, whilst the time history of the input F(t) can be easily determined by (7).


## 4. Numerical Examples

The suggested algorithm has been applied to several numerical examples including a truss structure, a 6-story shear frame, and a 12-story shear frame structure [25]. Since the observations for all numerical examples are similar, only the results for the 6-story shear frame structure are presented here in detail. The mass, stiffness, and damping coefficient of each story of the structure are shown in Table 1. Sinusoidal excitations are applied on the 4th and the 6th floor, which are f4(t)=6.74×105sin⁡(6t) and f6(t)=8.53×105sin⁡(4t), respectively. Therefore the set L in (18) is L={1,2,3,5} for this example. The dynamic responses of all six stories are first calculated in terms of displacement, velocity, and acceleration using the Newmark-β method [26]. Twenty ants are used in each iteration of the ACOR algorithm, and the convergence threshold in the objective function is set as 1.0×10-6 for all the cases. All the computation parameters used for ACOR are summarized in Table 2.Table 1
The basic parameters of the frame structure model.
Node number
1
2
3
4
5
6
Mass (kg)
610000
580000
530000
530000
530000
360000
Stiffness (kN/m)
271390
290360
282690
252480
233540
229940
Damping (kN·s/m)
2469.65
2874.56
2572.48
2297.57
2265.33
2092.45Table 2
Parameters used byACOR for HIM.
Parameters
Symbol
Value
Speed of convergence
ξ
0.85
Locality of the search process
q
0.2
Archive size
k
40
Number of dimensions
N
12


### 4.1. Noise-Free Measurements

The proposedACOR method is first applied to noise-free measurements. The parameter identification results of the three cases are shown in Table 3, where Cases 1 and 2 have the same initial parameters’ estimate but different measurement durations, and Cases 1 and 3 have the same measurement duration but different initial parameters. It is seen from Table 3 that, for all cases, the unknown parameters can be accurately identified by the proposed method with short duration of measurements and if the method is robust to obtain the initial estimated parameters. The input excitation can also be accurately identified by the ACOR method in the noise-free case.Table 3
Identified results with noise-free measurements.
Parameters and its real values
Case 1
Case 2
Case 3
θ
0
θ
⋄
θ
0
θ
⋄
θ
0
θ
⋄
k
1
271390
[1,1E6]
271390.00
[1,1E6]
271390.00
[1,1E7]
271389.99
k
2
290360
[1,1E6]
290360.00
[1,1E6]
290360.00
[1,1E7]
290360.00
k
3
282690
[1,1E6]
282690.00
[1,1E6]
282690.00
[1,1E7]
282689.99
k
4
252480
[1,1E6]
252490.00
[1,1E6]
252479.99
[1,1E7]
252479.99
k
5
233540
[1,1E6]
233540.00
[1,1E6]
233540.00
[1,1E7]
233540.00
k
6
229940
[1,1E6]
229939.99
[1,1E6]
229940.00
[1,1E7]
229940.00
c
1
2469.65
[1,1E6]
2469.6499
[1,1E6]
2469.6499
[1,1E7]
2469.6500
c
2
2874.56
[1,1E6]
2874.5600
[1,1E6]
2874.5600
[1,1E7]
2874.5599
c
3
2572.48
[1,1E6]
2572.4800
[1,1E6]
2572.4799
[1,1E7]
2572.4799
c
4
2297.57
[1,1E6]
2297.5699
[1,1E6]
2297.5700
[1,1E7]
2297.5700
c
5
2265.33
[1,1E6]
2265.3300
[1,1E6]
2265.3300
[1,1E7]
2265.3299
c
6
2092.45
[1,1E6]
2092.4499
[1,1E6]
2092.4499
[1,1E7]
2092.4500
L
100
20
100
Iteration number
410
417
968
θ
0: initial values range of the parameter vector, generate k solutions by uniform random sampling for archive T; θ⋄: identified parameter vector, L: number of sampling points used.


### 4.2. Noise-Pollution Measurements

White noise is numerically added to the calculated responses to simulate noisy measured data by the following equation:(21)x¨^=x¨+EpNnoiseMax⁡⁡[x¨],
where Ep is the noise level expressed as a percentage, Nnoise is a uniform distribution vector with interval [-1,1], and Max⁡⁡(x¨) is the maximum value of the calculated acceleration response. Four different noise levels at 1%, 5%, 10%, and 15% were considered in the calculation.The identification results are summarized in Table4, from which we can see that the proposed method can accurately identify the stiffness parameters for a noise-level up to 10%. The maximum parameter identification error of the stiffness parameter is lower than 0.1%, 0.6%, 1.6%, and 7.5% for noise levels 1%, 5%, 10%, and 15%, respectively. The identification accuracy for the damping ratio, however, is relatively low. This is not surprise since damping is very sensitive to measurement noise.Table 4
Identified results with noise-pollution measurements.
Case 1
Case 2
Case 3
Case 4
Parameters and its real values
(Noise = 1%)
(Noise = 5%)
(Noise = 10%)
(Noise = 15%)
θ
⋄
Error (%)
θ
⋄
Error (%)
θ
⋄
Error (%)
θ
⋄
Error (%)
k
1
k
1
271390
271212
−0.07
270655
−0.27
267110
−1.58
251259
−7.42
k
2
290360
290206
−0.05
288793
−0.54
288103
−0.78
273108
−5.94
k
3
282690
282877
0.07
281526
−0.41
278466
−1.49
276999
−2.01
k
4
252480
252271
−0.08
252238
−0.10
249497
−1.18
246445
−2.39
k
5
233540
233605
0.03
232816
−0.31
230917
−1.12
229261
−1.83
k
6
229940
229767
−0.08
229776
−0.07
230418
0.21
222248
−3.35
c
1
2469.65
2451.16
−0.75
2484.77
0.61
2061.82
−16.51
1591.87
−35.54
c
2
2874.56
2878.74
0.15
2865.31
−0.32
3047.32
6.01
3293.65
14.58
c
3
2572.48
2564.28
−0.32
2601.61
1.13
2124.07
−17.43
2192.95
−14.75
c
4
2297.57
2297.05
−0.02
2273.31
−1.06
2468.87
7.46
2490.99
8.42
c
5
2265.33
2266.12
0.03
2246.78
−0.82
2231.8
−1.48
2191.78
−3.25
c
6
2092.45
2095.25
0.13
2107.04
0.70
2134.34
2.00
2038.22
−2.59
L
200
300
200
300
For all cases, identical initial values range [1,1E7] is used for all stiffness and damping parameters; θ⋄: identified parameter vector; Error: relative error, L: number of sampling points used.It is interesting to compare the inverse input time history of the top floor with that of the second floor where no external force was applied. The identified input forces of the two floors for different noise levels are shown in Figure1, where the solid blue is the actual force curve and the dotted black line is the identified time history. Visually, even for noise level of 15% the identified force of the top floor matches well with the actual force. Amplitude of the identified curve of the second floor, on the other hand, is approximately zero compared with that of the top floor for noise level of 10%. This and several other numerical examples [25] have demonstrated the applicability and accuracy of the proposed identification method.Comparison of identified time histories of inputs with real input at different floor.
(a)
Noise level of 1%
(b)
Noise level of 5%
(c)
Noise level of 10%
(d)
Noise level of 15%


## 4.1. Noise-Free Measurements

The proposedACOR method is first applied to noise-free measurements. The parameter identification results of the three cases are shown in Table 3, where Cases 1 and 2 have the same initial parameters’ estimate but different measurement durations, and Cases 1 and 3 have the same measurement duration but different initial parameters. It is seen from Table 3 that, for all cases, the unknown parameters can be accurately identified by the proposed method with short duration of measurements and if the method is robust to obtain the initial estimated parameters. The input excitation can also be accurately identified by the ACOR method in the noise-free case.Table 3
Identified results with noise-free measurements.
Parameters and its real values
Case 1
Case 2
Case 3
θ
0
θ
⋄
θ
0
θ
⋄
θ
0
θ
⋄
k
1
271390
[1,1E6]
271390.00
[1,1E6]
271390.00
[1,1E7]
271389.99
k
2
290360
[1,1E6]
290360.00
[1,1E6]
290360.00
[1,1E7]
290360.00
k
3
282690
[1,1E6]
282690.00
[1,1E6]
282690.00
[1,1E7]
282689.99
k
4
252480
[1,1E6]
252490.00
[1,1E6]
252479.99
[1,1E7]
252479.99
k
5
233540
[1,1E6]
233540.00
[1,1E6]
233540.00
[1,1E7]
233540.00
k
6
229940
[1,1E6]
229939.99
[1,1E6]
229940.00
[1,1E7]
229940.00
c
1
2469.65
[1,1E6]
2469.6499
[1,1E6]
2469.6499
[1,1E7]
2469.6500
c
2
2874.56
[1,1E6]
2874.5600
[1,1E6]
2874.5600
[1,1E7]
2874.5599
c
3
2572.48
[1,1E6]
2572.4800
[1,1E6]
2572.4799
[1,1E7]
2572.4799
c
4
2297.57
[1,1E6]
2297.5699
[1,1E6]
2297.5700
[1,1E7]
2297.5700
c
5
2265.33
[1,1E6]
2265.3300
[1,1E6]
2265.3300
[1,1E7]
2265.3299
c
6
2092.45
[1,1E6]
2092.4499
[1,1E6]
2092.4499
[1,1E7]
2092.4500
L
100
20
100
Iteration number
410
417
968
θ
0: initial values range of the parameter vector, generate k solutions by uniform random sampling for archive T; θ⋄: identified parameter vector, L: number of sampling points used.


## 4.2. Noise-Pollution Measurements

White noise is numerically added to the calculated responses to simulate noisy measured data by the following equation:(21)x¨^=x¨+EpNnoiseMax⁡⁡[x¨],
where Ep is the noise level expressed as a percentage, Nnoise is a uniform distribution vector with interval [-1,1], and Max⁡⁡(x¨) is the maximum value of the calculated acceleration response. Four different noise levels at 1%, 5%, 10%, and 15% were considered in the calculation.The identification results are summarized in Table4, from which we can see that the proposed method can accurately identify the stiffness parameters for a noise-level up to 10%. The maximum parameter identification error of the stiffness parameter is lower than 0.1%, 0.6%, 1.6%, and 7.5% for noise levels 1%, 5%, 10%, and 15%, respectively. The identification accuracy for the damping ratio, however, is relatively low. This is not surprise since damping is very sensitive to measurement noise.Table 4
Identified results with noise-pollution measurements.
Case 1
Case 2
Case 3
Case 4
Parameters and its real values
(Noise = 1%)
(Noise = 5%)
(Noise = 10%)
(Noise = 15%)
θ
⋄
Error (%)
θ
⋄
Error (%)
θ
⋄
Error (%)
θ
⋄
Error (%)
k
1
k
1
271390
271212
−0.07
270655
−0.27
267110
−1.58
251259
−7.42
k
2
290360
290206
−0.05
288793
−0.54
288103
−0.78
273108
−5.94
k
3
282690
282877
0.07
281526
−0.41
278466
−1.49
276999
−2.01
k
4
252480
252271
−0.08
252238
−0.10
249497
−1.18
246445
−2.39
k
5
233540
233605
0.03
232816
−0.31
230917
−1.12
229261
−1.83
k
6
229940
229767
−0.08
229776
−0.07
230418
0.21
222248
−3.35
c
1
2469.65
2451.16
−0.75
2484.77
0.61
2061.82
−16.51
1591.87
−35.54
c
2
2874.56
2878.74
0.15
2865.31
−0.32
3047.32
6.01
3293.65
14.58
c
3
2572.48
2564.28
−0.32
2601.61
1.13
2124.07
−17.43
2192.95
−14.75
c
4
2297.57
2297.05
−0.02
2273.31
−1.06
2468.87
7.46
2490.99
8.42
c
5
2265.33
2266.12
0.03
2246.78
−0.82
2231.8
−1.48
2191.78
−3.25
c
6
2092.45
2095.25
0.13
2107.04
0.70
2134.34
2.00
2038.22
−2.59
L
200
300
200
300
For all cases, identical initial values range [1,1E7] is used for all stiffness and damping parameters; θ⋄: identified parameter vector; Error: relative error, L: number of sampling points used.It is interesting to compare the inverse input time history of the top floor with that of the second floor where no external force was applied. The identified input forces of the two floors for different noise levels are shown in Figure1, where the solid blue is the actual force curve and the dotted black line is the identified time history. Visually, even for noise level of 15% the identified force of the top floor matches well with the actual force. Amplitude of the identified curve of the second floor, on the other hand, is approximately zero compared with that of the top floor for noise level of 10%. This and several other numerical examples [25] have demonstrated the applicability and accuracy of the proposed identification method.Comparison of identified time histories of inputs with real input at different floor.
(a)
Noise level of 1%
(b)
Noise level of 5%
(c)
Noise level of 10%
(d)
Noise level of 15%


## 5. Experimental Example

The experimental data from a hammer test on a 3-story steel frame structure was adopted to further validate the applicability of the proposed identification method. The structural model is shown in Figure2 and dimensions of the model are given in Figure 3. The strength of the steel of the model was tested as 435 Mpa, and the elastic module was 200 Gpa. The overall dimension of the structure is 850*500*1450 mm (Figure 3). The structure’s column was made of a 9.5 × 75 mm steel plate and the floor was made by an 850 × 500 × 25 mm steel plate. Additional mass of 130 kg was installed on each floor leading a total floor mass to 230.3 Kg. The structural model can be treated as a shear-type building, and the shear stiffness of each floor was computed as 563651 N/m. Hammer tests, hand-shaking tests, and shaking table tests were conducted on the model. Only hammer tests and hand-shaking tests data are adopted in this paper. The hammer tests can be broadly divided into two categories. The first is an impact test where a certain floor of the building is hit once by a force hammer and then the building vibrates freely. The second is a continuous test wherein a certain floor of the building is hit continuously by the hammer.Figure 2
The 3-story building model in the test.Dimensions of the building model.
(a)
Elevation view
(b)
Plan view


### 5.1. Impact Test

In this test, the building model was hit by a hammer once on the 3rd floor and was then released for free decayed vibration in thex direction. Figure 4(a) shows the recorded accelerations for each floor. The accelerations were then integrated to obtain the corresponding velocity and displacement as shown in Figure 5. To reduce the measurement noise, the high-pass Butterworth signal filter was applied to the integration with the filter order N=4 where the lower cut-off frequency = 0.16 rad/sec.The recorded acceleration responses of each floor for impact test case.
(a)
The third floor
(b)
The second floor
(c)
The first floorFigure 5
The integrated velocity and displacement from accelerations shown in Figure4.The known information for this case isf1(t)≡0.0 and f2(t)≡0.0, and the structural parameters are to be identified by the proposed ACOR method with the following objective function:
(22)min⁡∑i=1Qz(ti)=∑i=1Q∥y1(ti)-y1⋄(ti)∥2+∑i=1Q∥y2(ti)-y2⋄(ti)∥2.Tables5(a), 5(b), and 5(c) show, respectively, the identified stiffness parameters from ten different response segments. Each has the same duration of 0.2, 0.33 and 1.0 second (i.e., 60, 100, and 300 data points for a sampling frequency 300 Hz). The same initial values for all the unknown parameters were used in the calculation. The maximum average identification error for all three cases is less than 2%, and the maximum identification error for a single parameter is less than 4%. Comparison between Tables 5(a), 5(b), and 5(c) also demonstrates that the identification accuracy is not sensitive to the duration of the response segment. Table 5(d) shows the identification results for response duration of 0.33 second (100 data points) and different initial values of unknown parameters. Comparison between Tables 5(b) and 5(d) shows that the proposed method is not sensitive to the initial estimates of the unknown parameters. All the above results indicate a very good identification accuracy achieved by the proposed method for a very short response duration, which has potential application for online structural health monitoring. It should be noted that the damping identification results are not satisfied in this case. The damping properties can be identified by other effective identification techniques such as the empirical mode decomposition plus Hilbert transform method suggested by the authors [27].Table 5
(a) Identified results for single hammer hit on the top floor (computation points = 60). (b) Identified results for single hammer hit on the top floor (computational points = 100). (c) Identified results for single hammer hit on the top floor (computation points = 300). (d) Identified results for single hammer hit on the top floor (computation points = 100 with different initial values for parameters).
(a)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~10.2 s
5.6027
−0.60
5.5243
−1.99
5.6330
−0.06
10.4~10.6 s
5.5882
−0.86
5.5804
−1.00
5.6430
0.12
10.8~11.0 s
5.5042
−2.35
5.7150
1.39
5.6406
0.07
11.2~11.4 s
5.6356
−0.02
5.4636
−3.07
5.6258
−0.19
11.6~11.8 s
5.4920
−2.56
5.6482
0.21
5.6706
0.60
12.0~12.2 s
5.5471
−1.59
5.6079
−0.51
5.6128
−0.42
12.4~12.6 s
5.6005
−0.64
5.5197
−2.07
5.6233
−0.23
12.8~13.0 s
5.4536
−3.25
5.6930
1.00
5.5977
−0.69
13.2~13.4 s
5.5688
−1.20
5.6138
−0.40
5.6182
−0.32
13.6~14.0 s
5.6279
−0.15
5.4499
−3.31
5.5719
−1.15
Average value
5.5621
−1.32
5.5816
−0.97
5.6237
−0.23
(b)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~10.33 s
5.5985
−0.67
5.6083
−0.50
5.6364
0.00
10.2~10.53 s
5.5308
−1.88
5.5969
−0.70
5.6462
0.17
10.4~10.73 s
5.5931
−0.77
5.5587
−1.38
5.6623
0.46
10.6~10.93 s
5.5365
−1.77
5.6322
−0.08
5.6652
0.51
10.8~11.13 s
5.5277
−1.93
5.6426
0.11
5.6669
0.54
11.0~11.33 s
5.5702
−1.18
5.5808
−0.99
5.6591
0.40
11.2~11.53 s
5.5274
−1.94
5.5934
−0.76
5.6655
0.51
11.4~11.73 s
5.5155
−2.15
5.6174
−0.34
5.6615
0.44
11.6~11.93 s
5.5655
−1.26
5.5664
−1.24
5.6284
−0.14
11.8~12.13 s
5.5815
−0.98
5.5713
−1.16
5.5944
−0.75
Average value
5.5547
−1.45
5.5968
−0.70
5.6486
0.21
(c)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~11.0 s
5.5588
−1.38
5.6011
−0.63
5.6436
0.13
10.4~11.4 s
5.5589
−1.38
5.6011
−0.63
5.6436
0.13
10.8~11.8 s
5.5202
−2.06
5.626
−0.19
5.6538
0.31
11.2~12.2 s
5.546
−1.61
5.5637
−1.29
5.6396
0.05
11.6~12.6 s
5.5457
−1.61
5.5805
−0.99
5.6326
−0.07
12.0~13.0 s
5.5195
−2.08
5.6151
−0.38
5.6089
−0.49
12.4~13.4 s
5.5409
−1.70
5.5773
−1.05
5.5951
−0.73
12.8~13.8 s
5.5409
−1.70
5.5773
−1.05
5.5951
−0.73
13.2~14.2 s
5.5392
−1.73
5.5713
−1.16
5.5805
−0.99
13.6~14.6 s
5.5578
−1.40
5.5061
−2.31
5.5519
−1.50
Average value
5.5428
−1.66
5.5820
−0.97
5.6145
−0.39
(d)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~10.33 s
5.5985
−0.67
5.6083
−0.50
5.6364
0.00
12.0~12.33 s
5.5665
−1.24
5.5997
−0.65
5.5795
−1.01
14.0~14.33 s
5.5667
−1.24
5.5997
−0.65
5.5795
−1.01
16.0~16.33 s
5.5005
−2.41
5.5325
−1.85
5.5381
−1.75
18.0~18.33 s
5.5362
−1.78
5.5068
−2.30
5.4753
−2.86
20.0~20.33 s
5.5362
−1.78
5.5068
−2.30
5.4752
−2.86
22.0~22.33 s
5.3287
−5.46
5.2234
−7.33
5.0971
−9.57
24.0~24.33 s
5.1685
−8.30
5.0521
−10.37
4.8238
−14.42
26.0~26.33 s
5.0674
−10.10
4.9758
−11.72
4.6406
−17.67
28.0~28.33 s
4.8131
−14.61
4.8464
−14.02
4.0776
−27.66Figures6 and 7 show the convergence procedure of all the stiffness parameters for the computational segments as 10.2~10.53 s and 11.0~11.33 s in Table 5(b). Note that the unknown stiffness parameters can be accurately identified using a short duration of response and the results are robust to the parameters’ initial guess.Parameter identification convergence procedure for Segments 10.2–10.53 sec in Table5(b).
(a)
Stiffnessk3
(b)
Stiffnessk2
(c)
Stiffnessk1Parameter identification convergence procedure for Segments 11.0–11.33 sec in Table5(b).
(a)
Stiffnessk3
(b)
Stiffnessk2
(c)
Stiffnessk1Figure8(a) shows the time history of the actual hammer force on the 3rd floor, and Figures 8(b) to 8(c) show the identified time histories of forces for the 3rd, 1st, and 2nd floor, respectively. Note that the identified input on the 3rd floor has a spike at the same time instant as the real input. The amplitudes of the identified inputs on the 1st and 2nd floor are almost zero compared to that of the 3rd floor.The real input on the top floor and the identified inputs of each floor.
(a)
Real input on the third floor
(b)
Identified input on the third floor
(c)
Identified input of the second floor
(d)
Identified input of the first floor


### 5.2. Continuous Hammer Test

In this test case, the hammer hit the 3rd floor continuously, and the resulting accelerations for each floor are shown in Figure9. The actual test duration was 60 seconds, but only the first 20 seconds were plotted for the sake of clarity. Using the same objective function as (22), the structural parameters were identified by the ACOR method, and the results are summarized in Table 6.Table 6
Identified results for continuous hammer hits on the top floor.
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
12.0~12.33 s
5.5252
−1.97
5.5454
−1.62
5.6421
0.10
13.0~13.33 s
5.5462
−1.60
5.6372
0.01
5.6888
0.93
27.0~27.33 s
5.5804
−1.00
5.6745
0.67
5.6446
0.14
28.0~28.33 s
5.5929
−0.77
5.4982
−2.45
5.7204
1.49
39.4~39.73 s
5.5038
−2.35
5.4916
−2.57
5.5697
−1.19
40.4~40.73 s
5.4955
−2.50
5.5077
−2.29
5.6238
−0.23
43.4~43.73 s
5.4719
−2.92
5.6114
−0.45
5.5332
−1.83
Average value
5.5308
−1.87
5.5666
−1.24
5.6318
−0.08The recorded acceleration responses of each floor for continuous hammer test.
(a)
Acceleration response of the third floor
(b)
Acceleration response of the second floor
(c)
Acceleration response of the first floorNote that the average identification error for each of the seven response segments is less than 2%, and the maximum identification error for each single parameter is less than 3%. Using the global average value of all the identified parameters, the inputs of each floor can be identified. Figures10(a) and 10(b) compare the identified input and real input on the 3rd floor. It is clear that the identified input has the continued impact spikes (peaks) at the same time as instants for the real input. The amplitudes of the identified inputs of the 2nd and 1st floor are almost zero compared to that of the 3rd floor. The result is consistent with the real test situation.The real input on the top floor and the identified inputs of each floor.
(a)
Real input on the third floor
(b)
Identified input on the third floor
(c)
Identified input of the second floor
(d)
Identified input of the first floor


### 5.3. Push the Building for a While and Release It to Free Vibration

In this test, the model was hand-shaked on the 2nd floor to vibrate for a while and was then released for free vibration. Figure11 shows the acceleration responses of each floor for this case. The known information is f1(t)≡0.0 and f3(t)≡0.0, where the objective function can then be established as follows:
(23)min⁡⁡∑i=1Qz(ti)=∑i=1Q∥y1(ti)-y1⋄(ti)∥2+∑i=1Q∥y3(ti)-y3⋄(ti)∥2.The recorded acceleration responses of each floor for hand-shaked test.
(a)
Acceleration response of the third floor
(b)
Acceleration response of the second floor
(c)
Acceleration response of the first floorTable7 shows the identification results for this case for seven different response segments. The maximum identification error is less than 1%. It is interesting to compare the identified results of input with the real input. Figures 12(a) to 12(c) show the identified inputs of 1st to 3rd floors, respectively. It is interesting to find that the input of the 2nd floor (Figure 12(b)) has a relatively large excitation force for the duration of 12 to 22 seconds. That is consistent with the test at that duration where the model was hand-shaked. The amplitudes of the inputs from the 1st and 3rd floor are relatively small compared to that of the 2nd floor. The above three cases demonstrate the applicability and feasibility of the proposed hybrid identification method using response-only measurements.Table 7
Identified results for hand-shaked tests on the second floor.
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
19.0~19.33 s
5.6022
−0.61
5.5341
−1.82
5.6282
−0.15
20.0~20.33 s
5.5308
−1.88
5.5969
−0.70
5.6462
0.17
24.0~24.33 s
5.6062
−0.54
5.6089
−0.49
5.6598
0.41
30.0~30.33 s
5.6013
−0.62
5.6068
−0.53
5.6424
0.10
40.0~40.33 s
5.6181
−0.33
5.6101
−0.47
5.6444
0.14
50.0~50.33 s
5.5968
−0.70
5.6174
−0.34
5.6403
0.07
60.0~60.33 s
5.5837
−0.94
5.6017
−0.62
5.6009
−0.63
Average value
5.5913
−0.80
5.5966
−0.71
5.6375
0.02The identified inputs of each floor.
(a)
Identified input on the third floor
(b)
Identified input on the second floor
(c)
Identified input on the third floor


## 5.1. Impact Test

In this test, the building model was hit by a hammer once on the 3rd floor and was then released for free decayed vibration in thex direction. Figure 4(a) shows the recorded accelerations for each floor. The accelerations were then integrated to obtain the corresponding velocity and displacement as shown in Figure 5. To reduce the measurement noise, the high-pass Butterworth signal filter was applied to the integration with the filter order N=4 where the lower cut-off frequency = 0.16 rad/sec.The recorded acceleration responses of each floor for impact test case.
(a)
The third floor
(b)
The second floor
(c)
The first floorFigure 5
The integrated velocity and displacement from accelerations shown in Figure4.The known information for this case isf1(t)≡0.0 and f2(t)≡0.0, and the structural parameters are to be identified by the proposed ACOR method with the following objective function:
(22)min⁡∑i=1Qz(ti)=∑i=1Q∥y1(ti)-y1⋄(ti)∥2+∑i=1Q∥y2(ti)-y2⋄(ti)∥2.Tables5(a), 5(b), and 5(c) show, respectively, the identified stiffness parameters from ten different response segments. Each has the same duration of 0.2, 0.33 and 1.0 second (i.e., 60, 100, and 300 data points for a sampling frequency 300 Hz). The same initial values for all the unknown parameters were used in the calculation. The maximum average identification error for all three cases is less than 2%, and the maximum identification error for a single parameter is less than 4%. Comparison between Tables 5(a), 5(b), and 5(c) also demonstrates that the identification accuracy is not sensitive to the duration of the response segment. Table 5(d) shows the identification results for response duration of 0.33 second (100 data points) and different initial values of unknown parameters. Comparison between Tables 5(b) and 5(d) shows that the proposed method is not sensitive to the initial estimates of the unknown parameters. All the above results indicate a very good identification accuracy achieved by the proposed method for a very short response duration, which has potential application for online structural health monitoring. It should be noted that the damping identification results are not satisfied in this case. The damping properties can be identified by other effective identification techniques such as the empirical mode decomposition plus Hilbert transform method suggested by the authors [27].Table 5
(a) Identified results for single hammer hit on the top floor (computation points = 60). (b) Identified results for single hammer hit on the top floor (computational points = 100). (c) Identified results for single hammer hit on the top floor (computation points = 300). (d) Identified results for single hammer hit on the top floor (computation points = 100 with different initial values for parameters).
(a)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~10.2 s
5.6027
−0.60
5.5243
−1.99
5.6330
−0.06
10.4~10.6 s
5.5882
−0.86
5.5804
−1.00
5.6430
0.12
10.8~11.0 s
5.5042
−2.35
5.7150
1.39
5.6406
0.07
11.2~11.4 s
5.6356
−0.02
5.4636
−3.07
5.6258
−0.19
11.6~11.8 s
5.4920
−2.56
5.6482
0.21
5.6706
0.60
12.0~12.2 s
5.5471
−1.59
5.6079
−0.51
5.6128
−0.42
12.4~12.6 s
5.6005
−0.64
5.5197
−2.07
5.6233
−0.23
12.8~13.0 s
5.4536
−3.25
5.6930
1.00
5.5977
−0.69
13.2~13.4 s
5.5688
−1.20
5.6138
−0.40
5.6182
−0.32
13.6~14.0 s
5.6279
−0.15
5.4499
−3.31
5.5719
−1.15
Average value
5.5621
−1.32
5.5816
−0.97
5.6237
−0.23
(b)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~10.33 s
5.5985
−0.67
5.6083
−0.50
5.6364
0.00
10.2~10.53 s
5.5308
−1.88
5.5969
−0.70
5.6462
0.17
10.4~10.73 s
5.5931
−0.77
5.5587
−1.38
5.6623
0.46
10.6~10.93 s
5.5365
−1.77
5.6322
−0.08
5.6652
0.51
10.8~11.13 s
5.5277
−1.93
5.6426
0.11
5.6669
0.54
11.0~11.33 s
5.5702
−1.18
5.5808
−0.99
5.6591
0.40
11.2~11.53 s
5.5274
−1.94
5.5934
−0.76
5.6655
0.51
11.4~11.73 s
5.5155
−2.15
5.6174
−0.34
5.6615
0.44
11.6~11.93 s
5.5655
−1.26
5.5664
−1.24
5.6284
−0.14
11.8~12.13 s
5.5815
−0.98
5.5713
−1.16
5.5944
−0.75
Average value
5.5547
−1.45
5.5968
−0.70
5.6486
0.21
(c)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~11.0 s
5.5588
−1.38
5.6011
−0.63
5.6436
0.13
10.4~11.4 s
5.5589
−1.38
5.6011
−0.63
5.6436
0.13
10.8~11.8 s
5.5202
−2.06
5.626
−0.19
5.6538
0.31
11.2~12.2 s
5.546
−1.61
5.5637
−1.29
5.6396
0.05
11.6~12.6 s
5.5457
−1.61
5.5805
−0.99
5.6326
−0.07
12.0~13.0 s
5.5195
−2.08
5.6151
−0.38
5.6089
−0.49
12.4~13.4 s
5.5409
−1.70
5.5773
−1.05
5.5951
−0.73
12.8~13.8 s
5.5409
−1.70
5.5773
−1.05
5.5951
−0.73
13.2~14.2 s
5.5392
−1.73
5.5713
−1.16
5.5805
−0.99
13.6~14.6 s
5.5578
−1.40
5.5061
−2.31
5.5519
−1.50
Average value
5.5428
−1.66
5.5820
−0.97
5.6145
−0.39
(d)
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
10.0~10.33 s
5.5985
−0.67
5.6083
−0.50
5.6364
0.00
12.0~12.33 s
5.5665
−1.24
5.5997
−0.65
5.5795
−1.01
14.0~14.33 s
5.5667
−1.24
5.5997
−0.65
5.5795
−1.01
16.0~16.33 s
5.5005
−2.41
5.5325
−1.85
5.5381
−1.75
18.0~18.33 s
5.5362
−1.78
5.5068
−2.30
5.4753
−2.86
20.0~20.33 s
5.5362
−1.78
5.5068
−2.30
5.4752
−2.86
22.0~22.33 s
5.3287
−5.46
5.2234
−7.33
5.0971
−9.57
24.0~24.33 s
5.1685
−8.30
5.0521
−10.37
4.8238
−14.42
26.0~26.33 s
5.0674
−10.10
4.9758
−11.72
4.6406
−17.67
28.0~28.33 s
4.8131
−14.61
4.8464
−14.02
4.0776
−27.66Figures6 and 7 show the convergence procedure of all the stiffness parameters for the computational segments as 10.2~10.53 s and 11.0~11.33 s in Table 5(b). Note that the unknown stiffness parameters can be accurately identified using a short duration of response and the results are robust to the parameters’ initial guess.Parameter identification convergence procedure for Segments 10.2–10.53 sec in Table5(b).
(a)
Stiffnessk3
(b)
Stiffnessk2
(c)
Stiffnessk1Parameter identification convergence procedure for Segments 11.0–11.33 sec in Table5(b).
(a)
Stiffnessk3
(b)
Stiffnessk2
(c)
Stiffnessk1Figure8(a) shows the time history of the actual hammer force on the 3rd floor, and Figures 8(b) to 8(c) show the identified time histories of forces for the 3rd, 1st, and 2nd floor, respectively. Note that the identified input on the 3rd floor has a spike at the same time instant as the real input. The amplitudes of the identified inputs on the 1st and 2nd floor are almost zero compared to that of the 3rd floor.The real input on the top floor and the identified inputs of each floor.
(a)
Real input on the third floor
(b)
Identified input on the third floor
(c)
Identified input of the second floor
(d)
Identified input of the first floor


## 5.2. Continuous Hammer Test

In this test case, the hammer hit the 3rd floor continuously, and the resulting accelerations for each floor are shown in Figure9. The actual test duration was 60 seconds, but only the first 20 seconds were plotted for the sake of clarity. Using the same objective function as (22), the structural parameters were identified by the ACOR method, and the results are summarized in Table 6.Table 6
Identified results for continuous hammer hits on the top floor.
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
12.0~12.33 s
5.5252
−1.97
5.5454
−1.62
5.6421
0.10
13.0~13.33 s
5.5462
−1.60
5.6372
0.01
5.6888
0.93
27.0~27.33 s
5.5804
−1.00
5.6745
0.67
5.6446
0.14
28.0~28.33 s
5.5929
−0.77
5.4982
−2.45
5.7204
1.49
39.4~39.73 s
5.5038
−2.35
5.4916
−2.57
5.5697
−1.19
40.4~40.73 s
5.4955
−2.50
5.5077
−2.29
5.6238
−0.23
43.4~43.73 s
5.4719
−2.92
5.6114
−0.45
5.5332
−1.83
Average value
5.5308
−1.87
5.5666
−1.24
5.6318
−0.08The recorded acceleration responses of each floor for continuous hammer test.
(a)
Acceleration response of the third floor
(b)
Acceleration response of the second floor
(c)
Acceleration response of the first floorNote that the average identification error for each of the seven response segments is less than 2%, and the maximum identification error for each single parameter is less than 3%. Using the global average value of all the identified parameters, the inputs of each floor can be identified. Figures10(a) and 10(b) compare the identified input and real input on the 3rd floor. It is clear that the identified input has the continued impact spikes (peaks) at the same time as instants for the real input. The amplitudes of the identified inputs of the 2nd and 1st floor are almost zero compared to that of the 3rd floor. The result is consistent with the real test situation.The real input on the top floor and the identified inputs of each floor.
(a)
Real input on the third floor
(b)
Identified input on the third floor
(c)
Identified input of the second floor
(d)
Identified input of the first floor


## 5.3. Push the Building for a While and Release It to Free Vibration

In this test, the model was hand-shaked on the 2nd floor to vibrate for a while and was then released for free vibration. Figure11 shows the acceleration responses of each floor for this case. The known information is f1(t)≡0.0 and f3(t)≡0.0, where the objective function can then be established as follows:
(23)min⁡⁡∑i=1Qz(ti)=∑i=1Q∥y1(ti)-y1⋄(ti)∥2+∑i=1Q∥y3(ti)-y3⋄(ti)∥2.The recorded acceleration responses of each floor for hand-shaked test.
(a)
Acceleration response of the third floor
(b)
Acceleration response of the second floor
(c)
Acceleration response of the first floorTable7 shows the identification results for this case for seven different response segments. The maximum identification error is less than 1%. It is interesting to compare the identified results of input with the real input. Figures 12(a) to 12(c) show the identified inputs of 1st to 3rd floors, respectively. It is interesting to find that the input of the 2nd floor (Figure 12(b)) has a relatively large excitation force for the duration of 12 to 22 seconds. That is consistent with the test at that duration where the model was hand-shaked. The amplitudes of the inputs from the 1st and 3rd floor are relatively small compared to that of the 2nd floor. The above three cases demonstrate the applicability and feasibility of the proposed hybrid identification method using response-only measurements.Table 7
Identified results for hand-shaked tests on the second floor.
Computation response segments
Identification results(×105) and identification errors
k
1
Error
k
2
Error
k
3
Error
19.0~19.33 s
5.6022
−0.61
5.5341
−1.82
5.6282
−0.15
20.0~20.33 s
5.5308
−1.88
5.5969
−0.70
5.6462
0.17
24.0~24.33 s
5.6062
−0.54
5.6089
−0.49
5.6598
0.41
30.0~30.33 s
5.6013
−0.62
5.6068
−0.53
5.6424
0.10
40.0~40.33 s
5.6181
−0.33
5.6101
−0.47
5.6444
0.14
50.0~50.33 s
5.5968
−0.70
5.6174
−0.34
5.6403
0.07
60.0~60.33 s
5.5837
−0.94
5.6017
−0.62
5.6009
−0.63
Average value
5.5913
−0.80
5.5966
−0.71
5.6375
0.02The identified inputs of each floor.
(a)
Identified input on the third floor
(b)
Identified input on the second floor
(c)
Identified input on the third floor


## 6. Conclusions

This paper combines the complete inverse method with the ant colony optimization algorithm to identify the input excitation and the structural parameters from the output-only measurements. The core idea of the complete inverse method is to convert the physical features of the excitation into a mathematical confinement condition. In theory, the proposed method has no limitation on the type of excitation. This paper applies the method to the situation where spatial locations of a limited number of excitations are known. Numerical examples were carried out to evaluate the feasibility of the proposed method. For the situation of noise-free output measurements, numerical studies show that the proposed method can reliably and efficiently identify both the structural parameters and the input time history using a short duration of measurements. Moreover, the accuracy and convergence of the method are robust to the initial values selected for the unknown parameters. Forf noise-pollution measurements, the stiffness parameters as well as the input excitation were identified satisfactorily even with high noise level. The identification accuracy of the damping coefficients is broadly acceptable at low noise level and become relatively poor at high noise level. The proposed method is then applied to experimental results of a three-story building model and the results also prove the feasibility of the method.



---

*Source: 185487-2014-03-04.xml*