王紅偉

(上海師范大學(xué) 建筑工程學(xué)院, 上海 201418)

基于NEXT/ERA時(shí)域聯(lián)合算法的模態(tài)參數(shù)識(shí)別

王紅偉

(上海師范大學(xué) 建筑工程學(xué)院, 上海 201418)

摘要:歸納總結(jié)了NEXT法和特征系統(tǒng)實(shí)現(xiàn)算法的基本思想,以及實(shí)施的基本步驟.在此基礎(chǔ)上,基于matlab編制了模態(tài)識(shí)別程序.以東海大橋監(jiān)測數(shù)據(jù)作為算例進(jìn)行數(shù)值分析,成功識(shí)別了大橋的豎向振動(dòng)動(dòng)力學(xué)參數(shù).識(shí)別結(jié)果可為橋梁動(dòng)力特性計(jì)算,抗震分析,健康監(jiān)測和損傷識(shí)別等提供基礎(chǔ)數(shù)據(jù).

關(guān)鍵詞:模態(tài)參數(shù); 環(huán)境激勵(lì); 結(jié)構(gòu)模態(tài)

近年來,橋梁的安全性、耐久性與正常使用日漸得到關(guān)注,促使橋梁健康監(jiān)測、損傷識(shí)別技術(shù)、模態(tài)參數(shù)識(shí)別技術(shù)迅速發(fā)展,并且日益成為國內(nèi)外學(xué)術(shù)界和工程界的研究熱點(diǎn).姜浩等研究了在環(huán)境激勵(lì)下對(duì)橋梁結(jié)構(gòu)進(jìn)行模態(tài)參數(shù)識(shí)別的一種時(shí)域聯(lián)合算法,通過對(duì)某一預(yù)應(yīng)力混凝土連續(xù)梁橋結(jié)構(gòu)的仿真,驗(yàn)證了方法的有效性[1].呂中亮等介紹了多點(diǎn)激勵(lì)模態(tài)參數(shù)識(shí)別方法的基礎(chǔ)上,總結(jié)了近十多年來多點(diǎn)激勵(lì)模態(tài)參數(shù)識(shí)別方法中的進(jìn)展,并分析了多點(diǎn)激勵(lì)參數(shù)識(shí)別方法在基于環(huán)境激勵(lì)下模態(tài)參數(shù)識(shí)別中的發(fā)展與應(yīng)用[2].很多大學(xué)的課題組都進(jìn)行了模態(tài)參數(shù)識(shí)別的研究[3-8],重點(diǎn)在理論方法上進(jìn)行研究.本研究在此基礎(chǔ)上,重點(diǎn)進(jìn)行程序編制和工程實(shí)例驗(yàn)證工作,可為橋梁抗震分析和健康監(jiān)測提供幫助.

1NExT法(響應(yīng)間的互相關(guān)函數(shù)理論)

NExT法的主要思想是受白噪聲環(huán)境激勵(lì)下結(jié)構(gòu)兩點(diǎn)之間響應(yīng)的互相關(guān)函數(shù)和脈沖響應(yīng)函數(shù)表達(dá)式相似,求出兩點(diǎn)之間響應(yīng)的互相關(guān)函數(shù)后,利用時(shí)域中模態(tài)識(shí)別方法進(jìn)行模態(tài)參數(shù)識(shí)別.

n個(gè)自由度線性時(shí)不變結(jié)構(gòu)的動(dòng)力特性方程可用振動(dòng)微分方程式(1)表示:

(1)

輸出信號(hào)x(t)可以利用標(biāo)準(zhǔn)變換為模態(tài)坐標(biāo)系下的形式:

(2)

假定是實(shí)模態(tài),式(2)兩邊同時(shí)乘以ΦT,則M、C、K是對(duì)角陣.則有:

(3)

由杜哈密積分,得出:

(4)

由式(2)和式(4)得到如下結(jié)果:

(5)

由式(3),當(dāng)系統(tǒng)k點(diǎn)受到外力f(t)激勵(lì)時(shí),系統(tǒng)l點(diǎn)得響應(yīng)Xlk(t)可表示為:

(6)

當(dāng)系統(tǒng)的k點(diǎn)受到單位脈沖力作用時(shí),由式(6)就得到系統(tǒng)l點(diǎn)的脈沖響應(yīng)hlk(t),可表示為:

(7)

在k點(diǎn)白噪聲激勵(lì)下,系統(tǒng)l,p兩測點(diǎn)的響應(yīng)分別為xlk(t)和xpk(t),這兩個(gè)響應(yīng)的互相關(guān)函數(shù)的表達(dá)式可以寫成:

(8)

將式(6)代入式(8)得到如下等式:

(9)

假設(shè)式中激勵(lì)f(t)為理想白噪聲,根據(jù)相關(guān)函數(shù)的定義,則有:

(10)

式中:αk為僅同激勵(lì)點(diǎn)k有關(guān)的常數(shù)項(xiàng).δ(t)是狄拉克(Dirac delta)函數(shù).

將式(10)代入式(9)中并積分,得到:

(11)

式中:β=t-τ.

利用式(4)的定義,考慮周期T,對(duì)式(11)進(jìn)行整理得:

(12)

將式(7)和式(12)進(jìn)行對(duì)照比較,可以看出,線性系統(tǒng)在白噪聲激勵(lì)下兩點(diǎn)響應(yīng)的互相關(guān)函數(shù)和脈沖激勵(lì)下的脈沖響應(yīng)的數(shù)學(xué)表達(dá)式在形式上是完全一致的.互相關(guān)函數(shù)確實(shí)可以表征為一系列復(fù)指數(shù)函數(shù)的疊加形式.在這點(diǎn)上,相關(guān)函數(shù)具有和系統(tǒng)的脈沖響應(yīng)函數(shù)同樣的性質(zhì).同時(shí)各測點(diǎn)的同階模態(tài)振型乘以同一因子時(shí),并不改變模態(tài)振型的特征.因此,互相關(guān)函數(shù)可以用來替代脈沖響應(yīng)函數(shù)進(jìn)行環(huán)境激勵(lì)下的模態(tài)參數(shù)識(shí)別.

為了使用該方法,將一個(gè)響應(yīng)信號(hào)作為參考信號(hào),該參考信號(hào)和其他響應(yīng)信號(hào)可以求得互譜密度函數(shù),對(duì)互譜密度函數(shù)進(jìn)行傅里葉逆變換就能夠得到互相關(guān)函數(shù).參考信號(hào)位置的選擇必須考慮到在那個(gè)位置所有的模態(tài)都能夠被觀測到.如果參考信號(hào)位置在某一模態(tài)的節(jié)線上,那么該模態(tài)就不能被識(shí)別出來.

2特征系統(tǒng)實(shí)現(xiàn)算法(ERA法)

特征系統(tǒng)實(shí)現(xiàn)算法的基本思想是根據(jù)系統(tǒng)的脈沖響應(yīng)函數(shù)(或自由響應(yīng))的離散時(shí)間序列,構(gòu)成廣義的Hankel矩陣,然后對(duì)它進(jìn)行奇異值分解,通過奇異值分解的結(jié)果得到系統(tǒng)的最小實(shí)現(xiàn),最后對(duì)最小實(shí)現(xiàn)的狀態(tài)矩陣進(jìn)行特征值分解,可得到系統(tǒng)動(dòng)力學(xué)參數(shù).

考慮m個(gè)輸入、n個(gè)輸出的N個(gè)自由度線性時(shí)不變系統(tǒng),在時(shí)間離散域上的狀態(tài)方程和觀測方程可分別寫成下列形式:

(13)

(14)

式中:x(k)為kΔt時(shí)刻系統(tǒng)的狀態(tài)矢量(2N維);y(k)為kΔt時(shí)刻系統(tǒng)響應(yīng)矢量(n維);p(k)為kΔt時(shí)刻系統(tǒng)激勵(lì)矢量(m維);A是系統(tǒng)矩陣(2N×2N);B是控制矩陣(2N×m);C是輸出矩陣(n×2N);其中kΔt為采樣時(shí)間間隔.

矩陣(A,B,C)反映了系統(tǒng)特征.特征系統(tǒng)實(shí)現(xiàn)算法是給定y(k)去構(gòu)成系統(tǒng)的(A,B,C),也就是要求辨識(shí)一組狀態(tài)方程的內(nèi)部描述去實(shí)現(xiàn)所要求的外部特性.但對(duì)同一輸入輸出關(guān)系,有無窮多個(gè)系統(tǒng)實(shí)現(xiàn),人們感興趣的是其中狀態(tài)空間模型階數(shù)最小的實(shí)現(xiàn),稱為最小系統(tǒng)實(shí)現(xiàn).對(duì)應(yīng)最小階系統(tǒng)實(shí)現(xiàn)的特征值集稱為固有特性參數(shù).

對(duì)式(13)和式(14)進(jìn)行z變換,得到:

由上面兩式可以得到:

(15)

系統(tǒng)傳遞函數(shù)矩陣(n×m):

(16)

或

(17)

因單位脈沖響應(yīng)函數(shù)矩陣H(k)與傳遞函數(shù)矩陣服從z變換對(duì),由式(17)可以得到下列脈沖響應(yīng)函數(shù)矩陣序列式(18),也稱Markov參數(shù).由Markov參數(shù)可組成下列(r+1)×(s+1)廣義Hankel矩陣,如式(19)所示.

(18)

(19)

或

(20)

如果系統(tǒng)是完全可控和可觀測的,則可控性矩陣Q和可觀測性矩陣P的秩都是2N.系統(tǒng)的階數(shù)為2N,則系統(tǒng)矩陣的秩為2N.理想情況下,Hankel矩陣的秩為2N.但由于噪聲干擾,由實(shí)測數(shù)據(jù)生成的H(k-1)會(huì)有秩的虧損.當(dāng)r與s的值足夠大后,矩陣的秩才保持不變.適當(dāng)?shù)剡x擇r與s的值,既保持矩陣的秩不變,又使H(k-1)尺寸最小.

由式(20),得到:

(21)

對(duì)上式進(jìn)行奇異值分解,得到:

(22)

矩陣H(0)的秩可由非零奇異值的個(gè)數(shù)確定.

(23)

式中:U2N、V2N分別由U和V前2N列構(gòu)成的矩陣.

假設(shè)存在一個(gè)矩陣H*,滿足下列關(guān)系:

(24)

式中:I2N是2N階單位矩陣.

由式(21)和式(24),得到:

上式說明:H*是H(0)的一種廣義逆矩陣,便可得到下式.

(25)

式中:0m、0n為m階、n階零矩陣;Im、In為階m、n階單位矩陣.

由式(19)和式(20)得到:

進(jìn)一步,上式可化為:

上式與式(18)比較,兩者由極其相似的形式,即h(k+1)=CAKB.

式中:

(26)

上列關(guān)系是特征系統(tǒng)實(shí)現(xiàn)算法的基本公式,系通過奇異值分解獲得最小階系統(tǒng)實(shí)現(xiàn).應(yīng)該指出,由于測量噪聲、結(jié)構(gòu)的非線性影響和計(jì)算截?cái)嗾`差影響,H(k)的秩不一定等于系統(tǒng)真實(shí)的階次,還需通過一定準(zhǔn)則剔除噪聲模態(tài).

解矩陣A的特征值問題,可得到特征值λi和相應(yīng)的特征矢量φi,有下式成立:

(27)

式中:λ=diag[λ1λ2…λ2N],Φ=diag[φ1φ2…φ2N].

考慮拉氏變換與z變換的關(guān)系:

(28)

可求下列模態(tài)參數(shù):

結(jié)構(gòu)模態(tài)與噪聲模態(tài)的鑒別.

引入模態(tài)幅值相關(guān)因子γ來評(píng)估模態(tài)置信度.在模態(tài)空間上,線性時(shí)不變系統(tǒng)的特征值和特征矢量是確定量,而噪聲模態(tài)具有隨機(jī)性.引入模態(tài)幅值相關(guān)因子γ,通過某階實(shí)測模態(tài)的幅值與真實(shí)模態(tài)的幅值在時(shí)間歷程上的相關(guān)程度鑒別結(jié)果模態(tài)和噪聲模態(tài),一般認(rèn)為,γ≥0.5即可認(rèn)為所識(shí)別的模態(tài)為結(jié)構(gòu)本身真實(shí)的模態(tài).

在t=kΔt(k=0,1,…,r-1)離散時(shí)刻,單位脈沖響應(yīng)函數(shù)矩陣序列的值為:

(29)

同樣

(30)

式中:β是由噪聲干擾產(chǎn)生的誤差矩陣.

將式(29)和式(30)進(jìn)行比較,得到:

(31)

考慮到第i階模態(tài),式(31)可寫成:

(32)

β越小,表明第i階模態(tài)具有理想線性特性,模態(tài)的純度越高.

實(shí)測的模態(tài)幅值時(shí)間歷程可根據(jù)Hankel矩陣分解得到:

(33)

利用辨識(shí)的特征值外推的模態(tài)幅值時(shí)間歷程為:

(34)

式中:[b1b2…bn]=Φ-1B.

這樣,第j階模態(tài)的幅值相關(guān)因子為:

(35)

γj在0和1之間取值.γj→1表明該階模態(tài)接近真實(shí)模態(tài);γj→0表明該階模態(tài)是噪聲模態(tài).

3NExT/ERA算法操作流程

NExT/ERA算法的應(yīng)用可以歸納為以下幾步:

(1) 選擇合適的參考點(diǎn);

(2) 計(jì)算各測點(diǎn)與參考點(diǎn)的互相關(guān)函數(shù);

(3) 選擇合適的互相關(guān)函數(shù)數(shù)據(jù)構(gòu)造Hankel矩陣H(0);

(4) 對(duì)Hankel矩陣H(0)進(jìn)行奇異值分解,計(jì)算出U2N、V2N及Σ2N;

(5) 構(gòu)造Hankel矩陣H(1),并且得到系統(tǒng)最小實(shí)現(xiàn)的矩陣;

(6) 計(jì)算矩陣A的特征值和特征向量;

(7) 剔除虛假模態(tài),最后得到結(jié)構(gòu)的模態(tài)參數(shù).

4NExT/ERA算法數(shù)值仿真算例

以東海大橋監(jiān)測數(shù)據(jù)作為算例進(jìn)行數(shù)值分析.

東海大橋是我國第一座真正意義上的跨海大橋,工程起點(diǎn)為蘆潮港客運(yùn)碼頭東側(cè)約5.6 km,靠北約1.4 km的海灘與現(xiàn)有老大堤的交接處,終點(diǎn)為浙江省嵊泗縣崎嶇列島的小城子山(進(jìn)入洋山深水港區(qū)的接線點(diǎn)),蘆潮港至小城子山線路總長度約32 km.其中大橋海上段約28 km.大橋設(shè)主通航孔橋一座,為雙塔單索面半漂浮體系疊合梁斜拉橋,主塔為倒Y型鋼筋混凝土結(jié)構(gòu),塔高150 m,主跨420 m,大橋南北走向,其從北向越東偏1.18°.

東海大橋作為上海國際航運(yùn)中心洋山深水港區(qū)重要的配套工程,為保證東海大橋交通暢通和提高大橋的維護(hù)管理水平,東海大橋上安裝了結(jié)構(gòu)健康監(jiān)測系統(tǒng).大橋健康監(jiān)測系統(tǒng)將大橋分8個(gè)區(qū)段,共計(jì)478個(gè)傳感器.其中主航道斜拉橋位于第5區(qū)段,安裝有169個(gè)傳感器,實(shí)時(shí)監(jiān)測大橋的加速度響應(yīng)、位移響應(yīng)、風(fēng)速/風(fēng)向、大氣溫度、結(jié)構(gòu)溫度、索力、結(jié)構(gòu)應(yīng)變,伸縮縫位移等(圖1).

主梁布設(shè)了14個(gè)豎向加速度傳感器和7個(gè)橫向加速度傳感器,分別布置在主跨的1/4、跨中和邊跨的跨中,共計(jì)7個(gè)斷面,每個(gè)斷面布設(shè)兩個(gè)豎向加速度傳感器和一個(gè)橫向加速度傳感器,其傳感器布設(shè)的橫向和立面圖如下:

用NexT/ERA算法識(shí)別系統(tǒng)的模態(tài)特性,獲得數(shù)據(jù)如表1所示:

表1中給出了所識(shí)別出的15階系統(tǒng)模態(tài)參數(shù),圖1至圖10是前10階模態(tài)振型圖.在計(jì)算中發(fā)現(xiàn),所取Hankel矩陣的階數(shù)越高,識(shí)別出的系統(tǒng)模態(tài)階次越多,但是,當(dāng)Hankel矩陣階數(shù)達(dá)到一定值時(shí),所識(shí)別出的模態(tài)階數(shù)不再變化,達(dá)到穩(wěn)定的狀態(tài).

本研究中所用的14個(gè)傳感器采集的都是東海大橋的豎向加速度,所以識(shí)別出的模態(tài)都是豎向彎曲模態(tài)或者扭轉(zhuǎn)模態(tài),橫向彎曲并未識(shí)別出來.

5總結(jié)

介紹了一種環(huán)境振動(dòng)模態(tài)分析的時(shí)域算法:自然環(huán)境激勵(lì)技術(shù)和特征系統(tǒng)實(shí)現(xiàn)算法結(jié)合,稱為NExT/ERA法.針對(duì)結(jié)構(gòu)振動(dòng)的離散狀態(tài)空間數(shù)學(xué)模型,逐步推演了NExT/ERA算法流程并給出了相應(yīng)的理論解釋.最后對(duì)一座雙塔單索面半漂浮體系疊合梁斜拉橋數(shù)值模型進(jìn)行仿真分析,通過模擬該結(jié)構(gòu)在平穩(wěn)隨機(jī)白噪聲激勵(lì)下的響應(yīng),應(yīng)用NExT/ERA算法識(shí)別出該斜拉橋的部分模態(tài)參數(shù),并對(duì)所得結(jié)果進(jìn)行了分析.

參考文獻(xiàn):

[1]JIANG H,GUO X D,YANG H L.Research on parameters identification of bridge structure ambient excitation[J].Journal of Vibration and Shock,2008(11):126-128.

[2]LV Z L,YANG C Q,AN P W,et al.Progress on modal parameter identification with Multipoint-excitation[J].Journal of Vibration and Shock,2011(1):197-202.

[3]YE X J.Modal parameter identification of Large-scale civil engineering structures based on ambient excitation[D].Guangzhou:South China University of Technology,2012.

[4]CHANG J.Modal parameters identification of bridge by stochastic subspace identification method[D].Shanghai:Tongji University,2006.

[5]GUO G H.Research on damage diagnosis of bridge structures based on the vibration parameters[D].Changsha:Hunan University,2001.

[6]GUO G W.Environment excitation,automatic identification of structural modal parameters and algorithm optimization[D].Chongqing:Chongqing University,2012.

[7]LUO K.Application of structural modal parameter identification based on natural excitation technique[D].Wuhan:Wuhan University of Technology,2010.

(責(zé)任編輯：顧浩然)

[8]HUI H H.Bridge structure under the incentive of modal parameter identification method research[D].Changsha:Central South University,2010.GUO X D,YANG H L.Research on parameters identification of bridge structure ambient excitation[J].Journal of Vibration and Shock,2008(11):126-128.

Modal parameter identification based on NEXT/ERA time-domain joint algorithmWANG Hongwei

(Architecture Engineering College, Shanghai Normal University, Shanghai 201418, China)

Abstract:This paper summarizes the basic idea of NEXT algorithm and features of the system, and the basic steps of implementation.On this basis,a modal identification program is written using Matlab.The monitoring data of East China Sea bridge are selected as an example of numerical analysis.Dynamic parameters of vertical vibration of the bridge are successfully identified.The recognition results can be served for the dynamic characteristics of bridge seismic analysis, and to provide the basis data for health monitoring and damage identification.

Key words:modal parameters; environmental incentive; structural modal

通信作者:王紅偉,中國上海市奉賢區(qū)海思路100號(hào),上海師范大學(xué)建筑工程學(xué)院結(jié)構(gòu)實(shí)驗(yàn)室,郵編:201418,E-mail:whw@shnu.edu.cn

基金項(xiàng)目:MTS土木結(jié)構(gòu)試驗(yàn)系統(tǒng)在橋墩低周反復(fù)試驗(yàn)中的準(zhǔn)確性和穩(wěn)定性研究(SK201420)

收稿日期:2014-02-21

中圖分類號(hào):TU 997

文獻(xiàn)標(biāo)志碼:A

文章編號(hào):1000-5137(2015)02-0146-08