基于靈敏度分析的結(jié)構(gòu)損傷識別方法研究

余潔歆,林 偉,陳 欣,謝本飛

(1.福建江夏學(xué)院 工程學(xué)院,福建 福州 350108;2.福州大學(xué) 土木工程學(xué)院,福建 福州 350108;3.福州桂武置業(yè)有限公司,福建 福州 350011)

隨著我國經(jīng)濟持續(xù)穩(wěn)定增長疊加城鎮(zhèn)化的高速推進(jìn),土木行業(yè)實現(xiàn)了快速發(fā)展。大量重大工程結(jié)構(gòu),如高層建筑、大跨橋梁、大型水利工程等不斷涌現(xiàn)。在其漫長的服役期內(nèi),由于環(huán)境侵蝕、日常超載、材料性能退化等原因可能導(dǎo)致結(jié)構(gòu)性能退化,或受到自然災(zāi)害如地震、臺風(fēng)等侵襲,導(dǎo)致不同程度的損傷、區(qū)域功能癱瘓,甚至倒塌。如2001年宜賓南門大橋因吊桿斷裂導(dǎo)致部分橋面垮塌[1];2004年法國巴黎戴高樂機場候機樓屋頂局部坍塌[2];2007年美國明尼蘇達(dá)州一座跨越密西西比河的I-35W大橋由于節(jié)點板破壞導(dǎo)致整橋坍塌,造成13人死亡和145人受傷[3]。

我國基礎(chǔ)設(shè)施體量龐大,隨著其服役時間的增加,逐漸進(jìn)入病害集中暴露期,迎來養(yǎng)護(hù)高峰,結(jié)構(gòu)維護(hù)管理費用的壓力愈發(fā)嚴(yán)重。為此,需利用健康監(jiān)測手段盡早發(fā)現(xiàn)隱患,保障工程結(jié)構(gòu)安全并實現(xiàn)預(yù)防性維護(hù)管理。

盡管不少結(jié)構(gòu)已經(jīng)安裝了健康監(jiān)測系統(tǒng),如日本明石海峽大橋、美國金門大橋,我國青馬大橋、虎門大橋、廣州電視塔等[4-10]。由于監(jiān)測數(shù)據(jù)與結(jié)構(gòu)損傷相關(guān)性弱,又缺乏對海量數(shù)據(jù)的有效處理方法,導(dǎo)致捕捉局部損傷猶如大海撈針。利用觀測數(shù)據(jù)對所監(jiān)測結(jié)構(gòu)進(jìn)行分析和反演,進(jìn)而定位與量化損傷,實則為一種逆分析理論。但工程結(jié)構(gòu)復(fù)雜、體型龐大,難以獲取完備信息,使得求解方程容易出現(xiàn)病態(tài),其準(zhǔn)確性受限于海量的冗余信息?，F(xiàn)有技術(shù)在實現(xiàn)結(jié)構(gòu)安全評價和預(yù)警方面仍存在眾多困難,如何從海量監(jiān)測數(shù)據(jù)分析結(jié)構(gòu)損傷性能成為重點要研究的問題[11-12]。

模態(tài)參數(shù)作為結(jié)構(gòu)的自身固有屬性,在結(jié)構(gòu)發(fā)生損傷時會隨之改變。為探求損傷對模態(tài)參數(shù)的影響,尋求能夠準(zhǔn)確評估結(jié)構(gòu)健康狀況的方法,本文從模態(tài)參數(shù)的靈敏度出發(fā),根據(jù)結(jié)構(gòu)模態(tài)參數(shù)與設(shè)計參數(shù)之間的近似關(guān)系,推導(dǎo)了結(jié)構(gòu)的一階靈敏度分析方程,以此建立基于靈敏度分析的損傷識別程序。對5層框架結(jié)構(gòu)進(jìn)行數(shù)值模擬,揭示模態(tài)參數(shù)的敏感性。通過選取合適的模態(tài)參數(shù),對結(jié)構(gòu)損傷位置和損傷程度進(jìn)行識別,驗證了該方法的有效性和準(zhǔn)確性。

1 基于靈敏度的有限元模型修正方法

1.1 靈敏度分析理論

有限元模型修正方法的基本思想是以結(jié)構(gòu)有限元模型和試驗?zāi)Ｐ椭g的差值為目標(biāo)函數(shù),通過修正設(shè)計參數(shù),使修正后的模型與試驗?zāi)Ｐ偷恼w信息趨于近似,廣泛應(yīng)用于各工程領(lǐng)域[13-14]。靈敏度理論是通過模態(tài)參數(shù)與損傷指標(biāo)之間的一階導(dǎo)數(shù),構(gòu)造兩者間的近似線性關(guān)系。

假定r為初始有限元模型設(shè)計參數(shù),fi為結(jié)構(gòu)模態(tài)參數(shù),且為r的函數(shù),當(dāng)參數(shù)r發(fā)生微小變動時,第i階模態(tài)參數(shù)f對r的一階泰勒展開式為:

(1)

式(1)可改寫為:

△f=S|m×n△r

(2)

式(2)中,△r為設(shè)計參數(shù)修正值向量;△f為模態(tài)參數(shù)的殘差向量,是受損結(jié)構(gòu)模態(tài)參數(shù)與初始有限元模型模態(tài)參數(shù)的差值。通過最小化殘差,實現(xiàn)損傷指標(biāo)反演。由于兩者間為近似關(guān)系,為確定損傷指標(biāo)的準(zhǔn)確幅值,需經(jīng)歷反復(fù)迭代。S為靈敏度矩陣,m表示所取特征參數(shù)的階數(shù),n表示待修正設(shè)計參數(shù)的個數(shù)。求解靈敏度矩陣S即是求特征參數(shù)對設(shè)計參數(shù)的導(dǎo)數(shù),如式3所示。

(3)

1.2 特征值的靈敏度

無阻尼結(jié)構(gòu)的第i階特征方程為:

(K-Mλi)Φi=0

(4)

式(4)中,K為結(jié)構(gòu)剛度矩陣,M為結(jié)構(gòu)的質(zhì)量矩陣,Φi為第i階特征向量,λi為第i階特征值。

對式(4)求關(guān)于設(shè)計參數(shù)rj的導(dǎo)數(shù):

(5)

(6)

結(jié)構(gòu)的剛度矩陣和質(zhì)量矩陣均為對稱矩陣,式(4)轉(zhuǎn)置可得:

(7)

因此式(6)中第三項為零,可簡化為:

(8)

采用正則振型對結(jié)構(gòu)質(zhì)量矩陣和剛度矩陣進(jìn)行歸一化,可得:

(9)

(10)

將式(9)代入式(8)可得第i階特征值關(guān)于設(shè)計參數(shù)r的靈敏度表達(dá)式:

(11)

1.3 特征向量的靈敏度

(12)

(13)

將式(9)和式(10)代入式(13)可得:

(14)

根據(jù)式(14)可得:

(15)

(16)

式(16)中的各項均為標(biāo)量,其中第一項表達(dá)式:

(17)

將式(17)代入式(16),可得:

(18)

將式(12)代入式(18),可得:

(19)

將式(9)和式(10)代入式(19),可以得到當(dāng)i=j時的系數(shù)為:

(20)

整理可得特征向量關(guān)于設(shè)計參數(shù)的靈敏度表達(dá)式:

(21)

2 損傷識別的程序?qū)崿F(xiàn)

使用ANSYS和MATLAB軟件編寫基于靈敏度分析的損傷識別程序,具體流程如圖1所示。

圖1 基于靈敏度分析的損傷識別流程圖

選取楊氏模量作為待修正參數(shù)r,通過降低某指定單元的楊氏模量模擬損傷單元。編寫特征值和特征向量對待修正參數(shù)r的一階偏導(dǎo)數(shù)程序,組成靈敏度矩陣Sk。利用結(jié)構(gòu)動力特性進(jìn)行待修正參數(shù)的修正,最小化待修正模型和損傷模型之間的差值,直至優(yōu)化過程收斂。采用剛度降低因子SRF作為識別輸出結(jié)果,表示修正參數(shù)的變化量△r與初始值r的比值,以確定損傷位置、量化損傷程度。

3 結(jié)構(gòu)損傷識別數(shù)值仿真

3.1 有限元模型

對5層空間框架結(jié)構(gòu)進(jìn)行數(shù)值研究,框架結(jié)構(gòu)如圖2所示。整體由24個節(jié)點組成,其中最底部4個節(jié)點為固定端,每個節(jié)點有6個自由度,共120個自由度。

圖2 5層框架結(jié)構(gòu)模型

圖2中編號為單元編號,共40個,單元特性:截面積為0.25m2,楊氏模量E為3.5×104MPa,質(zhì)量密度為2500kg/m3。利用MATLAB軟件編寫各節(jié)點的坐標(biāo)信息和單元信息,組成5層框架結(jié)構(gòu)的剛度矩陣和質(zhì)量矩陣,通過模態(tài)分析得到結(jié)構(gòu)頻率和振型。

3.2 損傷程度對特征參數(shù)的影響

為研究特征參數(shù)對損傷程度的敏感性,分別設(shè)置二層柱單元E25和二層梁單元E5(如圖2“D1”和“D2”)為損傷單元,將損傷程度由10%逐級遞增到90%來對單個單元進(jìn)行損傷模擬,觀察隨著損傷程度的提高,結(jié)構(gòu)低階特征參數(shù)的變化。變化情況采用特征值的相對變化率Diff和模態(tài)置信因子MAC(Modal Assurance Criterion)體現(xiàn),定義如式(22)和(23)所示。

(22)

(23)

當(dāng)MAC數(shù)值為1時表明兩組模態(tài)振型完全相關(guān);數(shù)值為0時表明兩組模態(tài)振型完全無關(guān)。

D1損傷前后的結(jié)構(gòu)特征參數(shù)的變化,限于篇幅,僅列舉損傷程度為50%和90%的工況下,結(jié)構(gòu)前7階模態(tài)和前20階均值分析結(jié)果如表1所示。

表1 不同損傷程度的結(jié)構(gòu)特征參數(shù)

由表1中可以看出,與未損傷狀態(tài)相比,損傷結(jié)構(gòu)的特征值和MAC明顯降低。當(dāng)損傷程度為50%時,引起的模態(tài)參數(shù)變化幅度不大,前20階Diff和MAC均值分別為2.405%和83.741%,說明局部小損傷對整體結(jié)構(gòu)的特征值和振型影響較小。當(dāng)損傷程度上升到90%時,特征參數(shù)變化幅度明顯增大。其中第2特征值受損傷影響較大,特征值相對變化率為14.580%,第6和7階特征向量受損傷影響較大,MAC值分別為88.339%和84.881%。結(jié)構(gòu)前20階Diff和MAC均值分別為9.425%和51.010%,振型變化較大。

損傷程度從10%遞增到90%的9種工況前7階特征值和特征向量的變化情況,D1損傷情況下,如圖3和圖4所示。

圖3 D1工況特征值隨損傷程度的變化

圖4 D1工況MAC隨損傷程度的變化

由圖3和圖4可以看出結(jié)構(gòu)特征參數(shù)的變化與損傷程度成正相關(guān)。當(dāng)損傷程度小于50%的工況下,前7階特征值相對變化率均小于4%,MAC最小值為99.393%,表明該程度損傷引起的特征參數(shù)變化幅度較小。隨著損傷程度的上升,特征參數(shù)變化幅度明顯增加。當(dāng)損傷程度為90%時,前6階特征值變化率均大于5%,各階特征向量產(chǎn)生較大變化,D2梁損傷單元的9種不同損傷程度工況下的特征參數(shù)變化,如圖5和圖6所示。

圖5 D2工況特征值隨損傷程度的變化

圖6 D2工況MAC隨損傷程度的變化

由圖5和圖6可以看出柱單元和梁單元受損的變化趨勢基本一致。從損傷引起的特征值變化幅度來看,當(dāng)D2損傷程度為50%的工況下,第1階特征值相對變化率最大為6.860%,其余階次特征值相對變化率均小于3%。當(dāng)損傷程度同為90%工況下,梁損傷單元最大特征值變化率為15.228%,與柱損傷引起的最大變化率相近。但其最小MAC值為98.276%,而柱單元最小MAC值為84.881%,說明梁損傷單元產(chǎn)生的結(jié)構(gòu)模態(tài)數(shù)據(jù)變化更為微小,容易被噪聲掩蓋。

3.3 損傷位置對特征參數(shù)的影響

為探究結(jié)構(gòu)特征參數(shù)對損傷位置的敏感性,保持單元受損程度不變,設(shè)置每種工況的損傷位置由結(jié)構(gòu)第一層依次上移到頂部第五層。受損單元類型分為柱單元和梁單元兩種情況討論,例如首層受損單元如圖2“D3”和“D4”,每個工況中構(gòu)件的損傷程度均為50%。特征參數(shù)隨損傷位置的變化情況,如圖7所示。

圖7 D3和D4工況特征值隨損傷位置的變化

當(dāng)損傷單元為柱時,結(jié)構(gòu)特征值的變化幅度隨著損傷樓層的上移而逐漸降低。損傷柱單元位于首層時,前20階特征值相對變化率均值為3.03%,而損傷柱單元位于5層時,相對變化率均值下降到1.72%。相較于底層柱,當(dāng)損傷單元位于高樓層時更難被發(fā)現(xiàn)。當(dāng)受損單元為梁時,結(jié)構(gòu)特征值沒有隨著損傷位置的移動發(fā)生明顯的改變,特征值變化幅度都較小,前20階特征值相對變化率均值在1.65%以下。由此可見梁單元對結(jié)構(gòu)低階動力特性的貢獻(xiàn)小于柱單元,在僅使用低階特征值時將難以被準(zhǔn)確識別。

通過損傷程度和位置對模態(tài)參數(shù)的影響可以發(fā)現(xiàn),損傷實質(zhì)是局部現(xiàn)象,對整體結(jié)構(gòu)模態(tài)數(shù)據(jù)產(chǎn)生的變化可能是微小的。而各階模態(tài)對損傷的敏感程度不同,在構(gòu)建目標(biāo)函數(shù)和靈敏度矩陣時宜選擇對損傷敏感的模態(tài),以提高計算效率。

3.4 損傷識別數(shù)值分析

以D2工況作為單一局部損傷識別,通過所提出的方法和程序識別損傷位置,量化損傷程度。首先僅以頻率作為損傷檢測指標(biāo),選取前8階特征值作為方案1進(jìn)行模型修正,方案1的損傷識別結(jié)果,如圖8所示。

圖8 單一損傷工況方案1的損傷識別結(jié)果

由圖8的識別結(jié)果可以看出,識別結(jié)果出現(xiàn)了誤判,將同層對稱梁單元E6識別為損傷單元,真實損傷單元E5的SRF值僅為25.93%,遠(yuǎn)小于實際損傷50%。

頻率是容易獲得且精度較高的全局量,但不同損傷工況可能導(dǎo)致相同的頻率改變,因此出現(xiàn)誤判。并且在損傷程度為50%的工況下,各階特征值變化幅度較小,前8階特征值變化率均值僅為1.858%。當(dāng)損傷引起特征值改變不明顯的情況下,僅利用特征值進(jìn)行損傷識別是不夠的。模態(tài)振型具有結(jié)構(gòu)振動的空間特性,可以提供結(jié)構(gòu)局部信息,因此利用頻率與振型共同作為損傷檢測指標(biāo)進(jìn)行識別。

在原特征值截取方案1的基礎(chǔ)上,增加前4階模態(tài)振型,形成方案2,構(gòu)建目標(biāo)函數(shù)和靈敏度矩陣?？紤]各階模態(tài)對損傷的敏感程度不同,通過對比損傷前后模態(tài)分析,發(fā)現(xiàn)節(jié)點扭轉(zhuǎn)自由度數(shù)值較小,為避免冗余信息干擾,產(chǎn)生病態(tài)反演及收斂性差等問題,在構(gòu)建靈敏度矩陣時從全部自由度中刪除貢獻(xiàn)較小的扭轉(zhuǎn)自由度。采用頻率與振型進(jìn)行模型修正,經(jīng)過77次迭代后的單一損傷工況方案2的損傷識別結(jié)果,如圖9所示。

圖9 單一損傷工況方案2的損傷識別結(jié)果

相比方案1,方案2的識別結(jié)果準(zhǔn)確定位了實際損傷單元,避免了誤判。E5單元相應(yīng)的SRF值為43.46%,識別精度高,收斂速度快,驗證了該方法在單一局部損傷工況下的可行性和準(zhǔn)確度。

以D3和D4作為損傷單元進(jìn)行多損傷工況識別,分別設(shè)置E21單元損傷程度為40%,E4單元損傷程度為50%,方案2采用頻率與振型進(jìn)行模型修正的損傷識別結(jié)果,如圖10所示。

圖10 多損傷工況方案2的損傷識別結(jié)果

通過200次迭代,準(zhǔn)確定位到了損傷單元,且E21和E4相應(yīng)的SRF值為38.97%與46.44%,與實際損傷值40%和50%相對誤差僅為2.58%和7.12%,其余未損傷單元SRF值均未達(dá)到10%,抗干擾強,損傷識別精度滿足多損傷工況下的要求。

4 結(jié)論

文中推導(dǎo)了結(jié)構(gòu)特征參數(shù)的靈敏度表達(dá)式,以此建立基于靈敏度分析的損傷識別程序。以5層框架結(jié)構(gòu)為數(shù)值模型,通過構(gòu)建不同損傷指標(biāo)方案反演模型,具體結(jié)論如下:

(1)結(jié)構(gòu)特征值和特征向量的變化與損傷程度成正相關(guān),若單元損傷程度較低,則引起的參數(shù)變化幅度較小。相較于底層柱,位于高樓層的柱損傷引起的特征值變化幅度更小。由于各階模態(tài)對損傷的敏感程度不同,識別過程涉及大量未知參數(shù),為減少數(shù)據(jù)量存儲,避免病態(tài)辨識問題,在構(gòu)建目標(biāo)函數(shù)和靈敏度矩陣時,需通過損傷前后模態(tài)差分析找出對損傷敏感的模態(tài),避免冗余信息干擾,產(chǎn)生病態(tài)反演及收斂性差等問題。

(2)從識別結(jié)果可以看出,本文提出的基于靈敏度分析的損傷識別方法對于結(jié)構(gòu)局部損傷識別是可行的,且滿足精度要求。損傷所導(dǎo)致自振頻率和振型的改變,可以反推結(jié)構(gòu)損傷位置和大小。在局部損傷引起特征值變化不明顯的情況下,僅利用低階特征值進(jìn)行損傷識別可能引起誤判。無論是單損傷或多損傷的工況下,同時利用頻率與振型進(jìn)行損傷識別,可以更好的反演結(jié)構(gòu)損傷信息,提升識別的準(zhǔn)確性。