# Identification of Dynamic Loads Based on Second-Order Taylor-Series Expansion Method

Authors: Xiaowang Li; Zhongmin Deng
Journal: Shock and Vibration
 (2016)

---

## Abstract

A new method based on the second-order Taylor-series expansion is presented to identify the structural dynamic loads in the time domain. This algorithm expresses the response vectors as Taylor-series approximation and then a series of formulas are deduced. As a result, an explicit discrete equation which associates system response, system characteristic, and input excitation together is set up. In a multi-input-multi-output (MIMO) numerical simulation study, sinusoidal excitation and white noise excitation are applied on a cantilever beam, respectively, to illustrate the effectiveness of this algorithm. One also makes a comparison between the new method and conventional state space method. The results show that the proposed method can obtain a more accurate identified force time history whether the responses are polluted by noise or not.


---

## Body

## 1. Introduction

With the continuous development of science and technology nowadays, the mechanical structure becomes more and more complex. To ensure the structure’s reliability and safety, it is essential to get the dynamic loads acting on the structure. Generally speaking, directly measuring the external load through a load transducer is preferred. However, in many cases of engineering, such as the missile subjected to the wind load and the tall building suffering the seismic load, direct measurement of dynamic loads exerted on the structure is difficult to achieve. To overcome this problem, force identification, as the second class of inverse problem in structural dynamics, is put forward. Along with the increase of complex engineering problems, force identification technology has become a crucial issue in structural dynamics. It is devoted to providing effective load reference for the design of the structure and plays a significant role in mechanical vibration analysis, structural reliability analysis, mechanical fault diagnosis, and other fields. Thanks to the rapid advancement of computer technology, the establishment of mathematical model becomes easier and the accuracy of the finite element simulation is improved obviously. Force identification technology, therefore, has made a great progress.There are mainly two categories of methods for dynamic load identification: the frequency domain method and the time domain method [1]. The frequency method has developed fast because of its simple principle. It is focused on searching for a system’s frequency response function (FRF) and response spectrums to calculate load spectrums [2]. Some classical frequency domain methods are introduced in [3–5]. Nevertheless, the frequency domain method has the drawback of low accuracy and being sensitive to the response noise [6], which makes its application limited in some cases.In recent years, a series of time domain methods have been studied. Compared with the frequency domain method, force identification technology in time domain is able to achieve relatively higher accuracy. Law and Fang [7] adopted the state space method in the control theory to perform force identification. This method is commonly utilized to identify stationary periodic loads. Yan [8] used Bayesian approach in statistics to calculate impact loads. The author first transformed the impact load identification problem to a parameter identification problem by representing the impact load using a set of parameters. Sun et al. [9] proposed a new improved regularization method for load identification, which can overcome the ill-condition of load reconstruction to some extent. In [10], the authors transformed the conventional implicit N
e
w
m
a
r
k
-
β algorithm into an explicit form for linear analysis of the structure. Liu et al. [11] presented a shape function method of moving least square fitting, by which the time domain of load is discretized and the local load is approximated by SFM under LSF. An analytical method was proposed in [12] to identify dynamic loads acting on stochastic structures, based on the Gegenbauer polynomial expansion theory and regularization method. In [13], the authors put forward an inverse method that combines the interval analysis with regularization. This algorithm is able to stably identify the bounds of dynamic load acting on the uncertain structures. Besides, the authors in [14] utilized support vector regression to identify nonlinear systems represented by Wiener models. Simulation results show that the method gives accurate models of systems. Polynomial interpolated Taylor series method was studied in [15]. It advances the technique in parameter identification of structures with significant nonlinear response dynamics.Taylor formula [16, 17], as an indispensable math tool in mathematical analysis, plays a key role on approximate calculation. It aims to transform a complex function into a concise polynomial function on the premise of maintaining a high approximation precision. Due to the distinguished advantages, Taylor-series expansion has been used in numerous areas to deal with sophisticated mathematical problems.In this paper, one proposes a new approach for the identification of dynamic loads, utilizing the formula of Taylor-series expansion. The proposed method expresses the response vectors as Taylor-series approximation in the neighborhood of timet and then deduces a series of formulas for force identification. Finally, a simple explicit equation which links system characteristic, system response, and input excitation together is established. The full use of outstanding feature of Taylor formula makes this approach reach a high identification accuracy. Complex iteration calculation is eliminated because of the establishment of explicit equation. In addition, this method is an implicit integration in essence; hence it has the merits of remaining unconditional stable. Multi-input-multi-output (MIMO) numerical simulations are carried out to illustrate this method. State space method is adopted to make a contrast with this method. The results indicate that the proposed method can obtain more satisfactory identified force time histories even in the case that noise is added into the responses.


## 2. Force Identification Based on Taylor-Series Expansion

For a general linear elastic structure with multiple dofs [18], the equation of motion is expressed as(1)
M
x
¨
t
+
C
x
˙
t
+
K
x
t
=
F
t
,where M, C, and K denote the mass, damping, and stiffness matrices. F
t denotes the vector of forces exerted on the structure. x
¨
t, x
˙
t, and x
t denote the vectors of acceleration, velocity, and displacement responses.Assuming Rayleigh damping [19], the damping matrix is(2)
C
=
α
M
+
β
K
,where α and β are the damping coefficients.


### 2.1. Method Deduction

The nature of this method is to express the acceleration vector as a Taylor-series approximation in the neighborhood of timet; namely,(3)
x
¨
t
+
τ
=
x
¨
t
+
P
·
G
τ
,where P
=
P
1
P
2
⋯
P
n and P
k
(
k
=
1,2
,
…
,
n
) are vectors to be determined. Also G
τ
=
τ
τ
2
/
2
⋯
τ
n
/
n
!
T. Then the vectors of velocity and displacement responses are written as(4)
x
˙
t
+
τ
=
x
˙
t
+
∫
0
τ
x
¨
t
+
τ
d
τ
=
x
˙
t
+
τ
x
¨
t
+
P
·
∫
0
τ
G
τ
d
τ
,
x
t
+
τ
=
x
t
+
∫
0
τ
x
˙
t
+
τ
d
τ
=
x
t
+
τ
x
˙
t
+
τ
2
2
x
¨
t
+
P
·
∫
0
τ
∫
0
τ
G
τ
d
τ
d
τ
.For the purpose of removing the integral symbol in (4) and carrying out follow-up deduction, one must set n as a specific number. Meanwhile, considering the process of force identification in this paper is completed by second-order Taylor-series expansion method, one lets n
=
2, and then(5)
x
¨
t
+
τ
=
x
¨
t
+
τ
P
1
+
τ
2
2
P
2
,
x
˙
t
+
τ
=
x
˙
t
+
τ
x
¨
t
+
τ
2
2
P
1
+
τ
3
6
P
2
,
x
t
+
τ
=
x
t
+
τ
x
˙
t
+
τ
2
2
x
¨
t
+
τ
3
6
P
1
+
τ
4
24
P
2
.The equation of motion at timet
+
Δ
t is(6)
M
x
¨
t
+
Δ
t
+
C
x
˙
t
+
Δ
t
+
K
x
t
+
Δ
t
=
F
t
+
Δ
t
.Replaceτ with Δ
t and substitute (5) into (6) to get(7)
χ
11
P
1
+
χ
12
P
2
=
Q
1
,where (8)
χ
11
=
Δ
t
M
+
Δ
t
2
2
C
+
Δ
t
3
6
K
,
χ
12
=
Δ
t
2
2
M
+
Δ
t
3
6
C
+
Δ
t
4
24
K
,
Q
1
=
F
t
+
Δ
t
-
M
+
Δ
t
C
+
Δ
t
2
2
K
x
¨
t
-
C
+
Δ
t
K
x
˙
t
-
K
x
t
.Replaceτ with ε
Δ
t and substitute (5) into (6) to get(9)
χ
21
P
1
+
χ
22
P
2
=
Q
2
,where(10)
χ
21
=
ε
Δ
t
M
+
ε
Δ
t
2
2
C
+
ε
Δ
t
3
6
K
,
χ
22
=
ε
Δ
t
2
2
M
+
ε
Δ
t
3
6
C
+
ε
Δ
t
4
24
K
,
Q
2
=
F
t
+
ε
Δ
t
-
M
+
ε
Δ
t
C
+
ε
Δ
t
2
2
K
x
¨
t
-
C
+
ε
Δ
t
K
x
˙
t
-
K
x
t
.F
t
+
ε
Δ
t is obtained by linear interpolation between F
t
+
Δ
t and F
t
+
2
Δ
t; namely,(11)
F
t
+
ε
Δ
t
=
2
-
ε
F
t
+
Δ
t
+
ε
-
1
F
t
+
2
Δ
t
1
<
ε
<
2
.Combining (7) and (9), one gets(12)
χ
11
χ
12
χ
21
χ
22
P
1
P
2
=
Q
1
Q
2
.Let(13)
η
=
χ
11
χ
12
χ
21
χ
22
-
1
=
η
11
η
12
η
21
η
22
;then P
1 and P
2 are written as(14)
P
1
=
η
11
Q
1
+
η
12
Q
2
=
L
0
F
i
+
1
+
L
1
F
i
+
2
+
L
2
x
i
+
L
3
x
˙
i
+
L
4
x
¨
i
,
P
2
=
η
21
Q
1
+
η
22
Q
2
=
S
0
F
i
+
1
+
S
1
F
i
+
2
+
S
2
x
i
+
S
3
x
˙
i
+
S
4
x
¨
i
,where(15)
L
0
=
η
11
+
2
-
ε
η
12
,
L
1
=
ε
-
1
η
12
,
L
2
=
-
η
11
+
η
12
K
,
L
3
=
-
η
11
+
η
12
C
-
Δ
t
η
11
+
ε
η
12
K
,
L
4
=
-
η
11
+
η
12
M
-
Δ
t
η
11
+
ε
η
12
C
-
Δ
t
2
2
η
11
+
ε
2
η
12
K
,
S
0
=
η
21
+
2
-
ε
η
22
,
S
1
=
ε
-
1
η
22
,
S
2
=
-
η
21
+
η
22
K
,
S
3
=
-
η
21
+
η
22
C
-
Δ
t
η
21
+
ε
η
22
K
,
S
4
=
-
η
21
+
η
22
M
-
Δ
t
η
21
+
ε
η
22
C
-
Δ
t
2
2
η
21
+
ε
2
η
22
K
.Substituting (14) into (3) and (4) one obtains(16)
x
i
+
1
=
A
1
F
i
+
1
+
A
2
F
i
+
2
+
A
d
x
i
+
A
v
x
˙
i
+
A
a
x
¨
i
,
(17)
x
˙
i
+
1
=
B
1
F
i
+
1
+
B
2
F
i
+
2
+
B
d
x
i
+
B
v
x
˙
i
+
B
a
x
¨
i
,
(18)
x
¨
i
+
1
=
C
1
F
i
+
1
+
C
2
F
i
+
2
+
C
d
x
i
+
C
v
x
˙
i
+
C
a
x
¨
i
,where(19)
A
1
=
Δ
t
3
6
L
0
+
Δ
t
4
24
S
0
,
A
2
=
Δ
t
3
6
L
1
+
Δ
t
4
24
S
1
,
A
d
=
I
+
Δ
t
3
6
L
2
+
Δ
t
4
24
S
2
,
A
v
=
Δ
t
I
+
Δ
t
3
6
L
3
+
Δ
t
4
24
S
3
,
A
a
=
Δ
t
2
2
I
+
Δ
t
3
6
L
4
+
Δ
t
4
24
S
4
,
B
1
=
Δ
t
2
2
L
0
+
Δ
t
3
6
S
0
,
B
2
=
Δ
t
2
2
L
1
+
Δ
t
3
6
S
1
,
B
d
=
Δ
t
2
2
L
2
+
Δ
t
3
6
S
2
,
B
v
=
I
+
Δ
t
2
2
L
3
+
Δ
t
3
6
S
3
,
B
a
=
Δ
t
I
+
Δ
t
2
2
L
4
+
Δ
t
3
6
S
4
,
C
1
=
Δ
t
L
0
+
Δ
t
2
2
S
0
,
C
2
=
Δ
t
L
1
+
Δ
t
2
2
S
1
,
C
d
=
Δ
t
L
2
+
Δ
t
2
2
S
2
,
C
v
=
Δ
t
L
3
+
Δ
t
2
2
S
3
,
C
a
=
I
+
Δ
t
L
4
+
Δ
t
2
2
S
4
.Combine (16), (17), and (18) to get(20)
x
i
+
1
x
˙
i
+
1
x
¨
i
+
1
=
A
1
A
2
B
1
B
2
C
1
C
2
F
i
+
1
F
i
+
2
+
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
x
i
x
˙
i
x
¨
i
.Then the responses at timet
i are expressed as(21)
x
i
x
˙
i
x
¨
i
=
∑
j
=
0
i
-
1
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
j
A
1
A
2
B
1
B
2
C
1
C
2
F
i
-
j
F
i
-
j
+
1
+
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
i
x
0
x
˙
0
x
¨
0
.Zero initial responses are assumed and (21) is simplified to(22)
y
i
=
x
i
x
˙
i
x
¨
i
=
∑
j
=
0
i
-
1
Z
j
W
j
F
i
-
j
F
i
-
j
+
1
,where(23)
Z
j
=
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
j
A
1
B
1
C
1
,
W
j
=
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
j
A
2
B
2
C
2
.Equation (22) can be rewritten as(24)
Y
=
H
F
,where the time duration is from t
1 to t
n
+
1 and (25)
Y
=
y
t
1
y
t
2
⋮
y
t
n
,
H
=
Z
0
W
0
0
0
⋯
0
Z
1
Z
0
+
W
1
W
0
0
⋯
0
Z
2
Z
1
+
W
2
Z
0
+
W
1
W
0
⋯
0
⋮
⋮
⋮
⋮
⋱
0
Z
n
-
1
Z
n
-
2
+
W
n
-
1
Z
n
-
3
+
W
n
-
2
Z
n
-
4
+
W
n
-
5
⋯
W
0
,
F
=
F
t
1
F
t
2
⋮
F
t
n
+
1
.Considering there may exist conditioning issues when using the Moore-Penrose inverse, one adopts the Tikhonov regularization method [20] to calculate F. Systematic error is written as(26)
e
=
Y
-
H
F
.To find out the minimum ofe, a penalty function J is introduced; namely,(27)
J
=
e
H
e
+
λ
F
H
F
.As the first-order derivative ofF is zero, e reaches a minimum. Then the solution of F is (28)
F
=
H
H
H
+
λ
I
-
1
H
H
Y
.λ is the regularization parameter. In this work, one applies the discrepancy principle [21] to determine λ.


### 2.2. Method Implementation Sequence

The Taylor-series expansion method for force identification is implemented in the sequence as follows:(i)
Determine parameterε and appropriate time step Δ
t. In principle the range of ε is 1~2. However, as τ is replaced with ε
Δ
t, then x
¨
t
+
τ, x
˙
t
+
τ, and x
t
+
τ are expressed as x
¨
t
+
Δ
t, x
˙
t
+
Δ
t, and x
t
+
Δ
t. This is to ensure the subsequent deduction carried out at the cost of a certain calculation error. In order to make the calculation error as small as possible, ε must be close to 1. But if ε
=
1, the matrix χ
11
χ
12
χ
21
χ
22 in (12) will be singular. This situation should be avoided. Weighing the above two aspects, one sets ε
=
1.2. Δ
t depends on the frequency of the excitation. In this work the excitation frequency is 20 Hz, so Δ
t
=
0.05
s.
(ii)
Compute matricesA
1, A
2, A
d, A
v, A
a, B
1, B
2, B
d, B
v, B
a, C
1, C
2, C
d, C
v, and C
a in (16), (17), and (18).
(iii)
Compute matricesZ
k and W
k.
(iv)
Assemble matrixH in (24).
(v)
Calculate the vectorF by (28). In addition, the relative error between identified forces and real forces can be computed as(29)
error
=
norm
2
⁡
F
id
-
F
real
norm
2
⁡
F
real
×
100
%
.


### 2.3. Selection of Expansion Order

In the process of method deduction in Section2.1, the acceleration vector is expanded to a second-order Taylor-series approximation. In order to verify the accuracy of second-order approximation, W
i
l
s
o
n
-
θ method [22] is compared with the proposed method. W
i
l
s
o
n
-
θ method is a classical stepwise integral method which expresses the acceleration, velocity, and displacement vectors as follows:(30)
x
¨
t
+
τ
=
x
¨
t
+
x
¨
t
+
θ
Δ
t
-
x
¨
t
θ
Δ
t
τ
,
(31)
x
˙
t
+
τ
=
x
˙
t
+
x
¨
t
τ
+
x
¨
t
+
θ
Δ
t
-
x
¨
t
θ
Δ
t
τ
2
2
,
(32)
x
t
+
τ
=
x
t
+
x
˙
t
τ
+
x
¨
t
τ
2
2
+
x
¨
t
+
θ
Δ
t
-
x
¨
t
θ
Δ
t
τ
3
6
.From (30), (31), and (32) it can be seen that the truncation errors of acceleration, velocity, and displacement are o
τ
2, o
τ
3, and o
τ
4, respectively. As for the new method, the truncation errors are o
τ
3, o
τ
4, and o
τ
5 in (5). The truncation error of the new method is one order higher than that of W
i
l
s
o
n
-
θ method, so two-order Taylor-series expansion method can achieve higher accuracy than W
i
l
s
o
n
-
θ method in the case of choosing the same time step.In addition, one makes a comparison between second-order approximation and first-order approximation. To implement a first-order approximation, letn
=
1 in (3). After a series of similar deductions, a formula as (24) is founded.A two-storey shear structure in Figure1 is set up. A sinusoidal force F
t is subjected to the structure, F
t
=
2
×
1
0
6
sin
⁡
10
π
t, m
1
=
m
2
=
1
0
5
k
g, k
1
=
3
×
1
0
4
k
N
/
m, and k
2
=
2
×
1
0
4
k
N
/
m. The displacement responses of two storeys are calculated by modal decomposition method. The dynamic load F
t is identified and the identification results are shown in Figure 2.Figure 1
A two-storey shear structure.Figure 2
Force identification results ofF
t.According to the identification results in Figure2, the method of expressing the acceleration vector as a first-order Taylor-series approximation will lead to prodigious identification error which is up to 158.37%, while in the case of adopting second-order Taylor-series, the relative error is 0.87% and the identification result is highly consistent with real force. So in this example adopting second-order Taylor-series is proper and higher order is unnecessary.According to above analysis, it can be considered that adopting second-order Taylor-series is appropriate.


## 2.1. Method Deduction

The nature of this method is to express the acceleration vector as a Taylor-series approximation in the neighborhood of timet; namely,(3)
x
¨
t
+
τ
=
x
¨
t
+
P
·
G
τ
,where P
=
P
1
P
2
⋯
P
n and P
k
(
k
=
1,2
,
…
,
n
) are vectors to be determined. Also G
τ
=
τ
τ
2
/
2
⋯
τ
n
/
n
!
T. Then the vectors of velocity and displacement responses are written as(4)
x
˙
t
+
τ
=
x
˙
t
+
∫
0
τ
x
¨
t
+
τ
d
τ
=
x
˙
t
+
τ
x
¨
t
+
P
·
∫
0
τ
G
τ
d
τ
,
x
t
+
τ
=
x
t
+
∫
0
τ
x
˙
t
+
τ
d
τ
=
x
t
+
τ
x
˙
t
+
τ
2
2
x
¨
t
+
P
·
∫
0
τ
∫
0
τ
G
τ
d
τ
d
τ
.For the purpose of removing the integral symbol in (4) and carrying out follow-up deduction, one must set n as a specific number. Meanwhile, considering the process of force identification in this paper is completed by second-order Taylor-series expansion method, one lets n
=
2, and then(5)
x
¨
t
+
τ
=
x
¨
t
+
τ
P
1
+
τ
2
2
P
2
,
x
˙
t
+
τ
=
x
˙
t
+
τ
x
¨
t
+
τ
2
2
P
1
+
τ
3
6
P
2
,
x
t
+
τ
=
x
t
+
τ
x
˙
t
+
τ
2
2
x
¨
t
+
τ
3
6
P
1
+
τ
4
24
P
2
.The equation of motion at timet
+
Δ
t is(6)
M
x
¨
t
+
Δ
t
+
C
x
˙
t
+
Δ
t
+
K
x
t
+
Δ
t
=
F
t
+
Δ
t
.Replaceτ with Δ
t and substitute (5) into (6) to get(7)
χ
11
P
1
+
χ
12
P
2
=
Q
1
,where (8)
χ
11
=
Δ
t
M
+
Δ
t
2
2
C
+
Δ
t
3
6
K
,
χ
12
=
Δ
t
2
2
M
+
Δ
t
3
6
C
+
Δ
t
4
24
K
,
Q
1
=
F
t
+
Δ
t
-
M
+
Δ
t
C
+
Δ
t
2
2
K
x
¨
t
-
C
+
Δ
t
K
x
˙
t
-
K
x
t
.Replaceτ with ε
Δ
t and substitute (5) into (6) to get(9)
χ
21
P
1
+
χ
22
P
2
=
Q
2
,where(10)
χ
21
=
ε
Δ
t
M
+
ε
Δ
t
2
2
C
+
ε
Δ
t
3
6
K
,
χ
22
=
ε
Δ
t
2
2
M
+
ε
Δ
t
3
6
C
+
ε
Δ
t
4
24
K
,
Q
2
=
F
t
+
ε
Δ
t
-
M
+
ε
Δ
t
C
+
ε
Δ
t
2
2
K
x
¨
t
-
C
+
ε
Δ
t
K
x
˙
t
-
K
x
t
.F
t
+
ε
Δ
t is obtained by linear interpolation between F
t
+
Δ
t and F
t
+
2
Δ
t; namely,(11)
F
t
+
ε
Δ
t
=
2
-
ε
F
t
+
Δ
t
+
ε
-
1
F
t
+
2
Δ
t
1
<
ε
<
2
.Combining (7) and (9), one gets(12)
χ
11
χ
12
χ
21
χ
22
P
1
P
2
=
Q
1
Q
2
.Let(13)
η
=
χ
11
χ
12
χ
21
χ
22
-
1
=
η
11
η
12
η
21
η
22
;then P
1 and P
2 are written as(14)
P
1
=
η
11
Q
1
+
η
12
Q
2
=
L
0
F
i
+
1
+
L
1
F
i
+
2
+
L
2
x
i
+
L
3
x
˙
i
+
L
4
x
¨
i
,
P
2
=
η
21
Q
1
+
η
22
Q
2
=
S
0
F
i
+
1
+
S
1
F
i
+
2
+
S
2
x
i
+
S
3
x
˙
i
+
S
4
x
¨
i
,where(15)
L
0
=
η
11
+
2
-
ε
η
12
,
L
1
=
ε
-
1
η
12
,
L
2
=
-
η
11
+
η
12
K
,
L
3
=
-
η
11
+
η
12
C
-
Δ
t
η
11
+
ε
η
12
K
,
L
4
=
-
η
11
+
η
12
M
-
Δ
t
η
11
+
ε
η
12
C
-
Δ
t
2
2
η
11
+
ε
2
η
12
K
,
S
0
=
η
21
+
2
-
ε
η
22
,
S
1
=
ε
-
1
η
22
,
S
2
=
-
η
21
+
η
22
K
,
S
3
=
-
η
21
+
η
22
C
-
Δ
t
η
21
+
ε
η
22
K
,
S
4
=
-
η
21
+
η
22
M
-
Δ
t
η
21
+
ε
η
22
C
-
Δ
t
2
2
η
21
+
ε
2
η
22
K
.Substituting (14) into (3) and (4) one obtains(16)
x
i
+
1
=
A
1
F
i
+
1
+
A
2
F
i
+
2
+
A
d
x
i
+
A
v
x
˙
i
+
A
a
x
¨
i
,
(17)
x
˙
i
+
1
=
B
1
F
i
+
1
+
B
2
F
i
+
2
+
B
d
x
i
+
B
v
x
˙
i
+
B
a
x
¨
i
,
(18)
x
¨
i
+
1
=
C
1
F
i
+
1
+
C
2
F
i
+
2
+
C
d
x
i
+
C
v
x
˙
i
+
C
a
x
¨
i
,where(19)
A
1
=
Δ
t
3
6
L
0
+
Δ
t
4
24
S
0
,
A
2
=
Δ
t
3
6
L
1
+
Δ
t
4
24
S
1
,
A
d
=
I
+
Δ
t
3
6
L
2
+
Δ
t
4
24
S
2
,
A
v
=
Δ
t
I
+
Δ
t
3
6
L
3
+
Δ
t
4
24
S
3
,
A
a
=
Δ
t
2
2
I
+
Δ
t
3
6
L
4
+
Δ
t
4
24
S
4
,
B
1
=
Δ
t
2
2
L
0
+
Δ
t
3
6
S
0
,
B
2
=
Δ
t
2
2
L
1
+
Δ
t
3
6
S
1
,
B
d
=
Δ
t
2
2
L
2
+
Δ
t
3
6
S
2
,
B
v
=
I
+
Δ
t
2
2
L
3
+
Δ
t
3
6
S
3
,
B
a
=
Δ
t
I
+
Δ
t
2
2
L
4
+
Δ
t
3
6
S
4
,
C
1
=
Δ
t
L
0
+
Δ
t
2
2
S
0
,
C
2
=
Δ
t
L
1
+
Δ
t
2
2
S
1
,
C
d
=
Δ
t
L
2
+
Δ
t
2
2
S
2
,
C
v
=
Δ
t
L
3
+
Δ
t
2
2
S
3
,
C
a
=
I
+
Δ
t
L
4
+
Δ
t
2
2
S
4
.Combine (16), (17), and (18) to get(20)
x
i
+
1
x
˙
i
+
1
x
¨
i
+
1
=
A
1
A
2
B
1
B
2
C
1
C
2
F
i
+
1
F
i
+
2
+
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
x
i
x
˙
i
x
¨
i
.Then the responses at timet
i are expressed as(21)
x
i
x
˙
i
x
¨
i
=
∑
j
=
0
i
-
1
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
j
A
1
A
2
B
1
B
2
C
1
C
2
F
i
-
j
F
i
-
j
+
1
+
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
i
x
0
x
˙
0
x
¨
0
.Zero initial responses are assumed and (21) is simplified to(22)
y
i
=
x
i
x
˙
i
x
¨
i
=
∑
j
=
0
i
-
1
Z
j
W
j
F
i
-
j
F
i
-
j
+
1
,where(23)
Z
j
=
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
j
A
1
B
1
C
1
,
W
j
=
A
d
A
v
A
a
B
d
B
v
B
a
C
d
C
v
C
a
j
A
2
B
2
C
2
.Equation (22) can be rewritten as(24)
Y
=
H
F
,where the time duration is from t
1 to t
n
+
1 and (25)
Y
=
y
t
1
y
t
2
⋮
y
t
n
,
H
=
Z
0
W
0
0
0
⋯
0
Z
1
Z
0
+
W
1
W
0
0
⋯
0
Z
2
Z
1
+
W
2
Z
0
+
W
1
W
0
⋯
0
⋮
⋮
⋮
⋮
⋱
0
Z
n
-
1
Z
n
-
2
+
W
n
-
1
Z
n
-
3
+
W
n
-
2
Z
n
-
4
+
W
n
-
5
⋯
W
0
,
F
=
F
t
1
F
t
2
⋮
F
t
n
+
1
.Considering there may exist conditioning issues when using the Moore-Penrose inverse, one adopts the Tikhonov regularization method [20] to calculate F. Systematic error is written as(26)
e
=
Y
-
H
F
.To find out the minimum ofe, a penalty function J is introduced; namely,(27)
J
=
e
H
e
+
λ
F
H
F
.As the first-order derivative ofF is zero, e reaches a minimum. Then the solution of F is (28)
F
=
H
H
H
+
λ
I
-
1
H
H
Y
.λ is the regularization parameter. In this work, one applies the discrepancy principle [21] to determine λ.


## 2.2. Method Implementation Sequence

The Taylor-series expansion method for force identification is implemented in the sequence as follows:(i)
Determine parameterε and appropriate time step Δ
t. In principle the range of ε is 1~2. However, as τ is replaced with ε
Δ
t, then x
¨
t
+
τ, x
˙
t
+
τ, and x
t
+
τ are expressed as x
¨
t
+
Δ
t, x
˙
t
+
Δ
t, and x
t
+
Δ
t. This is to ensure the subsequent deduction carried out at the cost of a certain calculation error. In order to make the calculation error as small as possible, ε must be close to 1. But if ε
=
1, the matrix χ
11
χ
12
χ
21
χ
22 in (12) will be singular. This situation should be avoided. Weighing the above two aspects, one sets ε
=
1.2. Δ
t depends on the frequency of the excitation. In this work the excitation frequency is 20 Hz, so Δ
t
=
0.05
s.
(ii)
Compute matricesA
1, A
2, A
d, A
v, A
a, B
1, B
2, B
d, B
v, B
a, C
1, C
2, C
d, C
v, and C
a in (16), (17), and (18).
(iii)
Compute matricesZ
k and W
k.
(iv)
Assemble matrixH in (24).
(v)
Calculate the vectorF by (28). In addition, the relative error between identified forces and real forces can be computed as(29)
error
=
norm
2
⁡
F
id
-
F
real
norm
2
⁡
F
real
×
100
%
.


## 2.3. Selection of Expansion Order

In the process of method deduction in Section2.1, the acceleration vector is expanded to a second-order Taylor-series approximation. In order to verify the accuracy of second-order approximation, W
i
l
s
o
n
-
θ method [22] is compared with the proposed method. W
i
l
s
o
n
-
θ method is a classical stepwise integral method which expresses the acceleration, velocity, and displacement vectors as follows:(30)
x
¨
t
+
τ
=
x
¨
t
+
x
¨
t
+
θ
Δ
t
-
x
¨
t
θ
Δ
t
τ
,
(31)
x
˙
t
+
τ
=
x
˙
t
+
x
¨
t
τ
+
x
¨
t
+
θ
Δ
t
-
x
¨
t
θ
Δ
t
τ
2
2
,
(32)
x
t
+
τ
=
x
t
+
x
˙
t
τ
+
x
¨
t
τ
2
2
+
x
¨
t
+
θ
Δ
t
-
x
¨
t
θ
Δ
t
τ
3
6
.From (30), (31), and (32) it can be seen that the truncation errors of acceleration, velocity, and displacement are o
τ
2, o
τ
3, and o
τ
4, respectively. As for the new method, the truncation errors are o
τ
3, o
τ
4, and o
τ
5 in (5). The truncation error of the new method is one order higher than that of W
i
l
s
o
n
-
θ method, so two-order Taylor-series expansion method can achieve higher accuracy than W
i
l
s
o
n
-
θ method in the case of choosing the same time step.In addition, one makes a comparison between second-order approximation and first-order approximation. To implement a first-order approximation, letn
=
1 in (3). After a series of similar deductions, a formula as (24) is founded.A two-storey shear structure in Figure1 is set up. A sinusoidal force F
t is subjected to the structure, F
t
=
2
×
1
0
6
sin
⁡
10
π
t, m
1
=
m
2
=
1
0
5
k
g, k
1
=
3
×
1
0
4
k
N
/
m, and k
2
=
2
×
1
0
4
k
N
/
m. The displacement responses of two storeys are calculated by modal decomposition method. The dynamic load F
t is identified and the identification results are shown in Figure 2.Figure 1
A two-storey shear structure.Figure 2
Force identification results ofF
t.According to the identification results in Figure2, the method of expressing the acceleration vector as a first-order Taylor-series approximation will lead to prodigious identification error which is up to 158.37%, while in the case of adopting second-order Taylor-series, the relative error is 0.87% and the identification result is highly consistent with real force. So in this example adopting second-order Taylor-series is proper and higher order is unnecessary.According to above analysis, it can be considered that adopting second-order Taylor-series is appropriate.


## 3. Numerical Examples

A linear elastic cantilever beam in Figure3 is built to perform the proposed method. The beam properties are shown in Table 1. Sinusoidal excitation and white noise excitation are applied on the beam, respectively; the system characteristics and responses are obtained by the way of finite element simulation.Table 1
Cantilever beam properties.
Young’s modulus,E
70 GPa
Density,ρ
2800k
g
/
m
3
Poisson’s ratio,μ
0.33
Length,l
0.6 m
Width,b
0.06 m
Height,h
0.03 mFigure 3
Cantilever beam model.


### 3.1. Sinusoidal Excitation without Noise

Two different sinusoidal excitations are applied on the beam together as shown in Figure3. The two external forces are (33)
F
5
t
=
600
sin
⁡
3
π
t
+
400
sin
⁡
1.5
π
t
,
F
9
t
=
1000
sin
⁡
2
π
t
+
800
sin
⁡
π
t
.The displacement responses of nodes 6, 8, and 10 are measured for the identification of loads. The identification results of F
5 and F
9 are shown in Figures 4 and 5, respectively. The relative identification errors are shown in Table 2.Table 2
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
2.81
4.36
F
9
5.93
7.12Figure 4
Force identification results ofF
5.Figure 5
Force identification results ofF
9.From above results in Figures4 and 5 and Table 2, it can be seen that for both methods the identified forces and real forces are highly consistent, and the relative errors are very small. With regard to the periodic sinusoidal load, it has the characteristics of simple form and no huge fluctuations. Hence the two algorithms are both able to precisely identify the force time histories of sinusoidal excitations.


### 3.2. Sinusoidal Excitation with Noise

The condition of this study is unchanged compared with last study, other than the fact that the displacement responses of nodes 6, 8, and 10 are mixed with a rand noise (SNR
=
100
dB). The identification results are shown in Figures 6 and 7, respectively. The relative identification errors are shown in Table 3.Table 3
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
12.26
23.15
F
9
6.51
8.37Figure 6
Force identification results ofF
5.Figure 7
Force identification results ofF
9.Figures6 and 7 and Table 3 show that the identification results obtained from the proposed method are satisfactory when the responses are mixed with noise, while the identification results of the state space method are not so exactly as that of the new method. In this example it can be seen that the new method has the merit of noise immunity.


### 3.3. White Noise Excitation without Noise

The condition of this study is same as Section3.1 except that the two sinusoidal excitations are replaced with two white noise excitations, which are generated separately, but both have the aptitude of 1000 N. The identification results are shown in Figures 8 and 9, respectively. The relative identification errors are shown in Table 4.Table 4
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
16.75
58.36
F
9
17.79
60.27Figure 8
Force identification results ofF
5.Figure 9
Force identification results ofF
9.Above results indicate that the identified forces calculated by the proposed method and real forces match well with each other, and the relative errors are acceptable, while the identified forces calculated by state space are unsatisfactory. Due to the complex form and huge fluctuation of the white noise load, it is normal that the identification error is larger than that of the periodic sinusoidal load. In this example, the new method has smaller identification error compared to the state space method. This is because the selection of the time stepΔ
t seriously impacts the precision of the state space method, which leads to a poor ability of identifying sophisticated loads. Nevertheless, the selection of Δ
t has little influence on the accuracy of the proposed method, which makes it superior to the state space method.


### 3.4. White Noise Excitation with Noise

The displacement responses in the last study are added into a rand noise (SNR
=
100
dB). The identification results are shown in Figures 10 and 11, respectively. The relative identification errors are shown in Table 5.Table 5
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
17.38
63.54
F
9
18.50
65.77Figure 10
Force identification results ofF
5.Figure 11
Force identification results ofF
9.From Figures10 and 11 and Table 5, the forces identified by the new method are close to the real forces in the case that responses are polluted by noise, and the relative errors are acceptable, while the identified forces obtained by state space method are inaccurate because of state space method’s weak noise resistance.


## 3.1. Sinusoidal Excitation without Noise

Two different sinusoidal excitations are applied on the beam together as shown in Figure3. The two external forces are (33)
F
5
t
=
600
sin
⁡
3
π
t
+
400
sin
⁡
1.5
π
t
,
F
9
t
=
1000
sin
⁡
2
π
t
+
800
sin
⁡
π
t
.The displacement responses of nodes 6, 8, and 10 are measured for the identification of loads. The identification results of F
5 and F
9 are shown in Figures 4 and 5, respectively. The relative identification errors are shown in Table 2.Table 2
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
2.81
4.36
F
9
5.93
7.12Figure 4
Force identification results ofF
5.Figure 5
Force identification results ofF
9.From above results in Figures4 and 5 and Table 2, it can be seen that for both methods the identified forces and real forces are highly consistent, and the relative errors are very small. With regard to the periodic sinusoidal load, it has the characteristics of simple form and no huge fluctuations. Hence the two algorithms are both able to precisely identify the force time histories of sinusoidal excitations.


## 3.2. Sinusoidal Excitation with Noise

The condition of this study is unchanged compared with last study, other than the fact that the displacement responses of nodes 6, 8, and 10 are mixed with a rand noise (SNR
=
100
dB). The identification results are shown in Figures 6 and 7, respectively. The relative identification errors are shown in Table 3.Table 3
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
12.26
23.15
F
9
6.51
8.37Figure 6
Force identification results ofF
5.Figure 7
Force identification results ofF
9.Figures6 and 7 and Table 3 show that the identification results obtained from the proposed method are satisfactory when the responses are mixed with noise, while the identification results of the state space method are not so exactly as that of the new method. In this example it can be seen that the new method has the merit of noise immunity.


## 3.3. White Noise Excitation without Noise

The condition of this study is same as Section3.1 except that the two sinusoidal excitations are replaced with two white noise excitations, which are generated separately, but both have the aptitude of 1000 N. The identification results are shown in Figures 8 and 9, respectively. The relative identification errors are shown in Table 4.Table 4
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
16.75
58.36
F
9
17.79
60.27Figure 8
Force identification results ofF
5.Figure 9
Force identification results ofF
9.Above results indicate that the identified forces calculated by the proposed method and real forces match well with each other, and the relative errors are acceptable, while the identified forces calculated by state space are unsatisfactory. Due to the complex form and huge fluctuation of the white noise load, it is normal that the identification error is larger than that of the periodic sinusoidal load. In this example, the new method has smaller identification error compared to the state space method. This is because the selection of the time stepΔ
t seriously impacts the precision of the state space method, which leads to a poor ability of identifying sophisticated loads. Nevertheless, the selection of Δ
t has little influence on the accuracy of the proposed method, which makes it superior to the state space method.


## 3.4. White Noise Excitation with Noise

The displacement responses in the last study are added into a rand noise (SNR
=
100
dB). The identification results are shown in Figures 10 and 11, respectively. The relative identification errors are shown in Table 5.Table 5
Relative errors of force identification.
Force
Proposed method (%)
State space method (%)
F
5
17.38
63.54
F
9
18.50
65.77Figure 10
Force identification results ofF
5.Figure 11
Force identification results ofF
9.From Figures10 and 11 and Table 5, the forces identified by the new method are close to the real forces in the case that responses are polluted by noise, and the relative errors are acceptable, while the identified forces obtained by state space method are inaccurate because of state space method’s weak noise resistance.


## 4. Conclusions

This paper has proposed a new approach for force identification. The main novelty of this algorithm is the use of Taylor-series expansion formula. Compared with conventional state space method, the new method has the following advantages and breakthroughs:(i)
The proposed method is implicit in nature but explicit in form, so it can eliminate complex iteration and remain unconditional convergent.
(ii)
The full use of Taylor formula’s prominent feature makes this approach achieve a higher accuracy for the identification of sinusoidal excitation and white noise excitation. The more complex the load form is, the more obvious superiority of the new method becomes.
(iii)
The identification results are satisfactory when the noise is taken into the responses, which indicates a lower noise sensibility.



---

*Source: 6461427-2016-01-14.xml*