雷 鷹，周 歡

（廈門大學 建筑與土木工程學院，福建 廈門 361005）

近年來，基于結(jié)構(gòu)振動數(shù)據(jù)的系統(tǒng)辨識和損傷識別方法的研究受到了廣泛的關(guān)注，產(chǎn)生了許多識別理論與方法［1－6］。當災害發(fā)生時，通過觀測的數(shù)據(jù)對結(jié)構(gòu)的損傷進行實時的診斷顯得越來越重要。目前，在時域范圍內(nèi)，對因結(jié)構(gòu)損傷產(chǎn)生的結(jié)構(gòu)參數(shù)改變的實時識別主要通過在最小二乘法或者卡爾曼濾波方法中引入一個常數(shù)遺忘因子［7－8］，但是用這種方法識別的結(jié)果跟遺忘因子大小有關(guān)，遺忘因子越小，追蹤能力比較好，但是對噪聲比較敏感，若遺忘因子較大，接近于1，識別結(jié)果雖然對噪聲不很敏感，但是追蹤能力又不好。Yang和Lin等在文獻［9－12］中提出了一個基于最小二乘算法的自適應追蹤方法來識別多層剪切框架的損傷情況，但是需要對結(jié)構(gòu)的速度，位移和加速度進行全觀測，在實際的工程中很難實現(xiàn)。本文作者提出了在部分觀測結(jié)構(gòu)響應的情況下，基于擴展卡爾曼濾波／預測對剪切框架和大型結(jié)構(gòu)的損傷診斷的方法［13－15］，但這些方法不能應用于對結(jié)構(gòu)損傷的實時診斷。Yang等［16］又提出了基于擴展卡爾曼濾波算法的自適應追蹤損傷的方法，只需對加速度部分觀測，用一個遺忘因子矩陣代替單個的遺忘因子來進行損傷識別。但自適應因子矩陣存在較繁瑣與多解等問題。

本文研究一種基于擴展卡爾曼濾波的結(jié)構(gòu)損傷實時診斷算法，在外荷載已知的情況下，該方法只需要對加速度進行部分的觀測，就能夠達到實時地追蹤結(jié)構(gòu)參數(shù)以及診斷其損傷。通過對比觀測加速度與預測加速度，可以判斷出結(jié)構(gòu)損傷的時間以及可能發(fā)生損傷的位置，再通過一個優(yōu)化函數(shù)，對可能損傷的剛度進行優(yōu)化，以減少兩個加速度之間的誤差，從而達到對損傷的位置以及程度的識別。

1 損傷診斷算法

對于m個自由度的結(jié)構(gòu)的動力方程可以表示為：

其中：x（t）是 m維的位移向量。M表示質(zhì)量矩陣，F(xiàn)c［，θ］和Fs［x，θ］分別為耗散力和恢復力向量。f表示外部激勵，B為力的位置矩陣。

本文提出的損傷實時診斷算法共分為三個階段。

1．1 第一階段

第一階段是識別結(jié)構(gòu)的未知參數(shù)。由于擴展卡爾曼濾波算法需要一定的時間才能收斂，因此第一階段是算法本身的收斂階段，可以認為結(jié)構(gòu)的參數(shù)是不變的，即損傷不會發(fā)生在此階段。

引入一個2m＋n維的向量：

其中：x是位移向量是速度向量，θT是n個未知的結(jié)構(gòu)參數(shù)。動力方程（1）可以被寫成以下形式，

通常情況下，結(jié)構(gòu)的加速度只能部分觀測到。因此擴展卡爾曼濾波中的觀測方程可以表示為

其中：h（Xk＋1，fk＋1）＝DM－1［－Fcx，θ）－Fs（x，θ）＋Bf］，yk＋1是觀測向量，D是和加速度傳感器放置位置有關(guān)的矩陣。vk＋1為測量噪聲。

通過擴展卡爾曼濾波，可以得到t＝（k＋1）×Δt時刻的擴展狀態(tài)向量：

其中＾k＋1k是 Xk＋1在 t＝k×Δt時刻的估計值，Kk是卡爾曼增益矩陣。

1．2 第二階段

在第一階段結(jié)構(gòu)的狀態(tài)和物理參數(shù)收斂后，就是判斷結(jié)構(gòu)剛度突變的時刻，以及可能發(fā)生損傷的位置。

結(jié)構(gòu)某一層的損傷，會對安置在相鄰兩層的加速度響應產(chǎn)生瞬時影響，但對較遠的加速度響應，影響則晚一點。如圖1所示，數(shù)值模擬n層剪切框架，受到一個頂部白噪聲的荷載作用結(jié)構(gòu)第i層剛度在第2．5 s突然發(fā)生了損傷。結(jié)構(gòu)的每層的剛度和質(zhì)量都相同，為ki＝12 kN／m，mi＝60 kg（i＝1，2，…，n），阻尼為黏滯阻尼 ci＝200 Ns／m。圖2（a）～2（d）是數(shù)值模擬得到的損傷附近層的樓層加速度響應在損傷前后時程變化與結(jié)構(gòu)不發(fā)生損傷情況的樓層加速度響應時程的對比。

圖1 n層剪切框架Fig．1 An n－story shear building

圖2 損傷位置附近樓層加速度變化Fig．2 Deviations of floor acceleration responses near locations of structural damage

可以看出，在第2．5 s的突變只影響到了和，而和則與無損時相比沒有明顯的變化，因為i層的剛度損傷需要一定的時間才能對其他層產(chǎn)生影響。對于多層剪切框架，加速度傳感器從第一層開始隔層放置，這樣就可以確保每一層的損傷，都能被觀測到。

第t＝（k＋1）×Δt時刻，加速度的誤差可以表示為：

yi，k＋1表示（k＋1）×Δt時刻的觀測加速度是t＝kΔt時刻預測的加速度，由下式算出：

結(jié)構(gòu)的損傷不能發(fā)生在結(jié)構(gòu)參數(shù)識別收斂之前，否則該算法不能對損傷情況進行診斷。定義時刻t＝（k＋1）Δt的加速度誤差平方和為：

其中：i表示加速度傳感器放置的樓層位置。m為放置的傳感器個數(shù)。

δ表示穩(wěn)定后得到的加速度差平方和的最大值。由圖2（a）～（d）可知，結(jié)構(gòu)的第 i和 i＋1層可能發(fā)生了損傷。

1．3 第三階段

第三階段通過優(yōu)化處理，來識別出結(jié)構(gòu)剛度損傷的具體位置和損傷的程度。

假設：

Δθi，Δθi＋1為待識別的結(jié)構(gòu)剛度參數(shù)的變化。利用MATLAB中的優(yōu)化函數(shù)FMINCON，該函數(shù)的作用是使目標函數(shù)取得最小值，將修正過的θ′＝θ－Δθ代入到擴展卡爾曼濾波算法中，循環(huán)之前的計算過程。

目標函數(shù)定義成以下的形式：

其中：j表示可能損傷的單元，p為可能損傷單元的總數(shù)。

約束函數(shù)定義為：

以上的優(yōu)化問題需要對每個可能損傷的單元進行單獨的優(yōu)化。優(yōu)化的過程是不斷代入可能的剛度損傷的位置和程度，使得目標函數(shù)最小，因為只有真實的剛度損傷值才能符合觀測值與預測值相等，這樣就可以得到準確的損傷的位置和損傷的程度。然后再實時的更新結(jié)構(gòu)的剛度參數(shù)，以便能對結(jié)構(gòu)的參數(shù)進行連續(xù)的識別追蹤。

地方自媒體已然成為縣級媒體建設中的活躍者，在兩微一端平臺發(fā)展迅速。以微信公眾號為例，某些縣甚至出現(xiàn)了百個城市的地方自媒體內(nèi)容矩陣公司，地方自媒體的高速發(fā)展獲取了相當一部分內(nèi)容流量和本地商機?；ヂ?lián)網(wǎng)公司最近幾年普遍嘗試下沉到四五線城市，下沉效應帶來用戶數(shù)量的快速上升，這將無可避免地對傳統(tǒng)媒體格局造成巨大沖擊。2014年以來的微博用戶月活超過4億，2018年快手用戶日活過億，這些用戶很大程度上來自下沉縣級市場的小鎮(zhèn)青年。鮮活的案例與龐大的數(shù)字無不昭示著一批蓄勢待發(fā)的潛力受眾，勢必在縣級融媒體建設進程中扮演著至關(guān)重要的角色。

2 數(shù)值算例

本文采用的是一個四層的剪切框架模型，如圖3。

圖3 四層剪切框架及部分加速度觀測Fig．3 A four-storey shear building with partial accleration observations

結(jié)構(gòu)的各個參數(shù)分別為：m1＝3 452 kg，m2＝m3＝2 652 kg，m4＝1 810 kg，k1＝k2＝k3＝k4＝67 900 kN／m，c1＝c2＝c3＝c4＝54．3 kN·s／m。結(jié)構(gòu)的質(zhì)量已知，加速度和外荷載可以通過觀測得到，加速度傳感器放置在結(jié)構(gòu)的第一層和第三層。所以需要識別的未知量為以下幾個量：X＝［x1，…，x4，x1，…，，k1，…，k4，c1，…，c4］T。通過上面介紹的擴展卡爾曼濾波的算法和損傷追蹤監(jiān)測的方法，就可以識別出結(jié)構(gòu)的剛度和阻尼，以及剛度因損傷產(chǎn)生的突變。

假設，結(jié)構(gòu)第一層剛度在第2．5 s發(fā)生突變，此時k1從67 900 kN／m降到了40 000 kN／m。用本文提出的方法對結(jié)構(gòu)損傷進行實時診斷，假設觀測加速度受到10%的白噪聲影響。數(shù)值模擬的過程及結(jié)果如下。

第一階段：用擴展卡爾曼濾波識別系統(tǒng)的時不變參數(shù)。假設算法在1 s之后趨于穩(wěn)定。

圖4 觀測與計算估計加速度誤差平方和Fig．4 Sum of quadratic errors between observed and predicted accelerations

第二階段：判斷損傷發(fā)生的時間以及可能的位置。

圖5 突變時各觀測層加速度誤差平方Fig．5 Quadratic errors of each observed and predicted acceleration at the time of abrupt change

通過對圖4的分析，可以發(fā)現(xiàn)，第2．5 s結(jié)構(gòu)的加速度誤差平方和發(fā)生突變，遠遠超過了1 s穩(wěn)定之后的最大值，說明結(jié)構(gòu)的某一位置發(fā)生了某種損傷。圖5則說明突變主要是由第一層加速度誤差產(chǎn)生的，所以可以判斷出可能發(fā)生損傷的位置為第一層和第二層。

第三階段：確定損傷的程度以及準確的位置。

根據(jù)目標函數(shù)可以確定損傷的位置和大小。圖6～圖9是各層剛度的收斂圖，據(jù)此可以得出本文中提出的基于擴展卡爾曼濾波的結(jié)構(gòu)損傷實時診斷的算法，可以準確的判斷出剛度突變的時間，并且能夠計算出突變的位置和突變的程度。

圖6 第一樓層剛度隨時間變化的識別Fig．6 Identification of time-variant k1

圖7 第二樓層剛度隨時間變化的識別Fig．7 Identification of time-variant k2

圖8 第三樓層剛度隨時間變化的識別Fig．8 Identification of time-variant k3

圖9 第四樓層剛度隨時間變化的識別Fig．9 Identification of time-variant k4

3 實驗驗證

實驗驗證采用的是由臺灣大學的Weng等［17］在臺灣地震研究中心做的三層鋼框架結(jié)構(gòu)損傷實驗數(shù)據(jù)。如圖10所示，該結(jié)構(gòu)底面積為3 m×2 m，層高3 m，梁柱節(jié)點為螺栓連接。一個鎖定系統(tǒng)用V形支撐安置在第一層框架上，在外界底部激勵作用下該鎖定結(jié)構(gòu)可以在任何瞬間釋放，用來模擬樓層剛度的突然變化。整個結(jié)構(gòu)放置在一個大型的振動臺上，其可以提供底部作用力。

地震作用下，結(jié)構(gòu)的運動方程可以寫成：

圖10 第一層帶有鎖定系統(tǒng)的剪切框架Fig．10 Shear frame with lock-up system in the 1st floor

圖11 各時刻加速度誤差平方和Fig．11 Sum of quadratic errors between observed and predicted accelerations

其中：m1＝m2＝m3＝6 000 kg，無鎖定系統(tǒng)的框架自振頻率為1．07 Hz，3．17 Hz和4．97 Hz，相對應的阻尼比為3．85%，0．01%和0．83%，無鎖定結(jié)構(gòu)的待識別樓層剛度為 1 568 kN／m（第一層），1 531 kN／m（第二層），1 616 kN／m（第三層）。對于加強剛度的結(jié)構(gòu)（在第一樓層安裝鎖定系統(tǒng)），樓層剛度為3 137 kN／m（第一層），1 515 kN／m（第二層），1 568 kN／m（第三層）。振動臺提供峰值為0．15 g的地面激勵。在t＝16．35 s將鎖定結(jié)構(gòu)釋放即該剛性元件和第一樓層斷開連接，用來模擬第一樓層突然出現(xiàn)剛度損傷。對結(jié)構(gòu)的第一層和第三層的加速度進行觀測，計算結(jié)果如圖11～圖15所示。

圖12 突變時各觀測層加速度誤差平方Fig．12 Quadratic errors of each observed and predicted acceleration at the time of abrupt change

圖13 第一樓層剛度隨時間變化的識別Fig．13 Identification of time-variant k1

圖14 第二樓層剛度隨時間變化的識別Fig．14 Identification of time-variant k2

圖15 第三樓層剛度隨時間變化的識別Fig．15 Identification of time-variant k3

通過以上這些圖可得，結(jié)構(gòu)的第一層剛度在16．35 s發(fā)生突變，突變程度與真實情況比較符合。綜上可得本文提出的基于擴展卡爾曼濾波的實時損傷診斷的算法，可以在部分觀測的基礎(chǔ)上較準確地判斷出突變時間，同時計算出突變的位置和突變的程度，比較好地滿足實驗結(jié)果。

4 結(jié) 論

本文提出了一個基于擴展卡爾曼濾波的結(jié)構(gòu)損傷實時診斷的方法用來識別結(jié)構(gòu)的參數(shù)以及診斷其損傷。結(jié)構(gòu)的損傷首先會對其相鄰兩層的加速度產(chǎn)生影響，基于該加速度的變化，文中提出一個優(yōu)化函數(shù)，通過對結(jié)構(gòu)的剛度參數(shù)進行優(yōu)化，來使得預測加速度與觀測加速度之間的誤差最小，進而實現(xiàn)對損傷的追蹤。通過對比已有文獻提出的算法，可以看出在對加速度部分觀測的基礎(chǔ)上，本文中提出的算法，對剛度進行實時追蹤，原理清晰明了。數(shù)值模擬和實驗的結(jié)果表明提出的方法具有較好的準確性和穩(wěn)定性來識別損傷。但該算法對于相鄰兩層的剛度同時損傷、對其他較復雜線性結(jié)構(gòu)、非線性結(jié)構(gòu)等情況的實時損傷識別，都還有待于后續(xù)更深入的研究。

［1］Bisht SS，Singh M P．Detecting sudden changes in stiffness using high-pass filters［J］．Structural Control＆ Health Monitoring，2012，19（3）：319－331．

［2］ Chu S Y，Lo S C．Application of real-time adaptive identification technique on damage detection and structural health monitoring［J］． Structural Control and Health Monitoring，2009，16（3）：154－177．

［3］Huang C S，Huang S L，Su W C，et al．Identification of time-variant modal parameters using Time-Varying Autoregressive with Exogenous Input and low-order polynomial function［J］． Computer-Aided Civil and Infrastructure Engineering，2009，24（7）：470－491．

［4］Loh C H，Huang S G．On-line physical parameter estimation and damage detection with adaptive Kalman filtering approach，Structural Health Monitoring 2007：Quantification，Validation，And Implementation，Vols 1 and 2：Structural Health Monitoring（SHM），2007：2023－2030．

［5］Smyth A W，Masri S F，Kosmatopoulos E B，et al．Development of adaptive modeling techniques for non-linear hysteretic systems［J］．International Journal of Non-Linear Mechanics，2007，37（8）：1435－1451．

［6］Shi Z Y，Law SS，Xu X．Identification of linear time-varying MDOF dynamic systems from Excitation using Hilbert transform and EMD method［J］．Journal of Sound and Vibration，2009，321（3－5）：572－589．

［7］Lin JW，Betti R．On-line identification and damage detection in non-linear structural systems using a variable forgetting factor approach［J］．Earthquake Engineering and Structural Dynamics，2004，33（4）：419－444．

［8］Smyth A W，Masri S F，Kosmatopoulos E B，et al．Development of adaptive modeling techniques for non-linear hysteretic systems［J］．International Journal of Non-Linear Mechanics，2002，37（8）：1435－1451．

［9］Yang JN，Lin S．Parametric identification of structures using an adaptive tracking technique［C］．Proceedings of the 16thASCE Engineering Mechanics Conference，Seattle，2003，CD ROM．

［10］Yang JN，Lin S．On-line identification of nonlinear hysteretic structures using an adaptive tracking technique［J］．International Journal of Non-linear Mechanics，2004，39：1481－1491．

［11］ Yang J N，Pan S，Lin S．Identification and tracking of tructural parameters with unknown excitations［J］．Proceedings of the 2004 American Control Conference，30 June－2 July，Boston，MA，2004，5：4189－4194．

［12］Yang J N，Lin S．Identification of parametric variations of structures based on least square estimation and adaptive tracking technique［J］．Journal of Engineering Mechanics（ASCE）2005，131（3）：290－298．

［13］雷鷹，倪萍禾，劉朝．有限觀測絕對加速度響應下剪切框架在未知地震作用下?lián)p傷診斷［J］．振動與沖擊，2012，31（6）：96－99．LEI Ying， NI Ping-hu， LIU Chao． Structural damage detection of a shear building under unknown earthquake excitation with limited observations of absolute acceleration reponses［J］．Vibration and Shock，2012，31（6）：96－99．

［14］雷鷹，毛亦可．部分觀測下基于子結(jié)構(gòu)的大型結(jié)構(gòu)損傷診斷法［J］．工程力學，2012（7）：180－185．LEI Ying， MAO Yi-ke．A distrubuted damage detection algorithm for large size structures under limited input and output measurements［J］．Engineering Mechanics，2012（7）：180－185．

［15］Lei Y，Jiang Y Q，Xu Z Q．Structural damage detection with limited input and output measurement signals［J］．Mechanical Systems and Signal Processing，2012（28）：229－243．

［16］Yang J N，Lin S，Huang H，et al．An adaptive extended Kalman filter for structural damage identification［J］．Structural Control and Health Monitoring，2006，13：849－867．

［17］Weng J H，Loh C H．Recursive subspace identification for on-line tracking of structural modal parameter［J］．Mechanical Systems and Signal Processing，2011，25（8）：2923－2937．