約束UKF初始參數(shù)對Bouc－Wen模型參數(shù)識別的影響

王 濤， 吳 斌， 孟麗巖， 許國山， 張 健， 尹曉黎

(1．黑龍江科技大學 建筑工程學院，哈爾濱150022;2．哈爾濱工業(yè)大學 土木工程學院，哈爾濱150090)

0 引 言

確定結構恢復力模型及其模型參數(shù)是進行結構非線性地震反應分析、健康監(jiān)測及結構混合實驗的關鍵技術。迄今為止，研究者們已提出一系列的恢復力模型［1］來描述力與變形之間的關系。從形式上可分為折線形或曲線形，從應用的對象可分為層間層次、構件層次、截面層次及材料層次的恢復力模型。其中，Bouc 等［2］建立了經(jīng)典的Bouc －Wen模型。之后，Baber 等［3］進一步發(fā)展了Bouc－Wen 模型，能夠考慮強度退化、剛度退化及滑移捏縮效應。

在過去的幾十年里，EKF 被廣泛地應用在大量土木工程應用中的結構損傷識別和非線性結構參數(shù)識別等領域［4］。由于EKF 需要在每一個遞推步中求解相應的雅克比矩陣，當非線性系統(tǒng)不可導時，則很難求解雅克比矩陣。為了解決EKF 存在的問題，Julier 等［5］提出了隱性卡爾曼濾波器(Unscented Kalman filter，UKF)。與EKF 不同，UKF 采用確定性采樣策略去逼近狀態(tài)量的統(tǒng)計特征值，即隱性變換(Unscented transformation，UT)，然后利用卡爾曼濾波的預測－更新遞推算法框架進行狀態(tài)估計。UKF 同時具有很強的非線性狀態(tài)估計能力和較高的計算效率。

2008年，Wu 等［6］應用UKF 在線識別考慮退化、滑移的非線性滯回模型參數(shù)，研究表明，UKF 算法具有較高的識別精度和計算效率，能夠實現(xiàn)實時狀態(tài)估計和參數(shù)識別。Xie 等［7］提出迭代UKF 算法，用來在線識別非線性結構參數(shù)，該算法提高了UKF 對觀測噪聲的魯棒性。2013年，王濤等［8－9］在UKF 算法基礎上提出了約束UKF(Constrained UKF，CUKF)算法，并將此算法應用于結構混合實驗中。

采用UKF 或CUKF 進行參數(shù)識別時，需要事先給出初始估計誤差協(xié)方差、過程噪聲協(xié)方差和觀測噪聲協(xié)方差等濾波器的初始參數(shù)，這些初始參數(shù)的設置會直接影響參數(shù)識別結果。對于Bouc －Wen模型這樣的非線性系統(tǒng)，采用解析的方式來分析濾波器初始參數(shù)對模型參數(shù)識別影響規(guī)律將變得異常困難。為此，筆者提出Bouc－Wen 模型在線參數(shù)識別方法，采用數(shù)值模擬方式分析約束UKF 濾波器的初始參數(shù)對Bouc －Wen 模型參數(shù)識別精度和收斂速度的影響規(guī)律。

1 約束UKF 方法

UKF 的基本思想就是將UT 變換和卡爾曼濾波框架相結合，進行非線性狀態(tài)估計。在卡爾曼濾波框架下，預測步中通過確定采樣的方法生成2n1+1個Sigma 點，將這些Sigma 點通過狀態(tài)方程，然后加權統(tǒng)計變換的Sigma 點，可以得到先驗狀態(tài)估計和協(xié)方差;在更新步，將狀態(tài)Sigma 點通過觀測方程，然后加權統(tǒng)計變換后的觀測Sigma 點，得到先驗的觀測先驗估計和協(xié)方差及狀態(tài)和觀測的互協(xié)方差，進而可得到卡爾曼增益矩陣，然后通過卡爾曼更新方程，得到狀態(tài)后驗估計和協(xié)方差矩陣。

UKF 方法具體算法［10］，將一個非線性動力系統(tǒng)表達為離散狀態(tài)方程:

觀測方程為

式(1)、(2)中，xk為系統(tǒng)狀態(tài)向量為狀態(tài)向量維數(shù);yk+1為系統(tǒng)觀測向量為觀測向量維數(shù);uk為系統(tǒng)輸入;vk為過程噪聲，假定為高斯白噪聲，其均值為零，協(xié)方差矩陣為Qk;wk+1為觀測噪聲，假定為高斯白噪聲，其均值為零，協(xié)方差矩陣為Rk+1。

狀態(tài)估計的目的就是根據(jù)當前第k 步及第k 步之前的觀測信息，及系統(tǒng)的狀態(tài)方程和觀測方程來計算第k+1 步系統(tǒng)狀態(tài)估計值。計算第k 步2n1+1 個確定性Sigma 點Xk，i和其相應的權重:

上面得到的第k +1 步系統(tǒng)狀態(tài)估計均值和協(xié)方差，并沒有考慮到第k +1 步系統(tǒng)觀測yk+1的影響。故可由觀測方程式(2)計算觀測向量Sigma點

應用卡爾曼濾波更新方程，計算第k+1 步系統(tǒng)狀態(tài)估計均值和協(xié)方差Pk+1:

式(15)、(16)中，Kk+1為卡爾曼濾波增益矩陣，

可以看出，標準的UKF 算法并未考慮有界限約束，或其他數(shù)學約束條件的狀態(tài)估計問題。若狀態(tài)估計值違反了約束條件，將會影響在線模型更新準確性，甚至會影響非線性模型本身的穩(wěn)定性。為了實現(xiàn)約束狀態(tài)估計，采用基于Sigma 點調整的約束UKF 算法［9］，該算法在標準UKF 算法基礎上進行以下兩個方面的改進:第一，在預測步中，將違反界限約束的Sigma 點沿著協(xié)方差的方向移至界限約束邊界上，同時對在約束界限內的Sigma 點作相應移動，以保持重新調整后的一組Sigma 點仍具有對稱性;第二，在校正步中，采用狀態(tài)更新方程產(chǎn)生變換后的Sigma 點，當更新的狀態(tài)估計違反約束時，將違反約束的Sigma 點垂直投影到約束邊界上，通過加權統(tǒng)計的方式計算后驗的狀態(tài)估計均值和協(xié)方差。

2 Bouc－Wen 模型參數(shù)識別方法

經(jīng)典Bouc－Wen 滯回模型的數(shù)學表達式［2］:

式(19)、(20)中，r 為恢復力;d、v、z 分別為位移、速度、滯變位移;k、α、β、γ、n 分別為控制滯回環(huán)形狀和尺寸的五個模型參數(shù)。當n =∞時，Bouc －Wen模型骨架曲線近似為雙折線模型，此時，k 和α 可以被理解為雙折線模型的初始剛度和第二剛度系數(shù)。然而，Bouc －Wen 模型參數(shù)的物理意義并不明確，這使得確定這些模型參數(shù)變得十分困難。我們通過在線參數(shù)識別來得到Bouc－Wen 模型的參數(shù)，即k、β、γ、n、α 和時變的滯變位移z。

采用約束UKF 方法在線識別模型參數(shù)，需要定義系統(tǒng)的狀態(tài)向量x:

系統(tǒng)觀測為恢復力，即

位移d 和速度v 為系統(tǒng)的已知輸入，則系統(tǒng)的連續(xù)狀態(tài)方程為

系統(tǒng)的觀測方程

設Bouc－Wen 模型參數(shù)真實值分別為k =135 kN/mm，β=0.2，γ =0.2，n =1，α =0.02。設真實系統(tǒng)恢復力觀測噪聲wk+1為離散的高斯白噪聲序列，噪聲標準差為真實系統(tǒng)輸出恢復力標準差的2%，其大小相當于真實恢復力峰值的0.84%。位移輸入為曾完成的防屈曲支撐混合實驗的位移數(shù)據(jù)，位移輸入時程如圖1a 所示;速度輸入采用位移向前差分方法計算，得到速度輸入時程如圖1b 所示。

圖1 Bouc－Wen 模型位移和速度輸入Fig．1 Displacement and velocity inputs of Bouc-Wen model

模型參數(shù)界限約束為:k≥0，β≥－γ，n≥1，0≤α≤1。過程噪聲協(xié)方差設為Q =10－6I6;觀測噪聲協(xié)方差設為R=7.839 6 kN2。

狀態(tài)初始估計值為

初始狀態(tài)估計誤差協(xié)方差為

建立了系統(tǒng)狀態(tài)空間模型和濾波器初始參數(shù)之后，就可以采用遞推形式的約束UKF 算法來在線識別Bouc－Wen 模型參數(shù)。

3 濾波器初始參數(shù)影響分析

3.1 初始狀態(tài)估計誤差協(xié)方差

為了分析初始協(xié)方差對參數(shù)識別的影響，設觀測噪聲協(xié)方差P0=kPP，分別對比kP為0.1、1．0、10．0 三種情況下模型參數(shù)識別結果，此時其余濾波器參數(shù)保持不變，即不同P0下參數(shù)識別值時程曲線對比如圖2 所示。參數(shù)識別終值及其相對誤差對比如表1 所示。

圖2 P0 對Bouc－Wen 模型參數(shù)識別值影響Fig．2 Estimated parameters of Bouc-Wen model with different P0

表1 不同P0 下參數(shù)識別終值及其相對誤差Table 1 Final identification values and relative errors with different P0

圖2 的結果表明:在無模型誤差情況下，隨著kP的增加參數(shù)識別值收斂速度明顯加快，識別終值的相對誤差絕對值變小。當kP=10．0 時，除了參數(shù)α 以外其余參數(shù)識別終值的相對誤差均小于2.5%。同時也可以看出，kP從0.1 變化到10．0，約束UKF 均可以正常工作，可見識別算法對P0的選擇具有較好的魯棒性。因此，當我們對參數(shù)初值信心不是很大時，選擇較大的P0可以提高參數(shù)識別收斂速度和精度。

3.2 觀測噪聲協(xié)方差

設觀測噪聲協(xié)方差Rk= kRR，分別對比kR為0.1、1．0、10．0 三種情況下模型參數(shù)識別結果，此時其余濾波器參數(shù)保持不變，即Q。參數(shù)識別終值及其相對誤差對比見表2。不同Rk下參數(shù)識別值時程曲線對比如圖3 所示。

表2 不同Rk 下參數(shù)識別終值及其相對誤差Table 2 Final identification values and relative errors with different Rk

圖3 的結果表明:在無模型誤差情況下，當kR=1．0時，即濾波器的觀測噪聲協(xié)方差與真實系統(tǒng)觀測噪聲協(xié)方差相等時，參數(shù)識別值收斂速度更快，識別終值的相對誤差絕對值變小。當kR=1．0 時，除了參數(shù)α 以外其余參數(shù)識別終值的相對誤差均小于7.5%。kR越小，說明我們更相信觀測，這樣前期參數(shù)識別受觀測噪聲的影響更大，識別值的波動性會加大;相反，kR越大，我們會更相信初值選擇，識別值的收斂速度會變慢。同時也可以看出，kR從0.1 變化到10．0，約束UKF 均可以正常工作，可見識別算法對Rk的選擇具有較好的魯棒性。因此，最好能夠事先統(tǒng)計真實的系統(tǒng)觀測噪聲協(xié)方差，將此作為濾波器的觀測噪聲協(xié)方差，這樣，可以提高參數(shù)識別收斂速度和精度。

圖3 Rk 對Bouc－Wen 模型參數(shù)識別值影響Fig．3 Estimated parameters of Bouc-Wen model with different Rk

3.3 過程噪聲協(xié)方差

設過程噪聲協(xié)方差Qk=kQQ，分別對比kQ為0、10．0、100．0 三種情況下模型參數(shù)識別結果，此時其余濾波器參數(shù)保持不變，即參數(shù)識別終值及其相對誤差對比見表3。不同Qk下參數(shù)識別值時程曲線對比如圖4 所示。

表3 不同Qk 下參數(shù)識別終值及其相對誤差Table 3 Final identification values and relative errors with different Qk

圖4 的結果表明:在無模型誤差情況下，當kQ=0 時，參數(shù)識別值收斂速度更快，識別終值的相對誤差絕對值變小。當kQ=0 時，除了參數(shù)α 以外其余參數(shù)識別終值的相對誤差均小于5.5%。kQ越大，參數(shù)識別值的波動性會加大，收斂速度會變慢。同時也可以看出，kQ從0 變化到10．0，約束UKF 均可以正常工作，而且對識別結果影響并不顯著，可見識別算法對Qk的選擇具有較好的魯棒性。因此，當真實系統(tǒng)無模型不確定性問題時，可以不用考慮過程噪聲，這樣可以提高參數(shù)識別收斂速度和精度。對于實際情況，系統(tǒng)模型一般不可能完全精確描述真實系統(tǒng)，可以通過設置較小的過程噪聲協(xié)方差來減小模型誤差影響。

3.4 初始狀態(tài)估計

圖4 Qk 對Bouc－Wen 模型參數(shù)識別值影響Fig．4 Estimated parameters of Bouc-Wen model with different Qk

圖5 對Bouc－Wen 模型參數(shù)識別值影響Fig．5 Estimated parameters of Bouc-Wen model with different

表4 不同 下參數(shù)識別終值及其相對誤差Table 4 Final identification values and relative errors with different

表4 不同 下參數(shù)識別終值及其相對誤差Table 4 Final identification values and relative errors with different

注:相對誤差(%)=100% × |識別值－真實值|/真實值。

參數(shù) 真實值 初始值kx =0.8識別值 誤差/%kx=1．0識別值 誤差/%kx =1.2識別值 誤差/%k 135.00 115．0 132.471 1.87 133.080 1.41 134.134 0.64 β 0.20 0.5 0.200 0 0.200 0 0.201 0.50 γ 0.20 0.5 0.178 11.00 0.185 －7.50 0.196 2.00 n 1.00 2.0 1.083 8.30 1.069 6.90 1.051 5.10 α 0.02 0.1 0.035 75.00 0.038 90.00 0.043 115.00

圖5 的結果表明:在無模型誤差情況下，kx從0.8 變化到1.2，約束UKF 均可以較好地識別模型參數(shù)，除了參數(shù)α 以外其余參數(shù)識別終值的相對誤差均不大于8.3%，可見識別算法對^x0的選擇具有較好的魯棒性。整體上，初始參數(shù)值與真實值越接近，識別值收斂會越快，識別終值的相對誤差絕對值也相應變小。因此，當對參數(shù)初值^x0的信心更有把握時，可以提高參數(shù)識別收斂速度和精度。

4 結 論

(1)與標準的UKF 算法相比，約束UKF 法在不增加計算負荷的基礎上，可以解決考慮有界限約束的非線性狀態(tài)估計問題。可以減小參數(shù)識別值的前期波動幅度，提高識別值的收斂速度，有助于提高在線模型更新準確性。文中給出基于約束UKF 的Bouc－Wen 模型參數(shù)識別方法。

(2)數(shù)值模擬研究表明，在無模型誤差的情況下，濾波器對于初始參數(shù)具有較好的魯棒性，初始參數(shù)可以具有較寬的取值范圍。

(3)給出在實際應用中濾波器參數(shù)選擇建議，通過適當?shù)脑龃蟪跏紶顟B(tài)估計協(xié)方差，減小過程噪聲，采用真實系統(tǒng)觀測噪聲協(xié)方差，以及減小初始參數(shù)值與真實值的偏差，可以有效提高參數(shù)識別收斂速度和精度。

(4)該研究可以為解決土木結構非線性模型在線參數(shù)識別和模型更新混合實驗提供參考。

［1］郭子雄，楊 勇．恢復力模型研究現(xiàn)狀及存在問題［J］．世界地震工程，2004，20(4):47 －51．

［2］WEN Y K．Method for random vibration of hysteretic systems［J］．Journal of the Engineering Mechanics Division ASCE，1976，102(2):249 －263．

［3］BABER T T，WEN Y K．Random vibration hysteretic，degrading systems［J］．Journal of the Engineering Mechanics Division ASCE，1981，107(6):1069 －1087．

［4］YANG J N，LIN S，HUANGL H，et al．An adaptive extended Kalman filter for structural damage identification［J］．Structural Control and Health Monitoring，2006，13(4):849 －867．

［5］JULIER S J，UHLMANN J K．Unscented filtering and nonlinear estimation［J］．Proceendings of the IEEE，2004，92(3):401－422．

［6］WU M，SMYTH A．Real-time parameter estimation for degrading and pinching hysteretic models［J］．International Journal of Non-Linear Mechanics，2008，43(9):822 －833．

［7］XIE Z B，F(xiàn)ENG J C．Real-time nonlinear structural system identification via iterated unscented Kalman filter［J］．Mechanical Systems and Signal Processing，2012，28:309 －322．

［8］王 濤，吳 斌．基于UKF 模型更新的混合試驗方法［J］．振動與沖擊，2013，32(5):138 －143．

［9］王 濤，吳 斌．基于約束UKF 模型更新的混合試驗方法［J］．地震工程與工程振動，2013，33(5):100 －109．

［10］WAN E，VAN DER MERWE R．Kalman filtering and neural networks［M］．Wiley:［s．n．］，2001:221 －280．