基于時序分析的銹蝕RC拱損傷識別：仿真與驗證

唐啟智,辛景舟,2,周建庭,付 雷,3,張 俊

(1.重慶交通大學 土木工程學院，重慶 400074; 2.廣西交通投資集團有限公司，廣西 南寧 530022; 3.貴州橋梁建設集團有限責任公司，貴州 貴陽 550000; 4.貴州畢節高速發展有限公司，貴州 畢節 551700)

0 引 言

拱橋因其結構形式多樣、造型美觀、剛度大等特點，在我國公路、鐵路系統中得到廣泛應用[1]。隨著服役年限的增加，拱橋性能劣化在所難免，識別拱橋的損傷狀態，掌握其運營狀況，對于保障國家、社會與人民群眾的生命財產安全，具有重要的科學研究意義與工程實用價值[2]。

國內外學者針對拱橋的損傷識別開展了大量的研究。B.ALEMDAR等[3]通過運營模態分析對圬工拱橋的固有頻率和阻尼比進行了研究；D.CAPECCHI等[4]提出基于模態頻率、振型與曲率的損傷識別方法，并在拋物線拱上進行了試驗驗證；M.M.ALAMDAN等[5]提出了基于譜矩的方法，以識別澳大利亞悉尼海港大橋子結構中可能存在的結構損傷或異常；周宇等[6]提出柔度曲率梯度的損傷特征指標并應用于剛架拱橋損傷識別；項貽強等[7]通過對鋼管混凝土中承式拱橋吊桿加速度響應信號進行離散小波變換，提出基于小波總能量相對變化的損傷識別方法；姜紹飛等[8]提出一種基于模態指標和數據融合的結構損傷識別方法，應用該法定位了鋼管混凝土拱橋的損傷。可見，目前拱橋的損傷識別主要集中在頻域方法且相關研究較少。

時間序列分析是最具代表性的時域方法[9]，由于其理論基礎成熟、無需進行繁瑣的頻譜分析，已逐漸成為損傷識別領域的熱點。國內外學者對基于時間序列的結構損傷識別進行了大量的研究。M.KRISHNAN等[10]提出了一種基于時變自回歸模型和主成分分析的損傷識別方法，實現了多自由度振動結構連續損傷在線檢測；D.ALESSIO等[11]對自回歸模型和馬氏距離平方法在結構健康監測領域的應用進行了全面的敏感性分析；H.SIMON等[12]使用遺傳算法對自回歸模型系數進行了優化選擇，構造了敏感性指標，實現了風力發電機葉片早期損傷的識別；杜永峰等[13]利用自回歸模型殘差的方差與待識別工況殘差的方差之比作為損傷指標，通過數值模擬的方法對一簡支梁進行損傷識別研究；劉綱等[14]基于向量自回歸模型，以模型系數的對角線元素的馬氏距離構造損傷指標，進行了損傷定位和程度的識別。由此可見，時間序列分析在結構損傷識別領域已取得長足發展，但研究對象僅僅是一些簡單的構件，如簡支梁等。

基于拱橋損傷識別的復雜性與既有方法的不足，筆者開展了基于時間序列的RC拱橋損傷識別研究。首先，介紹了時間序列模型的建模方法，將待識別工況與完好工況的模型殘差的方差之比作為損傷特征指標；其次，構建了鋼筋銹蝕與結構剛度退化的關系，建立了試驗拱肋的數值模型，利用靜力試驗數據對數值模型進行驗證之后，得到了不同損傷工況下的加速度響應；然后，對不同損傷工況下的響應數據建立時間序列模型，并提取損傷特征指標，獲取了銹蝕率與損傷特征指標之間的定量關系，提出了針對不同鋼筋直徑的識別模型；最后，將試驗數據送入所提識別模型進行預測，驗證了所提方法。圖1是全文的具體框架，圖中直接從銹蝕與剛度退化關系的建立開始，時間序列模型的建模方法會在第1節介紹。

圖1 全文框架

1 時序分析方法

1.1 時間序列模型

根據數據的平穩性，時間序列模型可以分為非平穩時間序列模型和平穩時間序列模型[15]。非平穩時間序列模型包括確定性趨勢、隨機性趨勢模型，由于其數據相關性不高，預測結果可靠性不強，一般可通過差分或消除趨勢來使其平穩化。平穩時間序列模型根據殘差項是否存在異方差性，可分為：線性時間序列模型，包括自回歸(AR)模型、滑動平均(MA)模型、自回歸滑動平均(ARMA)模型等；非線性時間序列模型，包括自回歸條件異方差(ARCH)模型、廣義自回歸條件異方差(GARCH)模型等。實際應用過程中，根據數據特點建立相應的時間序列模型即可。

如果一個隨機變量序列中的觀測量at可以表示為：

(1)

(2)

式中：c為常數；βi≥0，i=1,2,…,p；αj≥0，j=1,2,…,q分別為方差項和殘差項的系數。

ARCH(q)模型是GARCH模型在p=0時的一個特例，兩者可以認為是在ARMA模型的基礎上，對殘差項進行建模、修正的過程。

1.2 時序模型建模方法

由于平穩時間序列模型在實際應用中更為廣泛，這里主要介紹平穩時間序列模型的建模流程：

1)數據平穩性檢驗。平穩性檢驗即單位根檢驗，筆者采用ADF(augmented dickey-fuller)法進行平穩性檢驗。ADF檢驗有3個檢驗模型為：

(3)

模型中α、β、δ為回歸系數，原假設H0為3個模型中至少存在一個單位根(δ=0)，備選假設H1為序列不存在單位根。實際檢驗過程中從模型3至模型1依次檢驗，若拒絕H0，則認為序列不存在單位根，即序列平穩。若序列非平穩，可以通過差分等方式使數據平穩化。

2)模型識別，即判斷采用何種時間序列模型，一般運用序列的自相關函數(ACF)和偏自相關函數(PACF)圖，根據各個時間序列模型的截尾和拖尾特征來進行判別：MA模型表現為自相關函數截尾，偏自相關函數拖尾。AR模型表現為自相關函數拖尾，偏自相關函數截尾。ARMA模型表現為自相關函數拖尾，偏自相關函數拖尾。

3)模型定階，即確定模型滯后的階數。筆者采用AIC(Akaike’s information criterion)準則，選取使得該準則值較小的階數作為模型階數。

4)參數估計。模型參數估計的方法有矩估計、最大似然估計和最小二乘估計。筆者采用最小二乘估計，其基本思想是對于一般的ARMA(p,q)模型，尋找使得殘差平方和最小的一組參數值，即：

(4)

式中：s=max{p,q}+1。

5)模型檢驗。當參數估計完后，需要檢驗模型是否充分描述了數據，主要檢驗殘差是否存在異方差性。筆者采用ARCH LM檢驗法進行殘差異方差性檢驗，其基本思想是對時序模型殘差的平方運行如下的回歸：

(5)

ARCH LM檢驗的H0為殘差中直到p階都沒有異方差效應，H1為殘差存在p階異方差效應。當檢驗統計量小于假定顯著性水平的臨界值時則拒絕H0。

若殘差存在異方差性，即接受H1，高階時需要建立GARCH模型，低階時則建立ARCH模型來對殘差項進行修正。

時間序列模型建模流程如圖2。建立好時間序列模型后，構造損傷特征指標，通過特征指標與結構損傷間相聯關系的建立，實現結構的損傷診斷。

圖2 基于時間序列模型的損傷特征提取

1.3 損傷特征指標建立

綜合已有研究成果，發現國內外學者建立了各種各樣的損傷特征指標，大致包括基于殘差方差的損傷特征指標、基于自回歸系數的損傷特征指標等。文獻[13]采用基于殘差方差的損傷指標，取得了一定的效果，于是筆者采用基于殘差方差(這里指的是樣本方差，有別于前文所述殘差總體方差)的損傷特征指標。

(6)

2 模型試驗與數值仿真

2.1 模型試驗概況

模型試驗采用的是鋼筋混凝土矩形肋拱，為模擬結構真實的應力水平，根據相似理論，在拱肋相應位置施加配重，模型拱基本參數與測點布置如圖3。

圖3 模型拱基本參數

模型試驗共制作了兩片拱，一片未銹蝕，另一片銹蝕，分別記為0#、1#拱。采用通電加速銹蝕的方法對模型拱進行銹蝕，按10%的銹蝕率設定銹蝕時間，待達到預定時間之后，終止銹蝕過程。

動力試驗時，使用橡膠錘自由下落，在拱頂施加沖擊激勵，并采集加速度響應信號。靜力試驗時，使用油壓千斤頂，在拱肋3L/4(測點6)處按力加載的方式，加載至拱肋破壞。試驗過程中獲取了0#拱完好狀態、加載破壞和1#拱銹蝕損傷、加載破壞4種工況的加速度響應數據以及0#、1#拱加載破壞全過程的靜力數據，試驗過程見圖4。

圖4 試驗過程

試驗完成后，對1#拱破拆，分9段(圖3注釋)稱量鋼筋銹蝕后的重量，按式(7)計算各段的銹蝕率，各節段實際銹蝕率見表1。

表1 實際銹蝕率

(7)

式中：m、mc分別為銹蝕前后鋼筋的實際重量。

由于僅進行了1#拱的銹蝕試驗，未能獲取不同銹蝕程度下的加速度響應數據。為了獲取更多的損傷工況，以建立銹蝕與損傷特征指標之間的映射關系，故采用數值仿真獲取結構在不同銹蝕率下的響應數據，并與實測靜力數據對比來驗證數值仿真結果的正確性。

鋼筋銹蝕引起的銹脹裂縫是導致結構剛度退化的主要成因，從而造成動力響應的變化。然而在有限元模型中模擬大范圍的銹脹裂縫是較困難的，因此2.2節先給出銹蝕與剛度退化之間的定量關系，然后在2.3節進行有限元模型的建立。

2.2 銹蝕與剛度退化的關系

鋼筋銹蝕對截面剛度的影響主要從以下3個方面得以體現[16]：①鋼筋銹蝕會引起鋼筋自身截面的減少，從而影響鋼筋的力學性能；②鋼筋銹蝕產生的銹脹產物會在約束混凝土作用下形成銹脹力，使混凝土產生軟化效應，進而產生銹脹裂縫；③鋼筋銹蝕會降低鋼筋與混凝土之間的黏結作用。

需要指出的是，筆者僅考慮銹蝕作用下導致的截面剛度退化。

鋼筋銹蝕后，其屈服強度、彈性模量都會有所下降，梁巖等[17]給出了鋼筋銹蝕后材料性能計算公式為：

fyc=(1-1.727ηs)fy

(8)

Esc=(1-1.018ηs)Es

(9)

式中：fyc、Esc分別為銹蝕后的屈服強度和彈性模量；fy、Es分別為未銹蝕的屈服強度和彈性模量；ηs為鋼筋銹蝕率。

鋼筋銹蝕后，其截面半徑將減小，筆者采用均勻銹蝕下的理論計算公式[18]為：

(10)

則鋼筋銹蝕深度為：

(11)

式中：r為未銹蝕情況下的鋼筋半徑。

由于鋼筋銹蝕后銹脹力的作用，混凝土性能也將有所改變，筆者采用文獻[19]所提出的鋼筋銹蝕后混凝土性能計算公式為：

f′cu,k=(0.704 9+0.295 1e-4.131 7ηs)fcu,k

(12)

(13)

式中：f′cu,k為銹蝕后混凝土立方體抗壓強度；fcu,k為混凝土立方體抗壓強度；Ecc為銹蝕后混凝土彈性模量。

截面剛度按GB50010—2010《混凝土結構設計規范》所推薦的公式進行計算，并引入考慮鋼筋與混凝土黏結性能退化的剛度退化系數：

(14)

式中：B為截面剛度；β為考慮黏結退化而引入的剛度退化系數；其余符號與GB50010—2010《混凝土結構設計規范》中相同；對于銹蝕后的鋼筋，鋼筋截面尺寸、材料性能等取銹蝕后實際值。

β取值按孫彬等[20]通過試驗研究所提出的計算公式為：

(15)

式中：x為鋼筋銹蝕深度，按式(11)計算。

2.3 有限元模型建立

依據上述定量關系，筆者采用ANSYS建立了0#與1#拱的有限元模型，以彈性模量的折減來模擬截面剛度的退化，混凝土采用SOLID65單元，鋼筋采用LINK8單元，加載塊采用SOLID45單元，共劃分為2 888個單元，如圖5。混凝土的本構關系采用多線性等向強化模型(MISO)，并關閉壓碎，且考慮拉應力釋放，鋼筋的本構關系采用雙線性等向強化模型(BISO)。材料應力-應變關系的輸入采用GB50010—2010《混凝土結構設計規范》所推薦的公式，應力-應變關系曲線如圖6。

圖5 有限元模型

圖6 材料應力-應變關系

邊界條件以面約束的形式施加，約束兩拱腳所有平動及轉動位移；配重通過加載塊的自重來施加，并根據力值大小及加載塊體積換算等效密度。

靜力分析中加載位置與模型試驗一致，并以面荷載的形式施加。求解時，按力加載，并打開自動時間步和位移收斂準則，且將收斂條件放寬至5%。0#、1#拱極限荷載的有限元計算結果分別為64.2、58.1 kN，與試驗值61、55 kN相差5.2%和5.6%。圖7是0#和1#拱在3L/4截面的荷載-位移曲線仿真與試驗值的對比。從圖7中可以看出，有限元模型和實測結果吻合程度較高，即說明筆者所建立的有限元模型以及銹蝕率-剛度退化關系是合理的。

圖7 荷載-位移曲線對比

2.4 工況劃分及損傷識別模型建立

進一步地，為建立銹蝕與損傷特征指標之間的關系，在有限元模型中模擬了多種損傷工況，損傷為結構整體損傷，具體劃分情況見表2。在拱頂施加白噪聲激勵，獲取模型拱在不同損傷程度(銹蝕率)下的加速度時程響應數據。

表2 損傷工況劃分

獲取不同損傷工況下的加速度響應后，基于Matlab軟件平臺，編程建立了完好狀態下的時間序列基準模型，其余工況按此基準模型建立。以完好狀態下測點2的加速度響應數據為例，簡要介紹時序模型的建模過程：

1)平穩性檢驗。對樣本數據進行單位根檢驗，模型3、模型2、模型1的統計量都為-28.91，遠小于各自所對應的95%顯著性水平的臨界值-3.41、-2.86、-1.94，于是接受H1，即該數據序列不存在單位根，可進行平穩時間序列建模。

2)模型識別。由圖8可知，ACF拖尾，PACF截尾，故應采用AR模型。

圖8 ACF與PACF變化

3)模型定階。回歸不同滯后階數的AR模型，提取各模型的AIC值等。由圖9看出，模型的AIC值在滯后3階后逐漸趨于穩定，且擬合優度R2>0.99。綜合各節點響應數據的AIC值、系數顯著不為零以及滯后階數不宜過多的原則，選取AR(4)作為基準時間序列模型。

圖9 AIC與R2變化

4)參數估計。確定模型階數之后，利用最小二乘估計常數項和各系數的值。

5)異方差檢驗。對AR(4)模型的殘差進行異方差檢驗，由表3可知，滯后階數在1階后P值等于0.551，即不拒絕H0，認為異方差效應不明顯，故選擇AR(4)模型。

表3 異方差性檢驗

取每種損傷工況下所有測點DI值的均值作為該損傷工況所對應的DI值，結果如圖10。

圖10 DI與ηs關系

從圖10可以看出，隨著損傷程度的增加，各工況下的DI值逐漸減少，從而可以說明結構已經發生了損傷。將不同銹蝕率下的DI值進行回歸，得到銹蝕率與DI值之間的對應關系為：

(16)

利用式(16)，根據所計算的實際DI值即可推算出結構實際的銹蝕率，即實現了結構損傷程度的識別。從式(11)中可以知道，鋼筋銹蝕深度不僅與鋼筋銹蝕率有關，還與鋼筋直徑有關。

為此，筆者提出一種適用于不同鋼筋直徑的銹蝕率識別方法。仔細觀察式(16)后可以發現，式(16)所描述的分段函數的兩轉折點所對應的銹蝕率與式(15)是一致的，且轉折點之間呈線性關系，于是只要確定了各轉折點的銹蝕率ηs與損傷特征指標DI關系，就可以寫出整個分段函數的表達式。圖11是B1/B0與DI之間的關系，由非線性擬合得到。

圖11 DI與B1/B0關系

(17)

因此可按以下步驟識別不同直徑鋼筋的銹蝕率：

1)輸入鋼筋半徑r，根據式(11)、式(15)計算兩轉折點所對應的銹蝕率ηs1、ηs2(設ηs2>ηs1)。

2)取ηs3>ηs2，按式(14)、(17)計算ηs1、ηs2、ηs3所對應的DI值。

3)根據點(1,0)以及上面求得的3個點，即可分段寫出銹蝕率與DI之間的對應關系。

4)根據實測響應數據，按圖2流程計算DI，并用第3)步求得的函數關系計算銹蝕率。

2.5 銹蝕RC拱損傷識別

本小節基于1#拱實測加速度響應數據，按2.4節所建立的DI與ηs之間的定量關系，實現了銹蝕RC拱損傷識別，并與實測銹蝕率對比，以此驗證所提方法。

由于試驗條件的限制，未能在1#拱未銹前進行加速度響應數據的采集，故以0#拱的加速度數據作為1#拱在未銹狀態下的響應數據。

提取各測點(測點編號見圖3)響應數據中包含峰值點、衰減段的5段樣本數據，每個樣本段有500個數據點。選取AR(10)模型建立時序模型，并提取DI，送入式(16)進行預測，計算結果見表4。

表4 試驗計算結果

將表4的計算結果代入式(16)，得到計算的銹蝕率為9.73%。與表1實測的銹蝕率對比可知，存在2.86%的誤差，表明筆者提出的識別模型效果較好。分析造成此種誤差的原因，主要有以下幾點：拱圈并非均勻銹蝕，運用式(16)產生了相應的系統誤差；不能完全保證完好狀態和銹蝕狀態下的激勵一致，因為殘差方差不僅會因結構損傷發生變化，還會由于激勵大小、采樣頻率、模型階數、環境因素等不同而發生變化。這些也是后續研究中需要著手解決的問題。

3 結 論

1)針對鋼筋銹蝕導致的鋼筋混凝土拱結構損傷，提出了基于加速度時間序列分析的損傷識別方法，建立了鋼筋銹蝕與結構剛度退化的關系模型，并依此進行了數值仿真，在通過模型試驗靜力測試數據驗證數值模型的合理性的基礎上，得到了不同損傷工況下的結構動力響應，根據所提時序建模方法計算損傷特征指標，獲取了銹蝕率與損傷特征指標之間的定量關系，提出了針對不同鋼筋直徑的銹蝕損傷識別模型。

2)利用模型試驗所獲得的加速度響應數據，對提出的基于時序分析的RC拱橋損傷識別方法進行了驗證，計算結果與試驗結果差異僅為2.89%，證明了該識別方法的有效性。

筆者所采用的基于殘差方差的損傷特征指標對于激勵大小、環境因素的敏感性仍有待進一步的分析，這也是作者未來努力的方向。