# Constrained Unscented Kalman Filter for Structural Identification of Bouc–Wen Hysteretic System

Authors: Sung-Yong Kim; Jaemin Kim
Journal: Advances in Civil Engineering
 (2020)

---

## Abstract

The Bouc–Wen hysteresis model has been extensively used to identify earthquake damage on civil structures where hysteretic behavior is expected. One of the challenging problems in the identification of Bouc–Wen parameters is overfitting, which is caused by too many nonphysical model parameters. Limit on the search domain of parameters helps to solve the problem. The unscented Kalman filter (UKF) is widely developed to find Bouc–Wen parameters. In order to improve the accuracy and convergence speed of conventional UKF, this study developed a constrained UKF (CUKF) which can be used for the simultaneous identification of the Bouc–Wen hysteretic model. The proposed CUKF is very helpful when certain parameters can be constrained in a physical sense. We compare the proposed CUKF with conventional UKF to validate its robustness and efficacy in identifying Bouc–Wen parameters.


---

## Body

## 1. Introduction

Over the past few decades, several structural identification techniques have been successfully implemented for high-rise and long-span structures based on the development of computing and measurement technology. In the earthquake engineering field, structural identification techniques aim to quantify damage levels and guide potentially necessary immediate retrofitting efforts. Structural assessment in seismic engineering focuses on real-time hazard assessment using model updating techniques to find optimal parameters for structural models that can describe various complex nonlinear phenomena including yielding, postyielding, and strength and/or stiffness degradation. However, because it is difficult to describe complex inelastic behavior, monitoring of seismically affected structures is considered to be a challenging problem. Numerous models have been investigated for the description of nonlinear behaviors, and the Bouc–Wen model is one of the most popular models. It was originally introduced in [1] and extended in [2]. Owing to its versatility, the Bouc–Wen model has been successfully extended to create several variations, for assessing strength degradation, stiffness degradation [3], pinching effects [4–6], and asymmetric behavior [7–10]. It has also been widely adopted for the description of hysteretic structures [11–14].The Kalman filter (KF) is widely used to estimate system states in practical applications, which consist of two consecutive processes called prediction and estimation. In the prediction process, the KF calculates predicted states and covariances, and in the following estimation process, the KF calculates Kalman gain based on the predicted covariances. The KF is widely used due to its ability in facilitating with only current measured values, previous states, and covariances for each iteration without requiring any additional previous information. The KF, however, cannot be used to analyze structures that are subjected to a severe earthquake because such structures exhibit nonlinear hysteretic behavior while the KF cannot be applied to a nonlinear system. To estimate the states of the nonlinear structures accurately, various numerical techniques have been investigated including the extended Kalman filter (EKF) [15–18], unscented Kalman filter (UKF) [19–26], and particle filter [27]. Among them, the UKF has an advantage in comparison to the EKF and particle filter: the UKF does not require Jacobians because they transform a set of points using known nonlinear equations and combine the results to estimate the mean and covariance of a state. Compared with particle filters that do not need Jacobians, the number of points used in a particle filter generally needs to be much greater than the number of points in a UKF.However, as the demand for real-time applications has increased in many fields, the demand for technologies that provide accurate predictions in less computation has increased accordingly. Although the KF and its variations are powerful tools for state estimation on their own, they have the potential to facilitate faster and more accurate state estimation when using constraints on observed parameters identified by experts. Various studies have attempted to combine filtering techniques with nonlinear minimization algorithms [28–30] to facilitate the incorporation of prior knowledge regarding state parameters. In the health monitoring field, constrained dual EKFs are employed for the state estimation of linear structures. It has been proven numerically that constrained dual EKFs are efficient for state estimation. However, KFs with constraint functions have not yet been applied to nonlinear parameter estimation for the Bouc–Wen model.This study implemented constrained minimization scheme in the use of UKF in order to improve the performance of the conventional UKF for real-time parameter identification. This paper first reviews the original Bouc–Wen model and several variations that have been proposed to simulate strength and strength degradation. Next, we propose a constrained UKF (CUKF) that combines a constrained minimization technique with a conventional UKF to determine the appropriate constraints for state and parameter identification. The proposed method is then compared with conventional UKF-based parameter identification techniques to validate its robustness and efficacy in identifying Bouc–Wen parameters.


## 2. Parameter Identification via CUKF

### 2.1. Problem Formulation

Let us consider a nonlinear dynamic system modeled as(1)xk+1=Fkxk,θ+wk,(2)yk=Hkxk,θ+vk,where xk∈ℝN is the state variable, yk∈ℝM is the measured value, θ∈ℝP is the unknown parameter vector of interest, and wk∈ℝN and vk∈ℝM are additive process and sensor noises, respectively. We denote Fk:ℝN×ℝP↦ℝN as the nonlinear vector-valued transition dynamics of the system and Hk:ℝN×ℝP↦ℝM as a nonlinear vector-valued function, which transforms a state vector into the appropriate measurement vector. The subscript k∈ℕ stands for the k-th discrete time step.In order to estimate the unknown parameters, the unknown parameter vectorθ is modeled as a discretized vector θk, whose elements are the estimated parameters at the k-th step. We model the evolution of the estimated parameters as a Gaussian random walk process. An additional model for the parameter estimation can be considered as follows:(3)θk+1=θk+ξk,where ξk∈ℝP is a zero-mean Gaussian white noise process. Here, equation (5) is reformulated by defining augmented vectors as x˜k=xk;θkT∈ℝN+P, w˜k=wk;ξkT∈ℝN+P, F˜k:ℝN+P↦ℝN+P, and H˜k:ℝN+P↦ℝM. An augmented system is formulated as follows:(4)x˜k+1=F˜kx˜k+w˜k,(5)yk=H˜kx˜k+vk.We dispose the hat symbol for simplicity. We denote the total length of state vector asN in this section.


### 2.2. Stage 1: Conventional Unscented Kalman Filter

UKF is a conventional method that attempts to estimate the state of the next step of the nonlinear system, whose posterior distribution is calculated by deterministic sampling approximation. We suppose that a state variablexk∈ℝN follows a normal distribution with its estimate x^k and covariance Pk denoted by xk∼Nx^k,Pk. Under a conventional framework, a priori estimate and covariance for the next k+1 step are predicted by applying the following procedure.


#### 2.2.1. Step 1: Initialization

At thek-th step, a set of sampling points is declared as(6)χ1=x^k,χi+1=x^k+N+λPki,i=1,…,N,χi+N+1=x^k−N+λPki,i=1,…,N,where λ=σ2N+κ−N is a scaling parameter, σ stands for spread of the sampling points around x^k, κ is a secondary scaling factor, and ⋅i denotes the i-th row of the matrix square root. The corresponding weight values Wi for i=1,…,N are declared as(7)W1=λN+λ,Wi+1=12N+λ,i=1,…,N,Wi+N+1=12N+λ,i=1,…,N.


#### 2.2.2. Step 2: Prediction

The projected mean and covariance are subsequently evaluated by using the projected sigma points as follows:(8)x^k+1−=∑i=12N+1WiFkχi,Pk+1−=∑i=12N+1WiFkχi−x^k+1−Fkχi−x^k+1−T+Qk+1,where Qk+1 is the process noise covariance at the k+1-th step. The projected sigma points are further substituted into the observation equation to generate the predicted system responses as follows:(9)y^k+1−=∑i=12N+1WiHkχi.Pk+1yy=∑i=12N+1WiHkχi−y^k+1−Hkχi−y^k+1−T+Rk+1,where Rk+1 is the output noise covariance at the k+1-th step.


#### 2.2.3. Step 3: Correction

The Kalman gainKk+1 at the k+1-th step is calculated as(10)Pk+1xy=∑i=12N+1WiFkχi−x^k+1−Hkχi−y^k+1−T,Kk+1=Pk+1xyPk+1yy−1.The estimated statex^k+1+ and covariance Pk+1+ at the k+1-th step are finally calculated as(11)x^k+1+=x^k+1−+Kk+1yk+1−y^k+1−,(12)Pk+1+=Pk+1−−Kk+1Pk+1yyKk+1T.A conventional UKF completes the estimation of the mean and covariance for the next step by calculating equations (11) and (12) as x^k+1 and Pk+1, respectively. The estimated parameters, however, may be not in a feasible region such as negative values’ damping ratio or natural frequency. Such estimated parameters in the unacceptable range can result in lower accuracy in the estimation process or divergence problem. Hence, introducing an additional process that forces the corrected estimate x^k+1+ into the feasible range can improve the estimation performance. The next stage adopts constrained minimization in order to project the corrected estimate to improve the estimation accuracy.


### 2.3. Stage 2: Estimate Projection by Constrained Minimization

Having corrected an estimated statex^k+1+, a constrained estimate is obtained by projecting the unconstrained estimate onto a constraint surface. Interior-point methods are iterative schemes where the iterates approximate a local minimum from inside the feasible set. The constrained estimate can be obtained by solving the following minimization problem:(13)x^k+1=argminx∈Ωx−x^k+1+TSx−x^k+1+,where x^k and x˜k are the unconstrained estimate and constrained estimate of the state, respectively. We denote S∈S+N as a positive-definite weight matrix and Ω as the feasible region. A geometric illustration of estimate projection is presented in Figure 1, where the feasible set Ω is shaded and the contour lines of the convex quadratic objective function are depicted. It has been well established in previous studies [31] that setting S to IN×N, where IN×N is an identity matrix of size N×N, results in the maximum probability estimate of the state. Setting S to the covariance matrix Pk calculated in equation (12) results in the least square estimate of the state. The estimate projection described in equation (13) is a constrained nonlinear function minimization that can typically be solved using well-established algorithms, including interior-point methods [32] and the active set method [33].Figure 1
A geometric illustration of estimate projection diagram for constrained minimization of CUKF.The advantage of introducing the constrained minimization technique in the parameter estimation problem is that the predicted parameters are bounded by the reasonable constraints. The constraints would make the UKF estimate more accurate compared with the conventional UKF estimate. Many practical problems require parameter estimation for physically meaningful constraints, which include positive natural frequencies, damping ratios, and energy dissipation. In the next section, the above method is applied to the state and parameter estimation of the Bouc–Wen model, which includes a severe nonlinearity, and the results thereof are discussed.


## 2.1. Problem Formulation

Let us consider a nonlinear dynamic system modeled as(1)xk+1=Fkxk,θ+wk,(2)yk=Hkxk,θ+vk,where xk∈ℝN is the state variable, yk∈ℝM is the measured value, θ∈ℝP is the unknown parameter vector of interest, and wk∈ℝN and vk∈ℝM are additive process and sensor noises, respectively. We denote Fk:ℝN×ℝP↦ℝN as the nonlinear vector-valued transition dynamics of the system and Hk:ℝN×ℝP↦ℝM as a nonlinear vector-valued function, which transforms a state vector into the appropriate measurement vector. The subscript k∈ℕ stands for the k-th discrete time step.In order to estimate the unknown parameters, the unknown parameter vectorθ is modeled as a discretized vector θk, whose elements are the estimated parameters at the k-th step. We model the evolution of the estimated parameters as a Gaussian random walk process. An additional model for the parameter estimation can be considered as follows:(3)θk+1=θk+ξk,where ξk∈ℝP is a zero-mean Gaussian white noise process. Here, equation (5) is reformulated by defining augmented vectors as x˜k=xk;θkT∈ℝN+P, w˜k=wk;ξkT∈ℝN+P, F˜k:ℝN+P↦ℝN+P, and H˜k:ℝN+P↦ℝM. An augmented system is formulated as follows:(4)x˜k+1=F˜kx˜k+w˜k,(5)yk=H˜kx˜k+vk.We dispose the hat symbol for simplicity. We denote the total length of state vector asN in this section.


## 2.2. Stage 1: Conventional Unscented Kalman Filter

UKF is a conventional method that attempts to estimate the state of the next step of the nonlinear system, whose posterior distribution is calculated by deterministic sampling approximation. We suppose that a state variablexk∈ℝN follows a normal distribution with its estimate x^k and covariance Pk denoted by xk∼Nx^k,Pk. Under a conventional framework, a priori estimate and covariance for the next k+1 step are predicted by applying the following procedure.


### 2.2.1. Step 1: Initialization

At thek-th step, a set of sampling points is declared as(6)χ1=x^k,χi+1=x^k+N+λPki,i=1,…,N,χi+N+1=x^k−N+λPki,i=1,…,N,where λ=σ2N+κ−N is a scaling parameter, σ stands for spread of the sampling points around x^k, κ is a secondary scaling factor, and ⋅i denotes the i-th row of the matrix square root. The corresponding weight values Wi for i=1,…,N are declared as(7)W1=λN+λ,Wi+1=12N+λ,i=1,…,N,Wi+N+1=12N+λ,i=1,…,N.


### 2.2.2. Step 2: Prediction

The projected mean and covariance are subsequently evaluated by using the projected sigma points as follows:(8)x^k+1−=∑i=12N+1WiFkχi,Pk+1−=∑i=12N+1WiFkχi−x^k+1−Fkχi−x^k+1−T+Qk+1,where Qk+1 is the process noise covariance at the k+1-th step. The projected sigma points are further substituted into the observation equation to generate the predicted system responses as follows:(9)y^k+1−=∑i=12N+1WiHkχi.Pk+1yy=∑i=12N+1WiHkχi−y^k+1−Hkχi−y^k+1−T+Rk+1,where Rk+1 is the output noise covariance at the k+1-th step.


### 2.2.3. Step 3: Correction

The Kalman gainKk+1 at the k+1-th step is calculated as(10)Pk+1xy=∑i=12N+1WiFkχi−x^k+1−Hkχi−y^k+1−T,Kk+1=Pk+1xyPk+1yy−1.The estimated statex^k+1+ and covariance Pk+1+ at the k+1-th step are finally calculated as(11)x^k+1+=x^k+1−+Kk+1yk+1−y^k+1−,(12)Pk+1+=Pk+1−−Kk+1Pk+1yyKk+1T.A conventional UKF completes the estimation of the mean and covariance for the next step by calculating equations (11) and (12) as x^k+1 and Pk+1, respectively. The estimated parameters, however, may be not in a feasible region such as negative values’ damping ratio or natural frequency. Such estimated parameters in the unacceptable range can result in lower accuracy in the estimation process or divergence problem. Hence, introducing an additional process that forces the corrected estimate x^k+1+ into the feasible range can improve the estimation performance. The next stage adopts constrained minimization in order to project the corrected estimate to improve the estimation accuracy.


## 2.2.1. Step 1: Initialization

At thek-th step, a set of sampling points is declared as(6)χ1=x^k,χi+1=x^k+N+λPki,i=1,…,N,χi+N+1=x^k−N+λPki,i=1,…,N,where λ=σ2N+κ−N is a scaling parameter, σ stands for spread of the sampling points around x^k, κ is a secondary scaling factor, and ⋅i denotes the i-th row of the matrix square root. The corresponding weight values Wi for i=1,…,N are declared as(7)W1=λN+λ,Wi+1=12N+λ,i=1,…,N,Wi+N+1=12N+λ,i=1,…,N.


## 2.2.2. Step 2: Prediction

The projected mean and covariance are subsequently evaluated by using the projected sigma points as follows:(8)x^k+1−=∑i=12N+1WiFkχi,Pk+1−=∑i=12N+1WiFkχi−x^k+1−Fkχi−x^k+1−T+Qk+1,where Qk+1 is the process noise covariance at the k+1-th step. The projected sigma points are further substituted into the observation equation to generate the predicted system responses as follows:(9)y^k+1−=∑i=12N+1WiHkχi.Pk+1yy=∑i=12N+1WiHkχi−y^k+1−Hkχi−y^k+1−T+Rk+1,where Rk+1 is the output noise covariance at the k+1-th step.


## 2.2.3. Step 3: Correction

The Kalman gainKk+1 at the k+1-th step is calculated as(10)Pk+1xy=∑i=12N+1WiFkχi−x^k+1−Hkχi−y^k+1−T,Kk+1=Pk+1xyPk+1yy−1.The estimated statex^k+1+ and covariance Pk+1+ at the k+1-th step are finally calculated as(11)x^k+1+=x^k+1−+Kk+1yk+1−y^k+1−,(12)Pk+1+=Pk+1−−Kk+1Pk+1yyKk+1T.A conventional UKF completes the estimation of the mean and covariance for the next step by calculating equations (11) and (12) as x^k+1 and Pk+1, respectively. The estimated parameters, however, may be not in a feasible region such as negative values’ damping ratio or natural frequency. Such estimated parameters in the unacceptable range can result in lower accuracy in the estimation process or divergence problem. Hence, introducing an additional process that forces the corrected estimate x^k+1+ into the feasible range can improve the estimation performance. The next stage adopts constrained minimization in order to project the corrected estimate to improve the estimation accuracy.


## 2.3. Stage 2: Estimate Projection by Constrained Minimization

Having corrected an estimated statex^k+1+, a constrained estimate is obtained by projecting the unconstrained estimate onto a constraint surface. Interior-point methods are iterative schemes where the iterates approximate a local minimum from inside the feasible set. The constrained estimate can be obtained by solving the following minimization problem:(13)x^k+1=argminx∈Ωx−x^k+1+TSx−x^k+1+,where x^k and x˜k are the unconstrained estimate and constrained estimate of the state, respectively. We denote S∈S+N as a positive-definite weight matrix and Ω as the feasible region. A geometric illustration of estimate projection is presented in Figure 1, where the feasible set Ω is shaded and the contour lines of the convex quadratic objective function are depicted. It has been well established in previous studies [31] that setting S to IN×N, where IN×N is an identity matrix of size N×N, results in the maximum probability estimate of the state. Setting S to the covariance matrix Pk calculated in equation (12) results in the least square estimate of the state. The estimate projection described in equation (13) is a constrained nonlinear function minimization that can typically be solved using well-established algorithms, including interior-point methods [32] and the active set method [33].Figure 1
A geometric illustration of estimate projection diagram for constrained minimization of CUKF.The advantage of introducing the constrained minimization technique in the parameter estimation problem is that the predicted parameters are bounded by the reasonable constraints. The constraints would make the UKF estimate more accurate compared with the conventional UKF estimate. Many practical problems require parameter estimation for physically meaningful constraints, which include positive natural frequencies, damping ratios, and energy dissipation. In the next section, the above method is applied to the state and parameter estimation of the Bouc–Wen model, which includes a severe nonlinearity, and the results thereof are discussed.


## 3. Numerical Simulation

### 3.1. Example 1: Original Bouc–Wen Hysteresis

A single-degree-of-freedom (SDOF) system with original Bouc–Wen hysteresis is considered in this example. It should be note that for the identification of earthquake damage, one specific element where a hysteretic behavior is expected, such as a column or a slender brace, is selected for observation to reduce the computational complexity: the overestimation problem may occur if measurement and identification are performed at many points including the point where elastic behavior is likely to occur. From the practical perspective, this study examined a SDOF system. The equation of motion for the example system when subjected to ground acceleration is written as(14)u¨+2ζnωnu˙+αωn2u+1−αωn2z=−u¨g,where u is the displacement, ζn is the damping ratio, ωn is the natural frequency, ug is the ground displacement, α is the postyielding stiffness ratio, an overdot represents the time derivative, and z is an auxiliary variable for describing the Bouc–Wen hysteresis governed by following differential equation:(15)z˙=u˙1−znγ+βsgnu˙z.To perform parameter identification using the UKF, an augmented state vector is introduced as follows:(16)x=u,v,z,α,γ,β,nT,where v=u˙ is the velocity of the structure. Based on this augmented state vector, equations (14) and (15) can be reformulated into a state-space equation as follows:(17)x˙=fx=v−u¨g−2ζnωnv−αωn2u+1−αωn2zv1−znγ+βsgnvzO4×1,where ON×M∈ℝN×M is a zero matrix. In this example, the noise term for the state variable in equation (2) is neglected. Additionally, the damping ratio and natural frequency are discarded from the augmented vector because it is assumed that they were previously identified under ambient conditions. This point will be discussed in detail in the next example. We choose absolute acceleration as the measurement variable, which is common in real-world application. The output equation for the measurement can be written as follows:(18)y=u¨+u¨g+w=−2ζnωnv−αωn2u−1−αωn2z+u¨g+w,where w is the measurement noise. A Gaussian white noise process with root mean square (RMS) value of 1 m/s2 is considered based on its relatively large value in the context of the general measurement procedure. We choose the ground motion of the 1940 El Centro earthquake as the input excitation. It has a 50 Hz sampling rate and a 30-second duration. To perform system identification using the Kalman filter methods, we reformulate equations (17) and (18) into a discretized form of equation (2) using the fourth-order Runge–Kutta integration method. The initial values of the parameters u and v in the space vector x are set to zero. Table 1 lists the true model parameters and initial guessed model parameters for this example. The postyielding parameter α is set to 0.09, which is a practical value for steel. The shape parameters n, β, and γ are set to 1.2, 20, and 40, respectively. The initial guesses are determined to have the same or a similar order of magnitude compared to the corresponding true values. We configure the initial covariance matrix P in diagonal matrix form to ensure that the designed filter could be applied to a roughly estimated condition as follows:(19)P=diag10−3,10−3,10−3,10−1,103,103,10−1,where diag⋅ is a diagonal matrix with arguments along the main diagonal. To utilize the CUKF, the following constraint conditions are determined based on physical intuition. First, the parameters α, γ, β, and n are constrained to be positive. Second, the upper bound of the stiffness ratio α is set to 0.2. Finally, β+γ is bounded in the range of 35–70, meaning the limit value zm should be bounded in a range of 0.01–0.03 if the parameter n is near unity. The constraint conditions onto which the estimated state will be projected constitute the following feasible region:(20)Ω=x∈ℝ7:Ax<b,where(21)A=000−10000000−10000000−10000000−10000−1−100000110,b=0000−3570.Table 1
Model parameters used in Example 1.
ParameterTrueInitial guessEstimatedAbsolute errorUKFCUKFUKFCUKFζn† (%)2—————ωn† (rad/s)2π—————α (×10−2)91014.8010.345.801.34γ (×10m−1)22.52.772.240.770.24β (×10m−1)43.55.464.181.460.18n (—)1.211.191.150.010.05†These parameters are presumed to be determined.Therefore, the procedure for finding the constrained estimatex^k defined in equation (13) is a linearly constrained quadratic optimization problem. A variety of numerical tools can solve such problems. A MATLAB function called quadprog with an interior-point-convex algorithm is employed to solve the constrained quadratic problem in this example.Figure2 presents a comparison of the accelerations estimated by the UKF and CUKF. It is evident that the acceleration histories estimated by the UKF and CUKF correlate with the true simulated responses such that they cannot be distinguished from the true responses. Their absolute errors reach a maximum of 0.5 m/s2.Figure 2
Comparison of acceleration time histories identified by UKF and CUKF (a) and their absolute errors (b) for Example 1.
(a)(b)Figure3 presents the time histories of the parameters estimated by the UKF and CUKF. Table 1 summarizes the final estimated values for the parameters. By carefully inspecting the time histories for identifying parameter α, one can see that the UKF attempts to seek a value of α larger than 0.3 at approximately 1 s and 4 s, which is an infeasible region from a physical perspective. In contrast, the CUKF estimates parameters closer to the true values compared with the UKF by searching for the parameters within admissible physical constraints. The inaccuracy of parameters combined with good agreement with predicted measurements indicates overfitting for the UKF. The parameters of the Bouc–Wen model are empirical, meaning they are not derived from fundamental mechanics, and a given response may not determine parameters unambiguously [34]. Therefore, it is helpful for applying the CUKF to system parameter identification if certain parameters can be constrained in a physical sense. It can be seen that the proposed CUKF can provide accurate results for estimation of all parameters except for the parameter n, which shows a large error of about 4%. However, when we recall from the previous studies on the parameters of the Bouc model that the variable n only describes the sharpness of hysteretic behavior (for example, refer to [34]), it can be seen that the estimation results accurately predict the stiffness and strength of the hysteretic behavior.Figure 3
Results of parameter identification in Example 1.


### 3.2. Example 2: Overestimation Problem

To test the overestimation problem, an additional simulation is conducted in which the structure examined in the previous example is subjected to ambient ground motion. In the example, the governing equation is identical to that for Example 1, but the structure is assumed to be subjected to the ground motion of the 1940 El Centro earthquake at a scale of 1/10. We assume that the structure is elastic and the Bouc–Wen parameters should have no influence on structural behavior. The augmented vector includes the damping ratio and natural frequency. We define the augmented state vectorx as follows:(22)x=u,v,z,α,γ,β,n,ωn,ζnT.A Gaussian white noise process with an RMS value of 0.1 m/s2 is applied to the measurements. Just as in the previous example, the initial values for the original state vector x are set to zero. Table 2 lists the true and initially guessed model parameters for this example. The true natural frequency and damping ratio are 2π and 0.02, respectively, which are same as the corresponding values in Example 1. The postyielding parameter α and shape parameters n, β, and γ are set to 0.09, 1.2, 20, and 40, respectively.Table 2
Model parameters used in Example 2.
ParameterTrueInitial guessEstimatedAbsolute errorUKFCUKF1CUKF2UKFCUKF1CUKF2ζn (%)23−3.483.682.785.481.680.78ωn (rad/s)2π2.4π−4.206.286.3310.480.000.05α (×10−2)91033.945.25—†24.943.75—†γ (×10m−1)22.5−2.524.64—†4.522.64—†β (×10m−1)43.50.132.19—†3.871.81—†n (—)1.210.961.50—†0.300.24—††These values are undefined because they are constrained by equality constraints.Two CUKFs are utilized in the example. CUKF1 is identical to the CUKF used in Example 1, except that the natural frequency and damping ratio are included in the state vector and constrained to be positive. CUKF2 is a constrained UKF that includes several equality constraints based on physical intuition. In CUKF2, to suppress the Bouc–Wen parameters during the identification process, the parametersα, γ, β, and n are set to be constrained to 0.1, 1, 25, and 35, respectively. The constraint conditions onto which the estimated state is projected constitute the following feasible region:(23)Ω=x∈ℝ9:Ax<bandAeqx=beq,where(24)A=000−1000000000−1000000000−1000000000−1000000000−1000000000−10000−1−1000000011000,b=000000−3570.Aeq=000100000000010000000001000000000100,beq=1025351.Figure4 presents a comparison of the accelerations estimated by the UKF, CUKF1, and CUKF2. The acceleration histories estimated by the UKF and CUKFs correlate well with the true simulated responses. Figure 5 shows the time histories of the parameters estimated by the UKF and CUKFs. The UKF attempts to seek the values of ωn and ζn smaller than zero, which is infeasible from a physical perspective, resulting in overestimation. In contrast, both CUKFs estimate these parameters accurately by searching within a positive space. CUKF1 forcibly predicts the Bouc–Wen parameters, which are unnecessary for describing the elastic behavior of the examined structure. Therefore, as noted in the previous example, it is very helpful for applying the CUKF to system parameter identification if some parameters can be constrained in a physical sense. We summarize the final values estimated by the UKF and CUKFs in Table 2, which indicates that the CUKFs are superior for estimating accurate parameters.Figure 4
Comparison of acceleration time history results for Example 2 identified by the UKF, CUKF1, and CUKF2.Figure 5
Results of parameter identification for Example 2.


### 3.3. Example 3: Structure with Degrading Hysteresis

The final example uses the SDOF structure from Example 1, but degrading effects are included. In the model proposed by [4], the auxiliary variable z in equation (15) is modified in order to employ strength and stiffness degradation as follows:(25)z˙=u˙1−νznγ+βsgnu˙zη,where the parameters ν and η are expressed as linear functions of the dissipated energy and their slopes δν and δη such as(26)ν=1+δνε,η=1+δηε,with ε being the dissipated hysteretic energy as(27)ε=εt=1−αωn2∫0tzu˙dt.An augmented state vector is introduced to formulate the state-space equation as follows:(28)x=u,v,z,ε,α,γ,β,n,δν,δηT.Using this augmented state vector, equations (25), and (26) are reformulated as the following nonlinear state-space equation:(29)x˙=fx=v−u¨g−2ζnωnv−αωn2u−1−αωn2zvη1−νznγ+βsgnvzzvO6×1.Figure6 presents a comparison of the accelerations estimated by the UKF and CUKF. The acceleration histories estimated by the UKF and CUKF correlate well with the true simulated responses, and they cannot be distinguished from the true responses. Their absolute errors reach a maximum of only 0.5 m/s2. Figure 7 shows the time histories of the parameters estimated by the UKF and CUKF. The final estimated values for the parameters are summarized in Table 3. By carefully inspecting the time histories for identifying parameter α, the UKF attempts to seek a value greater than 0.3 at approximately 1 s and 4 s, which is unacceptable. In contrast, the CUKF estimates parameters closer to the true values compared to the UKF by searching for parameters within constraints. The estimated values in Table 3 are maximum likelihood estimates for the parameters of a normal distribution based on the estimated sample sequences from the simulations.Figure 6
Comparison of acceleration time history results for Example 3 identified by the UKF and CUKF (a) and the corresponding absolute errors (b).
(a)(b)Figure 7
Results of parameter identification for Example 3.Table 3
Model parameters used in Example 3.
ParameterTrueInitial guessEstimatedAbsolute errorUKFCUKFUKFCUKFζn† (%)2—————ωn† (rad/s)2π—————α (×10−2)121016.9611.334.961.33γ (×102m−1)40.10.795.030.770.24β (×102m−1)150.37.4411.947.563.06n (—)211.731.930.270.07δν (×10−2)0.501.20.50.70δη (×10−2)0.501.71.21.20.7†These parameters are presumed to be determined.


## 3.1. Example 1: Original Bouc–Wen Hysteresis

A single-degree-of-freedom (SDOF) system with original Bouc–Wen hysteresis is considered in this example. It should be note that for the identification of earthquake damage, one specific element where a hysteretic behavior is expected, such as a column or a slender brace, is selected for observation to reduce the computational complexity: the overestimation problem may occur if measurement and identification are performed at many points including the point where elastic behavior is likely to occur. From the practical perspective, this study examined a SDOF system. The equation of motion for the example system when subjected to ground acceleration is written as(14)u¨+2ζnωnu˙+αωn2u+1−αωn2z=−u¨g,where u is the displacement, ζn is the damping ratio, ωn is the natural frequency, ug is the ground displacement, α is the postyielding stiffness ratio, an overdot represents the time derivative, and z is an auxiliary variable for describing the Bouc–Wen hysteresis governed by following differential equation:(15)z˙=u˙1−znγ+βsgnu˙z.To perform parameter identification using the UKF, an augmented state vector is introduced as follows:(16)x=u,v,z,α,γ,β,nT,where v=u˙ is the velocity of the structure. Based on this augmented state vector, equations (14) and (15) can be reformulated into a state-space equation as follows:(17)x˙=fx=v−u¨g−2ζnωnv−αωn2u+1−αωn2zv1−znγ+βsgnvzO4×1,where ON×M∈ℝN×M is a zero matrix. In this example, the noise term for the state variable in equation (2) is neglected. Additionally, the damping ratio and natural frequency are discarded from the augmented vector because it is assumed that they were previously identified under ambient conditions. This point will be discussed in detail in the next example. We choose absolute acceleration as the measurement variable, which is common in real-world application. The output equation for the measurement can be written as follows:(18)y=u¨+u¨g+w=−2ζnωnv−αωn2u−1−αωn2z+u¨g+w,where w is the measurement noise. A Gaussian white noise process with root mean square (RMS) value of 1 m/s2 is considered based on its relatively large value in the context of the general measurement procedure. We choose the ground motion of the 1940 El Centro earthquake as the input excitation. It has a 50 Hz sampling rate and a 30-second duration. To perform system identification using the Kalman filter methods, we reformulate equations (17) and (18) into a discretized form of equation (2) using the fourth-order Runge–Kutta integration method. The initial values of the parameters u and v in the space vector x are set to zero. Table 1 lists the true model parameters and initial guessed model parameters for this example. The postyielding parameter α is set to 0.09, which is a practical value for steel. The shape parameters n, β, and γ are set to 1.2, 20, and 40, respectively. The initial guesses are determined to have the same or a similar order of magnitude compared to the corresponding true values. We configure the initial covariance matrix P in diagonal matrix form to ensure that the designed filter could be applied to a roughly estimated condition as follows:(19)P=diag10−3,10−3,10−3,10−1,103,103,10−1,where diag⋅ is a diagonal matrix with arguments along the main diagonal. To utilize the CUKF, the following constraint conditions are determined based on physical intuition. First, the parameters α, γ, β, and n are constrained to be positive. Second, the upper bound of the stiffness ratio α is set to 0.2. Finally, β+γ is bounded in the range of 35–70, meaning the limit value zm should be bounded in a range of 0.01–0.03 if the parameter n is near unity. The constraint conditions onto which the estimated state will be projected constitute the following feasible region:(20)Ω=x∈ℝ7:Ax<b,where(21)A=000−10000000−10000000−10000000−10000−1−100000110,b=0000−3570.Table 1
Model parameters used in Example 1.
ParameterTrueInitial guessEstimatedAbsolute errorUKFCUKFUKFCUKFζn† (%)2—————ωn† (rad/s)2π—————α (×10−2)91014.8010.345.801.34γ (×10m−1)22.52.772.240.770.24β (×10m−1)43.55.464.181.460.18n (—)1.211.191.150.010.05†These parameters are presumed to be determined.Therefore, the procedure for finding the constrained estimatex^k defined in equation (13) is a linearly constrained quadratic optimization problem. A variety of numerical tools can solve such problems. A MATLAB function called quadprog with an interior-point-convex algorithm is employed to solve the constrained quadratic problem in this example.Figure2 presents a comparison of the accelerations estimated by the UKF and CUKF. It is evident that the acceleration histories estimated by the UKF and CUKF correlate with the true simulated responses such that they cannot be distinguished from the true responses. Their absolute errors reach a maximum of 0.5 m/s2.Figure 2
Comparison of acceleration time histories identified by UKF and CUKF (a) and their absolute errors (b) for Example 1.
(a)(b)Figure3 presents the time histories of the parameters estimated by the UKF and CUKF. Table 1 summarizes the final estimated values for the parameters. By carefully inspecting the time histories for identifying parameter α, one can see that the UKF attempts to seek a value of α larger than 0.3 at approximately 1 s and 4 s, which is an infeasible region from a physical perspective. In contrast, the CUKF estimates parameters closer to the true values compared with the UKF by searching for the parameters within admissible physical constraints. The inaccuracy of parameters combined with good agreement with predicted measurements indicates overfitting for the UKF. The parameters of the Bouc–Wen model are empirical, meaning they are not derived from fundamental mechanics, and a given response may not determine parameters unambiguously [34]. Therefore, it is helpful for applying the CUKF to system parameter identification if certain parameters can be constrained in a physical sense. It can be seen that the proposed CUKF can provide accurate results for estimation of all parameters except for the parameter n, which shows a large error of about 4%. However, when we recall from the previous studies on the parameters of the Bouc model that the variable n only describes the sharpness of hysteretic behavior (for example, refer to [34]), it can be seen that the estimation results accurately predict the stiffness and strength of the hysteretic behavior.Figure 3
Results of parameter identification in Example 1.


## 3.2. Example 2: Overestimation Problem

To test the overestimation problem, an additional simulation is conducted in which the structure examined in the previous example is subjected to ambient ground motion. In the example, the governing equation is identical to that for Example 1, but the structure is assumed to be subjected to the ground motion of the 1940 El Centro earthquake at a scale of 1/10. We assume that the structure is elastic and the Bouc–Wen parameters should have no influence on structural behavior. The augmented vector includes the damping ratio and natural frequency. We define the augmented state vectorx as follows:(22)x=u,v,z,α,γ,β,n,ωn,ζnT.A Gaussian white noise process with an RMS value of 0.1 m/s2 is applied to the measurements. Just as in the previous example, the initial values for the original state vector x are set to zero. Table 2 lists the true and initially guessed model parameters for this example. The true natural frequency and damping ratio are 2π and 0.02, respectively, which are same as the corresponding values in Example 1. The postyielding parameter α and shape parameters n, β, and γ are set to 0.09, 1.2, 20, and 40, respectively.Table 2
Model parameters used in Example 2.
ParameterTrueInitial guessEstimatedAbsolute errorUKFCUKF1CUKF2UKFCUKF1CUKF2ζn (%)23−3.483.682.785.481.680.78ωn (rad/s)2π2.4π−4.206.286.3310.480.000.05α (×10−2)91033.945.25—†24.943.75—†γ (×10m−1)22.5−2.524.64—†4.522.64—†β (×10m−1)43.50.132.19—†3.871.81—†n (—)1.210.961.50—†0.300.24—††These values are undefined because they are constrained by equality constraints.Two CUKFs are utilized in the example. CUKF1 is identical to the CUKF used in Example 1, except that the natural frequency and damping ratio are included in the state vector and constrained to be positive. CUKF2 is a constrained UKF that includes several equality constraints based on physical intuition. In CUKF2, to suppress the Bouc–Wen parameters during the identification process, the parametersα, γ, β, and n are set to be constrained to 0.1, 1, 25, and 35, respectively. The constraint conditions onto which the estimated state is projected constitute the following feasible region:(23)Ω=x∈ℝ9:Ax<bandAeqx=beq,where(24)A=000−1000000000−1000000000−1000000000−1000000000−1000000000−10000−1−1000000011000,b=000000−3570.Aeq=000100000000010000000001000000000100,beq=1025351.Figure4 presents a comparison of the accelerations estimated by the UKF, CUKF1, and CUKF2. The acceleration histories estimated by the UKF and CUKFs correlate well with the true simulated responses. Figure 5 shows the time histories of the parameters estimated by the UKF and CUKFs. The UKF attempts to seek the values of ωn and ζn smaller than zero, which is infeasible from a physical perspective, resulting in overestimation. In contrast, both CUKFs estimate these parameters accurately by searching within a positive space. CUKF1 forcibly predicts the Bouc–Wen parameters, which are unnecessary for describing the elastic behavior of the examined structure. Therefore, as noted in the previous example, it is very helpful for applying the CUKF to system parameter identification if some parameters can be constrained in a physical sense. We summarize the final values estimated by the UKF and CUKFs in Table 2, which indicates that the CUKFs are superior for estimating accurate parameters.Figure 4
Comparison of acceleration time history results for Example 2 identified by the UKF, CUKF1, and CUKF2.Figure 5
Results of parameter identification for Example 2.


## 3.3. Example 3: Structure with Degrading Hysteresis

The final example uses the SDOF structure from Example 1, but degrading effects are included. In the model proposed by [4], the auxiliary variable z in equation (15) is modified in order to employ strength and stiffness degradation as follows:(25)z˙=u˙1−νznγ+βsgnu˙zη,where the parameters ν and η are expressed as linear functions of the dissipated energy and their slopes δν and δη such as(26)ν=1+δνε,η=1+δηε,with ε being the dissipated hysteretic energy as(27)ε=εt=1−αωn2∫0tzu˙dt.An augmented state vector is introduced to formulate the state-space equation as follows:(28)x=u,v,z,ε,α,γ,β,n,δν,δηT.Using this augmented state vector, equations (25), and (26) are reformulated as the following nonlinear state-space equation:(29)x˙=fx=v−u¨g−2ζnωnv−αωn2u−1−αωn2zvη1−νznγ+βsgnvzzvO6×1.Figure6 presents a comparison of the accelerations estimated by the UKF and CUKF. The acceleration histories estimated by the UKF and CUKF correlate well with the true simulated responses, and they cannot be distinguished from the true responses. Their absolute errors reach a maximum of only 0.5 m/s2. Figure 7 shows the time histories of the parameters estimated by the UKF and CUKF. The final estimated values for the parameters are summarized in Table 3. By carefully inspecting the time histories for identifying parameter α, the UKF attempts to seek a value greater than 0.3 at approximately 1 s and 4 s, which is unacceptable. In contrast, the CUKF estimates parameters closer to the true values compared to the UKF by searching for parameters within constraints. The estimated values in Table 3 are maximum likelihood estimates for the parameters of a normal distribution based on the estimated sample sequences from the simulations.Figure 6
Comparison of acceleration time history results for Example 3 identified by the UKF and CUKF (a) and the corresponding absolute errors (b).
(a)(b)Figure 7
Results of parameter identification for Example 3.Table 3
Model parameters used in Example 3.
ParameterTrueInitial guessEstimatedAbsolute errorUKFCUKFUKFCUKFζn† (%)2—————ωn† (rad/s)2π—————α (×10−2)121016.9611.334.961.33γ (×102m−1)40.10.795.030.770.24β (×102m−1)150.37.4411.947.563.06n (—)211.731.930.270.07δν (×10−2)0.501.20.50.70δη (×10−2)0.501.71.21.20.7†These parameters are presumed to be determined.


## 4. Summary and Conclusions

Parameter identification of the structure is challenging in the seismic engineering field due to its complex nonlinear behavior. The Bouc–Wen hysteresis model is the most widely used model to identify the parameter. However, it can suffer from overfitting problems when measurements are limited in practice. We adopt the constrained unscented Kalman filter (CUKF) for the state and parameter estimation problem of the Bouc–Wen model to overcome the overfitting problems. The applied CUKF combines a constrained minimization technique with a conventional UKF. The predicted parameters are bounded for the reasonable constraints, which make the UKF estimate more accurate compared with the conventional UKF. We validate the CUKF by exhibiting several numerical examples. We finally show that the proposed CUKF is very helpful, robust, and efficient for system parameter identification by the comparative result when certain parameters are constrained in a physical sense.



---

*Source: 8822239-2020-08-13.xml*