隨機橋梁地震可靠度分析的改進最大熵方法

陳志強， 鄭史雄, 2

(1.西南交通大學 土木工程學院，成都 610031； 2.陸地交通地質災害防治技術國家工程實驗室(西南交通大學)，成都 610031)

在基于性能的地震工程理論框架下，橋梁的地震可靠度分析對于橋梁的地震損傷評估和地震風險分析都具有重要意義[1]。一旦橋梁的地震可靠度確定，其地震風險可以通過失效概率對不同地震災害水平進行積分計算[2]。因此，從結構工程的觀點來看，橋梁地震風險分析中最為重要的部分就是橋梁的地震可靠度分析。

在實際工程中，地震動和橋梁結構本身通常都存在顯著的不確定性。原則上來說，在橋梁結構地震可靠度分析中最為合理的應該同時考慮這兩類不確定性，然而由于動力可靠度理論的不完善，現有的橋梁地震可靠度分析中通常采用隨機振動理論和隨機有限元理論兩類方法分別對地震動的不確定性和橋梁結構參數的不確定性單獨進行考慮[3]。當考慮地震動的不確定性時，橋梁的地震可靠度分析主要采用隨機振動理論進行[3-5]，然而現有的隨機振動理論主要是針對線彈性結構。雖然也有學者建立了針對非線性結構的隨機振動分析方法，如路徑積分法[6]、矩近似法[7]、等效線性化方法[8]以及概率密度演化理論[9-10]等，但是這些方法通常都只適用于簡單結構，大規模、強非線性結構的隨機地震響應分析仍然是一個巨大挑戰。當考慮橋梁結構參數的不確定性時，橋梁的地震可靠度分析通常采用隨機有限元方法[11-14]。然而，現有的方法不僅無法考慮地震動的不確定性，而且在復雜結構二階統計量的計算上都還存在較大難度[3]。另外，對于非線性橋梁結構而言，其地震響應具有顯著的非高斯性，僅采用一階和二階統計矩對其進行地震可靠度分析精度還遠遠不夠。當同時考慮橋梁結構和地震動的不確定性時，橋梁地震可靠度分析方法的唯一可行方法只有基于抽樣的數值模擬方法，如蒙特卡洛模擬(Monte Carlo simulation,MCS)[15]、子集模擬[16]以及重要性抽樣[17]等。這些方法雖然理論上可以獲得結構地震可靠度的精確結果，但是其所需要的樣本數量通常都與失效概率成反比。在小失效概率水平下，采用這些抽樣數值方法對橋梁結構進行地震可靠度分析時所需要的動力分析計算量非常大，同樣很難適用于大跨度斜拉橋、懸索橋以及高墩連續剛構橋等復雜橋梁結構的地震可靠度分析。

最近幾年，隨著分數階矩的出現，基于最大熵原理的結構動力可靠度分析獲得了研究者極大的關注[14,18-20]。文獻[18]分別采用單變量降維模型、雙變量降維模型以及稀疏網格積分對結構動力系統的分數矩進行估計，然后將其與最大熵模型相結合，建立了結構非線性地震可靠度分析的最大熵方法。文獻[14]采用乘法形式的降維模型對結構地震響應的矩進行估計，對一座典型高墩橋梁進行了地震可靠度分析。文獻[19]以最大熵原理為基礎，提出了結構動力可靠度分析中的非等權重擬蒙特卡洛方法。文獻[20]采用拉丁超立方抽樣對結構地震響應分數矩進行估計，并將其與最大熵原理相結合，建立了近斷層脈沖地震作用下高墩大跨連續剛構橋的地震可靠度分析方法。

由于分數階矩中含有高階中心矩信息，通過少數幾階分數矩就可以對結構地震響應極值分布進行準確模擬，因此基于分數矩的最大熵原理在橋梁地震可靠度分析中具有明顯的優勢。然而，由于現有的分數矩最大熵原理需要對多個變量進行迭代求解，在采用最大熵原理求解結構地震響應極值分布時求解過程受初值影響非常大，而且經常面臨著優化過程不穩定的問題[18]。為了克服這一問題，本文提出了一種基于似然函數的最大熵分布模型，并將該模型應用于隨機地震作用下，考慮非線性的隨機橋梁結構地震可靠度分析，建立了基于最大熵原理的橋梁非線性地震可靠度分析方法。最后，以一座典型的高墩大跨連續剛構橋為例，通過蒙特卡洛模擬對所提出的方法進行驗證。

1 橋梁地震可靠度問題闡述

對于一個考慮結構參數和地震動雙重不確定性的N自由度非線性動力系統，地震作用下，其動力方程可表示為

(1)

鑒于譜表示法具有簡單、高效、易于實現等優點，因此本文采用譜表示法對地震動進行模擬。根據非平穩隨機過程譜表示理論，一維地震動的表達式[21]通常為

(2)

式中：ωl=lΔω，Δω為頻率區間；{Xl,Yl}為2m個標準正交隨機變量；SU(ωl,t)為隨機地震動的時-頻演化功率譜，可以表示為

SU(ω,t)=|A(ω,t)|2S(ω)

(3)

其中A(ω,t)為地震動的時-頻調制函數。

根據譜表示理論，隨機變量{Xl,Yl}(l=1,2,…,m)應滿足[21]:

(4)

其中:E[·]表示期望算子，δrl為Kronecker-Delta函數。

將模擬的地震動代入式(1)中的動力學方程，并采用數值方法對式(1)進行求解，由此可以得到橋梁結構地震可靠度分析所關心的物理量(如墩頂位移、支座相對位移等)，其通常可以表示為

Z(t)=H(υ,χ,d,t)

(5)

式中:d為確定性參數向量，H(·)為確定性算子。

對于通常所關心的首次超越地震可靠度，其可以定義為

R=Pr{Z(t)

(6)

式中：Pr表示概率；T為地震動持時，ZT為橋梁結構抗震能力，其通常可以由規范規定或者根據Pushover分析進行確定。

由于式(6)中結構可靠度的求解非常困難，為了對其進行求解，這里通過等價極值原理[19]，將橋梁結構的地震可靠度轉化為

R=Pr{Zext

(7)

式中Zext為橋梁結構非線性地震響應的極值，可以定義為

(8)

(9)

橋梁結構的地震可靠度則定義為

(10)

式中pZext(z)為橋梁非線性地震響應極值Zext的概率密度，即橋梁結構的地震響應極值分布。同時，橋梁結構的失效概率則為

(11)

綜上所述，在橋梁地震可靠度分析中，最為重要的就是對其地震響應極值分布的估計。一旦橋梁地震響應的極值分布確定，其地震可靠度和失效概率就可以通過簡單的數值積分進行求解。此外，對于工程中的大跨度復雜橋梁結構，都是高可靠性結構，在設計地震作用下，失效的概率通常都非常小。為了對其地震可靠度進行估計，需要對橋梁結構地震響應極值分布的尾部進行精確計算。鑒于此，本文提出一種針對復雜非線性結構的地震響應極值分布估計的高效數值方法。

2 結構地震可靠度分析的最大熵方法

2.1 結構地震響應極值分布的最大熵模型

橋梁非線性地震響應的極值Zext在橋梁結構參數和地震動不確定性的影響下，通常是一個正的隨機變量，為簡單起見，這里將其表示為Z。根據最大熵原理，隨機變量Z的無偏概率分布可以在分數矩約束下，通過熵的最大化進行求解。因此，根據最大熵原理橋梁非線性地震響應極值分布pZext(z)可以通過以下約束非線性優化問題求解得到[18-19]：

(12)

式中μZαi為隨機變量Z的第αi階分數矩，J(z)為隨機變量Z的信息熵，定義為

(13)

通過采用拉格朗日乘數法對式(12)進行求解，可以得到橋梁地震響應極值分布的一般形式[18]為

(14)

式中λi(i=0,1,…,M)為M+1個拉格朗日乘子。

根據概率密度函數積分等于1的正則化條件可以得到

(15)

(16)

因此，pZext(z)的估計只需要對式(16)進行求解，得到分數指數向量α和拉格朗日乘子λ，并將其代入式(14)，即可得到橋梁非線性地震響應極值分布的無偏估計。

目前式(16)的求解普遍都是直接采用單純形算法進行迭代求解，所得到的都是局部解，求解過程對α和λ的初值具有很強的依賴性，而且求解過程中經常存在數值不穩定的問題。為了克服這一問題，對式(16)中的優化問題進行快速求解，并保證求解過程的穩定性，本文提出一種新的基于似然函數的求解方法。

2.2 基于似然函數的最大熵求解方法

根據分數矩的定義，對于隨機變量Z，其第αi階分數矩為

(17)

式中gi(z)=zαi。將式(14)代入式(17)，并采用分部積分法可以得到

(18)

(19)

根據式(18)、(19)，拉格朗日乘子向量可以通過以下線性方程組進行確定[18]：

λ=ω-1μZ

(20)

(21)

(22)

可以看到，通過步驟1～7，式(16)中分數矩最大熵模型轉化為了NQ,M個線性方程組的求解，因此其計算效率大大提高。另外，由于分數指數是在其整個潛在區間[αmin,αmax]離散，因此該方法得到的是全局最優解。

2.3 橋梁地震響應極值分數矩估計

Xl=Γ(cas(lχ1)),Yl=Γ(cas(lχ2))

(23)

式中：cas=sin(·)+cos(·)為Hartley正交基，Γ(·)為隨機映射算子。

根據式(23)，地震動模擬中的高維隨機變量就被減少到只有兩個基本隨機變量，因此地震動的不確定性將完全由隨機變量χ1和χ1所表征。此時，概率空間維度為d=d1+2 在結構地震響應極值分布分數矩的估計時，只需要在低維概率空間中進行即可。

對于結構非線性地震響應極值分布的第α階分數矩，其通常可以估計如下：

(24)

為了確定估計中的積分點，文獻[24]通過中心化的L2偏差對樣本點集進行篩選，提出了一種改進的偽相關性折減拉丁超立方抽樣方法(improved spurious correlation reduced LHS, ICLHS)。該方法在非線性結構動力響應的分數矩估計中取得了較好的效果，本文采用該方法對結構非線性地震響應極值分布的分數矩進行估計。

2.4 橋梁地震可靠度分析

采用2.3節中的ICLHS方法生成橋梁結構和地震隨機參數樣本點，并通過式(2)對地震動進行模擬，從而形成結構-地震動樣本對。然后采用非線性時程分析對結構地震響應進行求解，從而得到結構地震響應極值的分數矩。最后采用2.2節中所提出的最大熵方法對結構地震響應極值分布進行求解，即可得到橋梁結構非線性地震響應極值分布。再根據橋梁結構不同損傷極限狀態的定義，通過對極值分布進行數值積分即可得到橋梁結構的地震可靠度。橋梁結構地震可靠度分析流程如圖1所示。

圖1 橋梁地震可靠度分析流程圖

3 算例分析

3.1 橋梁概況及非線性有限元模型

以一座已經建成的高墩連續剛構橋為例，說明本文方法在復雜非線橋梁結構地震可靠度分析中的應用。該橋梁跨徑布置為90 m+170 m+90 m，橋梁總長350 m，左右兩個橋墩墩高分別為110.5 m和126.49 m，為典型的山區不規則高墩橋梁。橋梁上部主梁采用單箱單室預應力混凝土箱型梁，材料為C50預應力混凝土，主跨的跨中和橋墩的墩頂處截面梁高分別為3.7 m和10 m。為了適應主梁的變形，橋墩采用空心薄壁墩，墩身材料采用C40混凝土，橋梁結構詳細布置如圖2所示。

采用OpenSees作為分析平臺，建立橋梁非線性有限元分析模型。主梁本文采用基于位移的梁柱單元并結合彈性截面進行模擬，不考慮主梁的非線性行為。橋墩作為地震中的易損構件，強震作用下，墩頂和墩底有可能進入屈服，產生塑性鉸，本文采用基于力的非線性梁柱單元進行模擬，并通過纖維截面模擬橋墩的非線性行為。對于主梁和橋臺伸縮縫位置處，強震作用下其可能發生嚴重碰撞，本文采用接觸單元對梁-臺之間的碰撞行為進行模擬，碰撞單元材料采用Hertz-damp模型[25]，從而充分考慮碰撞過程中的剛度變化和能量耗散。在橋臺位置處，分別放置了兩個盆式橡膠支座對結構變形進行約束，支座是地震中的易損構件，其通常都會產生很強的非線性，這里采用零長度單元結合硬化材料對其非線性行為進行模擬。橋梁結構有限元模型如圖2所示。

圖2 橋梁結構有限元模型示意(m)

3.2 橋梁結構不確定性

橋梁結構的不確定性通常主要包括結構材料、幾何尺寸、邊界條件以及阻尼比等參數。既有研究表明，混凝土抗壓強度、鋼筋屈服強度以及橋梁結構阻尼比等參數對橋梁結構的地震需求具有顯著影響。本文在參考既有研究的基礎上[26-27]，選出7個影響橋梁結構地震既有最為顯著的參數作為考慮橋梁結構不確定性的隨機變量，其概率分布見表1。

表1 橋梁結構的隨機參數概率分布

3.3 隨機地震動功率譜模型

為了考慮地震動的時-頻非平穩性，這里采用文獻[28]中演化功率譜(evolutionary power spectral density, EPSD)模型生成地震記錄，從而對結構進行地震可靠度分析，其演化功率譜密度函數表示為

(25)

式中：A(t)為地震動非平穩調制函數，本文將其取為

(26)

其中c控制地震動峰值到達時間，d控制地震動的形狀，取值由橋位處所在場地確定。

此外，地震動演化功率譜模型中S0(t)為地震動譜強度因子，將其取為

(27)

(28)

圖3 地震演化功率譜

圖4 模擬地震動功率譜與目標值比較

3.4 橋梁結構地震可靠度分析方法驗證

根據橋梁結構隨機參數概率分布和地震動隨機參數，采用2.3節中的ICLHS方法生成400個積分點，通過非線性時程分析對結構地震響應進行求解即可得到400個橋梁結構地震響應極值。圖7以其中一個樣本為例，給出了地震作用下1號墩的墩底彎矩-曲率滯回曲線。從圖7中可以看到，在地震作用下橋梁結構產生了顯著的非線性變形。

為了對本文方法進行驗證，說明其在復雜結構非線性地震可靠度分析中的適用性。這里以1號墩的墩頂位移為例，圖8給出了1號墩的墩頂位移響應極值分數矩與MCS(104次)計算結果的對比。

圖5 典型地震記錄樣本

圖6 模擬地震記錄地震動均值和均方值與目標值比較

圖7 橋梁結構非線性地震響應典型樣本

從圖8中可以看到，采用ICLHS方法估計的分數矩與MCS的結果十分吻合，特別是在-2到1階內，二者相對誤差不到2%，由此說明了ICLHS方法在復雜非線性地震響應的分數階矩時具有較高的精度和效率。另外，從圖8中還可以看到當分數階數較小時，橋梁結構地震響應極值分布分數矩的相對誤差較小，而當分數階數較高時，分數矩相對誤差較大，由此說明了采用低階分數矩能夠更為準確地模擬結構地震響應極值分布。

圖8 橋梁非線性地震響應分數矩

若以橋梁結構非線性地震響應極值分數矩作為約束，采用2.2節中的方法對最大熵模型進行求解即可得到結構非線性地震響應極值分布。圖9給出了采用本文方法獲得的1號墩墩頂位移非線性地震響應極值分布的概率密度函數(probability density function,PDF)、累積密度分布(cumulative density function,CDF)以及超越概率(probability of exceedance,POE)與MCS結果的對比，同時，圖9中還給出了采用核密度估計(kernel density estimation,KDE)以及對數正態分布擬合得到的結果。

從圖9中可以看到，采用本文方法獲得的橋墩非線性地震響應極值分布與MCS結果具有較好的一致性，即使在小失效概率水平下(Pf<10-3)，二者也十分吻合。這說明了采用該方法能夠通過較少的幾百次非線性動力分析，對大規模復雜非線性結構地震響應極值分布進行精確估計。同時，從圖9中還可以看到，不論極值分布的主體還是尾部分布，采用核密度估計和對數正態分布擬合的結果與蒙特卡洛模擬結果相比都存在較大誤差，特別是對于尾部分布，二者間的誤差非常大，由此說明了采用核密度估計和對數正態分布擬合都無法對小失效概率水平下結構的地震可靠度進行有效估計。

(b)CDF

(c)POE

通常，對于非參數概率模型，其魯棒性通常較差，如核密度估計，其結果都嚴重依賴于帶寬的選擇。為了說明2.2節中基于似然函數的最大熵模型求解方法的魯棒性，圖10給出了采用不同分數指數增量Δα獲得的橋梁非線性地震響應極值分布的對比。從圖10中可以看到，當分數指數增量Δα=0.1時，該方法獲得的極值分布與蒙特卡洛模擬結果吻合得最好，二者基本完全重合。當Δα=0.2、0.25、0.3時，該方法獲得的結果與蒙特卡洛模擬結果仍然一致，這說明了本文提出的基于似然函數的最大熵模型求解方法具有較好的魯棒性。值得注意的是當Δα=0.2、0.25、0.3時，在尾部(10-4數量級)，本文方法獲得的超越概率曲線與MCS結果出現了非常細微的偏差。這主要是由于本文所采用的算例為實際的高墩連續剛構橋，模型較為復雜，動力分析計算量較大，因此在蒙特卡洛模擬時只進行了104次抽樣。在10-4數量級，蒙特卡洛模擬的結果自身就存在一定的誤差，因此由于Δα的不同帶來的這一點偏差基本是可以忽略不計的。同時，這也說明只要分數指數增量Δα的取值小于0.3，2.2節中基于似然函數的方法就能夠對最大熵模型進行精確求解。

(a) PDF

(b) CDF

4 結 論

1)建立了基于分數階矩最大熵原理的橋梁非線性地震可靠度分析方法，通過與蒙特卡洛模擬結果的對比說明了本文方法能夠通過少量的非線性動力分析實現小失效概率水平下結構非線性的地震響應極值分布的精確求解，該方法能夠為大規模復雜結構非線性地震可靠度分析提供一種有效途徑。

2)提出了一種基于似然函數的最大熵模型求解方法，能夠將最大熵模型求解過程中的無約束極值問題轉化為線性方程組的求解，從而克服了采用單純形方法求解最大熵分布迭代過程不穩定的問題。

3)核密度估計和對數正態分布對橋梁非線性地震響應極值分布的主體具有一定近似能力，但是二者都無法準確模擬結構極值分布的尾部分布，而最大熵分布不僅能夠模擬結構非線性地震響應極值分布的主體，而且對于尾部分布也具有非常好的近似能力。