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

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

0 引 言

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

在過去的幾十年里，EKF 被廣泛地應用在大量土木工程應用中的結(jié)構(gòu)損傷識別和非線性結(jié)構(gòu)參數(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 算法具有較高的識別精度和計算效率，能夠?qū)崿F(xiàn)實時狀態(tài)估計和參數(shù)識別。Xie 等［7］提出迭代UKF 算法，用來在線識別非線性結(jié)構(gòu)參數(shù)，該算法提高了UKF 對觀測噪聲的魯棒性。2013年，王濤等［8－9］在UKF 算法基礎上提出了約束UKF(Constrained UKF，CUKF)算法，并將此算法應用于結(jié)構(gòu)混合實驗中。

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

1 約束UKF 方法

UKF 的基本思想就是將UT 變換和卡爾曼濾波框架相結(jié)合，進行非線性狀態(tài)估計。在卡爾曼濾波框架下，預測步中通過確定采樣的方法生成2n1+1個Sigma 點，將這些Sigma 點通過狀態(tài)方程，然后加權(quán)統(tǒng)計變換的Sigma 點，可以得到先驗狀態(tài)估計和協(xié)方差;在更新步，將狀態(tài)Sigma 點通過觀測方程，然后加權(quán)統(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和其相應的權(quán)重:

上面得到的第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 點調(diào)整的約束UKF 算法［9］，該算法在標準UKF 算法基礎上進行以下兩個方面的改進:第一，在預測步中，將違反界限約束的Sigma 點沿著協(xié)方差的方向移至界限約束邊界上，同時對在約束界限內(nèi)的Sigma 點作相應移動，以保持重新調(diào)整后的一組Sigma 點仍具有對稱性;第二，在校正步中，采用狀態(tài)更新方程產(chǎn)生變換后的Sigma 點，當更新的狀態(tài)估計違反約束時，將違反約束的Sigma 點垂直投影到約束邊界上，通過加權(quán)統(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ù)識別結(jié)果，此時其余濾波器參數(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 的結(jié)果表明:在無模型誤差情況下，隨著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ù)識別結(jié)果，此時其余濾波器參數(shù)保持不變，即Q。參數(shù)識別終值及其相對誤差對比見表2。不同Rk下參數(shù)識別值時程曲線對比如圖3 所示。

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

圖3 的結(jié)果表明:在無模型誤差情況下，當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ù)識別結(jié)果，此時其余濾波器參數(shù)保持不變，即參數(shù)識別終值及其相對誤差對比見表3。不同Qk下參數(shù)識別值時程曲線對比如圖4 所示。

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

圖4 的結(jié)果表明:在無模型誤差情況下，當kQ=0 時，參數(shù)識別值收斂速度更快，識別終值的相對誤差絕對值變小。當kQ=0 時，除了參數(shù)α 以外其余參數(shù)識別終值的相對誤差均小于5.5%。kQ越大，參數(shù)識別值的波動性會加大，收斂速度會變慢。同時也可以看出，kQ從0 變化到10．0，約束UKF 均可以正常工作，而且對識別結(jié)果影響并不顯著，可見識別算法對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 的結(jié)果表明:在無模型誤差情況下，kx從0.8 變化到1.2，約束UKF 均可以較好地識別模型參數(shù)，除了參數(shù)α 以外其余參數(shù)識別終值的相對誤差均不大于8.3%，可見識別算法對^x0的選擇具有較好的魯棒性。整體上，初始參數(shù)值與真實值越接近，識別值收斂會越快，識別終值的相對誤差絕對值也相應變小。因此，當對參數(shù)初值^x0的信心更有把握時，可以提高參數(shù)識別收斂速度和精度。

4 結(jié) 論

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

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

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

(4)該研究可以為解決土木結(jié)構(gòu)非線性模型在線參數(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．