非一致地震激勵下大跨度橋梁彈塑性隨機(jī)響應(yīng)分析研究

劉小璐， 蘇 成， 聶 銘

(1. 廣東電網(wǎng)有限責(zé)任公司電力科學(xué)研究院, 廣州 510080； 2. 華南理工大學(xué) 土木與交通學(xué)院, 廣州 510640)

大跨度橋梁不同支承處的地震激勵存在顯著的非一致性，主要表現(xiàn)為行波效應(yīng)、失相干效應(yīng)和局部場地效應(yīng)，分別反映各支承處地震激勵間的相位差、相干性損失和頻譜差異[1]。同時(shí)，地震動具有本質(zhì)上的隨機(jī)性，應(yīng)該采用隨機(jī)振動方法對大跨度橋梁結(jié)構(gòu)開展非一致地震響應(yīng)分析研究。國內(nèi)外學(xué)者已采用不同隨機(jī)振動方法開展了大跨度橋梁非一致地震激勵問題研究。例如，Allam和Datta[2]與Dumanogluid和Soyluk[3]采用功率譜法計(jì)算了非一致地震激勵下大跨度斜拉橋的平穩(wěn)隨機(jī)響應(yīng)；Lin等[4]、焦常科和李愛群[5]與趙雷和劉寧國[6]均采用虛擬激勵法對非一致地震激勵下的大跨度斜拉橋開展隨機(jī)振動分析研究；劉小璐等[7]采用時(shí)域顯式法研究了非一致地震激勵下大跨度懸索橋的非平穩(wěn)隨機(jī)響應(yīng)。

上述研究主要解決非一致地震激勵下的線性隨機(jī)振動問題。強(qiáng)震作用下橋梁局部構(gòu)件會進(jìn)入彈塑性狀態(tài)，此時(shí)應(yīng)該采用非線性隨機(jī)振動方法進(jìn)行抗震分析。國內(nèi)外學(xué)者發(fā)展了多種非線性隨機(jī)振動方法，如FPK方程法[8]、隨機(jī)平均法[9]、矩方程法[10]、等效線性化法[11-12]、等效非線性系統(tǒng)法[13]、概率密度演化法[14]和Wiener路徑積分法[15]等。其中，大部分方法僅限于解決小規(guī)模非線性問題，難以應(yīng)用于大跨度橋梁結(jié)構(gòu)[16]。等效線性化法因其較強(qiáng)的適用性在實(shí)際工程中得到廣泛的應(yīng)用，例如Colangelo[17]利用等效線性化法計(jì)算了多層彈塑性框架在一致平穩(wěn)地震激勵下的隨機(jī)響應(yīng)；趙巖等[18]采用虛擬激勵法求解非線性結(jié)構(gòu)的等效線性化運(yùn)動方程，計(jì)算了連續(xù)梁橋在一致平穩(wěn)地震激勵下的彈塑性隨機(jī)響應(yīng)。盡管如此，該方法在求解非平穩(wěn)隨機(jī)問題時(shí)需要在不同時(shí)刻建立一系列等效線性系統(tǒng)，對于大規(guī)模工程結(jié)構(gòu)來說，計(jì)算量較大。該方法一般只能獲取結(jié)構(gòu)響應(yīng)的一、二統(tǒng)計(jì)矩，無法直接獲取非線性結(jié)構(gòu)響應(yīng)的尾部概率分布和峰值等更全面的統(tǒng)計(jì)信息。此外，隨機(jī)模擬法[19]是一種原理簡單且通用性強(qiáng)的方法，可以準(zhǔn)確獲取非線性結(jié)構(gòu)隨機(jī)響應(yīng)的任意階統(tǒng)計(jì)矩及概率分布，經(jīng)常用來檢驗(yàn)其它方法的計(jì)算精度。但隨機(jī)模擬法需涉及到大量反復(fù)的確定性樣本分析，若采用傳統(tǒng)彈塑性時(shí)程分析法進(jìn)行樣本分析將會非常耗時(shí)。

在作者前期研究[20]基礎(chǔ)上，針對非一致地震激勵下局部進(jìn)入彈塑性狀態(tài)的大跨度橋梁結(jié)構(gòu)，基于大質(zhì)量法[21]提出了一種時(shí)域顯式降維迭代算法，可以快速進(jìn)行非一致地震激勵下的彈塑性時(shí)程分析。將該方法與蒙特卡羅法結(jié)合，提出了一種時(shí)域顯式降維迭代-隨機(jī)模擬法，可高效獲取非一致地震激勵下大跨度橋梁結(jié)構(gòu)彈塑性響應(yīng)的統(tǒng)計(jì)矩和平均峰值。以主跨為1 200 m的某大跨度懸索橋?yàn)楣こ趟憷捎美w維梁柱單元模擬大橋彈塑性構(gòu)件，開展非一致地震激勵下彈塑性隨機(jī)振動分析研究，考察地震動空間效應(yīng)對彈塑性隨機(jī)響應(yīng)的影響。

1 非一致地震激勵

考慮地面運(yùn)動加速度為均勻調(diào)制的非平穩(wěn)零均值隨機(jī)過程，則結(jié)構(gòu)支承處的非一致地面運(yùn)動加速度向量可表示為

X(t)=Ψ(t)x(t)

(1)

式中

(2)

式中:m為非一致地震激勵個(gè)數(shù)；Xj(t)(j=1,2,…,m)為第j個(gè)非平穩(wěn)地面運(yùn)動加速度；ψj(t)為反映Xj(t)非平穩(wěn)特性的均勻調(diào)制函數(shù)；xj(t)為Xj(t)的平穩(wěn)部分，它們的互功率譜矩陣Sx(ω)可表達(dá)為

(3)

式中：Sjj(ω)為xj(t)的自功率譜，Sjk(ω)(j≠k)為xj(t)與xk(t)間的互功率譜。

根據(jù)式(3)所示的功率譜矩陣，采用諧波合成法生成非一致非平穩(wěn)地面運(yùn)動加速度為

(j=1,2,…,m)

(4)

式中:N為頻率分量個(gè)數(shù)；ωi為第i個(gè)頻率分量的圓頻率；φki為在(0,2π)范圍內(nèi)均勻分布的隨機(jī)相位角，且當(dāng)k≠r或i≠s時(shí)，φki和φrs相互獨(dú)立；ajk和θjk分別為第j與第k個(gè)地面運(yùn)動加速度間的相關(guān)幅值與相關(guān)相位角。

對平穩(wěn)地面運(yùn)動加速度x(t)的互功率譜矩陣Sx(ω)進(jìn)行Cholesky分解，可得

Sx(ω)=B(ω)B*T(ω)

(5)

式中

(6)

則式(4)中的ajk(ω)和θjk(ω)可分別表示為

(7)

(8)

式中:Δω為圓頻率間隔，函數(shù)Im、Re分別表示取復(fù)數(shù)的虛部、實(shí)部。

2 大質(zhì)量法

在一致地震激勵下，結(jié)構(gòu)內(nèi)力只與結(jié)構(gòu)相對地面的位移有關(guān)，因此一般基于結(jié)構(gòu)相對位移建立運(yùn)動方程，其中地震激勵為結(jié)構(gòu)各質(zhì)點(diǎn)上一致的地面運(yùn)動加速度。而在非一致地震激勵下，無法建立唯一的結(jié)構(gòu)-地面相對坐標(biāo)系，因此結(jié)構(gòu)運(yùn)動方程及其地震激勵施加方式與一致地震激勵情況大為不同。常用的處理方法有相對運(yùn)動法和大質(zhì)量法，其中相對運(yùn)動法由于需要進(jìn)行疊加運(yùn)算不適用于非線性情況，而大質(zhì)量法沒有這方面的限制，具有更強(qiáng)的適用性。

大質(zhì)量法是對結(jié)構(gòu)模型進(jìn)行動力等效的一種分析方法，該方法在處理非一致地震激勵問題時(shí)需要解除支承點(diǎn)沿地震作用方向的約束，并賦予該支承點(diǎn)大質(zhì)量，此時(shí)彈塑性結(jié)構(gòu)的運(yùn)動方程可寫為

(9)

式中

(10)

(11)

3 時(shí)域顯式降維迭代法

對于局部彈塑性結(jié)構(gòu)，將結(jié)構(gòu)非線性恢復(fù)力向量F(U)分解為線性部分和僅與彈塑性單元恢復(fù)力相關(guān)的部分為

F(U)=K0U+LnFn

(12)

式中:K0為結(jié)構(gòu)初始彈性剛度陣；Ln為彈塑性單元恢復(fù)力的指示矩陣；Fn為彈塑性單元恢復(fù)力中的非線性部分可表示為

Fn=R(Un)-Kn,0Un

(13)

式中:R為彈塑性單元的恢復(fù)力向量；Un為彈塑性單元的節(jié)點(diǎn)位移向量；Kn,0為彈塑性單元的初始彈性剛度矩陣，屬于結(jié)構(gòu)初始彈性剛度矩陣K0的一部分。將式(12)代入式(9)，并將非線性部分移至方程右端，可得

(14)

式中

(15)

L=[Lx-Ln]

(16)

Vi=Ai,0P0+Ai,1P1+…+Ai,iPi

(i=1,2,…,nt)

(17)

式中:Vi=V(ti)，Pi=P(ti)，i=1,2,…,nt，其中nt為總時(shí)間步數(shù)；Ai,j(j=0,1,…,i)為系數(shù)矩陣，它們只與結(jié)構(gòu)參數(shù)有關(guān)，可表達(dá)為如下閉合公式：

(18)

基于Newmark-β積分方法可獲取T,Q1和Q2為

(19)

式中:I為單位矩陣；取γ=0.5,β=0.25,則Newmark-β積分方法可以實(shí)現(xiàn)無條件穩(wěn)定。

對于當(dāng)前時(shí)刻ti，假設(shè)已逐步計(jì)算出結(jié)構(gòu)響應(yīng)V0,V1,…,Vi-1，則式(17)右端的為已知量P0,P1,…,Pi-1，然而P1中依然含有關(guān)于未知向量Un,i的非線性項(xiàng)，需要進(jìn)行迭代計(jì)算。式(17)共有2nd(nd為結(jié)構(gòu)總自由度數(shù))維，對于大規(guī)模結(jié)構(gòu)，直接迭代求解該式的計(jì)算量較大。為了提高計(jì)算效率，從式(17)中提取Un,i相關(guān)的行，可得

(20)

Ψ(Un,i)=Yi(i=1,2,…,nt)

(21)

式中

(22)

(i=1,2,…,nt)

(23)

采用Newton-Raphson方法求解式(21)，將所獲得的收斂結(jié)果Un,i代入式(17)的其余行中，則得到結(jié)構(gòu)當(dāng)前時(shí)刻其它單元的節(jié)點(diǎn)位移向量。例如對于結(jié)構(gòu)某一關(guān)鍵響應(yīng)r，可表示為

(24)

4 時(shí)域顯式降維迭代-隨機(jī)模擬法

將時(shí)域顯式降維迭代法與蒙特卡羅法相結(jié)合，利用時(shí)域顯式降維迭代法進(jìn)行樣本分析，高效獲取非一致地震激勵下大跨度橋梁結(jié)構(gòu)的彈塑性隨機(jī)響應(yīng)，具體計(jì)算步驟如下：

(2)建立結(jié)構(gòu)初始有限元模型，釋放地震激勵方向的支承約束，并在相應(yīng)自由度施加大質(zhì)量。

5 工程實(shí)例

5.1 有限元模型

虎門二橋大沙水道橋是一座主跨1 200 m的地錨式懸索結(jié)構(gòu)。該橋跨徑布置為(360+1 200+480)m，矢跨比1∶9.5，如圖 1所示。門式塔高191.1 m，設(shè)上、下兩道橫梁，如圖2所示。主塔為鋼筋混凝土結(jié)構(gòu)，塔底截面及其配筋情況如圖 3所示，其中縱向受力鋼筋強(qiáng)度等級為HRB400，混凝土強(qiáng)度等級為C50。利用ANSYS建立大橋的初始三維有限元模型。其中，主梁、主塔和邊墩均采用梁單元模擬，主纜和吊桿則通過桿單元模擬。主塔和邊墩底部固支，主纜通過大橋的兩側(cè)錨碇固定。主梁側(cè)向和豎向運(yùn)動受主塔和邊墩的約束，而縱向運(yùn)動不受主塔和邊墩約束。

圖1 懸索橋立面圖

圖2 主塔立面圖

圖3 主塔塔底截面及其配筋圖

歷史震害數(shù)據(jù)表明，大跨度懸索橋主塔相對其它部位更容易損傷[22]。虎門二橋大沙水道橋抗震專題研究表明，在重現(xiàn)期分別為945年和2450年的中震和大震作用下，整個(gè)結(jié)構(gòu)可以依然處于彈性狀態(tài)，然而主塔塔底截面承載能力的冗余度遠(yuǎn)低于其它關(guān)鍵截面。因此，當(dāng)發(fā)生強(qiáng)度高于大震的地震時(shí)，主塔塔底截面可能會率先進(jìn)入彈塑性狀態(tài)，應(yīng)該采用彈塑性單元模擬主塔底部，并預(yù)設(shè)塑性發(fā)展范圍為塔底10 m以內(nèi)，在后續(xù)分析中將對該塑性發(fā)展范圍進(jìn)行驗(yàn)證。在地震和恒載作用下，主塔底部主要受到軸力、剪力和彎矩作用，由于不允許發(fā)生剪切破壞，在此僅考慮軸力和彎矩引起的塑性行為，并采用可以模擬軸力-彎矩耦合作用的纖維彈塑性單元對預(yù)設(shè)的塔底塑性范圍進(jìn)行模擬。分別采用Menegotto-Pinto模型[23]和修正Kent-Park模型[24]模擬鋼筋纖維和混凝土纖維的軸向應(yīng)力應(yīng)變關(guān)系，鋼筋的屈服強(qiáng)度和混凝土的抗壓強(qiáng)度分別取為330 MPa和22.4 MPa。本文通過幾何非線性靜力分析計(jì)算恒載引起的大位移和應(yīng)力剛化效應(yīng)，在此基礎(chǔ)上開展地震激勵下的結(jié)構(gòu)動力分析，不考慮地震激勵新引起的幾何非線性效應(yīng)。

5.2 非一致地震激勵

考慮東、西塔底地震激勵的差異性，而同側(cè)塔底與錨碇相距較近，因此可以認(rèn)為同側(cè)塔底與錨碇的地震激勵具有一致性。為此，該結(jié)構(gòu)共有2個(gè)非一致地震激勵，式(1)中反映地面運(yùn)動加速度非平穩(wěn)的均勻調(diào)制函數(shù)可表示為

(25)

式中:ta=3.2 s，tb=17.0 s,c=0.13。地震激勵間的互功率譜可表示為

(26)

圖4 非一致地面運(yùn)動加速度時(shí)程樣本

5.3 彈塑性時(shí)程分析

針對圖 4所示的非一致地面運(yùn)動加速度時(shí)程樣本，分別采用本文的時(shí)域顯式降維迭代法和傳統(tǒng)非線性Newmark-β法對大橋進(jìn)行彈塑性時(shí)程分析。東塔底部截面Ⅰ和下橫梁上截面Ⅳ的彎矩時(shí)程分別如圖 5和圖 6所示。從圖中可以看出，兩種方法的結(jié)果完全重合，驗(yàn)證了時(shí)域顯式降維迭代法的準(zhǔn)確性。

圖5 東塔截面Ⅰ的彎矩時(shí)程

圖6 東塔截面Ⅳ的彎矩時(shí)程

傳統(tǒng)非線性Newmark-β法和時(shí)域顯式降維迭代法的計(jì)算時(shí)間分別為2 518.3 s和127.4 s。顯而易見，本文的時(shí)域顯式降維迭代法具有更高的計(jì)算效率，這是因?yàn)樵摲椒▋H僅需要對主塔底部4個(gè)彈塑性纖維單元進(jìn)行非線性降維迭代運(yùn)算，而傳統(tǒng)方法需要將橋梁所有自由度捆綁一起進(jìn)行迭代計(jì)算。需要指出的是，本文方法的計(jì)算時(shí)間主要包括兩部分：一部分用于構(gòu)建如式(17)所示動力響應(yīng)顯式表達(dá)式中的系數(shù)矩陣，為109.9 s；另一部分用于降維迭代運(yùn)算，僅為17.5 s。其中，系數(shù)矩陣只依賴于結(jié)構(gòu)本身信息，它們只需要計(jì)算一次，隨后可重復(fù)用于后續(xù)隨機(jī)模擬的所有樣本分析中。

為了考察強(qiáng)震作用下主塔塑性變形的發(fā)展情況，給出如圖 2所示東塔底部截面Ⅰ-Ⅲ的彎矩-曲率關(guān)系，分別如圖7(a)～(c)所示。從圖中可看出，非一致地震激勵下截面Ⅰ的彎矩-曲率曲線呈現(xiàn)出相對飽滿的旗形滯回環(huán)，表明截面Ⅰ表現(xiàn)出顯著的彈塑性行為；截面Ⅱ的滯回環(huán)面積明顯減小，非線性行為相對較弱；而截面Ⅲ的彎矩-曲率關(guān)系幾乎趨近于線彈性。以上結(jié)果表明，長度為10 m的彈塑性纖維單元可以完全捕捉主塔底部的彈塑性變形行為。截面Ⅰ-Ⅲ最外側(cè)鋼筋應(yīng)力-應(yīng)變關(guān)系如圖8(a)～(c)所示，同樣表明在所關(guān)注的地震水準(zhǔn)下，主塔底部彈塑性行為發(fā)展未超出預(yù)期范圍。

圖7 東塔底部截面Ⅰ-Ⅲ的彎矩-曲率關(guān)系

圖8 東塔底部截面Ⅰ-Ⅲ外側(cè)鋼筋應(yīng)力-應(yīng)變關(guān)系

5.4 彈塑性隨機(jī)振動分析

時(shí)域顯式降維迭代法用于彈塑性時(shí)程分析的計(jì)算精度和計(jì)算效率均得到了驗(yàn)證，可進(jìn)一步利用時(shí)域顯式降維迭代-隨機(jī)模擬法開展彈塑性抗震隨機(jī)振動分析。為了考察地震激勵非一致性的影響，同時(shí)計(jì)算懸索橋在一致和非一致順橋向地震激勵下的彈塑性隨機(jī)響應(yīng)，其中非一致地震激勵按5.2節(jié)中的參數(shù)確定，一致地震激勵取其東塔底處的激勵，隨機(jī)模擬樣本數(shù)取為500。

兩種地震激勵工況下西塔和東塔底部截面內(nèi)力標(biāo)準(zhǔn)差分別如圖9和圖10所示。從圖中可以觀察到，對于西塔底截面，非一致地震激勵下的軸力標(biāo)準(zhǔn)差偏小，而剪力和彎矩標(biāo)準(zhǔn)差偏大；對于東塔底截面，非一致地震激勵下的軸力和彎矩標(biāo)準(zhǔn)差偏小，剪力標(biāo)準(zhǔn)差與一致地震激勵下的結(jié)果非常接近。由此可知，地震激勵空間效應(yīng)會使得該大跨度懸索橋內(nèi)力結(jié)果既可能偏大也可能偏小。

圖9 西塔底截面內(nèi)力標(biāo)準(zhǔn)差

圖10 東塔底截面內(nèi)力標(biāo)準(zhǔn)差

此外，利用本文方法可進(jìn)一步獲取大橋關(guān)鍵響應(yīng)的平均峰值。主塔底截面軸力、剪力和彎矩的平均峰值如表1所示。為了考察地震激勵空間效應(yīng)對結(jié)構(gòu)響應(yīng)的影響，將兩種工況的相對差異也列于表1中。從表中可以看出，地震激勵空間效應(yīng)對響應(yīng)平均峰值的影響規(guī)律與對響應(yīng)標(biāo)準(zhǔn)差的影響規(guī)律基本相同。值得注意的是，當(dāng)考慮地震激勵空間效應(yīng)時(shí)，東塔底截面軸力的平均峰值比一致地震激勵下的結(jié)果小了9.4%；西塔底截面剪力和彎矩的平均峰值比一致地震激勵下的結(jié)果分別大了12.4%和11.8%。在該大橋彈塑性抗震設(shè)計(jì)中需要考慮地震激勵空間效應(yīng)對結(jié)構(gòu)內(nèi)力的影響。

表1 塔底截面內(nèi)力的平均峰值

時(shí)域顯式降維迭代-隨機(jī)模擬的計(jì)算時(shí)間如表2所示，取不同樣本數(shù)，平均每個(gè)樣本分析的時(shí)間基本相同，均是13 s左右。這表示只需要花13 s的時(shí)間就可以完成一次大跨度懸索橋確定性彈塑性時(shí)程分析。當(dāng)樣本數(shù)為500時(shí)，時(shí)域顯式降維迭代-隨機(jī)模擬法的計(jì)算時(shí)間為不到2個(gè)小時(shí)，而基于傳統(tǒng)非線性Newmark-β方法的隨機(jī)模擬法的計(jì)算時(shí)間則接近200個(gè)小時(shí)。本文的時(shí)域顯式降維迭代-隨機(jī)模擬法在計(jì)算效率方面具有明顯的優(yōu)越性。以上所有的計(jì)算均是在帶有Intel Core i5處理器和4 GB內(nèi)存的臺式機(jī)上完成的。

表2 時(shí)域顯式降維迭代-隨機(jī)模擬的計(jì)算時(shí)間

6 結(jié) 論

(1)基于大質(zhì)量法推導(dǎo)了非一致地震激勵下結(jié)構(gòu)動力響應(yīng)形式上的時(shí)域顯式表達(dá)式，利用該顯式表達(dá)式的降維列式優(yōu)勢，提出僅關(guān)于彈塑性單元節(jié)點(diǎn)自由度的時(shí)域顯式降維迭代算法，顯著提高非一致地震激勵下大跨度橋梁彈塑性時(shí)程分析的計(jì)算效率。并結(jié)合蒙特卡羅法，提出了非一致地震激勵下大跨度橋梁彈塑性隨機(jī)振動分析的時(shí)域顯式降維迭代-隨機(jī)模擬法。

(2)以虎門二橋大沙水道橋?yàn)楣こ虒?shí)例，開展順橋向非一致地震激勵下的彈塑性隨機(jī)振動分析，驗(yàn)證了本文方法的正確性和高效性。研究結(jié)果表明，對于該大跨度懸索橋，非一致地震激勵下主塔內(nèi)力的標(biāo)準(zhǔn)差和平均峰值比一致地震激勵下的結(jié)果既可能偏大也可能偏小，其中東塔底截面軸力的平均峰值偏小了9.4%，西塔底截面的剪力和彎矩平均峰值分別偏大了12.4%和11.8%。在該大橋彈塑性抗震設(shè)計(jì)中需要考慮地震激勵空間效應(yīng)對結(jié)構(gòu)內(nèi)力的影響。