王 云，許 斌

(1.湖南大學(xué) 土木工程學(xué)院，長沙 410082； 2.華僑大學(xué) 土木工程學(xué)院，福建 廈門 361021；3.華僑大學(xué) 福建省結(jié)構(gòu)工程與防災(zāi)重點實驗室，福建 廈門 361021)

部分觀測下結(jié)構(gòu)質(zhì)量識別迭代及滯回力識別方法

王 云1，許 斌2,3

(1.湖南大學(xué) 土木工程學(xué)院，長沙 410082； 2.華僑大學(xué) 土木工程學(xué)院，福建 廈門 361021；3.華僑大學(xué) 福建省結(jié)構(gòu)工程與防災(zāi)重點實驗室，福建 廈門 361021)

在強(qiáng)動力荷載作用下結(jié)構(gòu)滯回力可直接描述損傷的發(fā)生發(fā)展過程，并用于結(jié)構(gòu)耗能的定量評估。提出一種僅利用結(jié)構(gòu)部分自由度上加速度響應(yīng)時程且結(jié)構(gòu)剛度、阻尼、質(zhì)量均未知時的結(jié)構(gòu)質(zhì)量識別迭代和滯回力識別方法。首先基于質(zhì)量預(yù)估值利用擴(kuò)展卡爾曼濾波方法預(yù)測結(jié)構(gòu)的位移和速度響應(yīng)以及未知加速度響應(yīng)，利用二重切比雪夫多項式對質(zhì)量進(jìn)行更新并迭代直到收斂，最后基于收斂的質(zhì)量識別值描述結(jié)構(gòu)滯回性能，識別結(jié)構(gòu)滯回力。在一個多自由度集中質(zhì)量數(shù)值模型中引入磁流變阻尼器模擬非線性構(gòu)件，在結(jié)構(gòu)僅已知部分自由度上加速度響應(yīng)時程的條件下，數(shù)值模擬結(jié)果發(fā)現(xiàn)結(jié)構(gòu)質(zhì)量分布和阻尼器的滯回阻尼力均有良好的識別結(jié)果，同時考察了測量噪聲的影響。

振動與波：非線性恢復(fù)力識別；質(zhì)量識別；擴(kuò)展卡爾曼濾波；切比雪夫多項式模型；部分加速度測量

基于結(jié)構(gòu)動力響應(yīng)測量和結(jié)構(gòu)動力學(xué)基本原理對結(jié)構(gòu)物理參數(shù)進(jìn)行識別，借此描述結(jié)構(gòu)狀態(tài)和損傷，并對結(jié)構(gòu)進(jìn)行補(bǔ)強(qiáng)和加固進(jìn)而提升結(jié)構(gòu)的使用性能是土木工程領(lǐng)域緊迫的課題之一。

相對于非線性動力系統(tǒng)而言，對線性系統(tǒng)進(jìn)行參數(shù)化識別開展得較早，研究比較深入[1–3]。然而，結(jié)構(gòu)在地震等強(qiáng)動力荷載作用下的損傷發(fā)生發(fā)展以及破壞過程會呈現(xiàn)出典型的非線性特征。因此，研究適用于非線性動力系統(tǒng)的結(jié)構(gòu)識別方法具有重要意義。李宏男等對在結(jié)構(gòu)健康監(jiān)測和損傷識別中研究非線性系統(tǒng)的識別問題的重要性進(jìn)行了闡述[4–5]。Smyth等提出帶遺忘因子的修正最小二乘法來識別描述非線性滯回性能的Bouc-Wen模型中參數(shù)的識別，并通過在鏈?zhǔn)蕉嘧杂啥冉Y(jié)構(gòu)中引入Bouc-Wen非線性構(gòu)件驗證了方法的可行性[6]。此外實際土木工程結(jié)構(gòu)由于材料和結(jié)構(gòu)形式的差別其非線性行為差異明顯，難以用準(zhǔn)確的非線性數(shù)學(xué)模型建立。國外學(xué)者在結(jié)構(gòu)免模型識別上也得到了長足的發(fā)展。在免模型的結(jié)構(gòu)非線性狀態(tài)識別方面，Masri和Caughey最早提出了恢復(fù)力曲面法(Restoring Force Surface Method)，并利用該方法對單自由度非線性結(jié)構(gòu)的非線性恢復(fù)力進(jìn)行了有效識別[7]。之后該方法被推廣到多自由度系統(tǒng)，識別多自由度系統(tǒng)的非線性恢復(fù)力，即滯回性能[8]。Xu等利用已知激勵和結(jié)構(gòu)響應(yīng)，基于等價線性理論和最小二乘法對結(jié)構(gòu)非線性恢復(fù)力進(jìn)行識別，并通過裝有磁流變阻尼器的四層剪切型框架結(jié)構(gòu)模型的動力試驗測量數(shù)據(jù)，驗證了該方法的有效性，并考慮激勵部分未知時的實現(xiàn)方法[9–12]。許斌等利用數(shù)值模擬和試驗實測數(shù)據(jù)驗證了基于冪級數(shù)多項式和切比雪夫多項式模型的含磁流變阻尼器的非線性多自由度系統(tǒng)滯回力識別方法的有效性[13]。

本文提出基于擴(kuò)展卡爾曼濾波(Extended Kalman Filter，EKF)和切比雪夫多項式模型，在結(jié)構(gòu)僅部分自由度受激勵且僅部分自由度上加速度響應(yīng)已知時的質(zhì)量分布迭代識別和多自由度結(jié)構(gòu)滯回力識別方法。通過在一個多自由度線性系統(tǒng)中引入非線性構(gòu)件磁流變阻尼器，形成一個非線性結(jié)構(gòu)，且非線性結(jié)構(gòu)的非線性特征具有合適的數(shù)學(xué)物理模型來描述時，通過數(shù)值模擬驗證所提出的方法識別質(zhì)量分布和結(jié)構(gòu)滯回力的有效性，并考慮測量噪聲的影響。

1 非線性系統(tǒng)的擴(kuò)展卡爾曼濾波算法

為運用結(jié)構(gòu)部分加速度響應(yīng)測量進(jìn)行非線性系統(tǒng)的滯回性能識別，利用EKF算法預(yù)測結(jié)構(gòu)各個自由度位移和速度響應(yīng)以及未知加速度響應(yīng)，并識別非線性系統(tǒng)的物理參數(shù)。

對于一般的n自由度非線性動力系統(tǒng)，其運動方程可寫為

其中M為質(zhì)量矩陣，K為剛度矩陣，C為阻尼矩陣，fnon(θ)為非線性恢復(fù)力向量分別為結(jié)構(gòu)加速度、速度、位移向量，f(t)為外激勵向量。

定義以下擴(kuò)展?fàn)顟B(tài)向量

X(t)表示t時刻的狀態(tài)向量，、和分別為待識別的阻尼、剛度和非線性系統(tǒng)參數(shù)向量，不考慮參數(shù)的時變。于是得到如下的擴(kuò)展?fàn)顟B(tài)方程

式(3)改寫成如下形式

當(dāng)觀測部分加速度響應(yīng)時，對于離散采樣點，時刻tk=k t時的觀測響應(yīng)可以寫為如下形式

其中△t為采樣時間間隔，y(k)為tk時刻的觀測向量，V(k)為均值為零的觀測噪聲向量，滿足以下兩式要求

其中R為觀測噪聲的協(xié)方差矩陣，δij是Kronecker算子。

EKF算法通過在濾波估計值附近作線性化處理，逐步獲得下一時刻的最優(yōu)估計，其主要計算步驟如下。

1）狀態(tài)預(yù)測

2）誤差協(xié)方差預(yù)測方程

其中

3)增益矩陣

其中

4）狀態(tài)濾波方程

5）誤差協(xié)方差方程

在假定結(jié)構(gòu)質(zhì)量和已知部分結(jié)構(gòu)加速度的情況下，運用EKF，可以識別出非線性結(jié)構(gòu)的狀態(tài)向量以及參數(shù)向量，包括結(jié)構(gòu)位移、速度和位置加速度，以及結(jié)構(gòu)線性和非線性恢復(fù)力的模型參數(shù)。

以下通過引入二重切比雪夫多項式對結(jié)構(gòu)的質(zhì)量進(jìn)行識別與更新，并對多自由度結(jié)構(gòu)的滯回性能進(jìn)行描述。

2 基于二重切比雪夫多項式的多自由度結(jié)構(gòu)滯回力識別方法

2.1 二重切比雪夫多項式

余弦函數(shù)的多倍角公式如式（15）所示。

令θ=arccosα，則函數(shù)

其中n=0，1，2，…

可得到如下遞推關(guān)系

由式(17)遞推得到第一類切比雪夫多項式，其前幾階可以表示為

2.2 結(jié)構(gòu)質(zhì)量分布與滯回力的識別

將(1)式的運動方程改寫為下式

其中Rnon(t)為結(jié)構(gòu)總非線性恢復(fù)力，由線性系統(tǒng)的彈性恢復(fù)力、粘性阻尼力以及非線性構(gòu)件提供的非線性恢復(fù)力構(gòu)成，一般可用結(jié)構(gòu)相對速度和位移響應(yīng)來表示，如式(20)所示

式中x˙(t)，x(t)分別為結(jié)構(gòu)相對速度和位移響應(yīng)。

由于切比雪夫多項式的正交性，在非線性鏈?zhǔn)蕉嘧杂啥燃匈|(zhì)量模型中，兩集中質(zhì)量間的滯回力可近似表示為下式。

由于切比雪夫多項式中自變量的定義域為[-1，1]，需將相對位移和相對速度進(jìn)行標(biāo)準(zhǔn)化處理，v′i,i-1和S′i,i-1是經(jīng)過標(biāo)準(zhǔn)化層間相對速度和相對位移，通過下式確定。

于是，滯回力可用二維切比雪夫多項式模型來表示

考慮一般情況，即結(jié)構(gòu)僅有部分自由度上作用有激勵力。在這種非完整激勵條件下，以一個3自由度體系為例，假設(shè)外激勵僅作用在第3個自由度上，則各自由度的運動平衡方程可表示為式(25)－式(27)。

由式(25)可得

運用已知部分加速度響應(yīng)以及以上通過EKF預(yù)測的結(jié)構(gòu)位移、速度和加速度響應(yīng)，可以利用最小二乘優(yōu)化求解得到結(jié)構(gòu)的質(zhì)量以及滯回力切比雪夫多項式的參數(shù)。進(jìn)而根據(jù)牛頓第三定律可對其他自由度上質(zhì)量和滯回力進(jìn)行識別。

如根據(jù)自由度2與3之間的相互作用力的互等關(guān)系

由式(26)和式(29)得

同理識別出自由度2上的質(zhì)量以及自由度2與1之間的滯回力。同理，可得自由度1與地面之間的滯回力以及自由度1的質(zhì)量。于是得到結(jié)構(gòu)各自由度間的滯回力以及各質(zhì)點質(zhì)量識別結(jié)果。

3 基于EKF和二重切比雪夫多項式的質(zhì)量以及滯回力識別迭代算法

綜合運用EKF和二重切比雪夫多項式，通過迭代算法，可進(jìn)行結(jié)構(gòu)質(zhì)量和滯回力的識別。該迭代算法的框圖如圖1所示。

圖1 滯回力和質(zhì)量識別迭代算法流程

首先，將初始結(jié)構(gòu)質(zhì)量的預(yù)估值代入到EKF中預(yù)測出結(jié)構(gòu)位移和速度以及未知加速度響應(yīng)；其次，利用結(jié)構(gòu)位移和速度響應(yīng)，基于二重切比雪夫多項式，得到結(jié)構(gòu)的質(zhì)量的更新值；最后，引入學(xué)習(xí)因子，對結(jié)構(gòu)質(zhì)量近似值進(jìn)行修正，進(jìn)行迭代求解。在本算法中，學(xué)習(xí)因子α和β均為0.5。

4 數(shù)值算例

4.1 計算模型與參數(shù)

以一個具有4個集中質(zhì)量的結(jié)構(gòu)體系為例，為模擬結(jié)構(gòu)的非線性行為，在結(jié)構(gòu)第四層安裝一個非線性構(gòu)件磁流變阻尼器MR。線性結(jié)構(gòu)各層參數(shù)取值如下：質(zhì)量mi=200 kg，阻尼ci=100 N·s/m，剛度ki=3×104N/m ，(i=1，2，3，4)。不失一般性，假設(shè)外激勵僅水平方向作用在結(jié)構(gòu)的第3個自由度上，其他集中質(zhì)點上不施加荷載。

在數(shù)值模擬中MR阻尼器采用修正的Dahl模型[15]。其表達(dá)式為

式中參數(shù)σ用來控制阻尼器恢復(fù)力滯回曲線的形狀。在本文中，σ=50 s/m，K0=50 N/mC0=199 5 N·s/m，F(xiàn)d=34.85 N，f0=0 N。

圖2 非線性結(jié)構(gòu)計算模型

針對本問題，EKF的擴(kuò)展?fàn)顟B(tài)向量為

為了更好地模擬結(jié)構(gòu)非線性恢復(fù)力，本文采用的切比雪夫多項式模型中取k+q=3。觀測第1、第3、第4集中質(zhì)點上的加速度響應(yīng)作為已知加速度響應(yīng)，第2質(zhì)點上的加速度響應(yīng)未知。另外，作用在第3質(zhì)點上的外激勵已知。

本文選用的外激勵采用隨機(jī)激勵，其時域圖如圖3所示。

圖3 結(jié)構(gòu)隨機(jī)外激勵荷載時程

荷載作用時間為2 s。結(jié)構(gòu)在該動力荷載作用下的響應(yīng)采用4階龍格庫塔(Runge-Kutta)積分法計算，積分步長為0.000 5 s。其中，加速度混入了5%的量測噪聲，以考察該方法的抗噪性能。此外，為了驗證算法在不同初始條件下的收斂性能，本文分別考慮了初始質(zhì)量為理論值110%和90%兩種情況。

4.2 工況一：質(zhì)量大于真實值10%

(1)質(zhì)量、參數(shù)與響應(yīng)識別效果

基于擴(kuò)展卡爾曼濾波和切比雪夫多項式，識別出質(zhì)量的真實值收斂曲線如圖4所示。收斂準(zhǔn)則為后一次迭代與前一次的差的絕對值不大于0.1 kg。四個自由度對應(yīng)的質(zhì)量識別結(jié)果分別為m1=197 kg，m2=199 kg，m3=200 kg，m4=198 kg。與真實值吻合較好，本文方法可有效識別出結(jié)構(gòu)質(zhì)量分布。在迭代結(jié)束時，根據(jù)結(jié)構(gòu)第1、第3、第4集中質(zhì)點的加速度響應(yīng)觀測值，基于EKF，預(yù)測出結(jié)構(gòu)所有質(zhì)點的位移和速度響應(yīng)，并與積分計算所得結(jié)果進(jìn)行對比，結(jié)果如圖5和圖6所示。

對于未知的第2層加速度響應(yīng)，在得到速度預(yù)測值的情況下根據(jù)線加速度法的原理，可計算出第2層加速度響應(yīng)。第2層加速度響應(yīng)預(yù)測值與真實值的比較如圖7所示。

迭代結(jié)束時，結(jié)構(gòu)剛度、阻尼以及MR阻尼器的各參數(shù)的識別結(jié)果與真實值的比較見表1。比較發(fā)現(xiàn)，剛度與阻尼以及MR阻尼器的參數(shù)均有比較好的識別效果。

表1 結(jié)構(gòu)參數(shù)識別結(jié)果

圖4 質(zhì)量識別結(jié)果迭代過程

圖5 結(jié)構(gòu)位移響應(yīng)識別結(jié)果

(2)滯回力識別結(jié)構(gòu)

通過EKF方法識別的非線性參數(shù)和所有質(zhì)點的位移和速度響應(yīng)，MR阻尼器的非線性恢復(fù)力的識別值與真實值比較如圖8所示，對比發(fā)現(xiàn)，識別出的非線性恢復(fù)力有良好的識別精度。

圖6 結(jié)構(gòu)速度響應(yīng)識別結(jié)果

圖7 第2層加速度響應(yīng)識別值

圖8 MR阻尼器恢復(fù)力識別結(jié)果

工況二考慮了在質(zhì)量初始值小于理論值10%的情況。同樣在觀測的結(jié)構(gòu)加速度響應(yīng)中添加了5%的噪聲，仍采取圖3所示的隨機(jī)激勵，所得識別效果如下。

4.3 工況二：質(zhì)量小于真實值10%

(1)質(zhì)量、參數(shù)和響應(yīng)識別結(jié)果

基于EKF和二重切比雪夫多項式識別出的質(zhì)量為m1=202 kg，m2=200 kg，m3=200 kg，m4=202 kg，與真實值吻合較好，說明本文提供的方法在質(zhì)量初始值小于真實值10%的情況下可以有效識別出結(jié)構(gòu)的質(zhì)量分布，驗證了該方法的有效性，其質(zhì)量識別結(jié)果收斂曲線如圖9所示。

圖9 質(zhì)量識別迭代過程

迭代結(jié)束時通過擴(kuò)展卡爾曼濾波識別的參數(shù)結(jié)果如表2，識別出的參數(shù)與真實值吻合較好，在質(zhì)量小于10%的工況下，結(jié)構(gòu)的剛度、阻尼以及非線性參數(shù)均有良好的識別效果。

(2)滯回力識別結(jié)果

識別出的MR阻尼器滯回力結(jié)果如圖10所示。在質(zhì)量初始值低于真實值10%的情況下，MR阻尼器滯回力識別值與真實值吻合較好。

表2 結(jié)構(gòu)參數(shù)識別結(jié)果

圖10 MR阻尼器恢復(fù)力識別結(jié)果

5 結(jié)語

本文提出了一種在非完整激勵下僅利用結(jié)構(gòu)部分自由度上加速度響應(yīng)時程，在結(jié)構(gòu)質(zhì)量、剛度和阻尼信息未知時，基于EKF和二重切比雪夫多項式模型的質(zhì)量分布迭代識別以及結(jié)構(gòu)滯回力識別方法。

通過在一個4自由度線性模型中引入磁流變阻尼器作為非線性構(gòu)件，形成非線性結(jié)構(gòu)，運用其在一個集中動力荷載作用下3自由度上加速度響應(yīng)時程，未知結(jié)構(gòu)質(zhì)量分布，對結(jié)構(gòu)中質(zhì)量分布以及MR的滯回力進(jìn)行識別。考慮加速度測量噪聲以及不同質(zhì)量初始值下的計算結(jié)果。數(shù)值模擬結(jié)果表明，僅利用作用在部分自由度上的外激勵及其部分加速度響應(yīng)信息，所提出的算法在加速度測量噪聲以及不同質(zhì)量初始值下均可對結(jié)構(gòu)的質(zhì)量分布以及滯回力進(jìn)行識別，該算法具有較好的抗噪聲能力。

結(jié)構(gòu)滯回力的識別結(jié)果可以直接描述結(jié)構(gòu)構(gòu)件在強(qiáng)動力荷載作用下的非線性行為的發(fā)生發(fā)展過程，而且可以用于定量計算結(jié)構(gòu)構(gòu)件的耗能，對工程結(jié)構(gòu)在地震、爆炸等強(qiáng)動力荷載作用下的性能評估具有重要意義。所提出的方法可以運用于質(zhì)量未知情況下非線性系統(tǒng)的識別問題。

[1]WU Z S,XU B,HARADA T.Review on structural health monitoring forinfrastructure[J].Journalof Applied Mechanics,2003,70(6):1043-1054.

[2]孫鴻敏，李宏男.土木工程結(jié)構(gòu)健康監(jiān)測研究進(jìn)展[J].防災(zāi)減災(zāi)工程學(xué)報，2003，23(3)：92-98.

[3]宗周紅，任偉新，阮毅.土木工程結(jié)構(gòu)損傷診斷研究進(jìn)展[J].土木工程學(xué)報，2003，36(5)：105-110.

[4]李宏男，高東偉，伊延華.土木工程結(jié)構(gòu)健康監(jiān)測系統(tǒng)的研究狀況與進(jìn)展[J].力學(xué)進(jìn)展，2008，38(2)：151-166.

[5]李宏男，李東升.土木工程結(jié)構(gòu)安全性評估、健康監(jiān)測及診斷評述[J].地震工程與工程振動，2002，22(3)：82-90.

[6]A W SMYTH.On-line parametric identification of MDOF nonlinear hysteretic systems[J].Journal of Engineering Mechanics,1999,125(2):133-142.

[7]MASRISF,CAUGHEYTK.Anonparametric identification technique for nonlinear dynamic problem[J].Journal ofApplied Mechanics,1979,46(3):433-447.

[8]MASRI S F,SASSI H,CAUGHEY T K.Nonparametric identification of nearly arbitrary nonlinear systems[J].Journal ofApplied Mechanics,1982,49(5):619-628

[9]XU B,HE J,MASRI S F.Time domain data-based model free structural nonlinear performance identification[C].Proceeding of International Symposium on Innovation&Sustainability of Structures in Civil Engineering,2009:59-69.

[10]XU B,MASRI S F.Nonlinearity identification for a frame modelstructure with MR dampers underlimited excitations[J].Asia-pacific Young Researchers&Graduates Symposium,2010,159:167.

[11]XU B,HE J,AND MASRI S F.Data-based model-free hysteretic restoring force and mass identification for dynamic systems[J].Computer-Aided Civiland Infrastructural Engineering,30:2-18,2015

[12]XU B,HE J,MASRI S F.Data-based identification of nonlinearrestoring forceunderspatially incomplete excitationswith powerseriespolynomialmodel[J].Nonlinear Dynamics,2011,67(3):2063-2080.

[13]許斌，辛璐璐，賀佳.基于二重切比雪夫多項式的多自由度系統(tǒng)SMA非線性恢復(fù)力識別[J].振動與沖擊，2014，33(16)：6-13.

[14]周強(qiáng)，瞿偉廉.磁流變阻尼器的兩種力學(xué)模型和試驗驗證[J].地震工程與工程振動，2002，22(4)：144-150.

Mass and Hysteretic Force Identification for Structures Based on Limited Measurement Data

WANG Yun1,XU Bin2,3
(1.College of Civil Engineering,Hunan University,Changsha 410082,China;2.College of Civil Engineering,Huaqiao University,Xiamen 361021,Fujian China;3.Key Laboratory for Structural Engineering and Disaster Prevention of Fujian Province,Huaqiao University,Xiamen 361021,Fujian China)

Structural hysteretic force can provide a direct description for the initiation and development of structural damage and can be used for the quantitative structural energy consumption evaluation under dynamic loadings.In this study,a mass distribution and iterative structural hysteretic force identification approach is proposed.Acceleration responses at part of the degrees of freedoms(DOFs)of a nonlinear structure with unknown structure mass,stiffness and damping coefficients are employed for the hysteretic force and mass identification.First of all,the extended Kalman filter(EKF)is employed to forecast structure displacement,velocity in all DOFs and the unmeasured acceleration responses with the assessed mass matrix and the partially measured acceleration response.Then,a double Chebyshev polynomial model is employed to identify the hysteretic force and to update the mass.With the updated mass,the identification approach is repeated until the convergence criterion is met.Finally,a multi-degree-of-freedom(MDOF)structure equipped with a magnetorheological(MR)damper is employed to simulate the non-linear components.Based on the approach,the hysteretic force and mass distribution are identified successfully with partially known noise-polluted acceleration measurements.

vibration and wave;nonlinear restoring force identification;mass identification;extended Kalman filter;Chebyshev polynomial model;partial acceleration measurement

TH212；TH213.3 文獻(xiàn)標(biāo)示碼：A

10.3969/j.issn.1006-1355.2017.05.030

1006-1355(2017)05-0142-07

2017-03-14

國家自然科學(xué)基金資助項目(50978092)

王云(1992－)，男，安徽省舒城縣人，碩士生，主要研究方向為非線性恢復(fù)力識別。

E-mail:wangdasi4144@qq.com

許斌(1972－)，男，博士，教授，博士生導(dǎo)師。

E-mail:binxu@hqu.edu.cn