Next Article in Journal
Optimization of Heavy Reduction Process on Continuous-Casting Bloom
Previous Article in Journal
Fabrication and Characterization of Wire Arc Additively Manufactured AlSi5 Structures
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Dynamic Load Identification of Unspecified Metal Structures by Measuring Their Response

1
State Key Laboratory of Mechanics and Control for Mechanical Structures, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China
2
Institute for Infrastructure and Environment, Heriot-Watt University, Edinburgh EH14 4AS, UK
*
Author to whom correspondence should be addressed.
Metals 2022, 12(11), 1872; https://doi.org/10.3390/met12111872
Submission received: 12 September 2022 / Revised: 28 October 2022 / Accepted: 1 November 2022 / Published: 2 November 2022

Abstract

:
Many engineering structures are made of metal composite materials. External load information is a key issue for the design and condition monitoring of the structures. Due to the limitation of measurement technology and the external environment, it is difficult to directly measure dynamic loads on structures in many circumstances. This paper focuses on evaluating the external load applied on a structure with unknown dynamic properties. We proposed a novel dynamic load identification method that is based on the Bayesian principle coupled with the extended Kalman filter method. Firstly, the modal parameters are identified under ambient excitation using the Bayesian fast Fourier transform method (FFT). The posterior probability density function (PDF) and covariance of the modal parameters are obtained by the Fourier transform of the response data, and then the modal parameters of the structure are obtained based on unconstrained optimization. Next, the extended Kalman filter method in the modal space is used to update the modal parameters and identify the time-domain information of dynamic loads. The accuracy of the proposed theory was evaluated experimentally using a Bernoulli−Euler beam. The results showed that the method is feasible and efficient.

1. Introduction

Identifying the dynamic load acting on a metal structure is essential for a large variety of applications. For example, dynamic loads applied on a structure are responsible for the majority of failures during service [1,2,3,4]. Therefore, to improve the safety and the reliability of a structure, it is not enough to ensure that the static characteristics meet the design requirements. Dynamic loads should also be closely monitored in order to ensure that they do not exceed certain limits.
Current dynamic loads identification methods are mostly based on specified systems with known parameters, e.g., the stiffness and the damping [5,6,7]. These parameters are usually identified by applying artificial excitations. However, in many cases, this identification method may affect the normal operation of a structure (such as bridges and rotors in operation) [8,9]. Such artificial excitation may also cause structural damage. A practical approach to obtaining modal parameters is ambient vibration tests (AVTs) [10,11], which are usually economical and time-saving. Several modal parameter identification methods based on ambient excitation, such as the lbrahim time domain technique (ITD) [12], the autoregressive moving average model (ARMA) [13], and the random subspace method [14], are available in the literature. Among these approaches, the fast Fourier transform modal parameter identification method based on Bayesian probability [15,16] is the most efficient, as it makes full use of the response data and obtains the modal priori parameters. Furthermore, in the presence of noise, the identification results are often more accurate than other methods thanks to the Bayesian probability principle.
In the field of dynamic load identification, among Bayesian methods we have the Bayesian regularization method [17] and the Bayesian particle filter method [18]. Furthermore, the Kalman filter method [19] can also be considered as a special case of Bayesian methods where the Kalman filter is usually combined with the least square method [20]. Zhang et al. [21] used the Bayesian method to reconstruct a structural load with measurement noise and model uncertainty. Yan [22] identified the impact load on a composite plate by the Bayesian regularization method, and they also used the Kalman filter algorithm to identify the impact position. Faure et al. [23] used the Bayesian regularization method to identify the load position by assessing the displacement response. This method does not need the boundary conditions but is highly sensitive to the noise in the response measurement. The Kalman filter is developed in the modern control theory and has been successfully applied in the field of dynamic load identification. When considering a dynamic model that contains a Gaussian model error and the response signal is corrupted by a Gaussian noise, the Kalman filter load identification method [24,25,26,27,28] can accurately reconstruct the time-domain information of the load in both the physical space and the modal space. Several scholars have proposed the extended Kalman filter [20,29,30] and the unscented Kalman filter [27,31] to correct these errors and to identify loads. Lourens [25] proposed the extended Kalman filter to identify the load, which extends the unknown force vector to the state vector to form the augmented state vector, and uses the Kalman filter to identify the augmented state vector so as to identify the load. Naets et al. [26] studied the observability of this method and pointed out that the system is unobservable when only the acceleration response is used to identify the load. In order to solve this problem, they combined the virtual displacement technology with the acceleration response and improved the Kalman filter to identify the load. Yang et al. [27,28,29] composed the state vector and the model parameters into an augmented state vector. They used the extended Kalman filter to identify the augmented state vector and employed complex mathematical equations to identify the load.
The methods mentioned above are all based on the assumption of known system parameters. However, it remains to be studied how to identify the load when only the response data are known for an unspecified system under a random excitation. Moreover, these methods are only developed in the physical space and not in the modal space, which can lead to complicated computation when applied to large-scale structures. Developing a similar approach for the modal space is useful for large-scale engineering structures. In this paper, we firstly introduced the extended Kalman filter method into the modal space to identify the modal parameters and load of unspecified systems using only random excitation response data.
This paper is organized as follows: In Section 1, the introduction of the research background is described. In Section 2, the posterior PDF and covariance of modal parameters were obtained by the Fourier transform of the response data. Then, the modal parameters of the structure were obtained based on unconstrained optimization (Section 2.1). The extended Kalman filter load identification method in the modal space was derived to update the modal parameters with errors and identify the excitation (Section 2.2). In Section 3, a 3-DOF structure was simulated to verify the correctness of the extended Kalman filter. A simply supported beam under random loads was used to demonstrate the effectiveness and feasibility of the proposed methods in Section 2. Finally, the conclusions are presented in Section 5.

2. Dynamic Load Identification Method

The method proposed in this paper includes two stages:
Firstly, a Bayesian method is used to identify the modal parameters of the structure under ambient excitation. Secondly, the extended Kalman filter method is used to update the modal parameters and identify the time-domain information of the dynamic loads. Research manuscripts reporting large datasets that are deposited in a publicly available database should specify where the data have been deposited and provide the relevant accession numbers. If the accession numbers have not yet been obtained at the time of submission, please note that they will be provided during review. They must be provided prior to publication.

2.1. Bayesian FFT Modal Identification

In this section we introduced a multi-modal identification method based on Bayesian principle, which is a forward problem. The parameters optimal result was derived using the probability, while the rationality of the results was measured by the variance. It should be mentioned that the Bayesian method for multi-modal identification has the disadvantages of needing complex calculations and often involving large errors especially when considering high modal orders, i.e., orders exceed three. However, this is not the case in a single-mode identification. In this paper, the single-mode identification method was introduced.
We start by defining the vector { y j } , where { y j R n : j = 1 , N } is the acceleration time history measured at n DOFs of a structure. The number of sampling points per channel is N. The FFT of { y j } can then be defined as:
Y k = F k + i G k = ( 2 Δ t ) / N j = 1 N y j exp { 2 π i [ ( k 1 ) ( j 1 ) / N ] } ( k = 1 , N ) ,
where i 2 = 1 ; F k = Re Y k   and   G k = Im Y k are the real and imaginary parts of the FFT, respectively; and Δ t is the sampling interval. The FFT corresponds to the frequency f k = ( k 1 ) / N Δ t for k = 2 , 3 N d , where N d = int [ N / 2 ] + 1 is the number of frequencies used for the FFT. The scaling factor of the FFT makes the spectral density matrix unilateral invertible, and the frequency is measured in Hz. The frequency resolution after being converted to the frequency domain is Δ f = 1 / N Δ t .
The modal parameters θ to be identified include the modal frequency f , the modal damping ratio ζ , the modal mode shape Φ , the cross-spectral density of the modal excitation S i j , and the spectral density of the prediction error σ 2 . Let θ ^ be the most probable value (MPV). The posterior PDF of θ is proportional to the likelihood function p ( { Z k } | θ ) :
p ( θ | { Z k } ) p ( { Z k } | θ ) = ( 2 π ) ( N q 1 ) / 2 [ k = 2 N q det C k ( θ ) ] 1 / 2 × exp [ ( 1 / 2 ) k = 2 N q Z k T C k ( θ ) Z k ]
where
Z k = [ F k T , G k T ] T R 2 N s : k = 2 , , N d ,
C k is the covariance matrix of Z k [31] and described as:
C k = 1 2 [ Φ ( Re H k ) Φ T Φ ( Im H k ) T Φ T Φ ( Im H k ) Φ T Φ ( Re H k ) Φ T ] + σ 2 2 I 2 n ,
H k is the spectral density matrix of the model response with β i k = f i / f k and written as:
H k ( i , j ) = S i j [ ( β i k 2 1 ) + i ( 2 ζ i β i k ) ] 1 [ ( β j k 2 1 ) + i ( 2 ζ j β j k ) ] 1 ,
f i and ζ i are the natural frequency and the damping ratio of the i th mode, respectively; S i j is the cross spectral density between the i th- and j th-mode excitations. For convenience, Equation (2) is rewritten as:
p ( θ | { Z k } ) exp [ L ( θ ) ] ,
L ( θ ) = 1 2 k = 2 N q [ ln   det C k ( θ ) + Z k T C k ( θ ) 1 Z k ] .
When the sampling data are large enough, the posterior PDF can be represented as a Gauss distribution:
p ( θ | { Z k } ) exp [ 1 2 ( θ θ ^ ) T C ^ 1 ( θ ^ ) ( θ θ ^ ) ] ,
C ^ = H L ( θ ^ ) 1 ,
where H L ( θ ^ ) is the Hessian of L at the most optimal parameter value [27].
If all modes are completely separated and a selected frequency bandwidth contains only a single mode, a single-mode identification method can be used in this bandwidth. The FFT data of the selected bandwidth contain all the modal parameters of that mode. Modal parameter θ consists of the natural frequency f , the damping ratio ζ , the mode shape φ R n , the spectral density of the modal excitation S , and the spectral density of the prediction error σ 2 .
Next, L ( θ ) is given by:
L ( f , ζ , S , σ 2 , φ ) = n N f ln 2 + ( n 1 ) N f ln σ 2 + k ln ( S D k + σ 2 ) + σ 2 ( d Φ T A Φ )
where N f is the total number of the sampling points in the selected bandwidth,
D k = [ ( β k 2 1 ) 2 + ( 2 ζ β k ) 2 ] 1 ,
A = k [ 1 + ( σ 2 / S D k ) ] 1 D k ,
D k = F k F k T + G k G k T ,
d = k ( F k T F k + G k T G k ) = t r a c e ( A 0 ) ,
A 0 = k D k .
The sum in the equation adds N f terms of the bandwidth of the selected frequency. The most probable value (MPV) of the mode vector φ is the eigenvector corresponding to the maximum eigenvalue of A . Then, the unit normalization of the mode vector is carried out as φ / φ , φ 2 = φ T φ . Minimizing with respect to φ leads to L being only dependent on the remaining four parameters:
L ( f , ζ , S , σ 2 ) = n N f ln 2 + ( n 1 ) N f ln σ 2 + k ln ( S D k + σ 2 ) + σ 2 ( d λ ^ )
where λ ^ is the largest eigenvalue of A . The MPV of the four remaining parameters can be obtained using an unconstrained numerical optimization approach. To this end, we define the following factor as the signal-to-noise ratio at the frequency f k :
γ k = S D k / σ 2 .
In the selected bandwidth, if γ k > > 1 , the following relationship can be obtained:
( 1 + γ k 1 ) 1 = [ 1 + ( σ 2 / S D k ) ] 1 ~ 1 ( σ 2 / S D k ) ~ 1 .
Using the zero-order approximation, the matrix A becomes:
A ~ k D k = A 0 ,
which is a constant matrix. Then, the MPV of φ can be obtained from Equation (19), while the MPV of σ 2 and S can also be found as:
σ ^ 2 = ( d λ ^ 0 ) / ( n 1 ) N f ,
S ^ ( f , ς ) = N f 1 k D k ( f , ς ) 1 D ^ k ,
where
λ ^ 0 = φ ^ T A 0 φ ^ ,
D ^ k = φ ^ T D k φ ^ .
Thus, we see that L depends only on the remaining parameters f and ζ :
L ( f , ζ ) k ln D k ( f , ζ ) + N f ln [ N f 1 k D k ( f , ζ ) 1 D ^ k ] + const .
The steps to obtain modal parameters of a single mode can be summarized as:
  • Calculate the FFT of the sequence { F k , G k } and then evaluate the matrices D k using Equation (13) and the matrix A 0 using Equation (19);
  • Calculate the MPV as the eigenvector of A 0 using the largest eigenvalue and then the MPV of σ 2 using Equation (20);
  • Numerically optimize Equation (24) with respect to f and ζ . The initial guess for f can be obtained from the FFT spectrum of the response. The initial guess of ζ may be set to 1%. Then, the MPV of S can be obtained from Equation (21).

2.2. Dynamic Load Identification Method Based on Extended Kalman Filter in the Modal Space

The equation of motion for a dynamic system can be expressed as a set of ordinary differential equations:
M p ¨ ( t ) + C p ˙ ( t ) + K p ( t ) = B f f ( t ) ,
where f ( t ) is the excitation vector with B f R n × r being the matrix of the excitation positions composed of ones and zeros. Ones are used at the degrees of freedom corresponding to the excitation locations and zeros everywhere else. M , C , and K denote the mass, damping, and stiffness matrix, respectively. The state space model of the dynamic system, i.e., the equation of state and the equation of observation, can be written as Equations (26) and (27), respectively:
x ˙ ( t ) = [ p ˙ ( t ) p ¨ ( t ) θ ˙ ] + w ( t ) = [ p ˙ ( t ) M 1 [ B f f ( t ) C θ p ˙ ( t ) K θ p ( t ) ] [ 0 ] ] + w ( t ) = g ( x ( t ) , f ( t ) ) + w ( t ) ,
y ( t ) = H 0 p ¨ ( t ) + v ( t ) = h ( x ( t ) , f ( t ) ) + v ( t ) ,
where g and h are non-linear functions, x ( t ) = [ p ( t ) p ˙ ( t ) θ ] T is the augmented state vector, p ( t ) is the displacement vectors, θ = [ θ 1 θ 2 θ s ] is the modal parameter vectors, y ( t ) is the acceleration vector of the observation point, w ( t ) and v ( t ) are the model noise and the observation noise, respectively, the variances are G k and R k , and H 0 is the observation vectors.
Discretizing the above nonlinear equations in the time domain, we obtain:
x ˙ k = g ( x k , f k ) + w k ,
y k = h ( x k , f k ) + v k = H 0 p ¨ k + v k = H 0 M 1 [ B f f k C θ p ˙ k K θ p k ] + v k = h ¯ ( x k ) + H 0 M 1 B f f k + v k = h ¯ ( x k ) + D k f k + v k
Applying the first-order Taylor expansion to the non-linear functions g and h , we arrive at the following approximate linear functions:
g k ( x k , f k ) g k ( E ( x k ) , f k ) + g ( x , f ) x | x k | k , f k · ( x k E ( x k ) ) + g ( x , f ) f | x k | k , f k · ( f k E ( f k ) ) ,
h k ( x k , f k ) h k ( E ( x k ) , f k ) + h ( x , f ) x | x k | k , f k · ( x k E ( x k ) ) ,
with
g ( x , f ) x | x = x k | k , f = f k = [ g p g p ˙ g θ ] = [ [ 0 ] I [ 0 ] [ 0 ] M 1 K θ M 1 C θ M 1 K θ θ 1 p k | k M 1 C θ θ 1 p ˙ k | k M 1 K θ θ s p k | k M 1 C θ θ s p ˙ k | k [ 0 ] [ 0 ] [ 0 ] [ 0 ] ] = A
g ( x , f ) f | x = x k | k , f = f k = [ [ 0 ] M 1 B f [ 0 ] ] = B ,
h ( x , f ) x | x = x k | k 1 , f k = [ h p h p ˙ h θ ] = [ H 0 M 1 K θ H 0 M 1 C θ H 0 M 1 K θ θ 1 p k | k 1 H 0 M 1 C θ θ 1 p ˙ k | k 1       H 0 M 1 K θ θ s p k | k 1 H 0 M 1 C θ θ s p ˙ k | k 1 ] = H k
where K θ is a stiffness matrix with unknown stiffness coefficients and C θ is a damping matrix with unknown damping coefficients. The steps of the extended Kalman filter dynamic load identification method in the physical space are given as:
(1)
Start by choosing initial values x 0 | 0 and P 0 | 0 x ;
(2)
Set the model error variance matrix G k and the observation noise variance matrix R k ;
(3)
Locally linearize the observation equations H k and D k ;
(4)
Identify the load
K k = P k | k 1 x H k T ( H k P k | k 1 x H k T + R k ) 1
M k = [ D k T R k 1 ( I H k K k ) D k ] 1
f k = M k D k T R k 1 ( I H k K k ) ( y k h ¯ ( x k | k 1 ) ) ; )
(5)
Update the status x k | k = x k | k 1 + K k [ y k h ¯ ( x k | k 1 ) D k f k ] ;
(6)
Update the state covariance P k | k x = { I + K k D k M k D k T R k 1 H k } ( I K k H k ) P k | k 1 x ;
(7)
Update time x k + 1 | k = x k | k + t k t k + 1 g ( x ( t ) , f ( t ) ) d t ;
(8)
Locally linearize the state equation
A k = I + A Δ t , B k = B Δ t ;
(9)
Update the covariance of the state prediction
P k | k 1 x = A k P k 1 | k 1 A k T + G k .
Next, we aim to introduce the method into the modal space. The transformation equation of the physical space to the modal space is written as:
p ( t ) = Φ q ( t ) ,
where Φ is the mass normalized mode matrix and q ( t ) is the generalized coordinate. The dynamic equation is described as:
q ¨ ( t ) + Γ q ˙ ( t ) + Λ q ( t ) = Φ T B f f ( t ) ,
where with Λ being a generalized stiffness matrix, Λ = d i a g ( λ 1   λ 2     λ n ) , where λ j = ω j 2 and ω j is the natural frequency of order j, Γ is a generalized damping matrix, Γ = d i a g ( γ 1   γ 2     γ n ) , where γ j = 2 ζ j ω j is the j-order modal damping with ζ j being the modal damping ratio.
In the modal space, Equation (26) is transformed into the following equation:
x ˙ ( t ) = [ q ˙ ( t ) q ¨ ( t ) θ ˙ ] + w ( t ) = [ q ˙ ( t ) Φ T B f f ( t ) Γ θ q ˙ ( t ) Λ θ q ( t ) [ 0 ] ] + w ( t ) = g [ x ( t ) , f ( t ) ] + w ( t ) .
Splitting the function g into two parts, we obtain:
g [ x ( t ) , f ( t ) ] = [ q ˙ ( t ) Γ θ q ˙ ( t ) Λ θ q ( t ) [ 0 ] ] + [ [ 0 ] Φ T B f [ 0 ] ] f ( t ) = g ¯ [ x ( t ) ] + B f ( t ) .
Local linearization and time domain discretization leads to the following equations:
A = g ( x , f ) x | x = x k | k , f = f k = [ g q g q ˙ g θ ] = [ [ 0 ] I [ 0 ] [ 0 ] Λ θ Γ θ Λ θ θ 1 q k | k Γ θ θ 1 q ˙ k | k + ( Φ θ 1 ) T B f f k Λ θ θ s q k | k Γ θ θ s q ˙ k | k + ( Φ θ s ) T B f f k [ 0 ] [ 0 ] [ 0 ] [ 0 ] ]
A k = I + A Δ t ,
B = g ( x , f ) f | x = x k | k , f = f k = [ [ 0 ] Φ T B f [ 0 ] ] .
The observation equation is then obtained as:
y ( t ) = H 0 Φ q ¨ ( t ) + v ( t ) = H 0 Φ [ Φ T B f Λ θ q ( t ) Γ θ q ˙ ( t ) ] + v ( t ) = h [ x ( t ) , f ( t ) ] + v ( t ) .
Again, splitting the function h into two parts:
h [ x ( t ) , f ( t ) ] = [ H 0 Φ Λ θ H 0 Φ Γ θ 0 ] [ q ( t ) q ˙ ( t ) θ ] + H 0 Φ Φ T B f f ( t ) = h ¯ [ x ( t ) ] + D f ( t )
D = H 0 Φ Φ T B f .
The discrete form of the observation equation in the time domain is expressed as:
y k = h ¯ ( x k | k ) + D k f k .
The linearization of the observation function produces the following equations:
H k = h ( x ) x | x = x k | k 1 , f = f k 1 = [ h q h q ˙ h θ ] = [ H 0 Φ Λ θ H 0 Φ Γ θ H 0 Φ θ 1 Λ θ q k | k 1 H 0 Φ Λ θ θ 1 q k | k 1           H 0 Φ θ 1 Γ θ q ˙ k | k 1 H 0 Φ Γ θ θ 1 q ˙ k | k 1 + 2 H 0 Φ θ 1 Φ T B f f k 1           H 0 Φ θ s Λ θ q k | k 1 H 0 Φ Λ θ θ s q k | k 1 H 0 Φ θ s Γ θ q ˙ k | k 1 H 0 Φ Γ θ θ s q ˙ k | k 1 + 2 H 0 Φ θ s Φ T B f f k 1 ] ,
D k = h [ x k , f k ] f k = H 0 Φ Φ T B f
The error parameters in θ are eigenvalues and modal damping coefficients. Λ / θ s and Γ / θ s which are relatively easy to get while Φ / θ s can be used for the eigenvector sensitivity analysis.
The mode matrix Φ is differentiated by the variable structure parameter p i j :
Φ p i j = [ φ 1 p i j φ 2 p i j φ n p i j ] ,
where φ i is the mode vector of the ith order. If r s , the sensitivity of the eigenvectors to the stiffness coefficient k i j can be given as:
φ r k i j = s = 1 n β s φ s ,
β s = { 1 λ s λ r ( φ i s φ j r + φ i r φ j s ) ( i j ) 1 λ s λ r φ i s φ j r   ( i = j ) .
If r = s ,
φ r k i j = s = 1 n β s φ s
β s = 0
If r s , the sensitivity of eigenvectors to damping coefficient c i j is shown as:
φ r c i j = s = 1 n γ s φ s
γ s = { 1 λ s λ r ( φ i s φ j r + φ i r φ j s )                 ( i j ) 1 λ s λ r φ i s φ j r                                 ( i = j )
If r = s ,
γ s = { φ i r φ j r ( i j ) 1 2 φ i r 2 ( i = j )
where φ j s is the jth element of the mode vector of the order s.
In practical engineering applications, it is not possible to obtain all modes, and only the first m modes are usually used. If the sensitivity matrix is to be calculated by the above method, there will be a severe modal truncation error, which may not be calculated in this algorithm. For this reason, Wang [32] has proposed an incomplete modal method, in which the static modes are used to approximate the contribution of higher modes such that:
φ r p i j = s = 1 m β s ( γ s ) φ s + S R ,
where β s ( γ s ) means either β s or γ s are selected, with β s representing the derivation of the stiffness coefficient and γ s representing the derivation of the damping coefficient. In the next step, we write:
S R = K 1 F r s = 1 m 1 λ s φ j T F r φ j ,
F r = K p i j + λ r p i j C + λ r C p i j .
It should be noted that the matrices H k and A k obtained in the modal space are different from those obtained in the physical space. The load identification values in the previous step are used to calculate the matrices H k and A k in the modal space. Therefore, the order of the load identification step, the state-updating step, and the time-updating step of the load identification algorithm are changed in the modal space. Furthermore, the initial value needs to be set to include the augmented state vector and its covariance. The algorithm can then be summarized as:
(1)
Start by initializing x 0 | 0 , P 0 | 0 x , and f 0 ;
(2)
Set the model error variance matrix G k and the observation noise variance matrix R k ;
(3)
Update time x k + 1 | k = x k | k + t k t k + 1 g ( x , f k ) d τ ;
(4)
Locally linearize the state equation to obtain A k and B k ;
(5)
Update the covariance for the state prediction P k | k 1 x = A k P k 1 | k 1 A k T + G k ;
(6)
Locally linearize the observation equations H k and D k ;
(7)
Identify the load
K k = P k | k 1 x H k T ( H k P k | k 1 x H k T + R k ) 1 ,
M k = [ D k T R k 1 ( I H k K k ) D k ] 1 ,
f k = M k D k T R k 1 ( I H k K k ) ( y k h ¯ ( x k | k 1 ) ) ;
(8)
Update the status x k | k = x k | k 1 + K k [ y k h ¯ ( x k | k 1 ) D k f k ] ;
(9)
Update the state covariance P k | k x = { I + K k D k M k D k T R k 1 H k } ( I K k H k ) P k | k 1 x .

3. Numerical Study

A numerical example is given in this section to study the reliability of the proposed algorithm in the modal space. Using the proposed method, the error parameters of the model are corrected, and the load is identified. The considered dynamic system is shown in Figure 1. The system parameters are given in Table 1 and Table 2.
Computational acceleration response data were collected for 10 s, and the sampling interval was 0.01 s. The results showed that the first three natural frequencies of the multi-degree-of-freedom system were f 1 = 0.712   H z , f 2 = 1.513   H z , and f 3 = 1.823   H z . The corresponding generalized stiffness matrix was Λ = d i a g ( ω 1 2   ω 2 2   ω 3 2 ) = d i a g ( 20   90.4   131.3 ) . The normalized mode shapes of the mass were Φ = [ φ 1   φ 2   φ 3 ] , φ 1 = [ 0.309   0.413   0.206 ] T , φ 2 = [ 0.349   0.056   0.498 ] T , and φ 3 = [ 0.531   0.277   0.207 ] T . The excitation position matrix was Β f = [ 0   1   0 ] T . We considered 25%, 11%, and 5% errors to be applied to the eigenvalues of the generalized stiffness matrix. Thus, the generalized stiffness matrix was given as Λ = d i a g ( 25   80   128 ) . The acceleration responses of the three degrees of freedom were measured. To further test the robustness of the algorithm, the acceleration responses were contaminated with a zero-mean white Gaussian noise with the noise-to-signal ratio being 5%. The parameters to be identified are as follows:
(1)
Augmented state vector x = [ q 1   q 2   q 3   q ˙ 1   q ˙ 2   q ˙ 3   λ 1   λ 2   λ 3 ] T ;
(2)
Unknown excitation signal.
The initial value of the augmented state vector was set to x 0 | 0 = [ 0 0 0 0 0 0 25 80 128 ] T . The initial state prediction covariance was set to P x 0 | 0 = d i a g [ 1   1   1   1   1   1   10 8   10 8   10 8 ] T . The initial value of the excitation was zero. The measured noise covariance and the model noise covariance were set to R k = 10 2 Ι 2 , G k = 10 1 Ι 9 , respectively, where Ι m is the unit matrix. The identification results of the modal parameters and excitation are shown in Figure 2, Figure 3, Figure 4, Figure 5, Figure 6 and Figure 7.
The results in Figure 2, Figure 3 and Figure 4 showed that the modal parameters converged to stable values after several iterations and the recognition accuracy was high. The identification results of the modal parameters and their errors are shown in Table 3, where the identified eigenvalues were converted into natural frequencies and compared to the exact values.
Figure 5, Figure 6 and Figure 7 show the load identification results. It can be seen that the error of the load identification was relatively large in Figure 6. However, the accuracy of the results improved significantly, as the time progresses after the model parameters converged as can be seen in Figure 7. Clearly, the accuracy of the model is dependent on the modal parameters, and once they converge, the load can be identified with good accuracy.
Moreover, to evaluate the accuracy of the identified results, the mean square error (MSE) and the peak relative error method (PREM) were used. The equations to evaluate these errors are shown as follows:
M S E ( X , Y ) = i = 1 N [ Y ( i ) X ( i ) ] 2 / N ,
P R E M ( X , Y ) = | max Y ( i ) max X ( i ) max X ( i ) | × 100 % .
Load identification errors for different periods are shown in Table 4. The identified load matched the applied load with good accuracy, despite the noise in the measurements and the initially large errors in the identified modal parameters. Furthermore, the converged modal parameters remained to have some errors and the extended Kalman filter algorithm itself involves a truncation error resulting from considering a finite number of terms in the Taylor expansion. However, the proposed algorithm shows good stability and reliability for modal parameter identification as well as load identification.

4. Experimental Study

In this section, we aimed to evaluate the accuracy of the proposed algorithm using experimental validation. Many engineering metal structures can be abstracted into beams and plates. To this end, a simply supported beam made of steel was used. The acceleration responses of the beam under random loads were collected, and the modal parameters were obtained using the Bayesian modal parameter identification method. Furthermore, the extended Kalman filter in the modal space was used to correct the error parameters and to identify random loads.
The known parameters included the acceleration response of the simply supported beam under a random excitation (generated by a modal shaker) and the mass matrix of the beam. The geometric dimensions and the material properties are shown in Table 5.
The simply supported beam was meshed into 10 elements, each of which had a length of 0.07 m and a total of 11 nodes. The beam was supported in the vertical direction at the nodes of both ends. The acceleration sensors were placed at the other nine sections. The test setup is shown in Figure 8.
A random excitation with a frequency bandwidth from 0 to 1500 Hz was then applied on the beam. The sampling frequency of the sensor was 4096 Hz. The Fourier transformation of the response signal is shown in Figure 9. Figure 10 shows the Fourier spectra of the response signals.
It can be clearly seen that the spectrum had six peaks and the adjacent peaks were far apart. Hence, a single mode analysis can be performed. Furthermore, the seventh mode frequency was 1821 Hz, which was not considered in the simulation. We then considered the frequency corresponding to each peak as an initial value of the modal frequency in the unconstrained numerical optimization equation, which is also called a priori parameter in the Bayesian theory. The initial assumptions and the frequency bands are shown in Table 6. The identified modes are shown in Figure 11.
Four frequency bands were selected for each order modal analysis. The final modal frequency and damping ratio were the average values of the four calculations. Table 7 shows the identification values of the modal frequencies and the damping ratios for each order of the modal analysis compared with the results of the modal tests, as well as the MAC values of each mode.
The results shown in Table 7 indicated that the Bayesian modal parameter identification method has relatively large errors in the identification of high-order modal frequencies. Hence, the first two modal frequencies can be considered to be accurate. The third to sixth order modal frequencies were used as errors parameters to form an extended vector with the state vector. The extended Kalman filter load identification method in the modal space was then used to correct these errors parameters and to identify the time-domain load.
The modal error parameter vector was θ = [ λ 3 λ 4 λ 5 λ 6 ] . The extended state vector was then x = [ q 1 q 2 q 6 q ˙ 1 q ˙ 2   q ˙ 6 λ 3 λ 4 λ 5 λ 6 ] . The unknown quantities to be identified are as follows:
  • The extended vector;
  • The random excitation.
The initial values for the state variables were zero, and the initial assumptions for λ 3 , λ 4 , λ 5 , λ 6 were the identified values: λ i = ( 2 π f i ) 2 . The initial error covariance matrix of the extended state vector was taken to be P x 0 | 0 = d i a g [ 1   1   1   1   1   1   1   1   1   1   1   1   10 8   10 8   10 10 ] . The initial value of the excitation was zero. The covariance matrices of both the measurement noise and the model noise were chosen to be R k = 10 2 I 6 and G k = I 18 , respectively. The identified parameters are presented in Figure 12 and Figure 13. The identified eigenvalues and the natural frequencies and the error comparisons are shown in the Table 8. In the first few iterations of filtering, the identified value had a very large jump, but it converged to the exact value after several iterations.
Figure 14 shows the identified load for the first two seconds. Since the identification result was extremely abnormal in several iterative steps of the initial calculation, the statistics started from 0.02 s. Figure 15 shows the load from 0.02 to 0.1 s when the parameters were still being identified. Figure 16 shows the excitation identification results from 0.6 to 0.65 s, during which the modal parameters converged.
The MSE and the PREM were evaluated for the time span from 0.02 to 0.1 s, where the identification algorithm is still in the process of the model parameters updating. The errors measurements were then evaluated again for the time span from 0.6 to 0.65 s, when the load identification started to use converged parameters. The errors for both time spans are displayed in Table 9. The errors showed a significant improvement in the error once the algorithm has converged. The small errors for the time span from 0.6 to 0.65 s showed a good accuracy of the proposed algorithm in capturing the dynamic load.

5. Conclusions

In this paper, we proposed an algorithm for identifying the modal parameters and the dynamic load of an unspecified metal system using the system measured response. We used the Bayesian FFT to identify the modal parameters under the ambient excitation. Hence, the prior parameters of the modal frequency were obtained from the FFT spectrum of the measured response. Then, the modal parameters were refined using numerical unconstrained optimization. It should be noted that the resonance frequency band needs to be selected manually. The dynamic load identification method based on the extended Kalman filter was firstly derived in the modal space.
The results obtained from the numerical test of a 3-DOFs system and then an experiment of a random load applied on a simply supported steel beam demonstrated the advantages of the proposed method in identifying the applied dynamic loads. The main conclusions are as follows:
(1)
The modal parameters identified by the single-mode identification method had good accuracy. They were the initial values of the subsequent algorithm and can be further improved by iteration.
(2)
At the initial stage of load identification, the modal parameters were still in the progress of revision and did not converged. As a result, the error of the load was large. After many iterations, the error modal parameters converged to a stable value, and the recognition accuracy of the load was high.
(3)
The algorithm is suitable for large-scale structures with discrete model errors. It can reduce the amount of calculation while ensuring the calculation accuracy in the modal space compared with the methods in the time domain. It is simple, reliable and easy to be coded.

Author Contributions

Conceptualization and methodology, J.J.; investigation and formal analysis, N.S.; writing—original draft preparation, J.J. and N.S.; writing—review and editing, N.S. and M.S.M.; validation and supervision, F.Z. All authors have read and agreed to the published version of the manuscript.

Funding

Our research is supported by Supported by the Foundation of National Key Laboratory of Science and Technology on Rotorcraft Aeromechanics (No. 61422202105), Qing Lan Project and National Natural Science Foundation of China (No. 51775270). Financial support provided by the project of Qatar National Research Fund under the contract NPRP11S-1220- 170112 is gratefully acknowledged.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Bai, X.; Wei, X.; Ma, Q.; An, Z. Failure Rate Model of Materials under Uncertain Constant Amplitude Cyclic Load. Metals 2022, 12, 1181. [Google Scholar] [CrossRef]
  2. Wang, L.; Peng, Y.; Xie, Y.; Chen, B.; Du, Y. A new iteration regularization method for dynamic load identification of stochastic structures. Mech. Syst. Signal Process. 2021, 156, 107586. [Google Scholar] [CrossRef]
  3. Liu, R.; Dobriban, E.; Hou, Z.; Qian, K. Dynamic Load Identification for Mechanical Systems: A Review. Arch. Comput. Methods Eng. 2022, 29, 831–863. [Google Scholar] [CrossRef]
  4. Li, H.; Jiang, J.; Mohamed, M. Online Dynamic Load Identification Based on Extended Kalman Filter for Structures with Varying Parameters. Symmetry 2021, 13, 1372. [Google Scholar] [CrossRef]
  5. Jiang, J.; Ding, M.; Li, J. A novel time-domain dynamic load identification numerical algorithm for continuous systems. Mech. Syst. Signal Process. 2021, 160, 107881. [Google Scholar] [CrossRef]
  6. Yang, H.; Jiang, J.; Chen, G.; Mohamed, M.S.; Lu, F. A Recurrent Neural Network-Based Method for Dynamic Load Identification of Beam Structures. Materials 2021, 14, 7846. [Google Scholar] [CrossRef]
  7. Yang, H.; Jiang, J.; Chen, G.; Zhao, J. Dynamic load identification based on deep convolution neural network. Mech. Syst. Signal Process. 2023, 185, 109757. [Google Scholar] [CrossRef]
  8. Ronasi, H.; Johansson, H.; Larsson, F. A numerical framework for load identification and regularization with application to rolling disc problem. Comput. Struct. 2011, 89, 38–47. [Google Scholar] [CrossRef]
  9. Pan, C.D.; Yu, L.; Liu, H.L. Identification of moving vehicle forces on bridge structures via moving average Tikhonov regularization. Smart Mater. Struct. 2017, 26, 085041. [Google Scholar] [CrossRef]
  10. Huang, F.-L.; Wang, X.-M.; Chen, Z.-Q.; He, X.-H.; Ni, Y.-Q. A new approach to identification of structural damping ratios. J. Sound Vib. 2007, 303, 144–153. [Google Scholar] [CrossRef]
  11. Mirtaheri, M.; Salehi, F. Ambient vibration testing of existing buildings: Experimental, numerical and code provisions. Adv. Mech. Eng. 2018, 10, 1687814018772718. [Google Scholar] [CrossRef] [Green Version]
  12. Ibrahim, S.R. An upper hessenberg sparse matrix algorithm for modal identification on minicomputers. J. Sound Vib. 1987, 113, 47–57. [Google Scholar] [CrossRef]
  13. Lardies, J.; Larbi, N. A new method for model order selection and modal parameter estimation in time domain. J. Sound Vib. 2001, 245, 187–203. [Google Scholar] [CrossRef]
  14. Liu, K. Modal parameter estimation using the state space method. J. Sound Vib. 1996, 197, 387–402. [Google Scholar] [CrossRef]
  15. Yuen, K.-V.; Katafygiotis, L.S.; Beck, J.L. Spectral density estimation of stochastic vector processes. Probabilistic Eng. Mech. 2002, 17, 265–272. [Google Scholar] [CrossRef]
  16. Yuen, K.-V.; Katafygiotis, L.S. Bayesian Fast Fourier Transform Approach for Modal Updating Using Ambient Data. Adv. Struct. Eng. 2003, 6, 81–95. [Google Scholar] [CrossRef]
  17. Feng, W.; Li, Q.; Lu, Q.; Li, C.; Wang, B. Element-wise Bayesian regularization for fast and adaptive force reconstruction. J. Sound Vib. 2020, 490, 115713. [Google Scholar] [CrossRef]
  18. De Freitas, J.F.G.; Niranjan, M.; Gee, A.H.; Doucet, A. Sequential Monte Carlo Methods to Train Neural Network Models. Neural Comput. 2014, 12, 955–993. [Google Scholar] [CrossRef]
  19. Azam, S.E.; Chatzi, E.; Papadimitriou, C. A dual Kalman filter approach for state estimation via output-only acceleration measurements. Mech. Syst. Signal Process. 2015, 60–61, 866–886. [Google Scholar] [CrossRef]
  20. Lei, Y.; Jiang, Y.; Xu, Z. Structural damage detection with limited input and output measurement signals. Mech. Syst. Signal Process. 2012, 28, 229–243. [Google Scholar] [CrossRef]
  21. Zhang, E.; Antoni, J.; Feissel, P. Bayesian force reconstruction with an uncertain modal. J. Sound Vib. 2012, 331, 798–814. [Google Scholar] [CrossRef]
  22. Yan, G. A Bayesian approach for impact load identification of stiffened composite panel. Inverse Probl. Sci. Eng. 2014, 22, 940–965. [Google Scholar] [CrossRef]
  23. Faure, C.; Ablitzer, F.; Antoni, J.; Pézerat, C. Empirical and fully Bayesian approaches for the identification of vibration sources from transverse displacement measurements. Mech. Syst. Signal Process. 2017, 94, 180–201. [Google Scholar] [CrossRef]
  24. Hwang, J.S.; Lee, S.G.; Ji-hoon, P.; Eun-Jong, Y. Force identification from structural responses using kalman filter. Inst. Mater. Eng. 2008, 33, 257–266. [Google Scholar]
  25. Lourens, E.; Reynders, E.; De Roeck, G.; Degrande, G.; Lombaert, G. An augmented Kalman filter for force identification in structural dynamics. Mech. Syst. Signal Process. 2012, 27, 446–460. [Google Scholar] [CrossRef]
  26. Naets, F.; Cuadrado, J.; Desmet, W. Stable force identification in structural dynamics using Kalman filtering and dummy-measurements. Mech. Syst. Signal Process. 2015, 50–51, 235–248. [Google Scholar] [CrossRef]
  27. Yang, J.N.; Pan, S.; Lin, S. Least-Squares Estimation with Unknown Excitations for Damage Identification of Structures. J. Eng. Mech. 2007, 133, 12–21. [Google Scholar] [CrossRef]
  28. Yang, J.N.; Lin, S.; Huang, H.W.; Zhou, L. An adaptive extended Kalman filter for structural damage identification. Struct. Control Health Monit. 2006, 13, 849–867. [Google Scholar] [CrossRef]
  29. Yang, J.N.; Pan, S.; Huang, H.-W. An adaptive extended Kalman filter for structural damage identifications II: Unknown inputs. Struct. Control Health Monit. 2007, 14, 497–521. [Google Scholar] [CrossRef]
  30. Débarbouillé, A.; Renaud, F.; Dimitrijevic, Z.; Chojnacki, D.; Rota, L.; Dion, J.-L. Wheel forces estimation with an Augmented and Constrained Extended Kalman Filter applied on a nonlinear multi-body model of a half vehicle. Procedia Struct. Integr. 2022, 38, 342–351. [Google Scholar] [CrossRef]
  31. Ding, Y.; Zhao, B.; Wu, B.; Zhang, X.; Guo, L. Simultaneous identification of structural parameter and external excitation with an improved unscented kalman filter. Adv. Struct. Eng. 2015, 11, 5–18. [Google Scholar] [CrossRef]
  32. Wang, B.P. Improved Approximate Methods for Computing Eigenvector Derivatives in Structural Dynamics. AIAA J. 1991, 29, 1018–1020. [Google Scholar] [CrossRef]
Figure 1. DOF structure.
Figure 1. DOF structure.
Metals 12 01872 g001
Figure 2. λ 1 Identification result diagram.
Figure 2. λ 1 Identification result diagram.
Metals 12 01872 g002
Figure 3. λ 2 Identification result diagram.
Figure 3. λ 2 Identification result diagram.
Metals 12 01872 g003
Figure 4. λ 3 Identification result diagram.
Figure 4. λ 3 Identification result diagram.
Metals 12 01872 g004
Figure 5. Load identification results.
Figure 5. Load identification results.
Metals 12 01872 g005
Figure 6. Load identification results from 0–1 s.
Figure 6. Load identification results from 0–1 s.
Metals 12 01872 g006
Figure 7. Load identification results from 5–6 s.
Figure 7. Load identification results from 5–6 s.
Metals 12 01872 g007
Figure 8. Test setup.
Figure 8. Test setup.
Metals 12 01872 g008
Figure 9. Acceleration response signals of all channels.
Figure 9. Acceleration response signals of all channels.
Metals 12 01872 g009
Figure 10. Fourier spectra of the response signals.
Figure 10. Fourier spectra of the response signals.
Metals 12 01872 g010
Figure 11. Identification of the first six modes.
Figure 11. Identification of the first six modes.
Metals 12 01872 g011
Figure 12. Identified parameters of λ 3 and λ 4 .
Figure 12. Identified parameters of λ 3 and λ 4 .
Metals 12 01872 g012
Figure 13. Identified parameters of λ 5 and λ 6 .
Figure 13. Identified parameters of λ 5 and λ 6 .
Metals 12 01872 g013
Figure 14. Identified excitation from 0 to 2 s.
Figure 14. Identified excitation from 0 to 2 s.
Metals 12 01872 g014
Figure 15. Identified excitation from 0.02 to 0.1 s.
Figure 15. Identified excitation from 0.02 to 0.1 s.
Metals 12 01872 g015
Figure 16. Identified excitation from 0.6 to 0.65 s.
Figure 16. Identified excitation from 0.6 to 0.65 s.
Metals 12 01872 g016
Table 1. Stiffness values of the springs.
Table 1. Stiffness values of the springs.
k1k2k3k4
80   N / m 120   N / m 100   N / m 160   N / m
Table 2. Masses of the loads.
Table 2. Masses of the loads.
m1m2m3
2   k g 4   k g 3   k g
Table 3. Modification results of the modal parameters.
Table 3. Modification results of the modal parameters.
Modal ParameterConvergent ValueAccurate ValueError (%)
Characteristic value λ 1 20.6203
Characteristic value λ 2 90.690.40.2
Characteristic value λ 3 130.6131.26−0.7
Modal frequency f 1 (Hz)0.7220.7121.4
Modal frequency f 2 (Hz)1.5141.5130.06
Modal frequency f 3 (Hz)1.8181.823-0.2
Table 4. Error assessment of load identification results in different periods.
Table 4. Error assessment of load identification results in different periods.
TimeMSEPREM
0–1 s8.793.4
5–6 s0.05618.57
Table 5. Geometric dimensions and material properties of the simply supported beam.
Table 5. Geometric dimensions and material properties of the simply supported beam.
ParametersValue
length a = 0.7   m
Section width b = 0.04   m
Section height h = 0.008   m
Density ρ = 7900   k g / m 3
Modulus of elasticity E = 206   G P a
Poisson ratio μ = 0.3
Table 6. Initial assumptions and initial frequency bands.
Table 6. Initial assumptions and initial frequency bands.
ModeInitial Assumption (Hz)Frequency Band (Hz)
LowerUpper
1393840
2150.1145155
3339.2330350
4602.6590615
5931.8900960
61322.612901360
Table 7. Identification results of the modal parameters.
Table 7. Identification results of the modal parameters.
ModalModal Frequency MPV(Hz)Damping Ratio MPV (%)MAC
Identification ValueTest ValueIdentification ValueTest Value
138.1238.140.881.081
2149.1149.130.420.310.997
3336.4336.360.160.211
4597.6594.730.230.241
5928.8920.231.041.160.995
61314.21300.290.80.480.99
Table 8. Eigenvalue and modal frequency identification results.
Table 8. Eigenvalue and modal frequency identification results.
Modal Parameters λ 3 λ 4 λ 5 λ 6 f 3 ( Hz ) f 4 ( Hz ) f 5 ( Hz ) f 6 ( Hz )
Identification value4.48 × 1061.39 × 1073.29 × 1076.67 × 107336.96595.07914.521299.91
Test value4.46 × 1061.39 × 1073.34 × 1076.67 × 107336.36594.73920.231300.29
Error (%)2.570.14−2.66−0.060.170.057−0.65−0.15
Table 9. Statistical errors of load identification results in different periods.
Table 9. Statistical errors of load identification results in different periods.
Periods (s)MSEPREM
0.02–0.53.6174.4
0.6–0.651.3426.7
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Jiang, J.; Shen, N.; Mohamed, M.S.; Zhang, F. Dynamic Load Identification of Unspecified Metal Structures by Measuring Their Response. Metals 2022, 12, 1872. https://doi.org/10.3390/met12111872

AMA Style

Jiang J, Shen N, Mohamed MS, Zhang F. Dynamic Load Identification of Unspecified Metal Structures by Measuring Their Response. Metals. 2022; 12(11):1872. https://doi.org/10.3390/met12111872

Chicago/Turabian Style

Jiang, Jinhui, Nansun Shen, M. Shadi Mohamed, and Fang Zhang. 2022. "Dynamic Load Identification of Unspecified Metal Structures by Measuring Their Response" Metals 12, no. 11: 1872. https://doi.org/10.3390/met12111872

APA Style

Jiang, J., Shen, N., Mohamed, M. S., & Zhang, F. (2022). Dynamic Load Identification of Unspecified Metal Structures by Measuring Their Response. Metals, 12(11), 1872. https://doi.org/10.3390/met12111872

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop