基于高斯過(guò)程回歸的橋梁多變量地震易損性分析

閆業(yè)祥， 孫利民,2,3

(1.同濟(jì)大學(xué) 橋梁工程系，上海 200092；2.同濟(jì)大學(xué) 土木工程防災(zāi)國(guó)家重點(diǎn)實(shí)驗(yàn)室，上海 200092；3.上海期智研究院，上海 200232)

地震易損性分析表達(dá)了在給定地震動(dòng)強(qiáng)度的條件下工程結(jié)構(gòu)達(dá)到或超過(guò)某個(gè)性能極限狀態(tài)的概率，已廣泛應(yīng)用于土木基礎(chǔ)設(shè)施的震前風(fēng)險(xiǎn)預(yù)測(cè)、震后快速評(píng)估以及為抗震加固提供依據(jù)。在由太平洋地震工程研究中心發(fā)展的新一代基于性能的地震工程框架中，地震易損性分析是極其重要的一環(huán)，涵蓋了需求分析與損傷分析兩大模塊[1]。地震易損性分析的一個(gè)重要結(jié)果是地震易損性曲線，定義為結(jié)構(gòu)在不同的地震強(qiáng)度指標(biāo)(intensity measure, IM)下對(duì)各性能極限狀態(tài)的超越概率函數(shù)。但是，當(dāng)前廣泛使用的地震易損性曲線建立在單一的參數(shù)(即一個(gè)標(biāo)量IM，譬如峰值地面加速度PGA，或結(jié)構(gòu)第一階周期處的譜加速度)之上，其不足主要體現(xiàn)在兩個(gè)方面：① 單一的IM不足以表征地震動(dòng)對(duì)結(jié)構(gòu)的破壞能力，而矢量IM可能更有效[2]；② 結(jié)構(gòu)相關(guān)參數(shù)中存在的不確定性對(duì)易損性有不可小覷的影響[3]。因此，傳統(tǒng)的基于單參數(shù)的地震易損性曲線無(wú)法彌補(bǔ)上述不足，亟需一個(gè)新的工具來(lái)考慮多重變量對(duì)易損性的影響，這正是本文發(fā)展多變量地震易損性分析方法的背景與目的。

當(dāng)前已有一些學(xué)者對(duì)多變量地震易損性分析作了有益的探索[4-6]，其中基于代理模型的Logistic回歸是一種較為有效的方法。該方法主要利用代理模型(也稱元模型或近似模型，譬如響應(yīng)面模型、神經(jīng)網(wǎng)絡(luò)回歸等)來(lái)代替有限元“黑箱”分析，通過(guò)代理模型預(yù)測(cè)的大樣本地震需求來(lái)構(gòu)造生存-失效二值矢量(0代表安全，1代表超越了極限狀態(tài))，并在此基礎(chǔ)上擬合Logistic回歸獲得構(gòu)件的多變量地震易損性函數(shù)。上述方法的主要缺點(diǎn)在于流程上的繁瑣以及單一的IM仍被使用，而且易損性函數(shù)的精度極大地依賴于代理模型的預(yù)測(cè)精度以及Logistic回歸的擬合精度。不同于上述方法，本文提議了一個(gè)新的基于地震需求均值和需求方差模型的多變量地震易損性分析方法，并可以直接包含矢量IM與多重結(jié)構(gòu)參數(shù)。其中均值與方差采用高斯過(guò)程回歸來(lái)擬合，高斯過(guò)程回歸是一個(gè)非參數(shù)機(jī)器學(xué)習(xí)方法，最近一些研究已將其應(yīng)用于結(jié)構(gòu)的抗震性能評(píng)估中[7-9]。作為多變量易損性函數(shù)的直接擴(kuò)展，本文發(fā)展了均值與包絡(luò)易損性曲線的概念以充分考慮各種不確定參數(shù)的影響。最后，識(shí)別了不同的IM對(duì)易損性曲線離散性的影響，并確定了最優(yōu)的IM形式。

1 基于高斯過(guò)程回歸的多變量地震易損性模型

1.1 多變量概率地震需求模型

當(dāng)前的地震易損性評(píng)估方法主要采用非線性動(dòng)力分析進(jìn)行，即在一系列結(jié)構(gòu)參數(shù)與地震動(dòng)時(shí)程的隨機(jī)組合下進(jìn)行彈塑性動(dòng)力分析，然后基于結(jié)構(gòu)地震需求(或損傷指標(biāo))和輸入?yún)?shù)進(jìn)行后續(xù)的易損性分析。當(dāng)采用多元回歸的方式構(gòu)造地震需求與輸入?yún)?shù)之間的數(shù)學(xué)關(guān)系時(shí)，可建立如下表達(dá)式

lnyD=g(lnX)+ε

(1)

式(1)實(shí)質(zhì)上是使用回歸模型來(lái)代理隱式的非線性有限元模型，它將地震需求分為了兩部分：均值部分與誤差部分。均值部分呈現(xiàn)了對(duì)平均地震需求的預(yù)測(cè)；而誤差部分則主要來(lái)源于如下兩個(gè)方面：① 有限的IM難以全面呈現(xiàn)地震動(dòng)之間的變異性；② 回歸模型本身的擬合能力難以完全代替真實(shí)的有限元模型。因此，本文建議對(duì)均值與誤差分別建立回歸模型來(lái)預(yù)測(cè)結(jié)構(gòu)的地震需求

mD=gm(X)

(2)

(3)

r2=(yD-gm(X))2

(4)

當(dāng)然，一個(gè)常數(shù)方差估計(jì)也可利用下式得到

(5)

比較式(3)和式(5)不難發(fā)現(xiàn)，常數(shù)方差僅僅平均了殘差平方，而式(3)則是一個(gè)更加靈活的回歸模型，可以捕捉到可能存在的異方差模式，即方差隨著輸入?yún)?shù)的變化而出現(xiàn)一定的趨勢(shì)。

1.2 高斯過(guò)程回歸

式(1)提供了一個(gè)多變量需求模型的框架，其中回歸模型可以是任意的。本文使用高斯過(guò)程回歸作為特定的代理模型考慮到其強(qiáng)大的非線性擬合能力以及靈活的超參數(shù)求解技術(shù)。

所謂的高斯過(guò)程回歸即是假定地震需求在所有的N組樣本上服從一個(gè)均值矢量為mN和協(xié)方差矩陣為CN的高斯過(guò)程先驗(yàn)

yD～GP(mN,CN)

(6)

在實(shí)際使用時(shí)，mN可以設(shè)置為零元素矢量(即響應(yīng)數(shù)據(jù)減去其均值)，而協(xié)方差矩陣CN∈RN×N中的每個(gè)元素可以表達(dá)為

(7)

(8)

(9)

式中,κ*為核函數(shù)矢量，具有如下形式

(10)

(11)

結(jié)合式 (9)～(11) 可以看出在新樣本上的預(yù)測(cè)即為利用既有的N組地震響應(yīng)樣本進(jìn)行加權(quán)求和，距離新樣本點(diǎn)越近，則對(duì)預(yù)測(cè)的貢獻(xiàn)越大。

1.3 多變量地震易損性函數(shù)

(12)

式中：X*是由x*按行組成的矩陣；K*則是κ*按行構(gòu)造的矩陣；下標(biāo)m和β分別用來(lái)指代均值與方差的高斯過(guò)程回歸預(yù)測(cè)。式(12)之所以是多參數(shù)的，是因?yàn)樾碌臉颖军c(diǎn)中可以包括任意的IM和結(jié)構(gòu)參數(shù)矢量。之所以稱作易損性函數(shù)，是因?yàn)镮M可以設(shè)置為逐漸增加的形式。

式(12)可直接用于橋梁構(gòu)件層次的易損性模型，此時(shí)每個(gè)構(gòu)件對(duì)應(yīng)一個(gè)地震需求矢量yD。為了獲得一個(gè)參數(shù)化的系統(tǒng)層次易損性函數(shù)，本文提議一個(gè)基于多元高斯Copula函數(shù)的方法。首先根據(jù)結(jié)構(gòu)可靠度理論，多重構(gòu)件的串聯(lián)系統(tǒng)的易損性函數(shù)[11]可以表示為

(13)

式中：m為構(gòu)件的數(shù)目；Fi,Fj,…Fm分別為各構(gòu)件的易損性函數(shù)；P(FiFj…Fm)為所有構(gòu)件之間的聯(lián)合易損性函數(shù)，其余項(xiàng)類(lèi)似。式(13)的關(guān)鍵點(diǎn)在于聯(lián)合易損性函數(shù)的求解，用高斯Copula函數(shù)可以近似的代替為

P(FiFj…Fm)=CGauss(F)=

ΦR(Φ-1(F1),…,Φ-1(Fm))

(14)

式中：CGauss(·)為高斯Copula函數(shù)；ΦR(·)為協(xié)方差矩陣等于其相關(guān)系數(shù)矩陣R的多元高斯累計(jì)分布函數(shù)；Φ-1(·)為單變量標(biāo)準(zhǔn)正態(tài)分布的逆累積分布函數(shù)。式(14)的前提假定是將各構(gòu)件的易損性函數(shù)看作某個(gè)變量的累積分布函數(shù)，這與經(jīng)典易損性模型中的假定是一致的[12]。將式(12)和式(14)代入式(13)，即可獲得系統(tǒng)層次的易損性函數(shù)，且保留了多變量和參數(shù)化的特性。應(yīng)當(dāng)注意本文僅分析了基于高斯Copula函數(shù)的串聯(lián)系統(tǒng)多元易損性模型，至于其他的系統(tǒng)形式(如并聯(lián)和混合串并聯(lián)等)，可直接擴(kuò)展式(13)；對(duì)于其他Copula函數(shù)類(lèi)型(如Clayton、Frank、Gumbel等形式)，可直接擴(kuò)展式(14)。更多的發(fā)現(xiàn)將有賴于詳細(xì)的比較研究，限于篇幅本文暫未考慮。

1.4 均值與包絡(luò)易損性曲線

將任意給定的參數(shù)組合代入式(12)～(14)即可獲得一個(gè)特定極限狀態(tài)的超越概率，因此在連續(xù)變化的參數(shù)組合下其表現(xiàn)為高維空間中的易損性曲面。為了對(duì)應(yīng)經(jīng)典的單變量易損性曲線，可通過(guò)下式的高維積分得到

Fm(im)=FIM1[FIM(X)]=

(15)

式中：EIM1[·]為向標(biāo)量IM1求取期望的算子；IM*為除去IM1以外的地震強(qiáng)度指標(biāo)；fΘ和fIM*分別為結(jié)構(gòu)參數(shù)Θ與IM*的聯(lián)合概率密度函數(shù)?？梢?jiàn)式(15)為一個(gè)以積分求期望的過(guò)程(先積分結(jié)構(gòu)參數(shù)再積分強(qiáng)度指標(biāo))，但是當(dāng)變量維度很高或者變量間有相關(guān)性時(shí)變得十分復(fù)雜，甚至不可解。因此本文提議使用如下的非線性曲線擬合方式以獲得易損性曲線

(16)

圖1顯示了式(16)的擬合過(guò)程，首先將易損性投影到F-IM1平面內(nèi)，然后通過(guò)非線性最小二乘擬合出一個(gè)參數(shù)為λim和βim的對(duì)數(shù)正態(tài)分布的累積分布函數(shù)。式(16)所求的易損性曲線實(shí)質(zhì)上是考慮了所有參數(shù)隨機(jī)性的一個(gè)平均，因此將其稱作均值易損性曲線。

圖1 均值易損性曲線與包絡(luò)易損性曲線

既然投影到F-IM1平面內(nèi)的易損性是離散的，且式(16)代表了其均值，那么其包絡(luò)易損性曲線也可以容易獲得。首先提取出離散易損性點(diǎn)的包絡(luò)點(diǎn)，對(duì)上界和下界包絡(luò)點(diǎn)分別應(yīng)用式(16)即可獲得包絡(luò)易損性曲線，如圖1所示。

2 多變量地震易損性分析步驟

本文提議的多變量地震易損性分析的具體步驟,如圖2所示。

步驟1：選擇數(shù)量為N0的地面運(yùn)動(dòng)記錄，對(duì)于具有兩個(gè)正交水平分量的記錄可取其幾何均值作為新的記錄。確定用來(lái)縮放地面運(yùn)動(dòng)的地震強(qiáng)度指標(biāo)IMs(譬如PGA或Sa(T1))以及縮放目標(biāo)區(qū)間[IMsmin,IMsmax]。

圖2 多變量地震易損性分析流程圖

步驟2：選定隨機(jī)結(jié)構(gòu)參數(shù)Θ。

步驟3：利用拉丁超立方抽樣(LHS)技術(shù)抽取行數(shù)為N的結(jié)構(gòu)參數(shù)矩陣。注意將上述IMs設(shè)置成[IMsmin,IMsmax]內(nèi)的均勻分布，和結(jié)構(gòu)參數(shù)一起進(jìn)行抽樣。并以生成的IMs樣本為目標(biāo)對(duì)N0條地震波進(jìn)行縮放。如N0

步驟4：利用生成的N個(gè)結(jié)構(gòu)參數(shù)組合與N條縮放后的地震波代入到有限元模型進(jìn)行非線性動(dòng)力分析，得到每個(gè)構(gòu)件的地震需求矢量yD。

步驟5：將yD和X進(jìn)行對(duì)數(shù)變換后進(jìn)行高斯過(guò)程回歸，分別建立均值需求模型和方差需求模型。

步驟7：假定系統(tǒng)為串聯(lián)，利用高斯Copula函數(shù)建立系統(tǒng)的多變量地震易損性模型Fsys(X*)。

步驟8：分別擬合構(gòu)件與系統(tǒng)層次的均值易損性曲線與包絡(luò)易損性曲線。

需要強(qiáng)調(diào)的是步驟3中將縮放目標(biāo)IMs與結(jié)構(gòu)參數(shù)一起進(jìn)行拉丁超立方抽樣的目的是使生成的結(jié)構(gòu)參數(shù)樣本和目標(biāo)IM保持空間填充特征，這類(lèi)似于增量動(dòng)力分析(IDA)，區(qū)別在于IDA的縮放IM是離散的，而本文是在區(qū)間[IMsmin,IMsmax]內(nèi)均勻連續(xù)分布的。因此，這種均勻縮放的方法能夠以盡量少的動(dòng)力分析次數(shù)獲得盡可能高的易損性精度，而這也得到了后文的蒙特卡羅方法的驗(yàn)證。

3 橋梁結(jié)構(gòu)與有限元模型

3.1 有限元模型與本構(gòu)關(guān)系

本文的算例為一座60 m+100 m+60 m三跨鋼筋混凝土連續(xù)梁橋，采用單柱式橋墩，兩墩等高16 m，空心矩形截面，尺寸為縱橋向3 m，橫橋向8 m，壁厚0.8 m。上部結(jié)構(gòu)為變截面連續(xù)箱梁，寬17 m，高3 m。支座為盆式橡膠支座。由于本文目的是發(fā)展易損性分析方法，為簡(jiǎn)化問(wèn)題忽略了基礎(chǔ)對(duì)橋梁動(dòng)力響應(yīng)的貢獻(xiàn)。利用地震工程仿真平臺(tái)OpenSees建立了如圖3所示的橋梁有限元模型。上部結(jié)構(gòu)采用彈性梁?jiǎn)卧M；橋墩采用基于力的纖維截面梁柱單元模擬，其中混凝土纖維的本構(gòu)采用Concrete01材料模擬，核心區(qū)混凝土采用Mander模型來(lái)計(jì)算箍筋的約束效應(yīng)，鋼筋纖維采用Steel01材料模擬。橋臺(tái)與填土之間的被動(dòng)力-位移關(guān)系使用Hyperbolic Gap材料模型，橋臺(tái)與主梁之間的沖擊則利用Impact材料模擬。對(duì)于可滑動(dòng)支座，使用理想彈塑性模型來(lái)模擬其摩擦作用。圖3(b)列舉了這些本構(gòu)模型的參數(shù)，其具體解釋請(qǐng)參照OpenSees的用戶手冊(cè)，關(guān)于參數(shù)的計(jì)算方法可參見(jiàn)文獻(xiàn)[13]。利用該有限元模型，并將第3.2節(jié)中的隨機(jī)參數(shù)設(shè)置為其均值，可以計(jì)算出橋梁的動(dòng)力特征。由計(jì)算得橋梁的前兩階模態(tài)周期分別為0.71 s和0.42 s，其振型均包含橋墩的縱向彎曲。因此為節(jié)省篇幅，本文僅分析縱橋向上的易損性，并利用這兩階振型的頻率計(jì)算Rayleigh阻尼矩陣的質(zhì)量和剛度比例系數(shù)，用于隨后的動(dòng)力計(jì)算。

(a) 橋梁有限元模型

3.2 隨機(jī)參數(shù)

既然式(12)所示的易損性函數(shù)可以包含多個(gè)變量，因此可以將橋梁參數(shù)設(shè)置為隨機(jī)變量以考慮其不確定性對(duì)易損性的影響。本文共選取21個(gè)與橋梁地震響應(yīng)有關(guān)的隨機(jī)參數(shù)，如表1所示。這些參數(shù)包括4個(gè)無(wú)約束混凝土本構(gòu)參數(shù)(fcc-峰值受壓強(qiáng)度，εcc-與峰值強(qiáng)度對(duì)應(yīng)的應(yīng)變，fcu-極限受壓強(qiáng)度，εcu-極限受壓應(yīng)變)；1個(gè)與所有混凝土構(gòu)件相關(guān)的參數(shù)ρc-混凝土質(zhì)量密度；4個(gè)受約束混凝土的本構(gòu)參數(shù)，其參數(shù)符號(hào)意義與無(wú)約束下的相同；3個(gè)與鋼筋本構(gòu)相關(guān)的參數(shù)(fys-屈服應(yīng)力，Es-彈性模量，bs-屈服后剛度比)；3個(gè)與支座摩擦有關(guān)的參數(shù)(Fy1-橋墩上支座的屈服力，F(xiàn)y2-橋臺(tái)上支座的屈服力，δy-屈服位移)；3個(gè)與橋臺(tái)-主梁沖擊有關(guān)的參數(shù)(K1-初始剛度，K2-二次剛度，xy-屈服位移)；2個(gè)與橋臺(tái)-填土被動(dòng)擠壓有關(guān)的參數(shù)(Kmax-初始剛度，F(xiàn)ult-極限被動(dòng)土壓力)；1個(gè)橋梁阻尼比參數(shù)。

表1 各隨機(jī)結(jié)構(gòu)參數(shù)的概率分布特征

3.3 構(gòu)件的損傷指標(biāo)與極限狀態(tài)

本文共對(duì)五個(gè)構(gòu)件進(jìn)行易損性分析，包括橋臺(tái)后回填土，橋臺(tái)上支座以及墩上支座，縱向固定支座下的橋墩墩底和墩頂截面。橋臺(tái)和回填土之間的峰值壓縮位移，支座的縱向峰值滑動(dòng)位移以及橋墩截面的曲率Park-Ang指標(biāo)分別作為其損傷測(cè)度。曲率形式的Park-Ang指標(biāo)定義如下

(17)

式中：φy和φu分別為等效屈服和極限曲率，My為等效屈服彎矩，上述參數(shù)可由截面的彎矩-曲率分析得到；φm是地震產(chǎn)生的截面峰值曲率；ET是由滯回耗散的能量；因子βPA用來(lái)呈現(xiàn)循環(huán)荷載對(duì)結(jié)構(gòu)損傷的影響，本文取為0.294。

共定義輕微、中等、嚴(yán)重和完全破壞四種極限狀態(tài)(limit states, LS)，每個(gè)構(gòu)件的每種極限狀態(tài)所對(duì)應(yīng)的能力概率分布如表2所示。所有構(gòu)件的抗力被假定為對(duì)數(shù)正態(tài)分布，并為輕微和中等破壞采用較低的變異系數(shù)0.25，而為嚴(yán)重和完全破壞狀態(tài)指定更高的變異系數(shù)0.5，關(guān)于確定分布參數(shù)的細(xì)節(jié)參見(jiàn)文獻(xiàn)[12-14]。

表2 各構(gòu)件極限狀態(tài)與能力的均值

3.4 地震動(dòng)記錄

本文選取了報(bào)告FEMA P695[15]中整理的50條地震動(dòng)記錄，包括22條遠(yuǎn)場(chǎng)記錄，14條近場(chǎng)無(wú)脈沖記錄，以及14條近場(chǎng)含脈沖記錄。這50條地震波將根據(jù)第2章中的設(shè)計(jì)樣本總數(shù)N進(jìn)行縮放。譬如N=500，則先將50條地震波復(fù)制10次，然后根據(jù)拉丁超立方抽樣生成的500個(gè)IMs樣本進(jìn)行縮放，這樣既保證了地震波的隨機(jī)縮放特征，又保證了和結(jié)構(gòu)參數(shù)一起的空間填充特征，因此共進(jìn)行N次非線性動(dòng)力分析。關(guān)于地震波的細(xì)節(jié)參見(jiàn)FEMA P695，本文不再列出。

對(duì)于本文，考慮到PGA更加容易縮放，將其作為縮放IMs并在[0,2g]內(nèi)采樣。為了充分考慮到其余IM對(duì)易損性精度的貢獻(xiàn)，本文增加了另外3種地震動(dòng)強(qiáng)度指標(biāo)，分別為峰值地面加速度PGV，結(jié)構(gòu)第一階自振周期處的譜加速度Sa(T1)，以及周期區(qū)間[0.2T1,2T1]內(nèi)的幾何平均譜加速度Saavg。為了呈現(xiàn)縮放后的地震波的上述IM的分布特征，圖4顯示了N=500時(shí)各IM的分布，其中PGA服從[0,2g]內(nèi)的均勻分布，這正是本文所采用的均勻縮放所期望的結(jié)果，而其余三個(gè)IM擬合成了對(duì)數(shù)正態(tài)分布。對(duì)于算例橋梁T1=0.71 s，Sa(T1)與Saavg的均值分別為1.57和1.19g，而PGV與結(jié)構(gòu)無(wú)關(guān)，均值為1.44 m/s。

4 分析結(jié)果

4.1 地震需求均值與方差模型

按照第2章中的流程即可獲得各構(gòu)件的地震需求矢量yD與輸入?yún)?shù)矩陣X，對(duì)數(shù)變換后應(yīng)用高斯過(guò)程回歸即可獲得結(jié)構(gòu)的地震需求均值與方差模型。圖5顯示了構(gòu)件1和構(gòu)件5的地震需求以及殘差平方與PGA的關(guān)系，可以看出即使經(jīng)過(guò)了對(duì)數(shù)變換地震需求與PGA之間的關(guān)系仍然呈現(xiàn)出了非線性，意味著傳統(tǒng)的線性回歸是不適用的，而高斯過(guò)程回歸則可以較好地捕捉到非線性模式。此外，構(gòu)件5的殘差明顯表現(xiàn)出了異方差模式，即隨著PGA的增加方差逐漸增大。若使用常數(shù)方差估計(jì)，顯然會(huì)在較高的PGA處低估其方差，而再次使用高斯過(guò)程回歸進(jìn)行方差預(yù)測(cè)則可以靈活地適應(yīng)方差的變化。另外需要強(qiáng)調(diào)的是構(gòu)件1的異方差模式不顯著，則高斯過(guò)程回歸直接給出了常數(shù)方差估計(jì)，即如果從殘差中學(xué)習(xí)不到任何模式，其最終會(huì)給出一個(gè)均值估計(jì)。因此，高斯過(guò)程回歸可以靈活地適應(yīng)任何方差模式，這也是本文發(fā)展需求方差模型的原因。其余構(gòu)件與此類(lèi)似，不再贅述。

4.2 樣本尺寸與模型精度驗(yàn)證

實(shí)際上，樣本尺寸會(huì)影響高斯過(guò)程回歸的預(yù)測(cè)性能，進(jìn)而影響到由此生成的多變量易損性函數(shù)的精度。為了研究其影響以及驗(yàn)證所提方法的精度，使用蒙特卡洛仿真來(lái)直接生成各構(gòu)件的易損性。具體做法為將PGA在區(qū)間內(nèi)分成以0.1g為增量的20個(gè)離散的水平，在每一個(gè)PGA水平利用第2節(jié)中的抽樣技術(shù)抽取5 000個(gè)樣本進(jìn)行非線性動(dòng)力分析，各極限狀態(tài)下的超越概率可利用失效的樣本數(shù)直接除以總數(shù)(即5 000)得到。因此基于蒙特卡洛的易損性總共實(shí)施了10萬(wàn)次動(dòng)力分析。為了測(cè)量所提方法生成的均值易損性曲線與蒙特卡洛之間的誤差，定義如下的指標(biāo)

(18)

式中：Fi為所提議方法生成的均值易損性曲線在每個(gè)PGA水平處的值，F(xiàn)mc,i則對(duì)應(yīng)蒙特卡洛計(jì)算的失效概率，因此指標(biāo)Error測(cè)量了二者在所有PGA水平上的平均易損性偏差。

圖6為各構(gòu)件與系統(tǒng)層次的易損性誤差隨樣本尺寸N的變化情況，在N=50時(shí)所有構(gòu)件和系統(tǒng)有較大的誤差，這意味著較少數(shù)量的樣本無(wú)法為高斯過(guò)程回歸提供充足的訓(xùn)練。但是當(dāng)N=100時(shí)易損性誤差迅速減小，在N=300趨于穩(wěn)定，即不再隨著N的不斷增加而減小。事實(shí)上，這些不能被縮減的誤差是模型誤差，來(lái)源于高斯過(guò)程回歸的預(yù)測(cè)，以及在生成易損性函數(shù)時(shí)的正態(tài)性假定。若想繼續(xù)減小這些誤差，需要發(fā)展更加合理的代理模型以及易損性生成方法。

圖6 易損性精度與樣本尺寸的關(guān)系

對(duì)于本文算例來(lái)說(shuō)，N=300可以認(rèn)為是一個(gè)最優(yōu)的樣本尺寸，即在保證精度的同時(shí)最大程度地降低有限元分析的代價(jià)。圖7顯示了N=300時(shí)所有構(gòu)件與極限狀態(tài)的誤差，這些誤差都在0.02以下，這對(duì)于地震作用下的高失效概率來(lái)說(shuō)是可接受的。此外，觀察到系統(tǒng)層次的易損性誤差也是較低的，表明本文所提議的基于多元高斯Copula函數(shù)的系統(tǒng)易損性函數(shù)生成方法是可靠的。至于本文所提議方法的效率，由于蒙特卡洛法的總分析次數(shù)為10萬(wàn)，因此效率比為300/(1×105)=0.003，可見(jiàn)本文所提議的方法具有非常高的計(jì)算效率。

圖7 N=300時(shí)各構(gòu)件與系統(tǒng)的易損性誤差分布

4.3 IM的影響

4.2節(jié)驗(yàn)證了所提議方法的精度與效率，即生成的均值易損性曲線收斂于蒙特卡洛的解，不過(guò)使用的是PGA，使用PGA的原因是其作為縮放IM。但是其他的IM也可以用來(lái)生成均值易損性曲線，圖8顯示了由不同IM生成的均值易損性曲線的結(jié)果，每個(gè)IM的區(qū)間由圖4推斷得到。可以看出各圖最大的區(qū)別在于易損性投影點(diǎn)的離散性，其中Saavg的最小，其次Sa(T1)和PGV，離散性最大的是PGA。越高的離散性說(shuō)明還有更重要的參數(shù)影響著易損性，而其本身的影響則相對(duì)較小。因此，Saavg是一個(gè)更好的測(cè)度用來(lái)生成均值易損性曲線。盡管這樣，Saavg仍具有一定的離散性，說(shuō)明其它參數(shù)(包括其余IM以及結(jié)構(gòu)參數(shù))的影響也是不可忽略的，這也意味著任何單個(gè)參數(shù)都不足以完全呈現(xiàn)易損性，佐證了本文發(fā)展多變量易損性函數(shù)的必要性。此外，需要強(qiáng)調(diào)的是盡管PGA的離散程度較大，但是上節(jié)證明了其均值易損性曲線是精確的，因?yàn)槠淦骄怂胁淮_定參數(shù)的影響，而離散性反映的是其不確定程度，可以由包絡(luò)易損性曲線來(lái)呈現(xiàn)。

4.4 均值和包絡(luò)易損性曲線

根據(jù)4.3節(jié)的分析，圖9顯示了所有構(gòu)件以及系統(tǒng)層次以Saavg生成的均值與包絡(luò)易損性曲線。由圖可知幾乎所有構(gòu)件的所有極限狀態(tài)的均值易損性最終都達(dá)到了1，而且包絡(luò)易損性曲線的面積也幾乎縮減至0，意味著失效是必然的，沒(méi)有任何不確定性存在。這里Saavg的區(qū)間設(shè)置為[0,4g]，這由圖4中的均值和方差推斷而來(lái)，即μ+3σ=1.19+3×0.88=3.83g，這說(shuō)明雖然用來(lái)縮放的PGA限制在[0,2g]內(nèi)，但是Saavg卻可以高達(dá)4g左右，對(duì)于本文算例來(lái)說(shuō)是很高的地震強(qiáng)度。此外，上下包絡(luò)易損性之間的距離衡量了該極限狀態(tài)易損性的不確定性程度，可以發(fā)現(xiàn)這種不確定性呈現(xiàn)了兩端小(端部為0)，中間大的特點(diǎn)，而且不同的極限狀態(tài)之間會(huì)出現(xiàn)重疊現(xiàn)象，暗示著易損性曲線表征破壞狀態(tài)之間的模糊現(xiàn)象，需要未來(lái)更深入地探索。

表3列出了所有均值與包絡(luò)易損性曲線的中位值，即當(dāng)超越概率為0.5時(shí)所需的Saavg強(qiáng)度，其值越小表明更加易損。表3顯示橋臺(tái)填土和支座比橋墩更加容易遭受地震損傷；由于固定支座下橋墩與上部結(jié)構(gòu)表現(xiàn)為剛架形式，故而上端截面也承受較高的損傷概率，但稍低于底部截面；由于本文假定橋梁系統(tǒng)為串聯(lián)形式，所以其中位值Saavg低于所有構(gòu)件。

表3 均值與包絡(luò)易損性曲線的中位值

最后需要強(qiáng)調(diào)雖然Saavg被建議用來(lái)擬合均值與包絡(luò)易損性曲線，但是其依賴于結(jié)構(gòu)的自振周期，因此PGA等與結(jié)構(gòu)無(wú)關(guān)的指標(biāo)仍然被建議用來(lái)縮放地震波，考慮到縮放與生成均值、包絡(luò)易損性曲線是相互獨(dú)立的環(huán)節(jié)。

5 結(jié) 論

本文提議了一個(gè)基于高斯過(guò)程回歸的多變量地震易損性分析框架，可以將易損性表達(dá)為由地震動(dòng)強(qiáng)度指標(biāo)和結(jié)構(gòu)參數(shù)組成的矢量輸入的多元函數(shù)。生成的均值與包絡(luò)易損性曲線能夠全面考慮地震波與結(jié)構(gòu)參數(shù)中的不確定性。以一座三跨鋼筋混凝土連續(xù)箱梁橋?yàn)樗憷?，結(jié)果表明本文所提議的方法與直接蒙特卡洛仿真具有相一致的精度，并且具有極高的效率，大大減少了有限元分析數(shù)量。通過(guò)比較研究發(fā)現(xiàn)，平均譜加速度Saavg是一個(gè)生成均值與包絡(luò)易損性曲線更加合適的指標(biāo)，考慮到其具有最小的易損性離差，而且在區(qū)間端部不確定性縮減為零。此外，研究結(jié)果顯示任何單個(gè)參數(shù)均不能完全呈現(xiàn)結(jié)構(gòu)的易損性，需要在分析中包含多重變量的影響。