- Email: [email protected]
Parametric study on sequential deconvolution for force identification
Journal of Sound and Vibration 377 (2016) 76–89
Contents lists available at ScienceDirect
Journal of Sound and Vibration journal homepage: www.elsevier.com/locate/jsvi
Parametric study on sequential deconvolution for force identification Tao Lai a, Ting-Hua Yi a,n, Hong-Nan Li a,b a b
School of Civil Engineering, Dalian University of Technology, Dalian 116023, China School of Civil Engineering, Shenyang Jianzhu University, Shenyang 110168, China
a r t i c l e in f o
abstract
Article history: Received 20 November 2015 Received in revised form 5 February 2016 Accepted 8 May 2016 Handling Editor: L. G. Tham Available online 24 May 2016
The force identification can be mathematically viewed as the mapping from the observed responses to external forces through a matrix filled with system Markov parameters, which makes it difficult or even impossible for long time duration. A potentially efficient solution is to sequentially perform the identification processing. This paper presents a parametric study on the sequential deconvolution input reconstruction (SDR) method, which was proposed by Bernal. The behavior of the SDR method due to the effects of window parameters, noise levels and sensor configurations is investigated. In addition, a new normalized standard deviation of the reconstruction error over time is derived to cover the effect of only independent noise entries. The sinusoidal and band-limited white noise excitations are tested to be identified with good accuracy even with 10% noise. The simulation results yield various conclusions that may be helpful to engineering practitioners. & 2016 Elsevier Ltd. All rights reserved.
Keywords: Force identification Sequential deconvolution Parametric study
1. Introduction The force identification problem has recently attracted considerable attention in structural health monitoring, structural control and design communities [1–3]. Understanding the dynamic load exerted on structures is beneficial to engineering practitioners because structural integrity and safety are significantly involved with dynamic loading. However, in practice, directly measuring the dynamic force is generally difficult or even impossible, such as wind loads and earthquake excitations. The alternative is the recovery of inputs from vibration measurements. Nevertheless, the number of vibration sensors is usually limited for economic or technological reasons. Therefore, it is necessary to develop force identification techniques with incomplete measurements to evaluate the external loadings applied on the structures. Force identification from the observed responses is a typical inverse problem in structural dynamics. The reconstruction results usually suffer from numerical non-uniqueness and instability due to matrix ill-posedness and noise disturbance [4]. Regularization methods such as the widely used Tikhonov regularization and truncated singular value decomposition (TSVD) method are adopted to guarantee an accurate solution. A recent thorough review of the literature on force identification with observed measurements is found in Sanchez et al. [5]. The force identification problem is initially solved in frequency domain based on the fast Fourier transform (FFT) [6]. Ramkrishna et al. [7] proposed a method for vibration testing of flexible structures with electronic equipment in frequency n
Corresponding author at: School of Civil Engineering, Dalian University of Technology, Dalian 116023, China. Tel./fax: þ86 411 84706050. E-mail address: [email protected] (T.-H. Yi).
http://dx.doi.org/10.1016/j.jsv.2016.05.013 0022-460X/& 2016 Elsevier Ltd. All rights reserved.
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
77
domain. Chao et al. [8] presented a frequency-domain force identification technique considering the nonlinearity of the system. The FFT-based methods have two main adverse problems to address: these methods are limited to cases with zero initial conditions, and the hypotheses of periodic inputs result in a serious leakage problem particularly for short load durations or non-stationary signals [9–11]. In contrast, the more recent development of time-domain force identification techniques can avoid such problems. Mao et al. [12] estimated the inputs with high accuracy by introducing a precise timestep integration method in the state-space domain. Li et al. [13] decomposed the excitations and responses in wavelet domain and solved the force identification problem with Tikhonov method. Ding et al. [14] proposed a discrete average acceleration algorithm for force identification based on the Newmark-β algorithm. Liu et al. [15] extended this method and modified the conventional implicit Newmark-β algorithm into an explicit form for input estimation with a regularization method. Xu et al. [16] approximated the inputs and outputs locally by moving least-square method and transformed the differential equation of motion into an integral equation with the virtual work principle to iteratively identify inputs. Wang et al. [17] discretized the time element of inputs and outputs and recover the inputs in the state-space method based on weak Galerkin formulation. Qiao et al. [18] represented the desired force by the cubic B-spline scaling function based on the wavelet multi-resolution concept, and the force identification problem is equivalent to finding the coefficients of the scaling functions. This method was verified in an impact force identification experiment [11]. These studies of the time-domain techniques can be categorized into one group, which comprises two main steps: the first is to mathematically formulate the mapping from the sequences of inputs to outputs through an operation matrix; the second step is to solve the ill-posed inverse problem with the regularization method. However, the force identification problem is difficult or even impossible when the size of the operation matrix increases. A potentially effective and efficient method is to sequentially extract the inputs. Bernal et al. [19] proposed a sequential deconvolution input reconstruction (SDR) method, which can overcome this adverse problem and make online identification possible. In addition, Kalman filtering is another popular method that has been introduced for force identification [20–23]. In this paper, the parametric study on the SDR method is presented. The behavior of the SDR method due to the effects of window parameters, noise levels and sensor configurations is investigated. In addition, a new normalized standard deviation of the input identification error is derived to interpret the accuracy of the SDR method. The paper is organized as follows. Section 2 presents the mapping from the sequences of inputs to outputs. In Section 3, the first p values of the input vectors are uniquely extracted in the “window”. Section 4 is devoted to the stability and regularization problems. The parametric study of the SDR method is numerically illustrated in Section 5. Finally, some conclusions are summarized.
2. The mapping from the sequences of input to output 2.1. Fundamental formulation The equation of motion of a linear time-invariant system with n degrees of freedoms (DOFs) can be described by € þCuðtÞ _ þ KuðtÞ ¼ GfðtÞ ΜuðtÞ
(1)
where M A ℝnn ,C A ℝnn and K A ℝnn are the mass, damping and stiffness matrices, respectively; uðtÞ A ℝn1 is the displacement vector, fðtÞ A ℝr1 is the force vector applied to the system through the mapping matrix G A ℝnr , and r is the number of inputs. h iT _ T , Eq. (1) can be converted into the state-space representation as Defining the state vector x ¼ uT u _ ¼ Ac xðtÞ þBc fðtÞ xðtÞ
(2)
where Ac ¼
0 M
M
K
Bc ¼
I
1
M
C
A ℝ2n2n
(3-1)
0 1
1
G
A ℝ2nr
(3-2)
The general system output can be a linear combination of displacements, velocities and accelerations, which are grouped into an output vector yðtÞ A ℝm1 yðtÞ ¼ Cc xðtÞ þ Dc fðtÞ
(4)
where Cc A ℝ depends on the type of measurements, Dc A ℝ is the direct transmission term, and m is the number of outputs. The exact solution to Eq. (2) is given by Z t eAc ðt τÞ Bc fðτÞdτ (5) xðtÞ ¼ eAc t xð0Þ þ mn
mr
0
78
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Substituting Eq. (5) into Eq. (4) yields Z yðtÞ ¼ Cc eAc t xð0Þ þCc
t 0
eAc ðt τÞ Bc fðτÞdτ þ Dc fðtÞ
(6)
The continuous-time state-space model in Eqs. (5) and (6) can be discretized into the discrete-time state-space model, which is represented by Eqs. (7) and (8) with certain assumptions of the inter-sample behavior [24] (7) x k ¼ Ad x k 1 þBd f k 1 y k ¼ Cd x k þ Dd f k
(8)
2.2. Formulating the mapping from the sequences of inputs to outputs Based on the iterative equation in Eqs. (7) and (8), the state vector and observation at time kT, where T is the sampling time interval, can be formulated by kX 1 x k ¼ Akd x½0 þ Akd l 1 Bd f l
(9)
l¼0 kX 1 Cd Akd l 1 Bd f l þDd f k y k ¼ Cd Akd x½0 þ
(10)
l¼0
Stacking the sequences of inputs and outputs in columns yields h λ1 x½λ ¼ Aλd x½0 þ Ad Bd
2
Aλd 2 Bd
3 2 3 2 Dd Cd y ½ 0 7 6 6 7 6 6 Cd Bd 6 y ½ 1 7 6 Cd A d 7 7x½0 þ 6 6 7¼6 7 6 6 7 6 ⋮ ⋮ ⋮ 5 4 4 5 4 λ1 Cd A d Cd Aλd 2 Bd y½λ 1
3 f ½ 0 7 i6 6 f ½ 1 7 7 Bd 6 6 ⋮ 7 4 5 f ½λ 1
⋯
2
0
⋯
Dd
⋯
⋮
⋮
Cd Aλd 3 Bd
⋯
32
3 f ½0 76 7 0 76 f ½1 7 76 7 6 7 0 7 54 ⋮ 5 Dd f ½λ 1 0
(11)
(12)
Compactly write
where
x½λ ¼ Aλd x½0 þCoλ f~ 0jλ 1
(13)
y~ 0jλ 1 ¼ Obλ x½0 þHλ f~ 0jλ 1
(14)
h λ1 Coλ ¼ Ad Bd
Aλd 2 Bd 2
Cd
⋯
i Bd A ℝ2λnλr
3
6 7 6 Cd A d 7 7 A ℝλm2λn Obλ ¼ 6 6 7 ⋮ 4 5 λ1 Cd A d 2 6 6 Hλ ¼ 6 6 4
(15)
Dd
0
⋯
0
C d Bd ⋮
Dd ⋮
⋯ ⋮
0 0
Cd Aλd 2 Bd
Cd Aλd 3 Bd
⋯
Dd
(16)
3 7 7 7 A ℝλmλr 7 5
Eqs. (13) and (14) can be written in a general iterative formulation, which is used in later sections: x l þ λ ¼ Aλd x l þCoλ f~ ljl þ λ 1 y~ ljl þ λ 1 ¼ Obλ x l þ Hλ f~ ljl þ λ 1
(17)
(18) (19)
where Coλ is the controllability matrix (in reversed order), Obλ is the observability matrix; Ηλ is the lower block triangular Toeplitz matrix filled with system Markov parameters; l and λ are positive integers. From Eqs. (12) or (19), the sequences of inputs are mathematically mapped into the outputs through matrix Ηλ .
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
79
3. Force identification in a moving window From linear algebra theory, whether Eq. (19) has a solution fully depends on the characteristic of the matrix Hλ , to be more specific, the direct transmission term. In other words, if the sensors are collocated with the inputs, the full column rank of Ηλ ensures that the force identification problem has a unique solution in the least-square sense. On the contrary, matrix Ηλ is singular in the non-collocated situation, which results in an unstable solution. Therefore, the regularization method is used to impose an additional constraint and obtain a stabilized solution. Another concern is that the computation efficiency dramatically decreases when the size of matrix Ηλ increases. A potentially effective and efficient solution is to extract the inputs in a “window” that moves along the time axis, which makes online (or near online) force identification possible. The first step of the SDR method is to determine the length of the “window”. Herein, it is assumed that the window length is λ. Then, an identifiable matrix Q p is defined to recover the first p ðp r λÞ values of the input vector in the “window”. h i (20) Q p ¼ Ipr 0prðλ pÞr where Ipr represents the identity matrix of order p r. In essence, parameter p indicates how many time stations through which the “window” advances at one time. In this sense, p determines how “fast” the “window” can move. Simultaneously, the initial condition of the “window” in the next movement can be propagated through Eq. (18). From Eq. (19), the uniqueness of the solution within a window requires that the necessary and sufficient condition is fulfilled: the block diagonal element of Hλ in Eq. (17) must have full column rank, that is, m Zr, namely, the number of sensors must be no less than the number of inputs to be identified. This assumption is used throughout the paper. Then, the general solution of Eq. (19) is given by (21) f~ ljl þ λ 1 ¼ H†λ y~ ljl þ λ 1 Obλ x l þ Z U h where H†λ denotes the pseudo-inverse of Hλ , Z is the null space of Hλ , and h is an arbitrary coefficient vector with appropriate dimension. From a mathematical perspective, the solution has two parts: the optimal solution and an arbitrary combination of basis vector in the null space of Hλ . Then, the first p values of the input vector are selected by post-multiplying Eq. (21) by Eq. (20) f~ ljl þ p 1 ¼ Q p H†λ y~ ljl þ λ 1 Obλ x l þQ p Z U h (22) To guarantee the uniqueness of the identified force, the last term in Eq. (22) should be zero. Then, Q p Z should be zero because h is arbitrary. As a result, Zhas to be of a form given by " # 0prðÞ Z¼ ~ (23) Ζ ðλpÞrðÞ
4. Sequential deconvolution process From Eqs. (22) and (23), the input solution in the window can be uniquely identified by n o f~ ljl þ p 1 ¼ Q p H†λ y~ ljl þ p 1 Obλ x l
(24)
Once these inputs are identified in a time window, the initial condition in the next time window can be assessed through state propagation in Eq. (18). The identical process can be sequentially performed until the full duration of interest is covered. Unfortunately, the uncertainty of observations also makes the error propagate in the state vector and the identified input vector when the “window” slides. The stability of input reconstruction should be guaranteed during the sequential processing. Now, it is opportune to determine the window parameters λ and p. 4.1. Error propagation analysis It is necessary to investigate the error propagation of input reconstruction when the “window” advances along the time axis. Based on the stability analysis, the restrictive conditions on the parameters λ and p are expected. The stability of sequential deconvolution was originally derived in reference [19], and the derivation is repeated here for convenience in the subsequent discussion about conditioning. The measured responses are usually contaminated by noise, which leads to the uncertainty in the state vector and inputs. Using Eq. (24), one can write (25) f~ ljl þ p 1 þδf~ ljl þ p 1 ¼ Q p H†λ y~ ljl þ λ 1 þ v~ ljl þ λ 1 Obλ x l þ δx l
80
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Therefore δf~ ljl þ p 1 ¼ Q p H†λ v~ ljl þ λ 1 Obλ δx l
(26)
From Eq. (18), the error propagation in the state vector can be given by δx l þ p ¼ Apd δx l þ Cop δf~ ljl þ p 1
(27)
Substituting Eq. (26) into Eq. (27) produces
δx l þp ¼ Pp δx l þ Zp v~ ljl þ λ 1
(28)
Pp ¼ Apd Cop Q p H†λ Obλ
(29)
Zp ¼ Cop Q p H†λ
(30)
where
Assume that the “window” starts at time l ¼ 0 and after k ¼ jp time stations where j ¼ 1; 2; 3⋯, the propagation error in state vector can be tracked as j1 X δx k ¼ Pjp δx½0 þ Pjp l 1 Zp v~ lpjlp þ λ 1
(31)
l¼0
Eq. (31) expresses the uncertainty of the state vector as a combination of the uncertainty of initial state vector and noise sequence vectors. As shown in reference [19], the necessary condition of the stability requirement in Eq. (31) should be satisfied by
ρðPp Þ ¼ ηj max r 1 (32) where η is the eigenvalue of Pp , and ρ is the spectral radius of Pp . 4.2. Conditioning Based on Eq. (31), the expression for the accumulated input identification error remains a combination of the error in the initial condition and the uncertainty of noise sequences, as derived in reference [19]. However, for the reader's convenience, the noise sequence order is changed to be consistent with the time order. After substituting Eq. (31) into Eq. (26) and performing some rearrangements, one obtains δf~ jpjjp þ p 1 ¼ Q p H†λ Obλ Pjp δx½0 þ
j X
Dl v~ lpjlp þ λ 1
(33)
l¼0
where Dl ¼ Q p H†λ Obi Pjp l 1 Zp A ℝprλm
(34)
Dj ¼ Q p H†λ A ℝprλm
(35)
To simplify, it is assumed that the initial condition is exactly known, or the initial error can be neglected. Eq. (38) can be simplified as a combination of the measurement noise vectors δf~ jp pjjp 1 ¼ DV h
where D ¼ D1
D2
⋯
i
h
Dj A ℝprðj þ 1Þλm and V ¼ v~ T0jλ 1
T v~ pjp þ λ 1
(36) ⋯
T v~ jpjjp þ λ 1
iT
A ℝðj þ 1Þλm1 . If one interests a spe-
cific of inputs, like the qth, then δf~ jp pjjp 1 ¼ d V q
where
q f~ jp pjjp 1
q
(37)
q
is the qth input, and d is the qth row of D. Note that the vector V is not strictly random because its length is
ðj þ1Þλm while there are only ðjp þ λÞm independent noise entries. The noise entries of vector V are not independent leading to large input identification error if Eq. (37) is used as in reference [19]. Herein, a new normalized standard deviation of the reconstruction error over time is derived to cover the effect of only independent noise entries. Reshaping vector V to retain all independent noise entries and the corresponding q coefficients vector d , then Eq. (37) can be rewritten as q q δf~ jp pjjp 1 ¼ d~ V~
(38)
q where d~ A ℝprðjp þ λÞm is the recombination of coefficients in d , and V~ A ℝðjp þ λÞm1 is composed of independent noise entries. The reshaping processing is illustrated in Appendix A. q
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
81
Designating the covariance of the noise vector V~ as Σ v , the uncertain in the identified force can be assessed by the error propagation given by T
σ 2δf q ¼ dq Σ v dq
(39)
Note that Σ v is a diagonal matrix because the noise entries are independent. It is assumed that all channels have the same measurement noise σ v , then Eq. (39) can be evaluated as σ δf q ¼ σ v
qffiffiffiffiffiffiffiffiffiffiffi T dq dq
(40)
Define the normalized standard deviation of the input identification error as σ δf q σv
¼
qffiffiffiffiffiffiffiffiffiffiffi T dq dq
(41)
Eq. (41) can be considered a good indicator to determine if an acceptable accuracy is achievable before the input identification. If the standard deviation of the error is acceptably small and stable over time station, the force identification performance can be guaranteed. If not, a regularization method is needed. 4.3. TSVD regularization Herein, a simple scheme based on TSVD regularization method is used if regularization is needed [4]. From (Eqs. (18) and 19), one can see that the associated matrices Aλd ; Coλ ; Hλ and Obλ do not change as the “window” slides. This is helpful because TSVD only addresses the matrix Hλ once. For simplicity, the batch solution is approximately used to assess the
Fig. 1. Flowchart of the SDR method.
82
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
relative error of force identification. Taking the initial condition as zeros, Eq. (14) produces n o y~ 0jλ 1 þ v~ 0jλ 1 ¼ Hλ f~ 0jλ 1 þ δf~ 0jλ 1
(42)
From linear algebra theory, the relative error of force identification can be formulated as ‖δf~ 0jλ 1 ‖2 ‖v~ 0jλ 1 ‖2 r condðHλ Þ ‖y~ 0jλ 1 ‖2 ‖f~ ‖
(43)
0jλ 1 2
Define ε and ‖v~ 0jλ 1 ‖2 =‖y~ 0jλ 1 ‖2 as the bounded tolerance of the relative identification error and noise-to-signal ratio (NSR), respectively, then condðHλ Þ
‖v~ 0jλ 1 ‖2 σ max ¼ U NSR r ε ‖y~ 0jλ 1 ‖2 σ min
(44)
Hence σ min NSR Z ε σ max
(45)
Eq. (45) is considered as the criterion of TSVD regularization given the tolerance εand NSR. The bounded tolerance ε ¼ 1 is assumed throughout the paper. 4.4. Need for regularization From Eq. (41), one can plot a figure of the normalized standard deviation of input reconstruction error σ δf =σ v verse time stations and determine whether regularization is necessary. Because regularization controls the norm and distort the true result, it is important to maintain NSR=ε as small as possible [19]. Eq. (41) can be accounted for the input reconstruction performance when the regularization breaks the theoretical thread. Finally, a flowchart to implement the SDR method is schematically demonstrated in Fig. 1. Two schemes are available if the identification performance is not guaranteed. Scheme 1 is to adjust window parameters λ and p to reassess the identification performance. Alternatively, one can adopt the TSVD regularization in scheme 2.
5. Numerical illustration The parametric study of the SDR method is demonstrated with a numerical example given by Bernal [19]. For unit consistency, the physical parameters of the bar are set as follows: bar length L ¼ 10, cross-sectional area A ¼ 4, mass density ρd ¼ 1 and modulus of elasticity E ¼ 100. The bar is discretized into 40 elements with one end fixed and the other free shown in Fig. 2. The finite element model has a minimum and maximum frequency of 0:25 Hz and 12:73 Hz, respectively. The recommendation of time step size is T ¼ 0:05 s. Rayleigh damping with 2% in the first and last mode is used. The time histories of two decaying harmonic excitation applied at node 5 and 40, are given by f 5 ðt Þ ¼ 10e 0:3t ð1 cos 2πt Þ sin 6πt (46-1) f 40 ðt Þ ¼ 10e 0:2t ð1 cos πt Þ sin 4πt
(46-2)
The SDR accuracy is quantified by normalized root mean square (NRMS) error given by NRMS ¼
‖f e f r ‖2 100% ‖f r ‖2
(47)
where f e and f r are the estimated and real force vectors, respectively. Due to noise uncertainty, each simulation is independently repeated 50 times to account for the accuracy and robustness of the SDR method. 5.1. Sequential force identification results Herein, 4 accelerometers are assumed to be available located at coordinates {1, 10, 20, 30}. The window parameters λ ¼ 50; p ¼ 1 are used and the selection detail is discussed in a later section. To simulate the measurements, 2% RMS of white noises is added to the actual responses. As Section 4 mentioned, the accuracy of the force identification can be assessed by
Fig. 2. Cantilever bar with 40 DOFs.
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
83
Eq. (41) before input identification. Therefore, the plot of the normalized standard deviation of the input identification error versus time stations is first investigated, and the results are shown in Figs. 3 and 4. For the first input, the normalized standard deviation dramatically increases when the window slides, as shown in Fig. 3. In Fig. 4, the normalized standard deviation is stable after approximately 800 time stations but remains at a large value. The results indicate that the regularization method is necessary. Figs. 3 and 4 show that the normalized standard deviation after regularization is always maintained at a small value for both inputs. Finally, the input reconstruction results are illustrated in Figs. 5 and 6. It is observed that the identified inputs without regularization are deviating from the true inputs over time resulting in large error in force identification. However, after regularization, the identified inputs matches with the exact inputs quite well with NRMS error of 3:50% and 2:40%, respectively. 5.2. Effects analysis Herein, a detail study of the effects including window parameters, measurement noise levels and sensor configurations on the accuracy of the force identification is carried out. At the beginning, the selection of window parameters is investigated. Second, to account for the effect of measurement noise, four levels of Gaussian white noise, i.e.: 0%; 2%; 5% and 10% are studied. Finally, another three sets of sensor arrangements are considered except the sensor layout in the last subsection, as shown in Table 1. The sensors are approximately equal spaced at 1=3; 1=4 and 1=5 length of the bar for sensor cases S1, S2 and S3.
Fig. 3. Normalized standard deviation vs. time station for the first input.
Fig. 4. Normalized standard deviation versus time station for the second input.
84
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Fig. 5. SDR for the input reconstruction without regularization: (a) first input; (b) second input.
Fig. 6. SDR for the input reconstruction with regularization: (a) first input; (b) second input. Table 1 Sensor configuration schemes. Sensor configuration schemes
Sensor sets
S1
S2
S3
S4
{13,26}
{10,20,30}
{8,16,24,32}
{1,10,20,30}
(1) Window parameters λ; p . In Section 4, the stability requirement must be satisfied by Eq. (32), which is the only necessary condition. One can plot the 2D figure of the largest eigenvalue versus parameter pair λ; p , which is illustrated in Fig. 7. It is observed that there is not too much difference of the stabilized zoon, which is concentrated upon the bottom right in Fig. 7. However, the largest eigenvalue of Pp is quite different for regularization and
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
85
non-regularization situations. For illustration, one considers the specified situations p ¼ 1 and p ¼ 5 versus the varyingλ, as shown in Fig. 8. The zoom-in parts show that the largest eigenvalue with regularization is less than that without regularization after approximately λ ¼ 35. Mathematically, the smaller of the largest eigenvalue (less than 1), the faster of the convergence rate of matrix power of Pp . As a result, λ ¼ 30; 40; 50; 60; 70 and 80 with p ¼ 1 or 5 are considered. There are 12 numerical identification results in Figs. 9 and 10. One can see that the NRMS error decreases after λ ¼ 30 and remains at a low level for both p ¼ 1 and p ¼ 5. However, there are large errors for both cases without regularization. For computation efficiency, one takes λ ¼ 50 in the subsequent analysis. (2) Noise level. Herein, one takes λ ¼ 50; p ¼ 1 for the noise level effect analysis. The identification results are summarized in Table 2. The noise significantly affects the accuracy of the input reconstruction. The identification results with or without regularization are identical when there is absence of noise. When the noise levels increases, the performance with regularization slowly deteriorates in the beginning and dramatically increases with 10% noise. For the nonregularization situation, there is large identification error as noise level increases. It is demonstrated that the noise effect is an important factor to be considered for the ill-posed inverse problem. (3) Sensor configuration. Four schemes of sensor configuration are carried out. For each case, one takes λ ¼ 50; p ¼ 1, and 2% white noise is considered. The results are shown in Fig. 11. Compared with S1, scheme S2 has much better identification
Fig. 7. Largest eigenvalue of Pp versus λ; p : (a) non-regularization; (b) regularization.
Fig. 8. Largest eigenvalue Pp versus λ: (a) p ¼ 1; (b) p ¼ 5.
86
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Fig. 9. Input identification error versus λðp ¼ 1Þ: (a) regularization; (b) non-regularization.
Fig. 10. input identification error versus λðp ¼ 5Þ: (a) regularization (b) non-regularization. Table 2 Input identification error with different noise levels. NRMS error (%)
Non-regularization Regularization
Noise level
Input Input Input Input
I II I II
0%
2%
5%
10%
4.68 2.09 4.70 2.10
97.21 26.75 3.50 2.40
247.64 61.52 4.98 6.39
420.88 111.29 16.19 29.16
accuracy. However, S3 does not outperform the accuracy of S2, which uses one more sensor compared with S3. The scheme S4 has the best accuracy either with regularization or without regularization. 5.3. Band-limited white noise identification To further illustrate the applicability of the SDR method, band-limited white noise excitation is tested. The band-limited white noise that was applied at node 40 of the bar is generated in Simulink block and passes through a low-pass filter with 10 Hz cutoff frequency. Herein, it is assumed that λ ¼ 50; p ¼ 1 is used. 2% and 10% white noise levels are considered to mimic the true acceleration responses collected at coordinates {1, 10, 20, 30}. First, one should check the plot of the normalized standard deviation versus time stations. It is noteworthy that the normalized standard deviation is not related to the external excitation type. As a result, the trend of the normalized standard deviation over time stations is identical to that
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
87
Fig. 11. Input identification error with different sensor configuration schemes: (a) regularization; (b) non-regularization.
Fig. 12. SDR for band-limited white noise identification with regularization: (a) 2% noise; (b) 10% noise.
in Fig. 4, which indicates that the regularization is required. The identification results are shown in Fig. 12. It is observed that the identified band-limited white noise matches with the true input quite well with NRMS error of 4:61% under low noise level. In contrast, the identified input with 10% white noise has a large error with an NRMS of 18:45%.
6. Conclusion This paper presents a parametric study on the sequential deconvolution input reconstruction method. A new normalized standard deviation of the input identification error over time stations is defined to interpret the accuracy of reconstruction. The effects of the window parameters λ; p , noise levels and sensor configuration schemes are investigated. The sinusoidal and band-limited white noise excitations are tested. Numerical analysis of a simple bar is performed to demonstrate the accuracy and efficiency of the SDR method. The following conclusions are summarized as follows: (1) The selection of window parameters λ; p is crucial for the SDR method. Without the right choice, the input identification dramatically changed, which suggested a need for a larger λ and a smaller p on the premise of the necessary condition and computation efficiency consideration.
88
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
(2) The noise level significantly affects the accuracy of the input identification. One can achieve good reconstruction accuracy with 5% white noise and acceptable accuracy with 10% white noise after regularization. (3) The SDR accuracy dramatically varies with the sensor configurations. A good sensor configuration can reduce the illposedness of the input identification. The sensor layout remains an open problem for input reconstruction. (4) Single or multiple excitations of sinusoid and band-limited white noise are tested to be identified with good accuracy even with 10% noise.
Acknowledgments The authors would like to appreciate Professor Bernal for the innovation method and his help to clarify many related concepts and ideas. This research work was jointly supported by the 973 Program (Grant no. 2015CB060000), the National Natural Science Foundation of China (Grant nos. 51421064, 51478081), the Fok Ying Tong Education Foundation (Grant no. 141072), and the Science Fund for Distinguished Young Scholars of Dalian (Grant nos. 2014J11JH125, 2015J12JH209).
Appendix A. Illustration of the normalized standard deviation It is assumed that there are r ¼ 2 inputs and m ¼ 3 sensors available. Herein, one takes λ ¼ 5; p ¼ 1. Suppose that we want to characterize the reconstruction error at t ¼ 2T. As a result, we have j ¼ 2 in Eq. (33) written by δf~ 2pj3p 1 ¼ Q p H†λ Obλ Pjp δx½0 þ D0 v~ 0jλ 1 þD1 v~ pjp þ λ 1 þ D2 v~ 2pj2p þ λ 1
(A-1)
In matrix formulation, one has 2 3 ~1 δ f row#1 row#1 row#1 row#1 2pj3p 1 6 7 v~ 0jλ 1 þ v~ pjp þ λ 1 þ v~ 2pj2p þ λ 1 δx½0 þ 4 2 5¼ row#2 D1 row#2 D2 row#2 row#2 D0 δf~ 2pj3p 1
(A-2)
The rearrangement of Eq. (A-2) gives 1 δf~ 2pj3p 1
2 δf~ 2pj3p 1
h
¼ ½row#1δx½0 þ row#1D0
h
¼ ½row#2δx½0 þ row#2D0
row#1D1
row#2D1
i
2
v~ 0jλ 1
3
7 row#1D2 6 4 v~ pjp þ λ 1 5 v~ 2pj2p þ λ 1
row#2D2
i
(A-3)
2
3 v~ 0jλ 1 6 v~ 7 4 pjp þ λ 1 5 v~ 2pj2p þ λ 1
(A-4)
Vector V A ℝðj þ 1Þλm1 is composed of all noise sequencesv~ jpjjp þ λ 1 2 ~ 3 v 0jλ 1 6 7 V ¼ 4 v~ pjp þ λ 1 5 v~ 2pj2p þ λ 1 where
h ð1Þ v~ 0jλ 1 ¼ v0j4
2Þ vð0j4
3Þ vð0j4
1Þ vð1j5
2Þ vð1j5
3Þ vð1j5
1Þ vð2j6
2Þ vð2j6
(A-5)
3Þ vð2j6
1Þ vð3j7
h ð1Þ v~ pjp þ λ 1 ¼ v1j5
2Þ vð1j5
3Þ vð1j5
1Þ vð2j6
2Þ vð2j6
3Þ vð2j6
1Þ vð3j7
2Þ vð3j7
3Þ vð3j7
1Þ vð4j8
h ð1Þ v~ 2pj2p þ λ 1 ¼ v2j6
2Þ vð2j6
3Þ vð2j6
1Þ vð3j7
2Þ vð3j7
3Þ vð3j7
1Þ vð4j8
2Þ vð4j8
3Þ vð4j8
1Þ vð5j9
2Þ vð3j7 2Þ vð4j8 2Þ vð5j9
3Þ vð3j7 3Þ vð4j8 3Þ vð5j10
1Þ vð4j8 1Þ vð5j9 1Þ vð6j11
2Þ vð4j8
3Þ vð4j8
2Þ vð5j9 2Þ vð6j11
i
3Þ vð5j10
(A-6) i
3Þ vð6j11
(A-7) i
(A-8)
The corresponding coefficient vectors d are
h i 1 d ¼ row#1D0 row#1D1 row#1D2 h i ð1Þ ð1Þ ð1Þ d2 ⋯ d45 ¼ d1
(A-9)
h i 2 d ¼ row#1D0 row#1D1 row#1D2 h i ð2Þ ð2Þ ð2Þ d2 ⋯ d45 ¼ d1
(A-10)
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
89
Note that vector V is not strictly random because its length is ðjþ 1Þλm ¼ 45 while there are only ðjp þλÞm ¼ 21 independent noise entries. Now, the modified coefficient vector V~ A ℝ211 filled with independent noise entries is defined as V~ ¼
1Þ vð0j4
2Þ vð0j4
3Þ vð0j4
1Þ vð1j5
2Þ vð1j5
3Þ vð1j5
1Þ vð2j6
2Þ vð2j6
3Þ vð2j6
1Þ vð3j7
2Þ vð3j7
3Þ vð3j7
1Þ vð4j8
2Þ vð4j8
3Þ vð4j8
1Þ vð5j9
2Þ vð5j9
3Þ vð5j10
1Þ vð6j11
2Þ vð6j11
3Þ vð6j11
(A-11) The corresponding coefficient vectors d are reshaped as d~ 2 ð1Þ ð1Þ ð1Þ d1 d2 d3 6 1 ð1Þ ð1Þ ð1Þ ð1Þ ð1Þ d~ ¼ 6 4 d6 þ d18 d7 þ d19 þ d31 ⋯ ð1Þ ð1Þ ð1Þ ð1Þ ð1Þ d29 þ d41 d30 þd42 d43 2
ð2Þ
d1 6 ð2Þ ð2Þ ~ 6 d ¼ 4 d6 þ d18 ð2Þ ð2Þ d29 þ d41 2
ð2Þ
ð2Þ
d2 ð2Þ
ð2Þ
ð2Þ
d7 þ d19 þ d31 ð2Þ
ð1Þ
d30 þd42
ð1Þ
ð1Þ
ð1Þ
d4 þ d16 ð1Þ
ð1Þ
ð1Þ
d15 þ d27 þd39 ð1Þ
ð2Þ
ð2Þ
d4 þ d16
⋯
d15 þ d27 þd39
ð2Þ
d43
ð2Þ ð2Þ
d44
ð1Þ ð1Þ 7 d28 þ d40 7 5
(A-12)
d45
d3
ð2Þ
3
ð1Þ
d44 ð2Þ
ð1Þ
d5 þ d17
ð2Þ
d5 þ d17 ð2Þ
3
ð2Þ ð2Þ 7 d28 þ d40 7 5
(A-13)
ð2Þ
d45
1 2 Note that vectors d~ and d~ are row vectors, which are written as the form in Eqs. (A-12) and (A-13) for conciseness.
References [1] S.S. Law, J.Q. Bu, X.Q. Zhu, Time-varying wind load identification from structural responses, Engineering Structures 27 (2005) 1586–1598. [2] Y. Niu, C.P. Fritzen, H. Jung, I. Buethe, Y.Q. Ni, Y.W. Wang, Online simultaneous reconstruction of wind load and structural responses-theory and application to Canton Tower, Computer-Aided Civil and Infrastructure Engineering 30 (2015) 666–681. [3] A.L. Wu, C.H. Loh, J.N. Yang, J.H. Weng, C.H. Chen, T.S. Ueng, Input force identification: application to soil–pile interaction, Structural Control Health Monitoring 16 (2009) 223–240. [4] P.C. Hansen, Rank-Deficient and Discrete Ill-posed Problems: Numerical Aspects of Linear Inversion, SIAM, Philadelphia, 1998. [5] J. Sanchez, H. Benaroya, Review of force reconstruction techniques, Journal of Sound and Vibration 333 (2014) 2999–3018. [6] J.F. Doyle, Further developments in determining the dynamic contact law, Experimental Mechanics 24 (1984) 265–270. [7] D. Ramkrishna, A. Joshi, P.M. Mujumdar, Y. Krishna, Frequency domain input identification for vibration testing of flexible structures, Proceedings of the Institution of Mechanical Engineers Part G – Journal of Aerospace Engineering 227 (2013) 838–856. [8] M. Chao, H. Hongxing, X. Feng, The identification of external forces for a nonlinear vibration system in frequency domain, Proceedings of the Institution of Mechanical Engineers Part C – Journal of Mechanical Engineering Science 228 (2014) 1531–1539. [9] L.J.L. Nordstrom, T.P. Nordberg, A time delay method to solve non-collocated input estimation problems, Mechanical Systems and Signal Processing 18 (2004) 1469–1483. [10] M.S. Allen, T.G. Carne, Delayed, multi-step inverse structural filter for robust force identification, Mechanical Systems and Signal Processing 22 (2008) 1036–1054. [11] B.J. Qiao, X.F. Chen, X.F. Xue, X.J. Luo, R.N. Liu, The application of cubic B-spline collocation method in impact force identification, Mechanical Systems and Signal Processing 64-65 (2015) 413–427. [12] Y.M. Mao, X.L. Guo, Y. Zhao, A state space force identification method based on Markov parameters precise computation and regularization technique, Journal of Sound and Vibration 329 (2010) 3008–3019. [13] Z. Li, Z.P. Feng, F.L. Chu, A load identification method based on wavelet multi-resolution analysis, Journal of Sound and Vibration 333 (2014) 381–391. [14] Y. Ding, S.S. Law, B. Wu, G.S. Xu, Q. Lin, H.B. Jiang, Q.S. Miao, Average acceleration discrete algorithm for force identification in state space, Engineering Structures 56 (2013) 1880–1892. [15] K. Liu, S.S. Law, X.Q. Zhu, Y. Xia, Explicit form of an implicit method for inverse force identification, Journal of Sound and Vibration 333 (2014) 730–744. [16] X. Xu, J.P. Ou, Force identification of dynamic systems using virtual work principle, Journal of Sound and Vibration 337 (2015) 71–94. [17] T. Wang, Z.M. Wan, X.L. Wang, Y.J. Hu, A novel state space method for force identification based on the Galerkin weak formulation, Computers and Structures 157 (2015) 132–141. [18] B.J. Qiao, X.W. Zhang, X.J. Luo, X.F. Chen, A force identification method using cubic B-spline scaling functions, Journal of Sound and Vibration 337 (2015) 28–44. [19] D. Bernal, A. Ussia, Sequential deconvolution input reconstruction, Mechanical Systems and Signal Processing 50-51 (2015) 41–55. [20] E. Lourens, E. Reynders, G. De Roeck, G. Degrande, G. Lombaert, An augmented Kalman filter for force identification in structural dynamics, Mechanical Systems and Signal Processing 27 (2012) 446–460. [21] S.E. Azam, E. Chatzi, C. Papadimitriou, A dual Kalman filter approach for state estimation via output-only acceleration measurements, Mechanical Systems and Signal Processing 60-61 (2015) 866–886. [22] F. Naets, J. Cuadrado, W. Desmet, Stable force identification in structural dynamics using Kalman filtering and dummy-measurements, Mechanical Systems and Signal Processing 50-51 (2015) 235–248. [23] K. Maes, A.W. Smyth, G. De Roeck, G. Lombaert, Joint input-state estimation in structural dynamics, Mechanical Systems and Signal Processing 70-71 (2016) 445–466. [24] D. Bernal, Optimal discrete to nontinuous transfer for band limited inputs, Journal of Engineering Mechanics 133 (2007) 1370–1377.
Contents lists available at ScienceDirect
Journal of Sound and Vibration journal homepage: www.elsevier.com/locate/jsvi
Parametric study on sequential deconvolution for force identification Tao Lai a, Ting-Hua Yi a,n, Hong-Nan Li a,b a b
School of Civil Engineering, Dalian University of Technology, Dalian 116023, China School of Civil Engineering, Shenyang Jianzhu University, Shenyang 110168, China
a r t i c l e in f o
abstract
Article history: Received 20 November 2015 Received in revised form 5 February 2016 Accepted 8 May 2016 Handling Editor: L. G. Tham Available online 24 May 2016
The force identification can be mathematically viewed as the mapping from the observed responses to external forces through a matrix filled with system Markov parameters, which makes it difficult or even impossible for long time duration. A potentially efficient solution is to sequentially perform the identification processing. This paper presents a parametric study on the sequential deconvolution input reconstruction (SDR) method, which was proposed by Bernal. The behavior of the SDR method due to the effects of window parameters, noise levels and sensor configurations is investigated. In addition, a new normalized standard deviation of the reconstruction error over time is derived to cover the effect of only independent noise entries. The sinusoidal and band-limited white noise excitations are tested to be identified with good accuracy even with 10% noise. The simulation results yield various conclusions that may be helpful to engineering practitioners. & 2016 Elsevier Ltd. All rights reserved.
Keywords: Force identification Sequential deconvolution Parametric study
1. Introduction The force identification problem has recently attracted considerable attention in structural health monitoring, structural control and design communities [1–3]. Understanding the dynamic load exerted on structures is beneficial to engineering practitioners because structural integrity and safety are significantly involved with dynamic loading. However, in practice, directly measuring the dynamic force is generally difficult or even impossible, such as wind loads and earthquake excitations. The alternative is the recovery of inputs from vibration measurements. Nevertheless, the number of vibration sensors is usually limited for economic or technological reasons. Therefore, it is necessary to develop force identification techniques with incomplete measurements to evaluate the external loadings applied on the structures. Force identification from the observed responses is a typical inverse problem in structural dynamics. The reconstruction results usually suffer from numerical non-uniqueness and instability due to matrix ill-posedness and noise disturbance [4]. Regularization methods such as the widely used Tikhonov regularization and truncated singular value decomposition (TSVD) method are adopted to guarantee an accurate solution. A recent thorough review of the literature on force identification with observed measurements is found in Sanchez et al. [5]. The force identification problem is initially solved in frequency domain based on the fast Fourier transform (FFT) [6]. Ramkrishna et al. [7] proposed a method for vibration testing of flexible structures with electronic equipment in frequency n
Corresponding author at: School of Civil Engineering, Dalian University of Technology, Dalian 116023, China. Tel./fax: þ86 411 84706050. E-mail address: [email protected] (T.-H. Yi).
http://dx.doi.org/10.1016/j.jsv.2016.05.013 0022-460X/& 2016 Elsevier Ltd. All rights reserved.
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
77
domain. Chao et al. [8] presented a frequency-domain force identification technique considering the nonlinearity of the system. The FFT-based methods have two main adverse problems to address: these methods are limited to cases with zero initial conditions, and the hypotheses of periodic inputs result in a serious leakage problem particularly for short load durations or non-stationary signals [9–11]. In contrast, the more recent development of time-domain force identification techniques can avoid such problems. Mao et al. [12] estimated the inputs with high accuracy by introducing a precise timestep integration method in the state-space domain. Li et al. [13] decomposed the excitations and responses in wavelet domain and solved the force identification problem with Tikhonov method. Ding et al. [14] proposed a discrete average acceleration algorithm for force identification based on the Newmark-β algorithm. Liu et al. [15] extended this method and modified the conventional implicit Newmark-β algorithm into an explicit form for input estimation with a regularization method. Xu et al. [16] approximated the inputs and outputs locally by moving least-square method and transformed the differential equation of motion into an integral equation with the virtual work principle to iteratively identify inputs. Wang et al. [17] discretized the time element of inputs and outputs and recover the inputs in the state-space method based on weak Galerkin formulation. Qiao et al. [18] represented the desired force by the cubic B-spline scaling function based on the wavelet multi-resolution concept, and the force identification problem is equivalent to finding the coefficients of the scaling functions. This method was verified in an impact force identification experiment [11]. These studies of the time-domain techniques can be categorized into one group, which comprises two main steps: the first is to mathematically formulate the mapping from the sequences of inputs to outputs through an operation matrix; the second step is to solve the ill-posed inverse problem with the regularization method. However, the force identification problem is difficult or even impossible when the size of the operation matrix increases. A potentially effective and efficient method is to sequentially extract the inputs. Bernal et al. [19] proposed a sequential deconvolution input reconstruction (SDR) method, which can overcome this adverse problem and make online identification possible. In addition, Kalman filtering is another popular method that has been introduced for force identification [20–23]. In this paper, the parametric study on the SDR method is presented. The behavior of the SDR method due to the effects of window parameters, noise levels and sensor configurations is investigated. In addition, a new normalized standard deviation of the input identification error is derived to interpret the accuracy of the SDR method. The paper is organized as follows. Section 2 presents the mapping from the sequences of inputs to outputs. In Section 3, the first p values of the input vectors are uniquely extracted in the “window”. Section 4 is devoted to the stability and regularization problems. The parametric study of the SDR method is numerically illustrated in Section 5. Finally, some conclusions are summarized.
2. The mapping from the sequences of input to output 2.1. Fundamental formulation The equation of motion of a linear time-invariant system with n degrees of freedoms (DOFs) can be described by € þCuðtÞ _ þ KuðtÞ ¼ GfðtÞ ΜuðtÞ
(1)
where M A ℝnn ,C A ℝnn and K A ℝnn are the mass, damping and stiffness matrices, respectively; uðtÞ A ℝn1 is the displacement vector, fðtÞ A ℝr1 is the force vector applied to the system through the mapping matrix G A ℝnr , and r is the number of inputs. h iT _ T , Eq. (1) can be converted into the state-space representation as Defining the state vector x ¼ uT u _ ¼ Ac xðtÞ þBc fðtÞ xðtÞ
(2)
where Ac ¼
0 M
M
K
Bc ¼
I
1
M
C
A ℝ2n2n
(3-1)
0 1
1
G
A ℝ2nr
(3-2)
The general system output can be a linear combination of displacements, velocities and accelerations, which are grouped into an output vector yðtÞ A ℝm1 yðtÞ ¼ Cc xðtÞ þ Dc fðtÞ
(4)
where Cc A ℝ depends on the type of measurements, Dc A ℝ is the direct transmission term, and m is the number of outputs. The exact solution to Eq. (2) is given by Z t eAc ðt τÞ Bc fðτÞdτ (5) xðtÞ ¼ eAc t xð0Þ þ mn
mr
0
78
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Substituting Eq. (5) into Eq. (4) yields Z yðtÞ ¼ Cc eAc t xð0Þ þCc
t 0
eAc ðt τÞ Bc fðτÞdτ þ Dc fðtÞ
(6)
The continuous-time state-space model in Eqs. (5) and (6) can be discretized into the discrete-time state-space model, which is represented by Eqs. (7) and (8) with certain assumptions of the inter-sample behavior [24] (7) x k ¼ Ad x k 1 þBd f k 1 y k ¼ Cd x k þ Dd f k
(8)
2.2. Formulating the mapping from the sequences of inputs to outputs Based on the iterative equation in Eqs. (7) and (8), the state vector and observation at time kT, where T is the sampling time interval, can be formulated by kX 1 x k ¼ Akd x½0 þ Akd l 1 Bd f l
(9)
l¼0 kX 1 Cd Akd l 1 Bd f l þDd f k y k ¼ Cd Akd x½0 þ
(10)
l¼0
Stacking the sequences of inputs and outputs in columns yields h λ1 x½λ ¼ Aλd x½0 þ Ad Bd
2
Aλd 2 Bd
3 2 3 2 Dd Cd y ½ 0 7 6 6 7 6 6 Cd Bd 6 y ½ 1 7 6 Cd A d 7 7x½0 þ 6 6 7¼6 7 6 6 7 6 ⋮ ⋮ ⋮ 5 4 4 5 4 λ1 Cd A d Cd Aλd 2 Bd y½λ 1
3 f ½ 0 7 i6 6 f ½ 1 7 7 Bd 6 6 ⋮ 7 4 5 f ½λ 1
⋯
2
0
⋯
Dd
⋯
⋮
⋮
Cd Aλd 3 Bd
⋯
32
3 f ½0 76 7 0 76 f ½1 7 76 7 6 7 0 7 54 ⋮ 5 Dd f ½λ 1 0
(11)
(12)
Compactly write
where
x½λ ¼ Aλd x½0 þCoλ f~ 0jλ 1
(13)
y~ 0jλ 1 ¼ Obλ x½0 þHλ f~ 0jλ 1
(14)
h λ1 Coλ ¼ Ad Bd
Aλd 2 Bd 2
Cd
⋯
i Bd A ℝ2λnλr
3
6 7 6 Cd A d 7 7 A ℝλm2λn Obλ ¼ 6 6 7 ⋮ 4 5 λ1 Cd A d 2 6 6 Hλ ¼ 6 6 4
(15)
Dd
0
⋯
0
C d Bd ⋮
Dd ⋮
⋯ ⋮
0 0
Cd Aλd 2 Bd
Cd Aλd 3 Bd
⋯
Dd
(16)
3 7 7 7 A ℝλmλr 7 5
Eqs. (13) and (14) can be written in a general iterative formulation, which is used in later sections: x l þ λ ¼ Aλd x l þCoλ f~ ljl þ λ 1 y~ ljl þ λ 1 ¼ Obλ x l þ Hλ f~ ljl þ λ 1
(17)
(18) (19)
where Coλ is the controllability matrix (in reversed order), Obλ is the observability matrix; Ηλ is the lower block triangular Toeplitz matrix filled with system Markov parameters; l and λ are positive integers. From Eqs. (12) or (19), the sequences of inputs are mathematically mapped into the outputs through matrix Ηλ .
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
79
3. Force identification in a moving window From linear algebra theory, whether Eq. (19) has a solution fully depends on the characteristic of the matrix Hλ , to be more specific, the direct transmission term. In other words, if the sensors are collocated with the inputs, the full column rank of Ηλ ensures that the force identification problem has a unique solution in the least-square sense. On the contrary, matrix Ηλ is singular in the non-collocated situation, which results in an unstable solution. Therefore, the regularization method is used to impose an additional constraint and obtain a stabilized solution. Another concern is that the computation efficiency dramatically decreases when the size of matrix Ηλ increases. A potentially effective and efficient solution is to extract the inputs in a “window” that moves along the time axis, which makes online (or near online) force identification possible. The first step of the SDR method is to determine the length of the “window”. Herein, it is assumed that the window length is λ. Then, an identifiable matrix Q p is defined to recover the first p ðp r λÞ values of the input vector in the “window”. h i (20) Q p ¼ Ipr 0prðλ pÞr where Ipr represents the identity matrix of order p r. In essence, parameter p indicates how many time stations through which the “window” advances at one time. In this sense, p determines how “fast” the “window” can move. Simultaneously, the initial condition of the “window” in the next movement can be propagated through Eq. (18). From Eq. (19), the uniqueness of the solution within a window requires that the necessary and sufficient condition is fulfilled: the block diagonal element of Hλ in Eq. (17) must have full column rank, that is, m Zr, namely, the number of sensors must be no less than the number of inputs to be identified. This assumption is used throughout the paper. Then, the general solution of Eq. (19) is given by (21) f~ ljl þ λ 1 ¼ H†λ y~ ljl þ λ 1 Obλ x l þ Z U h where H†λ denotes the pseudo-inverse of Hλ , Z is the null space of Hλ , and h is an arbitrary coefficient vector with appropriate dimension. From a mathematical perspective, the solution has two parts: the optimal solution and an arbitrary combination of basis vector in the null space of Hλ . Then, the first p values of the input vector are selected by post-multiplying Eq. (21) by Eq. (20) f~ ljl þ p 1 ¼ Q p H†λ y~ ljl þ λ 1 Obλ x l þQ p Z U h (22) To guarantee the uniqueness of the identified force, the last term in Eq. (22) should be zero. Then, Q p Z should be zero because h is arbitrary. As a result, Zhas to be of a form given by " # 0prðÞ Z¼ ~ (23) Ζ ðλpÞrðÞ
4. Sequential deconvolution process From Eqs. (22) and (23), the input solution in the window can be uniquely identified by n o f~ ljl þ p 1 ¼ Q p H†λ y~ ljl þ p 1 Obλ x l
(24)
Once these inputs are identified in a time window, the initial condition in the next time window can be assessed through state propagation in Eq. (18). The identical process can be sequentially performed until the full duration of interest is covered. Unfortunately, the uncertainty of observations also makes the error propagate in the state vector and the identified input vector when the “window” slides. The stability of input reconstruction should be guaranteed during the sequential processing. Now, it is opportune to determine the window parameters λ and p. 4.1. Error propagation analysis It is necessary to investigate the error propagation of input reconstruction when the “window” advances along the time axis. Based on the stability analysis, the restrictive conditions on the parameters λ and p are expected. The stability of sequential deconvolution was originally derived in reference [19], and the derivation is repeated here for convenience in the subsequent discussion about conditioning. The measured responses are usually contaminated by noise, which leads to the uncertainty in the state vector and inputs. Using Eq. (24), one can write (25) f~ ljl þ p 1 þδf~ ljl þ p 1 ¼ Q p H†λ y~ ljl þ λ 1 þ v~ ljl þ λ 1 Obλ x l þ δx l
80
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Therefore δf~ ljl þ p 1 ¼ Q p H†λ v~ ljl þ λ 1 Obλ δx l
(26)
From Eq. (18), the error propagation in the state vector can be given by δx l þ p ¼ Apd δx l þ Cop δf~ ljl þ p 1
(27)
Substituting Eq. (26) into Eq. (27) produces
δx l þp ¼ Pp δx l þ Zp v~ ljl þ λ 1
(28)
Pp ¼ Apd Cop Q p H†λ Obλ
(29)
Zp ¼ Cop Q p H†λ
(30)
where
Assume that the “window” starts at time l ¼ 0 and after k ¼ jp time stations where j ¼ 1; 2; 3⋯, the propagation error in state vector can be tracked as j1 X δx k ¼ Pjp δx½0 þ Pjp l 1 Zp v~ lpjlp þ λ 1
(31)
l¼0
Eq. (31) expresses the uncertainty of the state vector as a combination of the uncertainty of initial state vector and noise sequence vectors. As shown in reference [19], the necessary condition of the stability requirement in Eq. (31) should be satisfied by
ρðPp Þ ¼ ηj max r 1 (32) where η is the eigenvalue of Pp , and ρ is the spectral radius of Pp . 4.2. Conditioning Based on Eq. (31), the expression for the accumulated input identification error remains a combination of the error in the initial condition and the uncertainty of noise sequences, as derived in reference [19]. However, for the reader's convenience, the noise sequence order is changed to be consistent with the time order. After substituting Eq. (31) into Eq. (26) and performing some rearrangements, one obtains δf~ jpjjp þ p 1 ¼ Q p H†λ Obλ Pjp δx½0 þ
j X
Dl v~ lpjlp þ λ 1
(33)
l¼0
where Dl ¼ Q p H†λ Obi Pjp l 1 Zp A ℝprλm
(34)
Dj ¼ Q p H†λ A ℝprλm
(35)
To simplify, it is assumed that the initial condition is exactly known, or the initial error can be neglected. Eq. (38) can be simplified as a combination of the measurement noise vectors δf~ jp pjjp 1 ¼ DV h
where D ¼ D1
D2
⋯
i
h
Dj A ℝprðj þ 1Þλm and V ¼ v~ T0jλ 1
T v~ pjp þ λ 1
(36) ⋯
T v~ jpjjp þ λ 1
iT
A ℝðj þ 1Þλm1 . If one interests a spe-
cific of inputs, like the qth, then δf~ jp pjjp 1 ¼ d V q
where
q f~ jp pjjp 1
q
(37)
q
is the qth input, and d is the qth row of D. Note that the vector V is not strictly random because its length is
ðj þ1Þλm while there are only ðjp þ λÞm independent noise entries. The noise entries of vector V are not independent leading to large input identification error if Eq. (37) is used as in reference [19]. Herein, a new normalized standard deviation of the reconstruction error over time is derived to cover the effect of only independent noise entries. Reshaping vector V to retain all independent noise entries and the corresponding q coefficients vector d , then Eq. (37) can be rewritten as q q δf~ jp pjjp 1 ¼ d~ V~
(38)
q where d~ A ℝprðjp þ λÞm is the recombination of coefficients in d , and V~ A ℝðjp þ λÞm1 is composed of independent noise entries. The reshaping processing is illustrated in Appendix A. q
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
81
Designating the covariance of the noise vector V~ as Σ v , the uncertain in the identified force can be assessed by the error propagation given by T
σ 2δf q ¼ dq Σ v dq
(39)
Note that Σ v is a diagonal matrix because the noise entries are independent. It is assumed that all channels have the same measurement noise σ v , then Eq. (39) can be evaluated as σ δf q ¼ σ v
qffiffiffiffiffiffiffiffiffiffiffi T dq dq
(40)
Define the normalized standard deviation of the input identification error as σ δf q σv
¼
qffiffiffiffiffiffiffiffiffiffiffi T dq dq
(41)
Eq. (41) can be considered a good indicator to determine if an acceptable accuracy is achievable before the input identification. If the standard deviation of the error is acceptably small and stable over time station, the force identification performance can be guaranteed. If not, a regularization method is needed. 4.3. TSVD regularization Herein, a simple scheme based on TSVD regularization method is used if regularization is needed [4]. From (Eqs. (18) and 19), one can see that the associated matrices Aλd ; Coλ ; Hλ and Obλ do not change as the “window” slides. This is helpful because TSVD only addresses the matrix Hλ once. For simplicity, the batch solution is approximately used to assess the
Fig. 1. Flowchart of the SDR method.
82
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
relative error of force identification. Taking the initial condition as zeros, Eq. (14) produces n o y~ 0jλ 1 þ v~ 0jλ 1 ¼ Hλ f~ 0jλ 1 þ δf~ 0jλ 1
(42)
From linear algebra theory, the relative error of force identification can be formulated as ‖δf~ 0jλ 1 ‖2 ‖v~ 0jλ 1 ‖2 r condðHλ Þ ‖y~ 0jλ 1 ‖2 ‖f~ ‖
(43)
0jλ 1 2
Define ε and ‖v~ 0jλ 1 ‖2 =‖y~ 0jλ 1 ‖2 as the bounded tolerance of the relative identification error and noise-to-signal ratio (NSR), respectively, then condðHλ Þ
‖v~ 0jλ 1 ‖2 σ max ¼ U NSR r ε ‖y~ 0jλ 1 ‖2 σ min
(44)
Hence σ min NSR Z ε σ max
(45)
Eq. (45) is considered as the criterion of TSVD regularization given the tolerance εand NSR. The bounded tolerance ε ¼ 1 is assumed throughout the paper. 4.4. Need for regularization From Eq. (41), one can plot a figure of the normalized standard deviation of input reconstruction error σ δf =σ v verse time stations and determine whether regularization is necessary. Because regularization controls the norm and distort the true result, it is important to maintain NSR=ε as small as possible [19]. Eq. (41) can be accounted for the input reconstruction performance when the regularization breaks the theoretical thread. Finally, a flowchart to implement the SDR method is schematically demonstrated in Fig. 1. Two schemes are available if the identification performance is not guaranteed. Scheme 1 is to adjust window parameters λ and p to reassess the identification performance. Alternatively, one can adopt the TSVD regularization in scheme 2.
5. Numerical illustration The parametric study of the SDR method is demonstrated with a numerical example given by Bernal [19]. For unit consistency, the physical parameters of the bar are set as follows: bar length L ¼ 10, cross-sectional area A ¼ 4, mass density ρd ¼ 1 and modulus of elasticity E ¼ 100. The bar is discretized into 40 elements with one end fixed and the other free shown in Fig. 2. The finite element model has a minimum and maximum frequency of 0:25 Hz and 12:73 Hz, respectively. The recommendation of time step size is T ¼ 0:05 s. Rayleigh damping with 2% in the first and last mode is used. The time histories of two decaying harmonic excitation applied at node 5 and 40, are given by f 5 ðt Þ ¼ 10e 0:3t ð1 cos 2πt Þ sin 6πt (46-1) f 40 ðt Þ ¼ 10e 0:2t ð1 cos πt Þ sin 4πt
(46-2)
The SDR accuracy is quantified by normalized root mean square (NRMS) error given by NRMS ¼
‖f e f r ‖2 100% ‖f r ‖2
(47)
where f e and f r are the estimated and real force vectors, respectively. Due to noise uncertainty, each simulation is independently repeated 50 times to account for the accuracy and robustness of the SDR method. 5.1. Sequential force identification results Herein, 4 accelerometers are assumed to be available located at coordinates {1, 10, 20, 30}. The window parameters λ ¼ 50; p ¼ 1 are used and the selection detail is discussed in a later section. To simulate the measurements, 2% RMS of white noises is added to the actual responses. As Section 4 mentioned, the accuracy of the force identification can be assessed by
Fig. 2. Cantilever bar with 40 DOFs.
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
83
Eq. (41) before input identification. Therefore, the plot of the normalized standard deviation of the input identification error versus time stations is first investigated, and the results are shown in Figs. 3 and 4. For the first input, the normalized standard deviation dramatically increases when the window slides, as shown in Fig. 3. In Fig. 4, the normalized standard deviation is stable after approximately 800 time stations but remains at a large value. The results indicate that the regularization method is necessary. Figs. 3 and 4 show that the normalized standard deviation after regularization is always maintained at a small value for both inputs. Finally, the input reconstruction results are illustrated in Figs. 5 and 6. It is observed that the identified inputs without regularization are deviating from the true inputs over time resulting in large error in force identification. However, after regularization, the identified inputs matches with the exact inputs quite well with NRMS error of 3:50% and 2:40%, respectively. 5.2. Effects analysis Herein, a detail study of the effects including window parameters, measurement noise levels and sensor configurations on the accuracy of the force identification is carried out. At the beginning, the selection of window parameters is investigated. Second, to account for the effect of measurement noise, four levels of Gaussian white noise, i.e.: 0%; 2%; 5% and 10% are studied. Finally, another three sets of sensor arrangements are considered except the sensor layout in the last subsection, as shown in Table 1. The sensors are approximately equal spaced at 1=3; 1=4 and 1=5 length of the bar for sensor cases S1, S2 and S3.
Fig. 3. Normalized standard deviation vs. time station for the first input.
Fig. 4. Normalized standard deviation versus time station for the second input.
84
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Fig. 5. SDR for the input reconstruction without regularization: (a) first input; (b) second input.
Fig. 6. SDR for the input reconstruction with regularization: (a) first input; (b) second input. Table 1 Sensor configuration schemes. Sensor configuration schemes
Sensor sets
S1
S2
S3
S4
{13,26}
{10,20,30}
{8,16,24,32}
{1,10,20,30}
(1) Window parameters λ; p . In Section 4, the stability requirement must be satisfied by Eq. (32), which is the only necessary condition. One can plot the 2D figure of the largest eigenvalue versus parameter pair λ; p , which is illustrated in Fig. 7. It is observed that there is not too much difference of the stabilized zoon, which is concentrated upon the bottom right in Fig. 7. However, the largest eigenvalue of Pp is quite different for regularization and
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
85
non-regularization situations. For illustration, one considers the specified situations p ¼ 1 and p ¼ 5 versus the varyingλ, as shown in Fig. 8. The zoom-in parts show that the largest eigenvalue with regularization is less than that without regularization after approximately λ ¼ 35. Mathematically, the smaller of the largest eigenvalue (less than 1), the faster of the convergence rate of matrix power of Pp . As a result, λ ¼ 30; 40; 50; 60; 70 and 80 with p ¼ 1 or 5 are considered. There are 12 numerical identification results in Figs. 9 and 10. One can see that the NRMS error decreases after λ ¼ 30 and remains at a low level for both p ¼ 1 and p ¼ 5. However, there are large errors for both cases without regularization. For computation efficiency, one takes λ ¼ 50 in the subsequent analysis. (2) Noise level. Herein, one takes λ ¼ 50; p ¼ 1 for the noise level effect analysis. The identification results are summarized in Table 2. The noise significantly affects the accuracy of the input reconstruction. The identification results with or without regularization are identical when there is absence of noise. When the noise levels increases, the performance with regularization slowly deteriorates in the beginning and dramatically increases with 10% noise. For the nonregularization situation, there is large identification error as noise level increases. It is demonstrated that the noise effect is an important factor to be considered for the ill-posed inverse problem. (3) Sensor configuration. Four schemes of sensor configuration are carried out. For each case, one takes λ ¼ 50; p ¼ 1, and 2% white noise is considered. The results are shown in Fig. 11. Compared with S1, scheme S2 has much better identification
Fig. 7. Largest eigenvalue of Pp versus λ; p : (a) non-regularization; (b) regularization.
Fig. 8. Largest eigenvalue Pp versus λ: (a) p ¼ 1; (b) p ¼ 5.
86
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
Fig. 9. Input identification error versus λðp ¼ 1Þ: (a) regularization; (b) non-regularization.
Fig. 10. input identification error versus λðp ¼ 5Þ: (a) regularization (b) non-regularization. Table 2 Input identification error with different noise levels. NRMS error (%)
Non-regularization Regularization
Noise level
Input Input Input Input
I II I II
0%
2%
5%
10%
4.68 2.09 4.70 2.10
97.21 26.75 3.50 2.40
247.64 61.52 4.98 6.39
420.88 111.29 16.19 29.16
accuracy. However, S3 does not outperform the accuracy of S2, which uses one more sensor compared with S3. The scheme S4 has the best accuracy either with regularization or without regularization. 5.3. Band-limited white noise identification To further illustrate the applicability of the SDR method, band-limited white noise excitation is tested. The band-limited white noise that was applied at node 40 of the bar is generated in Simulink block and passes through a low-pass filter with 10 Hz cutoff frequency. Herein, it is assumed that λ ¼ 50; p ¼ 1 is used. 2% and 10% white noise levels are considered to mimic the true acceleration responses collected at coordinates {1, 10, 20, 30}. First, one should check the plot of the normalized standard deviation versus time stations. It is noteworthy that the normalized standard deviation is not related to the external excitation type. As a result, the trend of the normalized standard deviation over time stations is identical to that
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
87
Fig. 11. Input identification error with different sensor configuration schemes: (a) regularization; (b) non-regularization.
Fig. 12. SDR for band-limited white noise identification with regularization: (a) 2% noise; (b) 10% noise.
in Fig. 4, which indicates that the regularization is required. The identification results are shown in Fig. 12. It is observed that the identified band-limited white noise matches with the true input quite well with NRMS error of 4:61% under low noise level. In contrast, the identified input with 10% white noise has a large error with an NRMS of 18:45%.
6. Conclusion This paper presents a parametric study on the sequential deconvolution input reconstruction method. A new normalized standard deviation of the input identification error over time stations is defined to interpret the accuracy of reconstruction. The effects of the window parameters λ; p , noise levels and sensor configuration schemes are investigated. The sinusoidal and band-limited white noise excitations are tested. Numerical analysis of a simple bar is performed to demonstrate the accuracy and efficiency of the SDR method. The following conclusions are summarized as follows: (1) The selection of window parameters λ; p is crucial for the SDR method. Without the right choice, the input identification dramatically changed, which suggested a need for a larger λ and a smaller p on the premise of the necessary condition and computation efficiency consideration.
88
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
(2) The noise level significantly affects the accuracy of the input identification. One can achieve good reconstruction accuracy with 5% white noise and acceptable accuracy with 10% white noise after regularization. (3) The SDR accuracy dramatically varies with the sensor configurations. A good sensor configuration can reduce the illposedness of the input identification. The sensor layout remains an open problem for input reconstruction. (4) Single or multiple excitations of sinusoid and band-limited white noise are tested to be identified with good accuracy even with 10% noise.
Acknowledgments The authors would like to appreciate Professor Bernal for the innovation method and his help to clarify many related concepts and ideas. This research work was jointly supported by the 973 Program (Grant no. 2015CB060000), the National Natural Science Foundation of China (Grant nos. 51421064, 51478081), the Fok Ying Tong Education Foundation (Grant no. 141072), and the Science Fund for Distinguished Young Scholars of Dalian (Grant nos. 2014J11JH125, 2015J12JH209).
Appendix A. Illustration of the normalized standard deviation It is assumed that there are r ¼ 2 inputs and m ¼ 3 sensors available. Herein, one takes λ ¼ 5; p ¼ 1. Suppose that we want to characterize the reconstruction error at t ¼ 2T. As a result, we have j ¼ 2 in Eq. (33) written by δf~ 2pj3p 1 ¼ Q p H†λ Obλ Pjp δx½0 þ D0 v~ 0jλ 1 þD1 v~ pjp þ λ 1 þ D2 v~ 2pj2p þ λ 1
(A-1)
In matrix formulation, one has 2 3 ~1 δ f row#1 row#1 row#1 row#1 2pj3p 1 6 7 v~ 0jλ 1 þ v~ pjp þ λ 1 þ v~ 2pj2p þ λ 1 δx½0 þ 4 2 5¼ row#2 D1 row#2 D2 row#2 row#2 D0 δf~ 2pj3p 1
(A-2)
The rearrangement of Eq. (A-2) gives 1 δf~ 2pj3p 1
2 δf~ 2pj3p 1
h
¼ ½row#1δx½0 þ row#1D0
h
¼ ½row#2δx½0 þ row#2D0
row#1D1
row#2D1
i
2
v~ 0jλ 1
3
7 row#1D2 6 4 v~ pjp þ λ 1 5 v~ 2pj2p þ λ 1
row#2D2
i
(A-3)
2
3 v~ 0jλ 1 6 v~ 7 4 pjp þ λ 1 5 v~ 2pj2p þ λ 1
(A-4)
Vector V A ℝðj þ 1Þλm1 is composed of all noise sequencesv~ jpjjp þ λ 1 2 ~ 3 v 0jλ 1 6 7 V ¼ 4 v~ pjp þ λ 1 5 v~ 2pj2p þ λ 1 where
h ð1Þ v~ 0jλ 1 ¼ v0j4
2Þ vð0j4
3Þ vð0j4
1Þ vð1j5
2Þ vð1j5
3Þ vð1j5
1Þ vð2j6
2Þ vð2j6
(A-5)
3Þ vð2j6
1Þ vð3j7
h ð1Þ v~ pjp þ λ 1 ¼ v1j5
2Þ vð1j5
3Þ vð1j5
1Þ vð2j6
2Þ vð2j6
3Þ vð2j6
1Þ vð3j7
2Þ vð3j7
3Þ vð3j7
1Þ vð4j8
h ð1Þ v~ 2pj2p þ λ 1 ¼ v2j6
2Þ vð2j6
3Þ vð2j6
1Þ vð3j7
2Þ vð3j7
3Þ vð3j7
1Þ vð4j8
2Þ vð4j8
3Þ vð4j8
1Þ vð5j9
2Þ vð3j7 2Þ vð4j8 2Þ vð5j9
3Þ vð3j7 3Þ vð4j8 3Þ vð5j10
1Þ vð4j8 1Þ vð5j9 1Þ vð6j11
2Þ vð4j8
3Þ vð4j8
2Þ vð5j9 2Þ vð6j11
i
3Þ vð5j10
(A-6) i
3Þ vð6j11
(A-7) i
(A-8)
The corresponding coefficient vectors d are
h i 1 d ¼ row#1D0 row#1D1 row#1D2 h i ð1Þ ð1Þ ð1Þ d2 ⋯ d45 ¼ d1
(A-9)
h i 2 d ¼ row#1D0 row#1D1 row#1D2 h i ð2Þ ð2Þ ð2Þ d2 ⋯ d45 ¼ d1
(A-10)
T. Lai et al. / Journal of Sound and Vibration 377 (2016) 76–89
89
Note that vector V is not strictly random because its length is ðjþ 1Þλm ¼ 45 while there are only ðjp þλÞm ¼ 21 independent noise entries. Now, the modified coefficient vector V~ A ℝ211 filled with independent noise entries is defined as V~ ¼
1Þ vð0j4
2Þ vð0j4
3Þ vð0j4
1Þ vð1j5
2Þ vð1j5
3Þ vð1j5
1Þ vð2j6
2Þ vð2j6
3Þ vð2j6
1Þ vð3j7
2Þ vð3j7
3Þ vð3j7
1Þ vð4j8
2Þ vð4j8
3Þ vð4j8
1Þ vð5j9
2Þ vð5j9
3Þ vð5j10
1Þ vð6j11
2Þ vð6j11
3Þ vð6j11
(A-11) The corresponding coefficient vectors d are reshaped as d~ 2 ð1Þ ð1Þ ð1Þ d1 d2 d3 6 1 ð1Þ ð1Þ ð1Þ ð1Þ ð1Þ d~ ¼ 6 4 d6 þ d18 d7 þ d19 þ d31 ⋯ ð1Þ ð1Þ ð1Þ ð1Þ ð1Þ d29 þ d41 d30 þd42 d43 2
ð2Þ
d1 6 ð2Þ ð2Þ ~ 6 d ¼ 4 d6 þ d18 ð2Þ ð2Þ d29 þ d41 2
ð2Þ
ð2Þ
d2 ð2Þ
ð2Þ
ð2Þ
d7 þ d19 þ d31 ð2Þ
ð1Þ
d30 þd42
ð1Þ
ð1Þ
ð1Þ
d4 þ d16 ð1Þ
ð1Þ
ð1Þ
d15 þ d27 þd39 ð1Þ
ð2Þ
ð2Þ
d4 þ d16
⋯
d15 þ d27 þd39
ð2Þ
d43
ð2Þ ð2Þ
d44
ð1Þ ð1Þ 7 d28 þ d40 7 5
(A-12)
d45
d3
ð2Þ
3
ð1Þ
d44 ð2Þ
ð1Þ
d5 þ d17
ð2Þ
d5 þ d17 ð2Þ
3
ð2Þ ð2Þ 7 d28 þ d40 7 5
(A-13)
ð2Þ
d45
1 2 Note that vectors d~ and d~ are row vectors, which are written as the form in Eqs. (A-12) and (A-13) for conciseness.
References [1] S.S. Law, J.Q. Bu, X.Q. Zhu, Time-varying wind load identification from structural responses, Engineering Structures 27 (2005) 1586–1598. [2] Y. Niu, C.P. Fritzen, H. Jung, I. Buethe, Y.Q. Ni, Y.W. Wang, Online simultaneous reconstruction of wind load and structural responses-theory and application to Canton Tower, Computer-Aided Civil and Infrastructure Engineering 30 (2015) 666–681. [3] A.L. Wu, C.H. Loh, J.N. Yang, J.H. Weng, C.H. Chen, T.S. Ueng, Input force identification: application to soil–pile interaction, Structural Control Health Monitoring 16 (2009) 223–240. [4] P.C. Hansen, Rank-Deficient and Discrete Ill-posed Problems: Numerical Aspects of Linear Inversion, SIAM, Philadelphia, 1998. [5] J. Sanchez, H. Benaroya, Review of force reconstruction techniques, Journal of Sound and Vibration 333 (2014) 2999–3018. [6] J.F. Doyle, Further developments in determining the dynamic contact law, Experimental Mechanics 24 (1984) 265–270. [7] D. Ramkrishna, A. Joshi, P.M. Mujumdar, Y. Krishna, Frequency domain input identification for vibration testing of flexible structures, Proceedings of the Institution of Mechanical Engineers Part G – Journal of Aerospace Engineering 227 (2013) 838–856. [8] M. Chao, H. Hongxing, X. Feng, The identification of external forces for a nonlinear vibration system in frequency domain, Proceedings of the Institution of Mechanical Engineers Part C – Journal of Mechanical Engineering Science 228 (2014) 1531–1539. [9] L.J.L. Nordstrom, T.P. Nordberg, A time delay method to solve non-collocated input estimation problems, Mechanical Systems and Signal Processing 18 (2004) 1469–1483. [10] M.S. Allen, T.G. Carne, Delayed, multi-step inverse structural filter for robust force identification, Mechanical Systems and Signal Processing 22 (2008) 1036–1054. [11] B.J. Qiao, X.F. Chen, X.F. Xue, X.J. Luo, R.N. Liu, The application of cubic B-spline collocation method in impact force identification, Mechanical Systems and Signal Processing 64-65 (2015) 413–427. [12] Y.M. Mao, X.L. Guo, Y. Zhao, A state space force identification method based on Markov parameters precise computation and regularization technique, Journal of Sound and Vibration 329 (2010) 3008–3019. [13] Z. Li, Z.P. Feng, F.L. Chu, A load identification method based on wavelet multi-resolution analysis, Journal of Sound and Vibration 333 (2014) 381–391. [14] Y. Ding, S.S. Law, B. Wu, G.S. Xu, Q. Lin, H.B. Jiang, Q.S. Miao, Average acceleration discrete algorithm for force identification in state space, Engineering Structures 56 (2013) 1880–1892. [15] K. Liu, S.S. Law, X.Q. Zhu, Y. Xia, Explicit form of an implicit method for inverse force identification, Journal of Sound and Vibration 333 (2014) 730–744. [16] X. Xu, J.P. Ou, Force identification of dynamic systems using virtual work principle, Journal of Sound and Vibration 337 (2015) 71–94. [17] T. Wang, Z.M. Wan, X.L. Wang, Y.J. Hu, A novel state space method for force identification based on the Galerkin weak formulation, Computers and Structures 157 (2015) 132–141. [18] B.J. Qiao, X.W. Zhang, X.J. Luo, X.F. Chen, A force identification method using cubic B-spline scaling functions, Journal of Sound and Vibration 337 (2015) 28–44. [19] D. Bernal, A. Ussia, Sequential deconvolution input reconstruction, Mechanical Systems and Signal Processing 50-51 (2015) 41–55. [20] E. Lourens, E. Reynders, G. De Roeck, G. Degrande, G. Lombaert, An augmented Kalman filter for force identification in structural dynamics, Mechanical Systems and Signal Processing 27 (2012) 446–460. [21] S.E. Azam, E. Chatzi, C. Papadimitriou, A dual Kalman filter approach for state estimation via output-only acceleration measurements, Mechanical Systems and Signal Processing 60-61 (2015) 866–886. [22] F. Naets, J. Cuadrado, W. Desmet, Stable force identification in structural dynamics using Kalman filtering and dummy-measurements, Mechanical Systems and Signal Processing 50-51 (2015) 235–248. [23] K. Maes, A.W. Smyth, G. De Roeck, G. Lombaert, Joint input-state estimation in structural dynamics, Mechanical Systems and Signal Processing 70-71 (2016) 445–466. [24] D. Bernal, Optimal discrete to nontinuous transfer for band limited inputs, Journal of Engineering Mechanics 133 (2007) 1370–1377.