# Indirect Reconstruction of Structural Responses Based on Transmissibility Concept and Matrix Regularization

Authors: Chudong Pan; Liwen Zhang; Zhuo Sun
Journal: Shock and Vibration
 (2021)

---

## Abstract

A novel method is proposed based on the transmissibility concept and matrix regularization for indirectly measuring the structural responses. The inputs are some measured responses that are obtained via physical sensors. The outputs are the structural responses corresponding to some critical locations where no physical sensors are installed. Firstly, the transmissibility concept is introduced for expressing the relationship between the measured responses and the indirectly measured ones. Herein, a transmissibility matrix is formulated according to the theory of force identification under unknown initial conditions. Then, in order to reduce the size of the transmissibility matrix, structural responses are reshaped in a form of a matrix by using the concept of moving time windows. According to the matrix form of input-output relationship, indirect reconstruction of responses is boiled down to an optimization equation. Since inverse problem may be ill-conditioned, matrix regularization such as F-norm regularization is then recommended for improving the optimization problem. Herein, the penalty function is defined by using a weighted sum of two F-norm values, which correspond to the estimated responses of physical sensors and the ones of the concerned critical locations, respectively. Numerical simulations and experimental studies are finally carried out for verifying the effectiveness and feasibility of the proposed method. Some results show that the proposed method can be applied for indirectly measuring the responses with good robustness.


---

## Body

## 1. Introduction

Measurement of structural responses is a fundamental task in structural health monitoring (SHM). The common measuring methods can be divided into two categories: contact and noncontact methods. Physical sensors are needed for the contact methods. The physical sensors are often limited due to some reasons such as the limited funds and the inaccessibility of locations for measurement. For the noncontact methods, structural responses can be usually measured by using the laser measuring technology, image/video-based measuring method, or even global positioning systems (GPS) [1–3]. Compared to the contact methods, noncontact technologies such as the image/video-based measuring method can be usually applied for monitoring multiple points at the same time. However, it is also hard for the noncontact methods to obtain the responses occurring in a covert area.To extrapolate the responses at some covert locations, structural response reconstruction has been drawing a lot of attention. In the earlier research, a representative research work is published by Kammer [4]. He proposed a method for estimating the structural responses by using a transformation matrix, which is formulated via system Markov parameters determined from a vibration test. Since then, technology for response reconstruction has developed rapidly, and many methods have been proposed [5]. In general, existing methods for response reconstruction can be divided into two major categories: data-driven methods and model-based methods. The data-driven methods are model-free. They often train black-box models based on the historical real-life measurement data and use them in the reconstructed process [6, 7]. The model-based methods often utilize an analysis model of the structure for describing the structural dynamic behavior. In view of this, the relationship between structural inputs and outputs can be expressed by using some mathematical models. Contrary to the data-driven methods, the model-based methods always require few real-life measured responses.The majority of existing model-based methods are proposed based on the transmissibility concept. Transmissibility is applied for reflecting the relationship between responses obtained from different locations of the structure [8–10]. The transmissibility concept can be expressed in the frequency domain [11]. Maia and his collaborators have published a lot of excellent papers for the theory and application of this topic [12–14]. A good review of the transmissibility-based research in the field of SHM has recently been done by Yan et al. [9]. For the application of transmissibility function, Law et al. [15], for example, have proposed a method for response reconstruction in a full structure or in a substructure. The transmissibility concept can be also developed based on the basic theory of force identification [16–18]. In these methods, the dynamic forces acting on the structure are firstly identified by using the measured responses. Then, the task of response reconstruction can be easily done as a forward problem [19, 20]. In fact, no matter what the specific form of the transmissibility concept is, the basic idea is using sensor measurements at a convenient location to extrapolate the responses at the desired locations via a transformation matrix or function. The transmissibility-based methods have a clear physical meaning. Most of them are applied for dealing with response reconstruction in a case when acting locations of the dynamic forces are known while time histories are unknown.Another type of response reconstruction is developed based on the theory of modal expansion. The basic idea is to estimate the model responses by using the measured structural responses. Then, according to the theory of modal expansion, structural responses occurring at the desired locations can be calculated [21, 22]. The number of real sensors applied to the modal expansion-based method should not be less than the number of considered modal orders.In addition to the above methods, advance signal processing technology has also been studied for response reconstruction. Empirical mode decomposition (EMD) and Kalman filter are two classical technologies. He et al. [23] proposed an empirical mode decomposition-based method for response reconstruction in time domain, wherein EMD is applied for estimating the model responses. This idea can be also seen in other research works [24, 25]. Kalman filter is an effective technology for estimating the vibration state of the structure; hence, it can be utilized for response reconstruction. Using moving windows, Zhang and Wu [26] proposed a Kalman filter-based method for predicting the structural responses under unknown measurement noises. Many articles relating to Kalman filter can be also found in the topic of structural response reconstruction [27–30]. One advantage of the Kalman filter-based strategy is that it can be used for dealing with different types of structural responses [31]. However, these methods solved the problem step by step in the time domain. Hence, some characteristics of signal cannot be effectively applied for analysis. In addition to the above methods, several model-based strategies have been also proposed for response reconstruction in recent decades [32–34].For the model-based methods, the model error of structure and the measurement noises are two important factors that could affect the reconstructed accuracy. The former one could affect the basic mathematical model that plays a foundation role. A refined model can often reduce the negative effect of model error to some extent. The latter ones can induce some estimated errors or even an unacceptable reconstructed result. This is because the problem of response reconstruction essentially belongs to the class of inverse problems. Many strategies can be applied for reducing the effects of measurement noises, for example, improving signal quality by optimizing sensor placement [35] and introducing the regularization method [36]. Herein, regularization method is a commonly used technology. It has been widely studied in the field of SHM [37, 38]. Furthermore, the transmissibility matrix-based methods, which are developed in the time domain, may not be an ideal technology when a long-time duration problem is considered. This is mainly because the size of transmissibility matrix will increase with increasing the number of sampling points. This characteristic will bring more memory cost and more computing time.This paper mainly focuses on the indirect reconstruction of structural responses in the time domain. Transmissibility matrix, which is formulated based on theory force identification under unknown initial condition, is applied for formulating a governing equation, so that it can use a transmissibility matrix with relatively small size to deal with a long-time duration problem. The problem of response reconstruction is boiled down to an optimization problem. Frobenius norm (F-norm) matrix regularization is introduced for reducing the effects of measurement noises. It is beneficial to note that the F-norm regularization applied in this study belongs to a matrix regularization method, but not the commonly used vector-based regularization method.The paper is organized as follows. Structural response reconstruction is reviewed in Section1. Theoretical background is introduced in Section 2. Numerical simulations and experimental verification are illustrated in Section 3 and Section 4, respectively. Finally, some conclusions are drawn in Section 5.


## 2. Theoretical Background

As shown in Figure1, a diagram of frame structure is taken as an example for expressing the main intention of this study. The real sensors are applied for directly measuring the structural responses. The virtual sensor does not exist in practice. However, obtaining the structural responses happening at the position of virtual sensor is the main purpose of this study. Because of this, this paper aims to propose a matrix regularization-based method for indirectly estimating the structural responses corresponding to the virtual sensor. Herein, the input signal of the proposed method is the structural responses obtained from the real sensors. To achieve this goal, the related theoretical background is firstly introduced in Section 2.Figure 1
Frame structure as example for expressing intention of the proposed method.


### 2.1. Structural Responses Induced by Excitation Sources

A relationship between structural responses and dynamic force is firstly given as follows. Herein, the structural initial conditions are assumed to be equal to zero. A dictionary consisted ofNf atoms is introduced for expanding the force as [38, 39](1)f=d1,d2,⋯,dNfα1α2⋮αNf=Dα,where f = [f(Δt), f(2Δt), …, f(nΔt)]T means a force vector. Δt is a time interval for the forward analysis. D is the introduced dictionary and di represents the i-th atom of dictionary D. Herein, it assume that nΔt is equal to tw, which means that the discussions of structural response mainly locate in a time range (0, tw]. α is a coefficient vector, whose element αi is a factor corresponding to di. Nf represents the number of atoms in the dictionary D. Some sensors installed on the structure are considered. The structural responses induced by dynamic force can be expressed as(2)bif=hi1,hi2,⋯,hiNfα1α2⋮αNf=Hifα,i=1,2,⋯,Ns,where bif = [b(Δt̅), b(2Δt̅), …, b(kΔt̅)]T is a response vector of the i-th sensor. Δt̅ (Δt̅ ≥ Δt) represents the sampling time interval. Hif is a system matrix of the i-th sensor corresponding to the force. hij is the j-th component of structural responses caused by atom dj. Herein, the subscript “i” means the i-th sensor. These response components can be calculated by applying the mode superposition method. Ns represents the number of sensors.Let us focus on the analysis of structural responses caused by the structural initial conditions. Modal space is used for expressing the initial conditions. Structural dynamic behavior is assumed that it can be approximately described by using the firstNm orders. Then the initial conditions can be organized as a vector like [40](3)y0=yN10,y˙N10,yN20,y˙N20,⋯,yNNm0,y˙NNm0T,where yNi(0) and ẏNi(0) are the i-th mode initial displacement and i-th mode initial velocity, respectively. The corresponding mode shape matrix Φ = [φ1, φ2, …, φNm] satisfies ΦTMΦ = I. M represents the mass matrix of the structure. According to equation (3), the structural responses caused by initial conditions can be expressed via linear equations like(4)biy=h¯i1,h¯i2,h¯i3,⋯,h¯i2Nm=Hiyy0,i=1,2,⋯,Ns,where biy represents the structural responses caused by the initial conditions. Hiy is the system matrix corresponding to the initial conditions. The column h¯ij means the structural responses caused by the unit j-th initial condition component, i.e., y(0) = [0(1), …, 0(j-1), 1(j), 0(j+1), …, 0(2Nm)]T. Hence, the matrix Hiy actually can be calculated by applying the classical linear vibration theory. In the above equation, subscript “i” represents the i-th sensor.Now, if both external force and initial conditions are considered as the excitation sources, then structural responses obtained from thei-th sensor can be expressed as(5)bi=bif+biy=Hif,Hiyαy0=Hiβ,i=1,2,⋯,Ns,where bi means structural responses of the i-th sensor. Hi = [Hif, Hiy] is a system matrix corresponding to two excitation sources. β = [αT, y(0)T]T represents the structural excitation sources. Furthermore, if more than one sensor is considered, then equation (5) can be rewritten into one equation such as(6)b1b2⋮bNs=H1H2⋮HNsβ⇒b=Hβ.In addition, in order to reduce the possible large difference of numerical values among columns of system matrixH, equation (5) can be further rewritten as [40](7)b=Hβ=AΛβ,where A = [A1, A2, …, AN] is a normalized system matrix, in which the i-th atom Ai is calculated as Ai = Hi/||Hi||2. Λ = diag(||H1||2, ||H2||2, …, ||HN||2) is a diagonal scaling matrix. Hi represents the i-th column in matrix H.


### 2.2. Theoretical Derivation of Transmissibility Matrix

Two sets of sensors are considered. As shown in Figure1, all real sensors are grouped into the first set, while the virtual sensor is considered as the second set. The virtual sensor means that its corresponding structural responses will be indirectly reconstructed from the responses corresponding to the real sensors. According to equation (7), the structural responses corresponding to two sets of sensors can be expressed as(8)br=ArΛrβ,(9)bv=AvΛvβ,where the subscripts “r” and “v” represent the real sensors and virtual sensor, respectively. According to equation (8), the structural excitation sources (force and initial conditions) can be approximately estimated as(10)βe=Λr−1Ar+br,where βe represents the structural inputs identified from the measured structural responses br. Ar+ is a pseudoinverse matrix of Ar. Substituting equation (10) into equation (9) yields(11)bve=AvΛvΛr−1Ar+br=Tvrbr,where Tvr=AvΛvΛr−1Ar+ is a transmissibility matrix calculated from real sensor set to virtual sensor set. bve represents the structural responses corresponding to virtual sensor estimated from structural responses br.According to the above analysis, it can be seen that the transmissibility matrix is developed based on the theory of force identification under unknown initial conditions. The first step is to estimate the excitation sources, and the second step is to reconstruct the responses by using the identified excitation sources. In view of this, the key issue for this definition is that the real sensors should provide enough information for effectively identifying the unknown excitation sources. The number of real sensors is better to be bigger than the number of considered forces. Herein, it is beneficial to note again that both dynamic force and initial conditions are applied for the calculation of the transmissibility matrix. This definition is different from the one proposed by Zhu et al. [16] due to the latter one is calculated without considering unknown initial conditions.From equation (11), it can be seen that structural responses of virtual sensor can be indirectly calculated by using the structural responses obtained from real sensors. However, equation (11) is not recommended to be directly applied for the response reconstruction. The reason comes from the influences of measurement noises. In real application, the measurement noises are inevitable. The reconstructed responses directly estimated from equation (11) can be expressed as(12)bve=Tvrbr+Tvrε,where ε is a vector used for representing the measured noises. Essentially speaking, applying the transmissibility matrix to calculate the structural responses in the covert area is actually an inverse problem. In view of this, the second term such as Tvrε shown in equation (12) may lead to a large fluctuation in the estimated results so that the reconstructed responses may be unacceptable. Furthermore, from equation (11) it is not hard to find that the size of the transmissibility matrix is related to the number of sampling points. It will consume more computing resources for the calculation if a problem with long-time duration is considered. To overcome these drawbacks, structural responses will be reshaped in a matrix form, and then the F-norm regularization method will be introduced to ensure that the reconstructed responses are robust to the measured noises.


### 2.3. Response Reconstruction Based on F-Norm Regularization

As shown in Figure2, moving time window with a length of k sampling points is applied for extracting local signal of the structural responses [41]. The length of overlap between any two adjacent windows contains k0 (k0 < k) sampling points. The structural responses extracted from each time window are stored as a column of a structural response matrix. The structural response matrix can be expressed as(13)B=b1b1k−ko+1b2k−ko+1⋯bnw−1k−ko+1b2b1k−ko+2b2k−ko+2⋯bnw−1k−ko+2b3b1k−ko+3b2k−ko+3⋯bnw−1k−ko+3⋮⋮⋮⋱⋮bkb1k−ko+kb2k−ko+k⋯bnw−1k−ko+k,where B means the structural response matrix corresponding to one measured point. nw represents the number of moving time windows. bi represents the i-th sampling point. Taking the matrix form of structural responses into consideration, equation (11) can be rewritten as(14)Bve=TvrBr,where Bve is the estimated structural responses corresponding to the virtual sensor. Br means the structural responses corresponding to the real sensors. It can be expressed as(15)Br=B1TB2T⋯BNsTT,where Bi represents the structural responses corresponding to the i-th real sensor. From equation (14), it is not hard to find that the size of the transmissibility matrix is only determined by the length of moving time window. This characteristic means that a transmissibility matrix with relatively small size can be used for dealing with a problem of response reconstruction even though a problem of long-time duration is considered.Figure 2
Extracting local signal of structural responses via moving time window.As we have briefly discussed above, it is not recommended that the relationship between two sets of sensors is directly used for the reconstruction of structural responses due to influences of the measured noises. Herein, least square method is introduced for estimating the structural responses corresponding to the real sensors. This step can be expressed as an optimization problem like(16)Bre=argminBreTrrBre−BrF2,where ||∗||F means an F-norm of the considered matrix. Bre is a denoised response matrix. Trr means a transmissibility matrix calculated from the real sensors to themselves. According to the physical meaning of Trr, this matrix actually can be selected as an identity matrix I.Equation (16) can be used for estimating the responses Bre. However, it may cause an overfitting. Regularization technology is introduced for avoiding the overfitting problem. Herein, we can find that the estimated responses are organized in a form of matrix. Hence, the classical vector-based regularization method, whose estimated variable is organized as a vector, cannot be directly used in this case. Matrix regularization is utilized in this study. Herein, matrix regularization generalizes notions of the vector-based regularization to cases where the object to be estimated is a matrix. In view of this, matrix regularization such as F-norm regularization can be applied as(17)Bre=argminBreBre−BrF2+λΓBreF2,where λ is a regularization parameter. ||ΓBre||F2 is a penalty function, in which Γ is a smooth matrix. The smooth matrix Γ can be used for representing some reasonable prior knowledge of the estimated result. It can be often selected as an identity matrix I. However, this selection is not recommended because it does not take any information of Bve = TvrBre into consideration. In fact, it should be noted that what we want is the estimation of structural responses Bve. The structural responses Bre actually play a role of intermediate variable. Hence, it is better to take some reasonable prior knowledge of Bve = TvrBre into the definition of the penalty function. To achieve this goal, however, it may be difficult to formulate a specific matrix of Γ in an explicit form.According to the above consideration, the smooth matrixΓ will be defined in an implicit form that is similar to the one proposed by Pan et al. [42]. In view of this, equation (17) can be rewritten as(18)Bre=argminBreBre−BrF2+λBreF2+ωTvrBreF2,where the penalty function is defined as a weighted sum of BreF2 and BveF2. ω is a weighted factor. Herein, it is not hard to find that the former part of penalty function can avoid the overfitting problem corresponding to real sensors, while the latter part can make sure that signal energy of the reconstructed responses corresponding to the virtual sensor will not be too large.Let us calculate the derivative of function shown in equation (18) with respect to Bre, and let the corresponding result be equal to 0. Then, the estimated response matrix solved by equation (18) can be obtained as(19)Bre=I+λI+ωTvrTTvr−1Br.From equation (19), we can find that the matrix (I+ωTvrTTvr) is an additional matrix for the original problem. The spectrum distribution of the matrix (I+ωTvrTTvr) should be relatively centralized; otherwise, the problem shown in equation (19) will be ill-conditioned if a big regularization parameter is considered. In view of this, the values of ω should be limited. As done by Pan et al. [42], a ratio factor defined between the maximal and minimal eigenvalues of (I+ωTvrTTvr) is used for reflecting the spectrum distribution. This factor can be calculated as(20)r=1+ωemax1+ωemin,where emax and emin are the maximal and minimal eigenvalues of TvrTTvr, respectively. Denoting the maximal acceptable value of r as rmax, then we can calculate the maximal acceptable value of ω as(21)ωmax=rmax−1emax−rmaxemin,emax>rmaxemin,∞,emax≤rmaxemin,where ωmax represents the maximal acceptable value of ω. In this study, the value of rmax will be set as rmax = 109. This value means that the requirement for the concentration of eigenvalues of matrix (I+ωTvrTTvr) is not very strict. In this way, the selection of weighted factor ω will be more flexible. Furthermore, equation (18) can be also rewritten as(22)Bre=argminBreBre−BrF2+λBreF2+ωTvrBreF2=argminBreBre−BrF2+λtrBreTI+ωTvrTTvrBre.The matrix (I+ωTvrTTvr) is a positive definite matrix so that it can be expressed as factorization of the matrix ΓTΓ. Therefore, equation (18) is equivalent to equation (17) even though the explicit form of matrix Γ is unknown to us. To sum up, the response reconstruction described by equation (18) is an F-norm regularization-based optimization problem. This is one of the main contributions of this study.Finally, according to equation (19), the matrix form of structural responses Bre can be firstly calculated. Then, the matrix form of responses corresponding to the virtual senor can be estimated as (23)Bve=TvrBre.The estimated structural responsesBve can be reorganized into a form of vector. Herein, some sampling points may have more than one estimated result due to overlapping areas of moving time windows. To obtain a single identified value, the measured result at one sampling point will be selected as an average value of all the estimated values corresponding to the considered sampling point. In view of this, we can indirectly measure the structural responses corresponding to the virtual senor by using the measured responses obtained from some real sensors.


### 2.4. Brief Summary of the Proposed Method

An F-norm regularization-based method is proposed for structural response reconstruction. The proposed method is a model-based method, and transmissibility matrix is applied for formulating a governing equation. One highlight is that the transmissibility matrix applied in this study is calculated based on the basic theory of force identification under unknown initial conditions. Furthermore, problem of response reconstruction is considered as an optimization problem, in which the optimization variable is considered as a response matrix. F-norm regularization is applied for ensuring that the estimated results are not very sensitive to the measurement noises. Hence, another highlight is the application of matrix regularization. According to the above analysis, the basic process of this method can be summarized as Figure3.Figure 3
Flowchart of the proposed method for indirect measurement of structural responses.


## 2.1. Structural Responses Induced by Excitation Sources

A relationship between structural responses and dynamic force is firstly given as follows. Herein, the structural initial conditions are assumed to be equal to zero. A dictionary consisted ofNf atoms is introduced for expanding the force as [38, 39](1)f=d1,d2,⋯,dNfα1α2⋮αNf=Dα,where f = [f(Δt), f(2Δt), …, f(nΔt)]T means a force vector. Δt is a time interval for the forward analysis. D is the introduced dictionary and di represents the i-th atom of dictionary D. Herein, it assume that nΔt is equal to tw, which means that the discussions of structural response mainly locate in a time range (0, tw]. α is a coefficient vector, whose element αi is a factor corresponding to di. Nf represents the number of atoms in the dictionary D. Some sensors installed on the structure are considered. The structural responses induced by dynamic force can be expressed as(2)bif=hi1,hi2,⋯,hiNfα1α2⋮αNf=Hifα,i=1,2,⋯,Ns,where bif = [b(Δt̅), b(2Δt̅), …, b(kΔt̅)]T is a response vector of the i-th sensor. Δt̅ (Δt̅ ≥ Δt) represents the sampling time interval. Hif is a system matrix of the i-th sensor corresponding to the force. hij is the j-th component of structural responses caused by atom dj. Herein, the subscript “i” means the i-th sensor. These response components can be calculated by applying the mode superposition method. Ns represents the number of sensors.Let us focus on the analysis of structural responses caused by the structural initial conditions. Modal space is used for expressing the initial conditions. Structural dynamic behavior is assumed that it can be approximately described by using the firstNm orders. Then the initial conditions can be organized as a vector like [40](3)y0=yN10,y˙N10,yN20,y˙N20,⋯,yNNm0,y˙NNm0T,where yNi(0) and ẏNi(0) are the i-th mode initial displacement and i-th mode initial velocity, respectively. The corresponding mode shape matrix Φ = [φ1, φ2, …, φNm] satisfies ΦTMΦ = I. M represents the mass matrix of the structure. According to equation (3), the structural responses caused by initial conditions can be expressed via linear equations like(4)biy=h¯i1,h¯i2,h¯i3,⋯,h¯i2Nm=Hiyy0,i=1,2,⋯,Ns,where biy represents the structural responses caused by the initial conditions. Hiy is the system matrix corresponding to the initial conditions. The column h¯ij means the structural responses caused by the unit j-th initial condition component, i.e., y(0) = [0(1), …, 0(j-1), 1(j), 0(j+1), …, 0(2Nm)]T. Hence, the matrix Hiy actually can be calculated by applying the classical linear vibration theory. In the above equation, subscript “i” represents the i-th sensor.Now, if both external force and initial conditions are considered as the excitation sources, then structural responses obtained from thei-th sensor can be expressed as(5)bi=bif+biy=Hif,Hiyαy0=Hiβ,i=1,2,⋯,Ns,where bi means structural responses of the i-th sensor. Hi = [Hif, Hiy] is a system matrix corresponding to two excitation sources. β = [αT, y(0)T]T represents the structural excitation sources. Furthermore, if more than one sensor is considered, then equation (5) can be rewritten into one equation such as(6)b1b2⋮bNs=H1H2⋮HNsβ⇒b=Hβ.In addition, in order to reduce the possible large difference of numerical values among columns of system matrixH, equation (5) can be further rewritten as [40](7)b=Hβ=AΛβ,where A = [A1, A2, …, AN] is a normalized system matrix, in which the i-th atom Ai is calculated as Ai = Hi/||Hi||2. Λ = diag(||H1||2, ||H2||2, …, ||HN||2) is a diagonal scaling matrix. Hi represents the i-th column in matrix H.


## 2.2. Theoretical Derivation of Transmissibility Matrix

Two sets of sensors are considered. As shown in Figure1, all real sensors are grouped into the first set, while the virtual sensor is considered as the second set. The virtual sensor means that its corresponding structural responses will be indirectly reconstructed from the responses corresponding to the real sensors. According to equation (7), the structural responses corresponding to two sets of sensors can be expressed as(8)br=ArΛrβ,(9)bv=AvΛvβ,where the subscripts “r” and “v” represent the real sensors and virtual sensor, respectively. According to equation (8), the structural excitation sources (force and initial conditions) can be approximately estimated as(10)βe=Λr−1Ar+br,where βe represents the structural inputs identified from the measured structural responses br. Ar+ is a pseudoinverse matrix of Ar. Substituting equation (10) into equation (9) yields(11)bve=AvΛvΛr−1Ar+br=Tvrbr,where Tvr=AvΛvΛr−1Ar+ is a transmissibility matrix calculated from real sensor set to virtual sensor set. bve represents the structural responses corresponding to virtual sensor estimated from structural responses br.According to the above analysis, it can be seen that the transmissibility matrix is developed based on the theory of force identification under unknown initial conditions. The first step is to estimate the excitation sources, and the second step is to reconstruct the responses by using the identified excitation sources. In view of this, the key issue for this definition is that the real sensors should provide enough information for effectively identifying the unknown excitation sources. The number of real sensors is better to be bigger than the number of considered forces. Herein, it is beneficial to note again that both dynamic force and initial conditions are applied for the calculation of the transmissibility matrix. This definition is different from the one proposed by Zhu et al. [16] due to the latter one is calculated without considering unknown initial conditions.From equation (11), it can be seen that structural responses of virtual sensor can be indirectly calculated by using the structural responses obtained from real sensors. However, equation (11) is not recommended to be directly applied for the response reconstruction. The reason comes from the influences of measurement noises. In real application, the measurement noises are inevitable. The reconstructed responses directly estimated from equation (11) can be expressed as(12)bve=Tvrbr+Tvrε,where ε is a vector used for representing the measured noises. Essentially speaking, applying the transmissibility matrix to calculate the structural responses in the covert area is actually an inverse problem. In view of this, the second term such as Tvrε shown in equation (12) may lead to a large fluctuation in the estimated results so that the reconstructed responses may be unacceptable. Furthermore, from equation (11) it is not hard to find that the size of the transmissibility matrix is related to the number of sampling points. It will consume more computing resources for the calculation if a problem with long-time duration is considered. To overcome these drawbacks, structural responses will be reshaped in a matrix form, and then the F-norm regularization method will be introduced to ensure that the reconstructed responses are robust to the measured noises.


## 2.3. Response Reconstruction Based on F-Norm Regularization

As shown in Figure2, moving time window with a length of k sampling points is applied for extracting local signal of the structural responses [41]. The length of overlap between any two adjacent windows contains k0 (k0 < k) sampling points. The structural responses extracted from each time window are stored as a column of a structural response matrix. The structural response matrix can be expressed as(13)B=b1b1k−ko+1b2k−ko+1⋯bnw−1k−ko+1b2b1k−ko+2b2k−ko+2⋯bnw−1k−ko+2b3b1k−ko+3b2k−ko+3⋯bnw−1k−ko+3⋮⋮⋮⋱⋮bkb1k−ko+kb2k−ko+k⋯bnw−1k−ko+k,where B means the structural response matrix corresponding to one measured point. nw represents the number of moving time windows. bi represents the i-th sampling point. Taking the matrix form of structural responses into consideration, equation (11) can be rewritten as(14)Bve=TvrBr,where Bve is the estimated structural responses corresponding to the virtual sensor. Br means the structural responses corresponding to the real sensors. It can be expressed as(15)Br=B1TB2T⋯BNsTT,where Bi represents the structural responses corresponding to the i-th real sensor. From equation (14), it is not hard to find that the size of the transmissibility matrix is only determined by the length of moving time window. This characteristic means that a transmissibility matrix with relatively small size can be used for dealing with a problem of response reconstruction even though a problem of long-time duration is considered.Figure 2
Extracting local signal of structural responses via moving time window.As we have briefly discussed above, it is not recommended that the relationship between two sets of sensors is directly used for the reconstruction of structural responses due to influences of the measured noises. Herein, least square method is introduced for estimating the structural responses corresponding to the real sensors. This step can be expressed as an optimization problem like(16)Bre=argminBreTrrBre−BrF2,where ||∗||F means an F-norm of the considered matrix. Bre is a denoised response matrix. Trr means a transmissibility matrix calculated from the real sensors to themselves. According to the physical meaning of Trr, this matrix actually can be selected as an identity matrix I.Equation (16) can be used for estimating the responses Bre. However, it may cause an overfitting. Regularization technology is introduced for avoiding the overfitting problem. Herein, we can find that the estimated responses are organized in a form of matrix. Hence, the classical vector-based regularization method, whose estimated variable is organized as a vector, cannot be directly used in this case. Matrix regularization is utilized in this study. Herein, matrix regularization generalizes notions of the vector-based regularization to cases where the object to be estimated is a matrix. In view of this, matrix regularization such as F-norm regularization can be applied as(17)Bre=argminBreBre−BrF2+λΓBreF2,where λ is a regularization parameter. ||ΓBre||F2 is a penalty function, in which Γ is a smooth matrix. The smooth matrix Γ can be used for representing some reasonable prior knowledge of the estimated result. It can be often selected as an identity matrix I. However, this selection is not recommended because it does not take any information of Bve = TvrBre into consideration. In fact, it should be noted that what we want is the estimation of structural responses Bve. The structural responses Bre actually play a role of intermediate variable. Hence, it is better to take some reasonable prior knowledge of Bve = TvrBre into the definition of the penalty function. To achieve this goal, however, it may be difficult to formulate a specific matrix of Γ in an explicit form.According to the above consideration, the smooth matrixΓ will be defined in an implicit form that is similar to the one proposed by Pan et al. [42]. In view of this, equation (17) can be rewritten as(18)Bre=argminBreBre−BrF2+λBreF2+ωTvrBreF2,where the penalty function is defined as a weighted sum of BreF2 and BveF2. ω is a weighted factor. Herein, it is not hard to find that the former part of penalty function can avoid the overfitting problem corresponding to real sensors, while the latter part can make sure that signal energy of the reconstructed responses corresponding to the virtual sensor will not be too large.Let us calculate the derivative of function shown in equation (18) with respect to Bre, and let the corresponding result be equal to 0. Then, the estimated response matrix solved by equation (18) can be obtained as(19)Bre=I+λI+ωTvrTTvr−1Br.From equation (19), we can find that the matrix (I+ωTvrTTvr) is an additional matrix for the original problem. The spectrum distribution of the matrix (I+ωTvrTTvr) should be relatively centralized; otherwise, the problem shown in equation (19) will be ill-conditioned if a big regularization parameter is considered. In view of this, the values of ω should be limited. As done by Pan et al. [42], a ratio factor defined between the maximal and minimal eigenvalues of (I+ωTvrTTvr) is used for reflecting the spectrum distribution. This factor can be calculated as(20)r=1+ωemax1+ωemin,where emax and emin are the maximal and minimal eigenvalues of TvrTTvr, respectively. Denoting the maximal acceptable value of r as rmax, then we can calculate the maximal acceptable value of ω as(21)ωmax=rmax−1emax−rmaxemin,emax>rmaxemin,∞,emax≤rmaxemin,where ωmax represents the maximal acceptable value of ω. In this study, the value of rmax will be set as rmax = 109. This value means that the requirement for the concentration of eigenvalues of matrix (I+ωTvrTTvr) is not very strict. In this way, the selection of weighted factor ω will be more flexible. Furthermore, equation (18) can be also rewritten as(22)Bre=argminBreBre−BrF2+λBreF2+ωTvrBreF2=argminBreBre−BrF2+λtrBreTI+ωTvrTTvrBre.The matrix (I+ωTvrTTvr) is a positive definite matrix so that it can be expressed as factorization of the matrix ΓTΓ. Therefore, equation (18) is equivalent to equation (17) even though the explicit form of matrix Γ is unknown to us. To sum up, the response reconstruction described by equation (18) is an F-norm regularization-based optimization problem. This is one of the main contributions of this study.Finally, according to equation (19), the matrix form of structural responses Bre can be firstly calculated. Then, the matrix form of responses corresponding to the virtual senor can be estimated as (23)Bve=TvrBre.The estimated structural responsesBve can be reorganized into a form of vector. Herein, some sampling points may have more than one estimated result due to overlapping areas of moving time windows. To obtain a single identified value, the measured result at one sampling point will be selected as an average value of all the estimated values corresponding to the considered sampling point. In view of this, we can indirectly measure the structural responses corresponding to the virtual senor by using the measured responses obtained from some real sensors.


## 2.4. Brief Summary of the Proposed Method

An F-norm regularization-based method is proposed for structural response reconstruction. The proposed method is a model-based method, and transmissibility matrix is applied for formulating a governing equation. One highlight is that the transmissibility matrix applied in this study is calculated based on the basic theory of force identification under unknown initial conditions. Furthermore, problem of response reconstruction is considered as an optimization problem, in which the optimization variable is considered as a response matrix. F-norm regularization is applied for ensuring that the estimated results are not very sensitive to the measurement noises. Hence, another highlight is the application of matrix regularization. According to the above analysis, the basic process of this method can be summarized as Figure3.Figure 3
Flowchart of the proposed method for indirect measurement of structural responses.


## 3. Numerical Simulations

### 3.1. Plane Frame Structure

Figure4(a) shows the finite element model (FEM) of a 7-story plane frame structure. The frame structure is modeled by 49 beam elements. The beam element has 2 nodes and 6 degrees of freedom. The properties of the material are simulated as follows: Young’s modulus E = 35 GPa, density of the material ρ = 2500 kg/m3. The cross section area is equal to 0.32 m2, and the moment of inertia is set as 0.017 m4. Each element has a length of 2 m [43]. The first five nature frequencies are 2.43 Hz, 7.81 Hz, 14.72 Hz, 22.56 Hz, and 23.94 Hz, respectively. The frame structure is subjected to a dynamic force, as shown in Figure 4(b). The time history of this force is simulated as a random process. Its signal energy is located in a frequency range from 2 Hz to 50 Hz. Two real sensors are simulated to be installed at nodes 38 and 34 for measuring the horizontal acceleration responses. In this paper, the structural responses are calculated by using the mode superposition method. Theoretical modal analysis is firstly carried out. Then, the modal force for each order is calculated according to the corresponding modal shape. For each order, the modal responses are calculated via using the simplified discrete convolution model. To ensure calculation accuracy, the first ten orders are used for analysis. The damping ratio for each order is set as 0.012. The analysis time interval is set as 0.001 s. Finally, the structural responses are calculated as a weighted combination of the modal responses. White noises are considered by using an equation shown as(24)bn=b+lev×1n∑i=1nbi×rand,where bn and b represent two response vectors for noise and noiseless, respectively. n is the element number in vector b. lev is the noise level. rand is a standard normal distribution vector. Figures 4(b)–4(d) show the dynamic force and the simulated acceleration responses obtained via two sensors. Herein, the noise level is set as lev = 10%. Figure 4(e) shows the simulated acceleration responses corresponding to the virtual sensor installed at node 30. The sampling frequency is considered as 1000 Hz, and the total time length is set as 60 s.Figure 4
Frame structure with corresponding force and acceleration responses. (a) Seven-story frame structure. (b) Simulated force. (c) Acceleration responses of real sensor 1 with 10% noise level. (d) Acceleration responses of real sensor 2 with 10% noise level. (e) Acceleration responses of virtual sensor installed at node 30.
(a)(b)(c)(d)(e)


#### 3.1.1. Response Reconstruction from Two Sensors

The acceleration responses obtained from two real sensors are taken as the input data. It is better for the number of real sensors to be more than the number of external force, because unknown initial conditions are considered in the proposed method. Hence, more measured data are needed to make sure that equation (8) is not an underdetermined system of equation. The virtual sensor is installed at node 30 for measuring the lateral acceleration responses. The structural responses shown in Figure 4 are reshaped in a form of matrix. The length of time window is selected as k = 1000 sampling points. The length of overlap between any two adjacent windows is set as k0 = 500 sampling points. In view of this, it is not hard to find that the number of moving time windows is 119. Time domain discretization is applied for expressing the force. Hence, the i-th atom of dictionary D, as shown in equation (1), can be selected as(25)di=d1d2⋯di−1didi+1⋯dNfT=00⋯010⋯0T,i=1,2,⋯,Nf,where the number of considered basis functions Nf is equal to the number of sampling points contained in one moving time window. Herein, a key point for the selection of force atom is that the space constructed by the selected atoms should be complete enough to reflect the main information of the force.According to equation (21), the maximal acceptable value of ω can be calculated as ωmax = 7.1174. In this case, all sensors including the virtual senor are used for measuring the acceleration responses. There may be no significant numerical difference between responses corresponding to real sensors and virtual senor. An intuitive example can be seen from Figures 4(c)–4(e). Hence, we will treat every sensor equally so that the weighted factor is selected as ω = 1 < ωmax. In addition, the idea of a posteriori criterion is used for the selection of regularization parameter [44]. A function is firstly defined for roughly representing a ratio value of the estimated noise energy to the signal energy of estimated responses such as(26)RVSλ=Breλ−BrF2/NreBreλF2+TvrBreλF2/Nre+Nve×100%,where Bre(λ) is a function of λ and it represents the estimated structural responses corresponding to the real sensors. Nre and Nve are the numbers of elements contained in matrices Bre and TvrBre, respectively. RVS(λ) is a function whose independent variable is the regularization parameter λ. According to equation (26), different values of λ located in a range from 10−6 to 101 are used for calculation, and then their corresponding values of RVS are plotted in Figure 5. Two additional horizontal curves, which are used for representing the acceptable levels of the RVS, are added to the figure. As shown in Figure 5, it is assumed that the signal energy of the estimated noises roughly is located in a range from 0.5% to 5% of the signal energy of the estimated structural responses. Herein, it should be noted that the values of lower limit (0.5%) and upper limit (5%) are not the only determined values. These two values actually play a role of some prior information, and they are usually selected according to experience and actual situations. From Figure 5, we can roughly determine an acceptable range of λ. On this basis, the regularization parameter is finally selected as λ = 0.05. Herein, it is beneficial to note again that λ = 0.05 is not the only option, but a reasonable selection.Figure 5
Diagram for selection of regularization parameter.Figure6 shows the reconstructed structural responses corresponding to node 30. From Figure 6, it can be seen that the proposed method can effectively reconstruct the lateral acceleration responses in this case. There are some numerical differences between the estimated responses and real ones. This phenomenon is reasonable because influences of the measured noises are considered. In spite of this, the responses shown in Figures 6(b) and 6(c) intuitively show that the reconstructed responses match well the real ones on the whole. This result, to some extent, shows the feasibility of the proposed method in this case.Figure 6
Reconstructed structural responses corresponding to node 30. (a) Total time range. (b) Local time range from 8 s to 10 s. (c) Local time range from 40 s to 42 s.


#### 3.1.2. Effect of Different Regularization Parameters

The former case is further considered in this section. However, the value of regularization parameter is no longer selected via equation (26). To intuitively understand the effect of different regularization parameters, six cases with different values of regularization parameter, namely, 0, 5 × 10−6, 5 × 10−3, 5 × 10−2, 5 × 10−1, and 5 × 100, are simulated for discussion. Figure 7 shows the local structural responses calculated from different cases.Figure 7
Local structural responses estimated considering different regularization parameters. (a)λ = 0. (b) λ = 5 × 10−6. (c) λ = 5 × 10−3. (d) λ = 5 × 10−2. (e) λ = 5 × 10−1. (f) λ = 5 × 10−0.
(a)(b)(c)(d)(e)(f)Figure7(a) clearly shows that the estimated responses fluctuate seriously when penalty function is not considered in this case; i.e., λ = 0. The reason may be due to the ill-conditioning of the considered problem. On the contrary, Figure 7(d) shows that the estimated responses, which are obtained considering an appropriate value of λ such as λ = 5 × 10−2 for the penalty term, are obviously close to the real ones. This is mainly because the penalty term can effectively reduce the influences of measurement noises. This phenomenon is actually the main reason why we do not recommend reconstructing the structural responses by directly using the transmissibility relationship.Figure7(b) shows that the estimated responses have a large fluctuation when the regularization parameter is selected as λ = 5 × 10−6 in the considered case. This result is much the same as the one shown in Figure 7(a). The reason behind this phenomenon is that a small value of regularization parameter indicates that the governing equation pays less attention to the added penalty term so that the matrix regularization cannot improve the estimated result effectively. This drawback can be overcome when a larger value of regularization parameter is considered, for example, the estimated result shown in Figures 7(c) and 7(d). However, this does not mean that a larger regularization parameter will naturally lead to a better estimated result. A too large value of regularization parameter may lead to an underfitting problem. As a result, the estimated structural responses will become smoother than the real structural responses, as shown in Figures 7(e) and 7(f).Furthermore, different noise levels such as 5%, 10%, 15%, 20%, 25%, 30%, and 50% are considered. Thirteen values of the regularization parameter are selected for formulating the matrix regularization-based governing equation. These values are located in a range from 1 × 10−3 to 1 × 10−0, as listed in Table 1. In order to evaluate the identified accuracy of each case, a relative percentage error (RPE) is defined as (27)RPE=br−be1br1×100%,where br and be represent the real structural responses and the estimated structural responses, respectively. According to equation (27), the RPEs are calculated and listed in Table 1.Table 1
Relative percentage errors calculated with different regularization parameters.
Regularization parametersRelative percentage errors (%)5% noise10% noise15% noise20% noise25% noise30% noise50% noiseλ = 1 × 10−347.0492.58138.39184.27230.21276.16460.04λ = 2.5 × 10−329.5656.6284.06111.64139.27166.94277.80λ = 5 × 10−320.2936.3252.9669.8186.77103.78172.11λ = 7.5 × 10−317.1128.1640.0452.2164.5376.93126.93λ = 1 × 10−215.8523.9833.1142.6152.3062.08101.72λ = 2.5 × 10−216.2818.5321.6725.3329.3133.5051.26λ = 5 × 10−220.1720.8121.8223.1424.7126.4834.88λ = 7.5 × 10−223.9424.2424.7125.3626.1527.0831.88λ = 1 × 10−127.3927.5627.8428.2228.6929.2532.29λ = 2.5 × 10−142.7742.8042.8642.9343.0243.1343.73λ = 5 × 10−157.3957.4057.4257.4457.4757.5057.68λ = 7.5 × 10−165.9465.9565.9665.9765.9866.0066.08λ = 1 × 10−071.6071.6071.6071.6171.6271.6371.68Bold values in each column represent the minimum value.From Table1, three obvious phenomena can be found. The first one is that the estimated accuracy decreases with increasing the noise level in the considered cases if a specific regularization parameter is applied. This result shows that the noise level could bring negative effects to the estimated responses. The second one is that the RPE first falls and then rises with increasing the value of regularization parameter if a specific noise level is considered. The reason behind these results can be intuitively seen from Figure 7. That is, a small value of regularization parameter cannot avoid the appearance of overfitting, while a large one could cause an underfitting result. The third phenomenon is that the optimal regularization parameters are not exactly equal to each other among different conditions of noise level. Herein, the optimal regularization parameter represents a value of λ corresponding to the smallest RPE in the considered cases whose noise levels are the same. The optimal regularization parameter, on the whole, increases with increasing the noise level in the considered cases. This is reasonable since the reduction of greater measured noises needs a higher weight for the penalty term in the matrix regularization model. In addition, the above results also indicate that it is difficult to select an optimal regularization parameter if the noise level is unknown to us. This is the reason why a preliminary estimation of noise energy is required in this paper for selecting the regularization parameter, as shown in Figure 5.To sum up, the above results clearly show the importance of the selection of regularization parameter. This is actually a common problem in the application of regularization-based method. It is beneficial to note that the optimal selection of regularization parameter is still a key issue in this paper, even though there is a posteriori method proposed in this paper.


#### 3.1.3. Effect of Different Positions of Virtual Sensors

Three cases with different positions of the virtual sensor are considered. As shown in Figure4(a), three virtual sensors are installed at nodes 30, 36, and 40, respectively. The acceleration responses corresponding to each sensor are calculated independently. Three noise levels, namely, 5%, 10%, and 15%, are considered for each case. Regularization parameter λ for each case is selected as that we have done in Section 3.1.1. Figure 8 shows a diagram for the selection of regularization parameter for the cases with 15% noise levels considered. From Figure 8, it can be seen that λ = 0.05 can ensure that its corresponding value of RVS is intuitively located in a range from RVS = 0.5% to RVS = 5% in the case of node 30, while the regularization parameter playing a similar role in the other cases can be selected as λ = 0.075. Because of this, the regularization parameters are roughly selected as λ = 0.05 for the case of node 30 and λ = 0.075 for the other two cases. Other parameters are set the same as the ones that have been done in the previous section.Figure 8
Diagram for selection of regularization parameter in cases with 15% noise level considered.Table2 lists the RPEs of the structural responses estimated considering different positions of virtual sensors. From Table 2, it can be seen that the estimated accuracies, on the whole, are close to each other in the cases considering the same virtual sensor, even though the RPE increases with increasing the noise level. This phenomenon indicates again that the measured noises may bring negative effects to the estimated responses. However, the estimated result is not very sensitive to the measured noises in the considered cases. In addition, the RPEs corresponding to virtual sensors at node 36 and node 40 are smaller than the ones corresponding to virtual sensor at node 30. It seems that the position of virtual sensor could affect the accuracy of the estimated result in the considered cases. However, the reason behind this may be due to many factors; for example, it may be related to the selection of regularization parameter. As listed in Table 1, the RPE value can decrease to 15.85% if the regularization parameter is selected as λ = 1 × 10−2 for the case that is denoted as node 30 and 5% noise level.Table 2
Relative percentage errors calculated with different virtual sensors.
Positions of virtual sensorRelative percentage errors (%)5% noise10% noise15% noiseNode 3020.1720.8121.82Node 3614.5915.0215.70Node 4014.0414.1714.39For more detailed study, Figure9 shows the local responses calculated from different virtual sensors and noise levels. From the first row of Figure 9, it can be seen that the acceleration responses calculated considering different noise levels are all close to the real structural responses. This phenomenon, on the whole, means that the estimated result is not very sensitive to the measured noises in the case of virtual sensor installed at node 30. Similar results can also be seen from the other rows of Figure 9. This advantage mainly benefits from the penalty function such as a term of (BreF2+ωTvrBreF2). Furthermore, Figure 9 clearly shows that the proposed method can be used for reconstructing the acceleration responses at nodes 30, 36, and 40 with a good accuracy in these cases. Combining the estimated accuracies shown in Table 2, we can see that the proposed method is not very sensitive to the positions of the considered virtual sensors, although there are some differences between the RPEs in the cases considering different positions of the virtual sensors. The reason behind this result may be that structural responses obtained from two simulated real sensors can provide enough information for reflecting the signal features of the simulated force. In such situation, the structural responses can be naturally reconstructed with a relatively better accuracy even though different positions of the virtual sensors are considered in the above cases.Figure 9
Local responses calculated considering different virtual sensors and noise levels. (a) Node 30 and 5% noise level. (b) Node 30 and 10% noise level. (c) Node 30 and 15% noise level. (d) Node 36 and 5% noise level. (e) Node 36 and 10% noise level. (f) Node 36 and 15% noise level. (g) Node 40 and 5% noise level. (h) Node 40 and 10% noise level. (i) Node 40 and 15% noise level.
(a)(b)(c)(d)(e)(f)(g)(h)(i)


#### 3.1.4. Effect of Different Lengths of Moving Time Window

The length of moving time window is one of the significant parameters of the proposed method. This parameter directly relates to the size of the transmissibility matrixTvr. To get an intuitive understanding, five cases with different lengths of moving time windows, as listed in Table 3, are carried out. The lengths of moving time windows considered in Case 1, Case 2, Case 3, Case 4, and Case 5 actually correspond to 1 s, 2 s, 3 s, 4 s, and 5 s, respectively. The length of overlapping part between any two adjacent windows is auto-selected as half the length of the moving window. The regularization parameter is selected as λ = 0.05. Other parameters are selected as the ones that have been done in Section 3.1.1.Table 3
Comparison of calculation performance of cases with different lengths of time windows.
CasesNumber of sampling pointsSize of matricesCalculation performancesEach windowOverlap partTvrBreCalculation time (s)Ratio value (tcase i/tcase 1)Case 110005001000 × 20002000 × 1192.0040911.00Case 2200010002000 × 40004000 × 599.8422164.91Case 3300015003000 × 60006000 × 3929.99130814.97Case 4400020004000 × 80008000 × 2966.06360232.96Case 5500025005000 × 1000010000 × 23132.70321366.22All the cases are calculated via the same program with different calculation parameters and a personal computer with information CPU (AMD Ryzen 7 2700X Eight-Core Processor 3.70 GHz), RAM 16.0 GB, Windows 10, Matlab R2015a. The analysis times for five cases are gathered and listed in Table3. From Table 3, it is not hard to find that the computational time increases rapidly with increasing the length of the moving time window in the considered cases. The direct cause is that a longer length of moving time window actually leads to a larger size of transmissibility matrix, as listed in Table 3. In addition, Case 1 and Case 5 are picked up as two representative cases. Comparison of the reconstructed responses obtained from these two cases is plotted in Figure 10. From Figure 10, it can be seen that there are no significant differences between two reconstructed responses on the whole. Overall, Case 1 can spend less calculation time on reconstructing the acceleration responses. Therefore, with a comprehensive consideration of calculation time and reconstructed accuracy, it is not hard to find that the case with window length of 1000 sampling points is the best solution way among the five cases listed in Table 3. To some extent, the advantage of Case 1 actually benefits from the consideration of unknown initial conditions. It is because this consideration is the theoretical basis for the formulation of equation (18), which is a governing equation expressed in a matrix form.Figure 10
Reconstructed responses at node 30 with different lengths of moving time windows. (a) Total time range. (b) Local time range from 19 s to 21 s.Furthermore, different noise levels such as 5%, 10%, and 15% are considered. The RPEs are calculated and listed in Table4. From Table 4, two phenomena can be obtained. The first one is that the measured noises could bring negative effects to the estimated results in the considered cases. This phenomenon is the same as the one discussed in the former sections. The second phenomenon is that the estimated accuracy increases with increasing the length of time windows in the considered cases. In fact, the estimated errors at the beginning of each moving time window may be larger than the ones at other zones due to the interference of structural initial conditions. Hence, if more moving time windows are applied for covering the total time range, the estimated responses may have a lower accuracy. To provide a more intuitive understanding, Case 1 and Case 5 are picked up as two representative cases, and their estimated errors between the estimated responses and the real responses are plotted in Figure 11. From Figure 11, it can be seen that more spikes appear in the estimated errors in Case 1. These spikes actually correspond to the beginning of each time window. That is why the coordinate values of the spikes shown in Figure 11(b) are close to the integral multiple of 0.5 seconds, while the coordinate values of the spikes shown in Figure 11(d) are near to the integral multiple of 2.5 seconds. The above phenomenon is actually the direct reason why estimated accuracy of Case 5 is better than that of Case 1.Table 4
Relative percentage errors calculated with different lengths of time windows and noise levels.
CasesNumber of sampling pointsRelative percentage errors (%)Each windowOverlap part5% noise10% noise15% noiseCase 1100050020.1720.8121.82Case 22000100015.7716.5517.76Case 33000150014.0214.9016.22Case 44000200013.2914.1815.54Case 55000250013.0013.9115.30Figure 11
Estimated errors between estimated responses and real responses with 10% noise level considered. (a) Total length in Case 1. (b) Local length in Case 1. (c) Total length in Case 5. (d) Local length in Case 5.


### 3.2. Further Verification on Plane Truss Structure

A plan truss structure is simulated for further verifying the proposed method. As shown in Figure12, the truss is simply supported with pin joint at one end and roller support at the other end. Both horizontal and vertical members have a length of 1.0 m. The cross section area is 0.0016 m2. The density and Young’s modulus are 7850 kg·m−3 and 2.06 GPa, respectively [45]. Lumped mass matrix is considered. The first five nature frequencies are 1.56 Hz, 5.17 Hz, 6.94 Hz, 11.98 Hz, and 17.72 Hz, respectively. One dynamic force is considered, and its time history is simulated as(28)ft=501−cos10πtsin30πt.Figure 12
Plan truss for numerical simulations.Four real sensors installed at nodes 6, 8, 10, and 12 are used for measuring the vertical acceleration responses. The sampling frequency is set as 1000 Hz. The structural responses are calculated using the mode superposition method and considering the first fifteen modal orders. A basic process for simulating the structural responses is the same as the one described in Section3.1. Three noise levels such as 5%, 10%, and 15% levels are considered. The discrete time interval is selected as 0.001 s, and the total sampling time length is set as 10 s. The damping ratio for each order is selected as 0.012. The initial conditions are set as zeros.The objective is to reconstruct the acceleration responses at nodes 15, 17, and 19, respectively. Effects of different numbers of real sensors are considered. For more details, structural responses corresponding to each virtual sensor will be estimated by using structural responses measured from two sensors at nodes 6 and 8; three sensors at nodes 6, 8, and 10; and four sensors at nodes 6, 8, 10, and 12, respectively. Hence, it is not hard to find that there are twenty-seven cases simulated in this section. Acceleration responses are reshaped in a matrix form by using parameters such as length of time window of 1000 sampling points and length of overlap part of 500 sampling points. The weighted factor is set asω = 1. It should be noted that ω = 1 satisfies the condition ω < ωmax in every considered case. The regularization parameter λ is carefully selected as the one that has been done in the previous section. Herein, the lower limit value of RVS is reduced appropriately, and finally λ will be selected as λ = 0.05.According to equation (27), the RPEs corresponding to nine cases are calculated and listed in Table 5. From Table 5, it can be seen that the effect of the measured noises is almost the same as the ones discussed in the previous sections. In addition, Table 5 also shows that the estimated accuracy of each virtual sensor increases with increasing the real sensors in the simulated cases. The reason may be that more sensors can provide more information for reflecting the external force so that the reconstructed acceleration responses are closer to the real ones. However, more sensors will lead to a larger size of transmissibility matrix. Along with this, more computing resources are needed. In view of this, a good balance between the calculation accuracy and efficiency should be carefully considered in practical applications. In addition, Figure 13 shows the time histories of reconstructed acceleration responses that are calculated from nodes 6, 8, 10, and 12 considering 15% noise level. From this figure, it can be seen that the reconstructed responses match pretty well the real ones for each virtual sensor. This result means that the proposed method is feasible and effective in the considered cases to some extent.Table 5
Relative percentage errors calculated from different sets of real sensors.
Real sensorsRPEs at node 15 (%)RPEs at node 17 (%)RPEs at node 19 (%)5% noise10% noise15% noise5% noise10% noise15% noise5% noise10% noise15% noiseNodes 6 and 810.1110.7511.7215.0015.2215.548.9111.6815.20Nodes 6, 8, and 106.917.117.439.689.8410.066.288.1910.63Nodes 6, 8, 10, and 125.996.146.386.766.867.005.526.918.77Figure 13
Reconstructed responses calculated from nodes 6, 8, 10, and 12 considering 15% noise level. (a) Total time range and node 15. (b) Local time range and node 15. (c) Total time range and node 17. (d) Local time range and node 17. (e) Total time range and node 19. (f) Local time range and node 19.Overall, the above discussions show that the proposed method can be used for indirectly measuring the acceleration responses corresponding to the considered virtual sensors. The proposed method can reduce influences of the measured noises to some extent. Hence, the estimated responses shown in this section also have good robustness.


## 3.1. Plane Frame Structure

Figure4(a) shows the finite element model (FEM) of a 7-story plane frame structure. The frame structure is modeled by 49 beam elements. The beam element has 2 nodes and 6 degrees of freedom. The properties of the material are simulated as follows: Young’s modulus E = 35 GPa, density of the material ρ = 2500 kg/m3. The cross section area is equal to 0.32 m2, and the moment of inertia is set as 0.017 m4. Each element has a length of 2 m [43]. The first five nature frequencies are 2.43 Hz, 7.81 Hz, 14.72 Hz, 22.56 Hz, and 23.94 Hz, respectively. The frame structure is subjected to a dynamic force, as shown in Figure 4(b). The time history of this force is simulated as a random process. Its signal energy is located in a frequency range from 2 Hz to 50 Hz. Two real sensors are simulated to be installed at nodes 38 and 34 for measuring the horizontal acceleration responses. In this paper, the structural responses are calculated by using the mode superposition method. Theoretical modal analysis is firstly carried out. Then, the modal force for each order is calculated according to the corresponding modal shape. For each order, the modal responses are calculated via using the simplified discrete convolution model. To ensure calculation accuracy, the first ten orders are used for analysis. The damping ratio for each order is set as 0.012. The analysis time interval is set as 0.001 s. Finally, the structural responses are calculated as a weighted combination of the modal responses. White noises are considered by using an equation shown as(24)bn=b+lev×1n∑i=1nbi×rand,where bn and b represent two response vectors for noise and noiseless, respectively. n is the element number in vector b. lev is the noise level. rand is a standard normal distribution vector. Figures 4(b)–4(d) show the dynamic force and the simulated acceleration responses obtained via two sensors. Herein, the noise level is set as lev = 10%. Figure 4(e) shows the simulated acceleration responses corresponding to the virtual sensor installed at node 30. The sampling frequency is considered as 1000 Hz, and the total time length is set as 60 s.Figure 4
Frame structure with corresponding force and acceleration responses. (a) Seven-story frame structure. (b) Simulated force. (c) Acceleration responses of real sensor 1 with 10% noise level. (d) Acceleration responses of real sensor 2 with 10% noise level. (e) Acceleration responses of virtual sensor installed at node 30.
(a)(b)(c)(d)(e)


### 3.1.1. Response Reconstruction from Two Sensors

The acceleration responses obtained from two real sensors are taken as the input data. It is better for the number of real sensors to be more than the number of external force, because unknown initial conditions are considered in the proposed method. Hence, more measured data are needed to make sure that equation (8) is not an underdetermined system of equation. The virtual sensor is installed at node 30 for measuring the lateral acceleration responses. The structural responses shown in Figure 4 are reshaped in a form of matrix. The length of time window is selected as k = 1000 sampling points. The length of overlap between any two adjacent windows is set as k0 = 500 sampling points. In view of this, it is not hard to find that the number of moving time windows is 119. Time domain discretization is applied for expressing the force. Hence, the i-th atom of dictionary D, as shown in equation (1), can be selected as(25)di=d1d2⋯di−1didi+1⋯dNfT=00⋯010⋯0T,i=1,2,⋯,Nf,where the number of considered basis functions Nf is equal to the number of sampling points contained in one moving time window. Herein, a key point for the selection of force atom is that the space constructed by the selected atoms should be complete enough to reflect the main information of the force.According to equation (21), the maximal acceptable value of ω can be calculated as ωmax = 7.1174. In this case, all sensors including the virtual senor are used for measuring the acceleration responses. There may be no significant numerical difference between responses corresponding to real sensors and virtual senor. An intuitive example can be seen from Figures 4(c)–4(e). Hence, we will treat every sensor equally so that the weighted factor is selected as ω = 1 < ωmax. In addition, the idea of a posteriori criterion is used for the selection of regularization parameter [44]. A function is firstly defined for roughly representing a ratio value of the estimated noise energy to the signal energy of estimated responses such as(26)RVSλ=Breλ−BrF2/NreBreλF2+TvrBreλF2/Nre+Nve×100%,where Bre(λ) is a function of λ and it represents the estimated structural responses corresponding to the real sensors. Nre and Nve are the numbers of elements contained in matrices Bre and TvrBre, respectively. RVS(λ) is a function whose independent variable is the regularization parameter λ. According to equation (26), different values of λ located in a range from 10−6 to 101 are used for calculation, and then their corresponding values of RVS are plotted in Figure 5. Two additional horizontal curves, which are used for representing the acceptable levels of the RVS, are added to the figure. As shown in Figure 5, it is assumed that the signal energy of the estimated noises roughly is located in a range from 0.5% to 5% of the signal energy of the estimated structural responses. Herein, it should be noted that the values of lower limit (0.5%) and upper limit (5%) are not the only determined values. These two values actually play a role of some prior information, and they are usually selected according to experience and actual situations. From Figure 5, we can roughly determine an acceptable range of λ. On this basis, the regularization parameter is finally selected as λ = 0.05. Herein, it is beneficial to note again that λ = 0.05 is not the only option, but a reasonable selection.Figure 5
Diagram for selection of regularization parameter.Figure6 shows the reconstructed structural responses corresponding to node 30. From Figure 6, it can be seen that the proposed method can effectively reconstruct the lateral acceleration responses in this case. There are some numerical differences between the estimated responses and real ones. This phenomenon is reasonable because influences of the measured noises are considered. In spite of this, the responses shown in Figures 6(b) and 6(c) intuitively show that the reconstructed responses match well the real ones on the whole. This result, to some extent, shows the feasibility of the proposed method in this case.Figure 6
Reconstructed structural responses corresponding to node 30. (a) Total time range. (b) Local time range from 8 s to 10 s. (c) Local time range from 40 s to 42 s.


### 3.1.2. Effect of Different Regularization Parameters

The former case is further considered in this section. However, the value of regularization parameter is no longer selected via equation (26). To intuitively understand the effect of different regularization parameters, six cases with different values of regularization parameter, namely, 0, 5 × 10−6, 5 × 10−3, 5 × 10−2, 5 × 10−1, and 5 × 100, are simulated for discussion. Figure 7 shows the local structural responses calculated from different cases.Figure 7
Local structural responses estimated considering different regularization parameters. (a)λ = 0. (b) λ = 5 × 10−6. (c) λ = 5 × 10−3. (d) λ = 5 × 10−2. (e) λ = 5 × 10−1. (f) λ = 5 × 10−0.
(a)(b)(c)(d)(e)(f)Figure7(a) clearly shows that the estimated responses fluctuate seriously when penalty function is not considered in this case; i.e., λ = 0. The reason may be due to the ill-conditioning of the considered problem. On the contrary, Figure 7(d) shows that the estimated responses, which are obtained considering an appropriate value of λ such as λ = 5 × 10−2 for the penalty term, are obviously close to the real ones. This is mainly because the penalty term can effectively reduce the influences of measurement noises. This phenomenon is actually the main reason why we do not recommend reconstructing the structural responses by directly using the transmissibility relationship.Figure7(b) shows that the estimated responses have a large fluctuation when the regularization parameter is selected as λ = 5 × 10−6 in the considered case. This result is much the same as the one shown in Figure 7(a). The reason behind this phenomenon is that a small value of regularization parameter indicates that the governing equation pays less attention to the added penalty term so that the matrix regularization cannot improve the estimated result effectively. This drawback can be overcome when a larger value of regularization parameter is considered, for example, the estimated result shown in Figures 7(c) and 7(d). However, this does not mean that a larger regularization parameter will naturally lead to a better estimated result. A too large value of regularization parameter may lead to an underfitting problem. As a result, the estimated structural responses will become smoother than the real structural responses, as shown in Figures 7(e) and 7(f).Furthermore, different noise levels such as 5%, 10%, 15%, 20%, 25%, 30%, and 50% are considered. Thirteen values of the regularization parameter are selected for formulating the matrix regularization-based governing equation. These values are located in a range from 1 × 10−3 to 1 × 10−0, as listed in Table 1. In order to evaluate the identified accuracy of each case, a relative percentage error (RPE) is defined as (27)RPE=br−be1br1×100%,where br and be represent the real structural responses and the estimated structural responses, respectively. According to equation (27), the RPEs are calculated and listed in Table 1.Table 1
Relative percentage errors calculated with different regularization parameters.
Regularization parametersRelative percentage errors (%)5% noise10% noise15% noise20% noise25% noise30% noise50% noiseλ = 1 × 10−347.0492.58138.39184.27230.21276.16460.04λ = 2.5 × 10−329.5656.6284.06111.64139.27166.94277.80λ = 5 × 10−320.2936.3252.9669.8186.77103.78172.11λ = 7.5 × 10−317.1128.1640.0452.2164.5376.93126.93λ = 1 × 10−215.8523.9833.1142.6152.3062.08101.72λ = 2.5 × 10−216.2818.5321.6725.3329.3133.5051.26λ = 5 × 10−220.1720.8121.8223.1424.7126.4834.88λ = 7.5 × 10−223.9424.2424.7125.3626.1527.0831.88λ = 1 × 10−127.3927.5627.8428.2228.6929.2532.29λ = 2.5 × 10−142.7742.8042.8642.9343.0243.1343.73λ = 5 × 10−157.3957.4057.4257.4457.4757.5057.68λ = 7.5 × 10−165.9465.9565.9665.9765.9866.0066.08λ = 1 × 10−071.6071.6071.6071.6171.6271.6371.68Bold values in each column represent the minimum value.From Table1, three obvious phenomena can be found. The first one is that the estimated accuracy decreases with increasing the noise level in the considered cases if a specific regularization parameter is applied. This result shows that the noise level could bring negative effects to the estimated responses. The second one is that the RPE first falls and then rises with increasing the value of regularization parameter if a specific noise level is considered. The reason behind these results can be intuitively seen from Figure 7. That is, a small value of regularization parameter cannot avoid the appearance of overfitting, while a large one could cause an underfitting result. The third phenomenon is that the optimal regularization parameters are not exactly equal to each other among different conditions of noise level. Herein, the optimal regularization parameter represents a value of λ corresponding to the smallest RPE in the considered cases whose noise levels are the same. The optimal regularization parameter, on the whole, increases with increasing the noise level in the considered cases. This is reasonable since the reduction of greater measured noises needs a higher weight for the penalty term in the matrix regularization model. In addition, the above results also indicate that it is difficult to select an optimal regularization parameter if the noise level is unknown to us. This is the reason why a preliminary estimation of noise energy is required in this paper for selecting the regularization parameter, as shown in Figure 5.To sum up, the above results clearly show the importance of the selection of regularization parameter. This is actually a common problem in the application of regularization-based method. It is beneficial to note that the optimal selection of regularization parameter is still a key issue in this paper, even though there is a posteriori method proposed in this paper.


### 3.1.3. Effect of Different Positions of Virtual Sensors

Three cases with different positions of the virtual sensor are considered. As shown in Figure4(a), three virtual sensors are installed at nodes 30, 36, and 40, respectively. The acceleration responses corresponding to each sensor are calculated independently. Three noise levels, namely, 5%, 10%, and 15%, are considered for each case. Regularization parameter λ for each case is selected as that we have done in Section 3.1.1. Figure 8 shows a diagram for the selection of regularization parameter for the cases with 15% noise levels considered. From Figure 8, it can be seen that λ = 0.05 can ensure that its corresponding value of RVS is intuitively located in a range from RVS = 0.5% to RVS = 5% in the case of node 30, while the regularization parameter playing a similar role in the other cases can be selected as λ = 0.075. Because of this, the regularization parameters are roughly selected as λ = 0.05 for the case of node 30 and λ = 0.075 for the other two cases. Other parameters are set the same as the ones that have been done in the previous section.Figure 8
Diagram for selection of regularization parameter in cases with 15% noise level considered.Table2 lists the RPEs of the structural responses estimated considering different positions of virtual sensors. From Table 2, it can be seen that the estimated accuracies, on the whole, are close to each other in the cases considering the same virtual sensor, even though the RPE increases with increasing the noise level. This phenomenon indicates again that the measured noises may bring negative effects to the estimated responses. However, the estimated result is not very sensitive to the measured noises in the considered cases. In addition, the RPEs corresponding to virtual sensors at node 36 and node 40 are smaller than the ones corresponding to virtual sensor at node 30. It seems that the position of virtual sensor could affect the accuracy of the estimated result in the considered cases. However, the reason behind this may be due to many factors; for example, it may be related to the selection of regularization parameter. As listed in Table 1, the RPE value can decrease to 15.85% if the regularization parameter is selected as λ = 1 × 10−2 for the case that is denoted as node 30 and 5% noise level.Table 2
Relative percentage errors calculated with different virtual sensors.
Positions of virtual sensorRelative percentage errors (%)5% noise10% noise15% noiseNode 3020.1720.8121.82Node 3614.5915.0215.70Node 4014.0414.1714.39For more detailed study, Figure9 shows the local responses calculated from different virtual sensors and noise levels. From the first row of Figure 9, it can be seen that the acceleration responses calculated considering different noise levels are all close to the real structural responses. This phenomenon, on the whole, means that the estimated result is not very sensitive to the measured noises in the case of virtual sensor installed at node 30. Similar results can also be seen from the other rows of Figure 9. This advantage mainly benefits from the penalty function such as a term of (BreF2+ωTvrBreF2). Furthermore, Figure 9 clearly shows that the proposed method can be used for reconstructing the acceleration responses at nodes 30, 36, and 40 with a good accuracy in these cases. Combining the estimated accuracies shown in Table 2, we can see that the proposed method is not very sensitive to the positions of the considered virtual sensors, although there are some differences between the RPEs in the cases considering different positions of the virtual sensors. The reason behind this result may be that structural responses obtained from two simulated real sensors can provide enough information for reflecting the signal features of the simulated force. In such situation, the structural responses can be naturally reconstructed with a relatively better accuracy even though different positions of the virtual sensors are considered in the above cases.Figure 9
Local responses calculated considering different virtual sensors and noise levels. (a) Node 30 and 5% noise level. (b) Node 30 and 10% noise level. (c) Node 30 and 15% noise level. (d) Node 36 and 5% noise level. (e) Node 36 and 10% noise level. (f) Node 36 and 15% noise level. (g) Node 40 and 5% noise level. (h) Node 40 and 10% noise level. (i) Node 40 and 15% noise level.
(a)(b)(c)(d)(e)(f)(g)(h)(i)


### 3.1.4. Effect of Different Lengths of Moving Time Window

The length of moving time window is one of the significant parameters of the proposed method. This parameter directly relates to the size of the transmissibility matrixTvr. To get an intuitive understanding, five cases with different lengths of moving time windows, as listed in Table 3, are carried out. The lengths of moving time windows considered in Case 1, Case 2, Case 3, Case 4, and Case 5 actually correspond to 1 s, 2 s, 3 s, 4 s, and 5 s, respectively. The length of overlapping part between any two adjacent windows is auto-selected as half the length of the moving window. The regularization parameter is selected as λ = 0.05. Other parameters are selected as the ones that have been done in Section 3.1.1.Table 3
Comparison of calculation performance of cases with different lengths of time windows.
CasesNumber of sampling pointsSize of matricesCalculation performancesEach windowOverlap partTvrBreCalculation time (s)Ratio value (tcase i/tcase 1)Case 110005001000 × 20002000 × 1192.0040911.00Case 2200010002000 × 40004000 × 599.8422164.91Case 3300015003000 × 60006000 × 3929.99130814.97Case 4400020004000 × 80008000 × 2966.06360232.96Case 5500025005000 × 1000010000 × 23132.70321366.22All the cases are calculated via the same program with different calculation parameters and a personal computer with information CPU (AMD Ryzen 7 2700X Eight-Core Processor 3.70 GHz), RAM 16.0 GB, Windows 10, Matlab R2015a. The analysis times for five cases are gathered and listed in Table3. From Table 3, it is not hard to find that the computational time increases rapidly with increasing the length of the moving time window in the considered cases. The direct cause is that a longer length of moving time window actually leads to a larger size of transmissibility matrix, as listed in Table 3. In addition, Case 1 and Case 5 are picked up as two representative cases. Comparison of the reconstructed responses obtained from these two cases is plotted in Figure 10. From Figure 10, it can be seen that there are no significant differences between two reconstructed responses on the whole. Overall, Case 1 can spend less calculation time on reconstructing the acceleration responses. Therefore, with a comprehensive consideration of calculation time and reconstructed accuracy, it is not hard to find that the case with window length of 1000 sampling points is the best solution way among the five cases listed in Table 3. To some extent, the advantage of Case 1 actually benefits from the consideration of unknown initial conditions. It is because this consideration is the theoretical basis for the formulation of equation (18), which is a governing equation expressed in a matrix form.Figure 10
Reconstructed responses at node 30 with different lengths of moving time windows. (a) Total time range. (b) Local time range from 19 s to 21 s.Furthermore, different noise levels such as 5%, 10%, and 15% are considered. The RPEs are calculated and listed in Table4. From Table 4, two phenomena can be obtained. The first one is that the measured noises could bring negative effects to the estimated results in the considered cases. This phenomenon is the same as the one discussed in the former sections. The second phenomenon is that the estimated accuracy increases with increasing the length of time windows in the considered cases. In fact, the estimated errors at the beginning of each moving time window may be larger than the ones at other zones due to the interference of structural initial conditions. Hence, if more moving time windows are applied for covering the total time range, the estimated responses may have a lower accuracy. To provide a more intuitive understanding, Case 1 and Case 5 are picked up as two representative cases, and their estimated errors between the estimated responses and the real responses are plotted in Figure 11. From Figure 11, it can be seen that more spikes appear in the estimated errors in Case 1. These spikes actually correspond to the beginning of each time window. That is why the coordinate values of the spikes shown in Figure 11(b) are close to the integral multiple of 0.5 seconds, while the coordinate values of the spikes shown in Figure 11(d) are near to the integral multiple of 2.5 seconds. The above phenomenon is actually the direct reason why estimated accuracy of Case 5 is better than that of Case 1.Table 4
Relative percentage errors calculated with different lengths of time windows and noise levels.
CasesNumber of sampling pointsRelative percentage errors (%)Each windowOverlap part5% noise10% noise15% noiseCase 1100050020.1720.8121.82Case 22000100015.7716.5517.76Case 33000150014.0214.9016.22Case 44000200013.2914.1815.54Case 55000250013.0013.9115.30Figure 11
Estimated errors between estimated responses and real responses with 10% noise level considered. (a) Total length in Case 1. (b) Local length in Case 1. (c) Total length in Case 5. (d) Local length in Case 5.


## 3.1.1. Response Reconstruction from Two Sensors

The acceleration responses obtained from two real sensors are taken as the input data. It is better for the number of real sensors to be more than the number of external force, because unknown initial conditions are considered in the proposed method. Hence, more measured data are needed to make sure that equation (8) is not an underdetermined system of equation. The virtual sensor is installed at node 30 for measuring the lateral acceleration responses. The structural responses shown in Figure 4 are reshaped in a form of matrix. The length of time window is selected as k = 1000 sampling points. The length of overlap between any two adjacent windows is set as k0 = 500 sampling points. In view of this, it is not hard to find that the number of moving time windows is 119. Time domain discretization is applied for expressing the force. Hence, the i-th atom of dictionary D, as shown in equation (1), can be selected as(25)di=d1d2⋯di−1didi+1⋯dNfT=00⋯010⋯0T,i=1,2,⋯,Nf,where the number of considered basis functions Nf is equal to the number of sampling points contained in one moving time window. Herein, a key point for the selection of force atom is that the space constructed by the selected atoms should be complete enough to reflect the main information of the force.According to equation (21), the maximal acceptable value of ω can be calculated as ωmax = 7.1174. In this case, all sensors including the virtual senor are used for measuring the acceleration responses. There may be no significant numerical difference between responses corresponding to real sensors and virtual senor. An intuitive example can be seen from Figures 4(c)–4(e). Hence, we will treat every sensor equally so that the weighted factor is selected as ω = 1 < ωmax. In addition, the idea of a posteriori criterion is used for the selection of regularization parameter [44]. A function is firstly defined for roughly representing a ratio value of the estimated noise energy to the signal energy of estimated responses such as(26)RVSλ=Breλ−BrF2/NreBreλF2+TvrBreλF2/Nre+Nve×100%,where Bre(λ) is a function of λ and it represents the estimated structural responses corresponding to the real sensors. Nre and Nve are the numbers of elements contained in matrices Bre and TvrBre, respectively. RVS(λ) is a function whose independent variable is the regularization parameter λ. According to equation (26), different values of λ located in a range from 10−6 to 101 are used for calculation, and then their corresponding values of RVS are plotted in Figure 5. Two additional horizontal curves, which are used for representing the acceptable levels of the RVS, are added to the figure. As shown in Figure 5, it is assumed that the signal energy of the estimated noises roughly is located in a range from 0.5% to 5% of the signal energy of the estimated structural responses. Herein, it should be noted that the values of lower limit (0.5%) and upper limit (5%) are not the only determined values. These two values actually play a role of some prior information, and they are usually selected according to experience and actual situations. From Figure 5, we can roughly determine an acceptable range of λ. On this basis, the regularization parameter is finally selected as λ = 0.05. Herein, it is beneficial to note again that λ = 0.05 is not the only option, but a reasonable selection.Figure 5
Diagram for selection of regularization parameter.Figure6 shows the reconstructed structural responses corresponding to node 30. From Figure 6, it can be seen that the proposed method can effectively reconstruct the lateral acceleration responses in this case. There are some numerical differences between the estimated responses and real ones. This phenomenon is reasonable because influences of the measured noises are considered. In spite of this, the responses shown in Figures 6(b) and 6(c) intuitively show that the reconstructed responses match well the real ones on the whole. This result, to some extent, shows the feasibility of the proposed method in this case.Figure 6
Reconstructed structural responses corresponding to node 30. (a) Total time range. (b) Local time range from 8 s to 10 s. (c) Local time range from 40 s to 42 s.


## 3.1.2. Effect of Different Regularization Parameters

The former case is further considered in this section. However, the value of regularization parameter is no longer selected via equation (26). To intuitively understand the effect of different regularization parameters, six cases with different values of regularization parameter, namely, 0, 5 × 10−6, 5 × 10−3, 5 × 10−2, 5 × 10−1, and 5 × 100, are simulated for discussion. Figure 7 shows the local structural responses calculated from different cases.Figure 7
Local structural responses estimated considering different regularization parameters. (a)λ = 0. (b) λ = 5 × 10−6. (c) λ = 5 × 10−3. (d) λ = 5 × 10−2. (e) λ = 5 × 10−1. (f) λ = 5 × 10−0.
(a)(b)(c)(d)(e)(f)Figure7(a) clearly shows that the estimated responses fluctuate seriously when penalty function is not considered in this case; i.e., λ = 0. The reason may be due to the ill-conditioning of the considered problem. On the contrary, Figure 7(d) shows that the estimated responses, which are obtained considering an appropriate value of λ such as λ = 5 × 10−2 for the penalty term, are obviously close to the real ones. This is mainly because the penalty term can effectively reduce the influences of measurement noises. This phenomenon is actually the main reason why we do not recommend reconstructing the structural responses by directly using the transmissibility relationship.Figure7(b) shows that the estimated responses have a large fluctuation when the regularization parameter is selected as λ = 5 × 10−6 in the considered case. This result is much the same as the one shown in Figure 7(a). The reason behind this phenomenon is that a small value of regularization parameter indicates that the governing equation pays less attention to the added penalty term so that the matrix regularization cannot improve the estimated result effectively. This drawback can be overcome when a larger value of regularization parameter is considered, for example, the estimated result shown in Figures 7(c) and 7(d). However, this does not mean that a larger regularization parameter will naturally lead to a better estimated result. A too large value of regularization parameter may lead to an underfitting problem. As a result, the estimated structural responses will become smoother than the real structural responses, as shown in Figures 7(e) and 7(f).Furthermore, different noise levels such as 5%, 10%, 15%, 20%, 25%, 30%, and 50% are considered. Thirteen values of the regularization parameter are selected for formulating the matrix regularization-based governing equation. These values are located in a range from 1 × 10−3 to 1 × 10−0, as listed in Table 1. In order to evaluate the identified accuracy of each case, a relative percentage error (RPE) is defined as (27)RPE=br−be1br1×100%,where br and be represent the real structural responses and the estimated structural responses, respectively. According to equation (27), the RPEs are calculated and listed in Table 1.Table 1
Relative percentage errors calculated with different regularization parameters.
Regularization parametersRelative percentage errors (%)5% noise10% noise15% noise20% noise25% noise30% noise50% noiseλ = 1 × 10−347.0492.58138.39184.27230.21276.16460.04λ = 2.5 × 10−329.5656.6284.06111.64139.27166.94277.80λ = 5 × 10−320.2936.3252.9669.8186.77103.78172.11λ = 7.5 × 10−317.1128.1640.0452.2164.5376.93126.93λ = 1 × 10−215.8523.9833.1142.6152.3062.08101.72λ = 2.5 × 10−216.2818.5321.6725.3329.3133.5051.26λ = 5 × 10−220.1720.8121.8223.1424.7126.4834.88λ = 7.5 × 10−223.9424.2424.7125.3626.1527.0831.88λ = 1 × 10−127.3927.5627.8428.2228.6929.2532.29λ = 2.5 × 10−142.7742.8042.8642.9343.0243.1343.73λ = 5 × 10−157.3957.4057.4257.4457.4757.5057.68λ = 7.5 × 10−165.9465.9565.9665.9765.9866.0066.08λ = 1 × 10−071.6071.6071.6071.6171.6271.6371.68Bold values in each column represent the minimum value.From Table1, three obvious phenomena can be found. The first one is that the estimated accuracy decreases with increasing the noise level in the considered cases if a specific regularization parameter is applied. This result shows that the noise level could bring negative effects to the estimated responses. The second one is that the RPE first falls and then rises with increasing the value of regularization parameter if a specific noise level is considered. The reason behind these results can be intuitively seen from Figure 7. That is, a small value of regularization parameter cannot avoid the appearance of overfitting, while a large one could cause an underfitting result. The third phenomenon is that the optimal regularization parameters are not exactly equal to each other among different conditions of noise level. Herein, the optimal regularization parameter represents a value of λ corresponding to the smallest RPE in the considered cases whose noise levels are the same. The optimal regularization parameter, on the whole, increases with increasing the noise level in the considered cases. This is reasonable since the reduction of greater measured noises needs a higher weight for the penalty term in the matrix regularization model. In addition, the above results also indicate that it is difficult to select an optimal regularization parameter if the noise level is unknown to us. This is the reason why a preliminary estimation of noise energy is required in this paper for selecting the regularization parameter, as shown in Figure 5.To sum up, the above results clearly show the importance of the selection of regularization parameter. This is actually a common problem in the application of regularization-based method. It is beneficial to note that the optimal selection of regularization parameter is still a key issue in this paper, even though there is a posteriori method proposed in this paper.


## 3.1.3. Effect of Different Positions of Virtual Sensors

Three cases with different positions of the virtual sensor are considered. As shown in Figure4(a), three virtual sensors are installed at nodes 30, 36, and 40, respectively. The acceleration responses corresponding to each sensor are calculated independently. Three noise levels, namely, 5%, 10%, and 15%, are considered for each case. Regularization parameter λ for each case is selected as that we have done in Section 3.1.1. Figure 8 shows a diagram for the selection of regularization parameter for the cases with 15% noise levels considered. From Figure 8, it can be seen that λ = 0.05 can ensure that its corresponding value of RVS is intuitively located in a range from RVS = 0.5% to RVS = 5% in the case of node 30, while the regularization parameter playing a similar role in the other cases can be selected as λ = 0.075. Because of this, the regularization parameters are roughly selected as λ = 0.05 for the case of node 30 and λ = 0.075 for the other two cases. Other parameters are set the same as the ones that have been done in the previous section.Figure 8
Diagram for selection of regularization parameter in cases with 15% noise level considered.Table2 lists the RPEs of the structural responses estimated considering different positions of virtual sensors. From Table 2, it can be seen that the estimated accuracies, on the whole, are close to each other in the cases considering the same virtual sensor, even though the RPE increases with increasing the noise level. This phenomenon indicates again that the measured noises may bring negative effects to the estimated responses. However, the estimated result is not very sensitive to the measured noises in the considered cases. In addition, the RPEs corresponding to virtual sensors at node 36 and node 40 are smaller than the ones corresponding to virtual sensor at node 30. It seems that the position of virtual sensor could affect the accuracy of the estimated result in the considered cases. However, the reason behind this may be due to many factors; for example, it may be related to the selection of regularization parameter. As listed in Table 1, the RPE value can decrease to 15.85% if the regularization parameter is selected as λ = 1 × 10−2 for the case that is denoted as node 30 and 5% noise level.Table 2
Relative percentage errors calculated with different virtual sensors.
Positions of virtual sensorRelative percentage errors (%)5% noise10% noise15% noiseNode 3020.1720.8121.82Node 3614.5915.0215.70Node 4014.0414.1714.39For more detailed study, Figure9 shows the local responses calculated from different virtual sensors and noise levels. From the first row of Figure 9, it can be seen that the acceleration responses calculated considering different noise levels are all close to the real structural responses. This phenomenon, on the whole, means that the estimated result is not very sensitive to the measured noises in the case of virtual sensor installed at node 30. Similar results can also be seen from the other rows of Figure 9. This advantage mainly benefits from the penalty function such as a term of (BreF2+ωTvrBreF2). Furthermore, Figure 9 clearly shows that the proposed method can be used for reconstructing the acceleration responses at nodes 30, 36, and 40 with a good accuracy in these cases. Combining the estimated accuracies shown in Table 2, we can see that the proposed method is not very sensitive to the positions of the considered virtual sensors, although there are some differences between the RPEs in the cases considering different positions of the virtual sensors. The reason behind this result may be that structural responses obtained from two simulated real sensors can provide enough information for reflecting the signal features of the simulated force. In such situation, the structural responses can be naturally reconstructed with a relatively better accuracy even though different positions of the virtual sensors are considered in the above cases.Figure 9
Local responses calculated considering different virtual sensors and noise levels. (a) Node 30 and 5% noise level. (b) Node 30 and 10% noise level. (c) Node 30 and 15% noise level. (d) Node 36 and 5% noise level. (e) Node 36 and 10% noise level. (f) Node 36 and 15% noise level. (g) Node 40 and 5% noise level. (h) Node 40 and 10% noise level. (i) Node 40 and 15% noise level.
(a)(b)(c)(d)(e)(f)(g)(h)(i)


## 3.1.4. Effect of Different Lengths of Moving Time Window

The length of moving time window is one of the significant parameters of the proposed method. This parameter directly relates to the size of the transmissibility matrixTvr. To get an intuitive understanding, five cases with different lengths of moving time windows, as listed in Table 3, are carried out. The lengths of moving time windows considered in Case 1, Case 2, Case 3, Case 4, and Case 5 actually correspond to 1 s, 2 s, 3 s, 4 s, and 5 s, respectively. The length of overlapping part between any two adjacent windows is auto-selected as half the length of the moving window. The regularization parameter is selected as λ = 0.05. Other parameters are selected as the ones that have been done in Section 3.1.1.Table 3
Comparison of calculation performance of cases with different lengths of time windows.
CasesNumber of sampling pointsSize of matricesCalculation performancesEach windowOverlap partTvrBreCalculation time (s)Ratio value (tcase i/tcase 1)Case 110005001000 × 20002000 × 1192.0040911.00Case 2200010002000 × 40004000 × 599.8422164.91Case 3300015003000 × 60006000 × 3929.99130814.97Case 4400020004000 × 80008000 × 2966.06360232.96Case 5500025005000 × 1000010000 × 23132.70321366.22All the cases are calculated via the same program with different calculation parameters and a personal computer with information CPU (AMD Ryzen 7 2700X Eight-Core Processor 3.70 GHz), RAM 16.0 GB, Windows 10, Matlab R2015a. The analysis times for five cases are gathered and listed in Table3. From Table 3, it is not hard to find that the computational time increases rapidly with increasing the length of the moving time window in the considered cases. The direct cause is that a longer length of moving time window actually leads to a larger size of transmissibility matrix, as listed in Table 3. In addition, Case 1 and Case 5 are picked up as two representative cases. Comparison of the reconstructed responses obtained from these two cases is plotted in Figure 10. From Figure 10, it can be seen that there are no significant differences between two reconstructed responses on the whole. Overall, Case 1 can spend less calculation time on reconstructing the acceleration responses. Therefore, with a comprehensive consideration of calculation time and reconstructed accuracy, it is not hard to find that the case with window length of 1000 sampling points is the best solution way among the five cases listed in Table 3. To some extent, the advantage of Case 1 actually benefits from the consideration of unknown initial conditions. It is because this consideration is the theoretical basis for the formulation of equation (18), which is a governing equation expressed in a matrix form.Figure 10
Reconstructed responses at node 30 with different lengths of moving time windows. (a) Total time range. (b) Local time range from 19 s to 21 s.Furthermore, different noise levels such as 5%, 10%, and 15% are considered. The RPEs are calculated and listed in Table4. From Table 4, two phenomena can be obtained. The first one is that the measured noises could bring negative effects to the estimated results in the considered cases. This phenomenon is the same as the one discussed in the former sections. The second phenomenon is that the estimated accuracy increases with increasing the length of time windows in the considered cases. In fact, the estimated errors at the beginning of each moving time window may be larger than the ones at other zones due to the interference of structural initial conditions. Hence, if more moving time windows are applied for covering the total time range, the estimated responses may have a lower accuracy. To provide a more intuitive understanding, Case 1 and Case 5 are picked up as two representative cases, and their estimated errors between the estimated responses and the real responses are plotted in Figure 11. From Figure 11, it can be seen that more spikes appear in the estimated errors in Case 1. These spikes actually correspond to the beginning of each time window. That is why the coordinate values of the spikes shown in Figure 11(b) are close to the integral multiple of 0.5 seconds, while the coordinate values of the spikes shown in Figure 11(d) are near to the integral multiple of 2.5 seconds. The above phenomenon is actually the direct reason why estimated accuracy of Case 5 is better than that of Case 1.Table 4
Relative percentage errors calculated with different lengths of time windows and noise levels.
CasesNumber of sampling pointsRelative percentage errors (%)Each windowOverlap part5% noise10% noise15% noiseCase 1100050020.1720.8121.82Case 22000100015.7716.5517.76Case 33000150014.0214.9016.22Case 44000200013.2914.1815.54Case 55000250013.0013.9115.30Figure 11
Estimated errors between estimated responses and real responses with 10% noise level considered. (a) Total length in Case 1. (b) Local length in Case 1. (c) Total length in Case 5. (d) Local length in Case 5.


## 3.2. Further Verification on Plane Truss Structure

A plan truss structure is simulated for further verifying the proposed method. As shown in Figure12, the truss is simply supported with pin joint at one end and roller support at the other end. Both horizontal and vertical members have a length of 1.0 m. The cross section area is 0.0016 m2. The density and Young’s modulus are 7850 kg·m−3 and 2.06 GPa, respectively [45]. Lumped mass matrix is considered. The first five nature frequencies are 1.56 Hz, 5.17 Hz, 6.94 Hz, 11.98 Hz, and 17.72 Hz, respectively. One dynamic force is considered, and its time history is simulated as(28)ft=501−cos10πtsin30πt.Figure 12
Plan truss for numerical simulations.Four real sensors installed at nodes 6, 8, 10, and 12 are used for measuring the vertical acceleration responses. The sampling frequency is set as 1000 Hz. The structural responses are calculated using the mode superposition method and considering the first fifteen modal orders. A basic process for simulating the structural responses is the same as the one described in Section3.1. Three noise levels such as 5%, 10%, and 15% levels are considered. The discrete time interval is selected as 0.001 s, and the total sampling time length is set as 10 s. The damping ratio for each order is selected as 0.012. The initial conditions are set as zeros.The objective is to reconstruct the acceleration responses at nodes 15, 17, and 19, respectively. Effects of different numbers of real sensors are considered. For more details, structural responses corresponding to each virtual sensor will be estimated by using structural responses measured from two sensors at nodes 6 and 8; three sensors at nodes 6, 8, and 10; and four sensors at nodes 6, 8, 10, and 12, respectively. Hence, it is not hard to find that there are twenty-seven cases simulated in this section. Acceleration responses are reshaped in a matrix form by using parameters such as length of time window of 1000 sampling points and length of overlap part of 500 sampling points. The weighted factor is set asω = 1. It should be noted that ω = 1 satisfies the condition ω < ωmax in every considered case. The regularization parameter λ is carefully selected as the one that has been done in the previous section. Herein, the lower limit value of RVS is reduced appropriately, and finally λ will be selected as λ = 0.05.According to equation (27), the RPEs corresponding to nine cases are calculated and listed in Table 5. From Table 5, it can be seen that the effect of the measured noises is almost the same as the ones discussed in the previous sections. In addition, Table 5 also shows that the estimated accuracy of each virtual sensor increases with increasing the real sensors in the simulated cases. The reason may be that more sensors can provide more information for reflecting the external force so that the reconstructed acceleration responses are closer to the real ones. However, more sensors will lead to a larger size of transmissibility matrix. Along with this, more computing resources are needed. In view of this, a good balance between the calculation accuracy and efficiency should be carefully considered in practical applications. In addition, Figure 13 shows the time histories of reconstructed acceleration responses that are calculated from nodes 6, 8, 10, and 12 considering 15% noise level. From this figure, it can be seen that the reconstructed responses match pretty well the real ones for each virtual sensor. This result means that the proposed method is feasible and effective in the considered cases to some extent.Table 5
Relative percentage errors calculated from different sets of real sensors.
Real sensorsRPEs at node 15 (%)RPEs at node 17 (%)RPEs at node 19 (%)5% noise10% noise15% noise5% noise10% noise15% noise5% noise10% noise15% noiseNodes 6 and 810.1110.7511.7215.0015.2215.548.9111.6815.20Nodes 6, 8, and 106.917.117.439.689.8410.066.288.1910.63Nodes 6, 8, 10, and 125.996.146.386.766.867.005.526.918.77Figure 13
Reconstructed responses calculated from nodes 6, 8, 10, and 12 considering 15% noise level. (a) Total time range and node 15. (b) Local time range and node 15. (c) Total time range and node 17. (d) Local time range and node 17. (e) Total time range and node 19. (f) Local time range and node 19.Overall, the above discussions show that the proposed method can be used for indirectly measuring the acceleration responses corresponding to the considered virtual sensors. The proposed method can reduce influences of the measured noises to some extent. Hence, the estimated responses shown in this section also have good robustness.


## 4. Experimental Verifications

Experimental studies of an aluminum beam are conducted in a laboratory, as shown in Figure14. The beam has a length of 3 m. The cross section is a rectangle with 150 mm width and 10 mm height. Four accelerometers (B&K) installed at 0.6 m, 1.2 m, 1.8 m, and 2.4 m are stuck on the lower surface of the beam in order to measure the structural responses. Sensitivity values of the accelerometers are set as the original factory defaults. A data acquisition system made by China Orient Institute of Noise & Vibration is used for measuring the structural responses. The sampling frequency is set as 1024 Hz. Furthermore, there are two lightweight foam guides installed on the upper surface of the beam. Some strain gauges are installed on the lower surface of the beam. However, these guides and strain gauges are actually used for another experiment that is related to moving force identification.Figure 14
Experimental setup. (a) Aluminum beam. (b) Data acquisition system. (c) Experimental support. (d) Sensor placement in experimental test. (e) Finite element model.
(a)(b)(c)(d)(e)


### 4.1. Modal Analysis

A modal experiment analysis (MEA) is firstly carried out. In this step, the beam is impacted by using a hammer with a hollow soft plastic hammerhead. The impact position is set asx = 1.65 m. Acceleration responses are measured by using the installed sensors. Herein, the total sampling time is selected as 20 s. The vibration responses are then applied for estimating the modal parameters by using software of modal analysis. The measured frequencies are listed in Table 6.Table 6
Comparison of the first four natural frequencies.
OrdersFrequencies (Hz)Relative errors|fe − fa|/fe × 100%Damping ratios (experimental)Experimental (fe)Analytical (fa)1st3.0633.0051.89%1.022nd10.68210.8221.31%0.173rd23.89423.8070.36%0.224th42.27141.9780.69%0.16An FEM is built for simulating the experimental beam, as shown in Figure14(e). The FEM contains 20 beam elements. Each element has 2 nodes and 4 degrees of freedom. Torsion spring is added to the end of the beam for restraining the bending deformation. The reason is mainly due to the installation precision of the beam. In fact, our installation process is not always meticulous. The end of the real beam and the screws installed on the support are in contact. Hence, the beam is actually constrained by the screws to some extent. The weight of the experiment beam is obtained by using an electronic scale. Then, the density of unite length can be obtained as ρA = 4.15 kg·m−1. The bending stiffness of the beam and the stiffness of torsion springs are carefully selected via adjustment many times, and finally these two values are set as EI = 9.17 × 102N·m2 and kr = 2.8 × 102N·m·rad−1. Comparison of the first four natural frequencies is presented in Table 6. From Table 6, it can be seen that the calculated frequencies are close to the measured ones. These results actually show the applicability of the proposed FEM to some extent.


### 4.2. Vibration Experiment and Response Reconstruction

The aluminum beam is continuously impacted by using the hollow soft plastic hammer. The impact location is set asx = 1.65 m. The acceleration responses are recorded by using the installed sensors. All the measured responses are filtered by using a low-pass filter with stop frequency of 200 Hz. In the following analysis, the responses obtained from two sensors installed at x = 0.6 m and x = 1.8 m are applied as the input of algorithm. The responses obtained from the remaining sensors are used for the evaluation of the reconstructed results. Figure 15 shows the details of the sensor groups and the measured responses obtained from sensors S1 and S3.Figure 15
Sensor groups and measured responses. (a) Diagrammatic sketch of sensor groups. (b) Acceleration responses of sensor S1. (c) Acceleration responses of sensor S3.
(a)(b)(c)Two experimental cases are considered here for accessing the effectiveness and feasibility of the proposed method. The first one is used to reconstruct the responses at sensor S2 (Case 1), and the second case is used to estimate responses at sensor S4 (Case 2). The length of moving window is selected ask = 1024 sampling points. The length of overlap between any two adjacent windows is set as k0 = 512 sampling points. Time domain discretization is applied for expressing the force as we have selected in Section 3.1.1. The first eight orders of the beam are applied for the calculation of the proposed transmissibility matrix. Herein, the damping ratios of the first four orders are set as the ones listed in Table 6, while the ones for other orders are set as 0.16%. The discrete time interval is set as 1/1024 s.Two maximal acceptable values of the weighted factorω are firstly calculated. They are ωmax = 2.0571 × 108 and ωmax = 2.0330 × 108 for Case 1 and Case 2, respectively. As the consideration discussed in previous section, the weighted factor ω will be set as ω = 1 (ω < ωmax) in two experimental cases. The regularization parameter is selected as we have done in previous section. Herein, it should be noted that the measurement noises in the laboratory are relatively small. Hence, the values of lower limit and upper limit of the RVS, shown in equation (26), can be adjusted appropriately. The acceleration responses obtained via S1, as shown in Figure 15(b), are taken for roughly determining the above two limit values. The responses corresponding to (0 s, 2 s] and the ones corresponding to (10 s, 12 s] are firstly extracted. Herein, acceleration responses in the former segment are regarded as the measurement noises. A ratio of signal energies corresponding to (0 s, 2 s] and (10 s, 12 s] is calculated. It is equal to 6.0274 × 10−6. In view of this, two values of lower limit and upper limit of the RVS can be set as 10−7 and 10−5, respectively. Figure 16 shows the RVS curves for selecting the regularization parameters. According to Figure 16, two regularization parameters are roughly selected as 5 × 10−4 for two experimental cases.Figure 16
Diagram for selection of regularization parameter in test cases. (a) Case 1. (b) Case 2.
(a)(b)Comparisons of the measured responses and the reconstructed responses are shown in Figure17. From Figures 17(a) and 17(d), it can been seen that the reconstructed responses match well the measured responses throughout the time range in the considered cases. More details about the reconstructed responses can further be seen from the other subfigures of Figure 17. These four subfigures actually show us a clearer result; that is, the reconstructed responses are close to the measured ones. This phenomenon means that the proposed method is feasible and effective in solving the problem of response reconstruction in the experimental cases. Furthermore, Figure 18 illustrates the power spectral densities of the responses. These curves are estimated by using a Matlab function named pwelch( ). From Figure 18, it can be seen that two curves match each other well in some ranges, where they are near to the peaks of the curve. This phenomenon shows us that the proposed method can effectively reconstruct the main signal components of the structural responses. There is a relatively large difference between two curves in some ranges, where they are near to the valley values. This is mainly because signal energies corresponding to these frequency ranges are relatively small. As a result, the response components located in these ranges are easier to be weakened, when the regularization-based method is applied for reducing the effects of the measurement noises.Figure 17
Comparison of measured responses and reconstructed responses. (a) Total time range and responses at S2. (b, c) Local time range and responses at S2. (d) Total time range and responses at S4. (e, f) Local time range and responses at S4.Figure 18
Power spectral density (PSD). (a) PSD corresponding to S2. (b) PSD corresponding to S4.
(a)(b)To sum up, the reconstructed results shown in Figures17 and 18 actually show that good correspondence is found in both the time and frequency domains for demonstrating the effectiveness of the proposed F-norm regularization-based method.


## 4.1. Modal Analysis

A modal experiment analysis (MEA) is firstly carried out. In this step, the beam is impacted by using a hammer with a hollow soft plastic hammerhead. The impact position is set asx = 1.65 m. Acceleration responses are measured by using the installed sensors. Herein, the total sampling time is selected as 20 s. The vibration responses are then applied for estimating the modal parameters by using software of modal analysis. The measured frequencies are listed in Table 6.Table 6
Comparison of the first four natural frequencies.
OrdersFrequencies (Hz)Relative errors|fe − fa|/fe × 100%Damping ratios (experimental)Experimental (fe)Analytical (fa)1st3.0633.0051.89%1.022nd10.68210.8221.31%0.173rd23.89423.8070.36%0.224th42.27141.9780.69%0.16An FEM is built for simulating the experimental beam, as shown in Figure14(e). The FEM contains 20 beam elements. Each element has 2 nodes and 4 degrees of freedom. Torsion spring is added to the end of the beam for restraining the bending deformation. The reason is mainly due to the installation precision of the beam. In fact, our installation process is not always meticulous. The end of the real beam and the screws installed on the support are in contact. Hence, the beam is actually constrained by the screws to some extent. The weight of the experiment beam is obtained by using an electronic scale. Then, the density of unite length can be obtained as ρA = 4.15 kg·m−1. The bending stiffness of the beam and the stiffness of torsion springs are carefully selected via adjustment many times, and finally these two values are set as EI = 9.17 × 102N·m2 and kr = 2.8 × 102N·m·rad−1. Comparison of the first four natural frequencies is presented in Table 6. From Table 6, it can be seen that the calculated frequencies are close to the measured ones. These results actually show the applicability of the proposed FEM to some extent.


## 4.2. Vibration Experiment and Response Reconstruction

The aluminum beam is continuously impacted by using the hollow soft plastic hammer. The impact location is set asx = 1.65 m. The acceleration responses are recorded by using the installed sensors. All the measured responses are filtered by using a low-pass filter with stop frequency of 200 Hz. In the following analysis, the responses obtained from two sensors installed at x = 0.6 m and x = 1.8 m are applied as the input of algorithm. The responses obtained from the remaining sensors are used for the evaluation of the reconstructed results. Figure 15 shows the details of the sensor groups and the measured responses obtained from sensors S1 and S3.Figure 15
Sensor groups and measured responses. (a) Diagrammatic sketch of sensor groups. (b) Acceleration responses of sensor S1. (c) Acceleration responses of sensor S3.
(a)(b)(c)Two experimental cases are considered here for accessing the effectiveness and feasibility of the proposed method. The first one is used to reconstruct the responses at sensor S2 (Case 1), and the second case is used to estimate responses at sensor S4 (Case 2). The length of moving window is selected ask = 1024 sampling points. The length of overlap between any two adjacent windows is set as k0 = 512 sampling points. Time domain discretization is applied for expressing the force as we have selected in Section 3.1.1. The first eight orders of the beam are applied for the calculation of the proposed transmissibility matrix. Herein, the damping ratios of the first four orders are set as the ones listed in Table 6, while the ones for other orders are set as 0.16%. The discrete time interval is set as 1/1024 s.Two maximal acceptable values of the weighted factorω are firstly calculated. They are ωmax = 2.0571 × 108 and ωmax = 2.0330 × 108 for Case 1 and Case 2, respectively. As the consideration discussed in previous section, the weighted factor ω will be set as ω = 1 (ω < ωmax) in two experimental cases. The regularization parameter is selected as we have done in previous section. Herein, it should be noted that the measurement noises in the laboratory are relatively small. Hence, the values of lower limit and upper limit of the RVS, shown in equation (26), can be adjusted appropriately. The acceleration responses obtained via S1, as shown in Figure 15(b), are taken for roughly determining the above two limit values. The responses corresponding to (0 s, 2 s] and the ones corresponding to (10 s, 12 s] are firstly extracted. Herein, acceleration responses in the former segment are regarded as the measurement noises. A ratio of signal energies corresponding to (0 s, 2 s] and (10 s, 12 s] is calculated. It is equal to 6.0274 × 10−6. In view of this, two values of lower limit and upper limit of the RVS can be set as 10−7 and 10−5, respectively. Figure 16 shows the RVS curves for selecting the regularization parameters. According to Figure 16, two regularization parameters are roughly selected as 5 × 10−4 for two experimental cases.Figure 16
Diagram for selection of regularization parameter in test cases. (a) Case 1. (b) Case 2.
(a)(b)Comparisons of the measured responses and the reconstructed responses are shown in Figure17. From Figures 17(a) and 17(d), it can been seen that the reconstructed responses match well the measured responses throughout the time range in the considered cases. More details about the reconstructed responses can further be seen from the other subfigures of Figure 17. These four subfigures actually show us a clearer result; that is, the reconstructed responses are close to the measured ones. This phenomenon means that the proposed method is feasible and effective in solving the problem of response reconstruction in the experimental cases. Furthermore, Figure 18 illustrates the power spectral densities of the responses. These curves are estimated by using a Matlab function named pwelch( ). From Figure 18, it can be seen that two curves match each other well in some ranges, where they are near to the peaks of the curve. This phenomenon shows us that the proposed method can effectively reconstruct the main signal components of the structural responses. There is a relatively large difference between two curves in some ranges, where they are near to the valley values. This is mainly because signal energies corresponding to these frequency ranges are relatively small. As a result, the response components located in these ranges are easier to be weakened, when the regularization-based method is applied for reducing the effects of the measurement noises.Figure 17
Comparison of measured responses and reconstructed responses. (a) Total time range and responses at S2. (b, c) Local time range and responses at S2. (d) Total time range and responses at S4. (e, f) Local time range and responses at S4.Figure 18
Power spectral density (PSD). (a) PSD corresponding to S2. (b) PSD corresponding to S4.
(a)(b)To sum up, the reconstructed results shown in Figures17 and 18 actually show that good correspondence is found in both the time and frequency domains for demonstrating the effectiveness of the proposed F-norm regularization-based method.


## 5. Conclusions

This paper investigates the indirect reconstruction of structural responses. A novel method based on the concept of transmissibility and matrix regularization is proposed. Firstly, a linear relationship between different measured points is formulated based on the theory of force identification. In this step, a highlight is that structural responses are reshaped in a form of matrix so that the proposed method has the ability to deal with the long-time duration problem. Frobenius norm matrix regularization is then innovatively introduced for ensuring that the estimated results are robust to the measured noises. Herein, the penalty term is defined as a weighted combination of the measured structural responses and the estimated ones. This is different from the direct application of the regularization mathematical model. Finally, the problem of response reconstruction is boiled down to an optimization problem in mathematics. The optimization problem aims at searching an optimal matrix that represents the structural responses corresponding to the real sensors. To verify the effectiveness and feasibility of the proposed method, numerical simulations and experimental studies are carried out. Some contributions and conclusions can be summarized as follows:(1)
Time histories of dynamic force and initial vibration state of the structure are not required for the proposed method.(2)
The proposed method can use a transmissibility matrix with a small size to formulate a governing equation so that it can spend less cost on solving the long-time duration problem. Herein, the length of moving time window is a key determinant of the size of transmissibility matrix.(3)
Numerical and experimental studies show that the proposed method is correct, reasonable, and suitable for indirectly measuring the structural responses with a good accuracy.(4)
Case studies indicate that measured noises could bring negative effects to the estimated structural responses. However, the negative effects are limited, and the proposed method is not very sensitive to the measured noises.(5)
Illustrated results also show that selection of regularization parameter is very important for the application of the proposed method. Regularization parameter with too small value or too large value could cause a poor estimated accuracy.(6)
In the current investigation, the number of dynamic forces is equal to one. In addition, the structure is also relatively simple. These two may be different from the complex engineering practice. Furthermore, the selection of regularization parameter still depends on experience. Hence, future works are expected to conduct in-depth research on the above problems.



---

*Source: 5176613-2021-10-14.xml*