基于BOX-COX變換的橋梁結(jié)構(gòu)地震易損性分析

張鵬輝， 郭軍軍， 周連緒， 袁萬城

(土木工程防災(zāi)國家重點實驗室(同濟大學(xué))， 上海 200092)

理論易損性曲線作為確定橋梁結(jié)構(gòu)潛在地震風(fēng)險的重要工具，建立了地震動強度指標(biāo)(Intensity Measure, IM)與結(jié)構(gòu)在某一損傷極限狀態(tài)(Limit State, LS)下失效概率之間的聯(lián)系[1]，而云圖法由于其計算效率上的優(yōu)越性已經(jīng)成為獲取結(jié)構(gòu)理論易損性曲線的主要途徑之一[2-3]。然而以往的研究發(fā)現(xiàn)云圖法在建立概率地震需求模型(Probabilistic Seismic Demand Model, PSDM)時，為簡化計算所作出的三個基本假設(shè)常常與現(xiàn)實相違背，存在以下三個問題，即：(1) 特定地震動強度水平下的工程需求參數(shù)(Engineering Demand Parameter, EDP)不服從對數(shù)正態(tài)分布；(2) ln(IM)與ln(EDP)間不滿足線性相關(guān)；(3) ln(EDP)的方差隨IM取值變化[4-6]。Karamlou等[5]考慮結(jié)構(gòu)不確定性，通過大量數(shù)值模擬研究了上述三個問題對理論易損性曲線計算結(jié)果的影響，研究表明與問題(1)相比后兩者對易損性曲線計算結(jié)果的影響更嚴重。為改善問題(2)，Padgett等[7]考慮高效性、實用性、充分性以及災(zāi)害的可計算性對IM的選取進行優(yōu)化，Pan等[8]采用二次函數(shù)對ln(IM)與ln(EDP)進行回歸，Mackie等[9]采用雙線性函數(shù)對ln(IM)與ln(EDP)進行回歸，均取得了較好的效果。此外，為了避免云圖法三個基本假設(shè)對易損性曲線計算結(jié)果的影響，研究人員通過增量動力分析(Incremental Dynamic Analysis, IDA)方法獲取易損性曲線，但該方法存在對地震動時程的放縮且需要的時程分析次數(shù)劇增等問題，故應(yīng)用受到限制[10]。如何改進理論易損性分析方法，使其能夠不受上述三個基本假設(shè)的限制，同時又不增加時程分析次數(shù)是當(dāng)前理論易損性分析亟待解決的問題。

本文旨在不顯著增加計算量的前提下提高易損性分析結(jié)果的可靠性，通過引入BOX-COX變換使概率地震需求模型滿足對數(shù)正態(tài)分布、線性和同方差性假設(shè)，結(jié)合蒙特卡洛抽樣得到結(jié)構(gòu)易損性曲線，并以一座三跨連續(xù)梁橋為例，驗證該方法的有效性。

1 基于BOX-COX變換的地震易損性分析方法

1.1 云圖法的基本假設(shè)和BOX-COX變換

Cornell等[11]假設(shè)EDP的均值(SD)與IM之間滿足對數(shù)線性關(guān)系，并建立普通最小二乘(Ordinary Least Square, OLS)回歸模型

(1)

式中，α1、α2均為回歸系數(shù)。因此式(1)需服從OLS回歸中的對e作出的正態(tài)性和同方差性假設(shè)，由此得到概率地震需求模型

(2)

式中：D為構(gòu)件的需求;Φ(·)為標(biāo)準(zhǔn)正態(tài)分布函數(shù)；βD為構(gòu)件需求的對數(shù)標(biāo)準(zhǔn)差。

(3)

若進一步假定構(gòu)件在某一損傷極限狀態(tài)(LSj)的抗震能力C服從均值為SC對數(shù)標(biāo)準(zhǔn)差為βC的對數(shù)正態(tài)分布，則云圖法得到的構(gòu)件在LSj下的失效概率為

(4)

由此可見，云圖法中所作的三個基本假設(shè)，本質(zhì)是為OLS回歸模型式(1)服務(wù)的。當(dāng)ln(SD)與ln(IM)之間關(guān)系非線性，或回歸的殘差e不再服從正態(tài)性和同方差假設(shè)時，回歸系數(shù)α1、α2的估計不再是最小方差線性無偏估計，式(1)的應(yīng)用效果將不理想，進而影響構(gòu)件易損性計算的可靠性[12]。

BOX-COX變換由Box等[13]提出，在回歸分析領(lǐng)域應(yīng)用廣泛，通過對因變量y=ln(EDP)進行適當(dāng)?shù)淖儞Q，可在不丟失任何信息的前提下，使回歸模型滿足線性、正態(tài)性和同方差假設(shè)

(5)

其中a為使y+a的各個分量都大于0的系數(shù)，當(dāng)min(y)>0時a=0；當(dāng)min(y)≤0時a=roundup(|min(y)|)，roundup(·)為向上取整函數(shù)[14]。引入BOX-COX變換后，式(1)更新為BOX-COX回歸模型

(6)

(7)

記x=ln(IM)，那么y=(y1,y2,…,yn)的聯(lián)合概率密度函數(shù)可表示為

f(y)=

(8)

對式(8)兩邊取對數(shù)可得極大似然函數(shù)

(9)

考慮到式(9)求偏導(dǎo)較為困難，無法得到λ的解析解，因此本文以0.01作為搜索步長，找到使L(λ)達到最大值的λ作為其數(shù)值最優(yōu)解。

1.2 基于BOX-COX變換和蒙特卡洛抽樣的地震易損性分析

(10)

則構(gòu)件在LSj下的失效概率為

(11)

值得指出的是，由于蒙特卡洛抽樣所增加的計算時間是有限的，相較于時程分析來說可忽略不計。

2 算例分析

2.1 有限元模型

本文以一座預(yù)應(yīng)力混凝土連續(xù)梁橋為例，如圖1所示?？鐝讲贾脼?0 m+100 m+60 m，主梁采用C50混凝土的預(yù)應(yīng)力箱梁；支座為球鋼支座，支座布置如圖1所示?；瑒又ёΣ烈驍?shù)取0.02，屈服位移取2 mm[17]，屈服力為最大靜摩擦力；橋墩采用C40混凝土，鋼筋種類為HRB335，縱向鋼筋配筋率0.8%，體積配箍率1%；以梅花形布置的直徑2 m鉆孔灌注樁作為基礎(chǔ)；兩側(cè)橋臺均為樁接蓋梁橋臺。使用OpenSees程序建立全橋非線性有限元模型(如圖2)。主梁采用彈性梁柱單元模擬；橋墩采用彈塑性纖維單元模擬，混凝土本構(gòu)采用Kent-Park-Scott模型[18]，不考慮混凝土的抗拉性能，鋼筋采用考慮強化的雙線性本構(gòu)。支座、樁基、碰撞均采用零長單元模擬，采用M法計算土彈簧六個方向的剛度，碰撞單元采用Muthukumar提出的本構(gòu)模型[19]，伸縮縫采用碰撞單元模擬，碰撞單元的初始縫隙根據(jù)橋臺處伸縮縫縫寬取為15 cm。由于橋臺對構(gòu)件的地震響應(yīng)有顯著影響[20]，本文參考Nielson給出的橋臺模擬方法[21]，被動土壓力方向的剛度由樁基和橋臺共同提供，主動土壓力方向剛度僅由樁基提供，以零長單元進行模擬。

圖1 橋梁結(jié)構(gòu)總體布置圖(m)

圖2 有限元模型示意圖

2.2 地震動時程

地震動時程的選取是易損性分析的關(guān)鍵之一，Baker等[22]提出的非條件譜選波方法能夠根據(jù)地震動預(yù)測模型(Ground Motion Prediction Model, GMPM)得到的場地均值譜和方差選取符合目標(biāo)反應(yīng)譜概率分布的一組波。采用該選波方法，首先按焦馳宇[23]給出的中美剪切波速換算公式得到Ⅲ類場地30 m深度范圍內(nèi)的土層平均剪切波速V30；然后根據(jù)Boore等[24]提出的地震動預(yù)測模型生成適用于中國區(qū)域震級M=7，斷層距Rjb=15 km，V30=205 m/s的反應(yīng)譜概率分布；最后從NGA-West2強震數(shù)據(jù)庫中按M=6～8，Rjb=0～100 km，V30=150～260 m/s選出100條符合目標(biāo)反應(yīng)譜分布的原場單向地震動分量。圖3畫出了它們的加速度反應(yīng)譜，圖4給出了這些地震動時程的震級和震源距分布直方圖。將這100條地震動時程按縱橋向輸入，得到結(jié)構(gòu)各構(gòu)件的響應(yīng)。

圖3 地震波的加速度反應(yīng)譜(阻尼比5%)

圖4 震級和震源距分布直方圖

2.3 概率地震需求模型建立與比較

Riddell[25]指出對于短周期結(jié)構(gòu)選用加速度相關(guān)型的地震動強度指標(biāo)(IM)是合適的，Moschonas等[26-27]采用PGA(Peak Ground Acceleration)作為IM對一系列中小跨徑梁橋和大跨度鋼構(gòu)橋進行易損性分析，并取得了較好的效果。本文模型橋一階周期0.95 s，結(jié)構(gòu)周期較短；同時PGA的地震危險性曲線目前已經(jīng)能夠得到，即具有良好的災(zāi)害可計算性，而結(jié)構(gòu)相關(guān)型的地震動強度指標(biāo)Sa(T1)的地震危險性曲線目前還未給出，故選用PGA(單位：g)作為地震動強度指標(biāo)。連續(xù)梁橋的易損部位為橋臺、支座和墩柱，故選用被動土壓力方向橋臺位移δp、主動土壓力方向橋臺位移δa、滑動支座位移δb、固定墩墩頂漂移率Dr作為工程需求參數(shù)[28]。將時程分析結(jié)果，分別按式(1)和式(6)進行OLS回歸和BOX-COX回歸得到概率地震需求模型。

為評價殘差e的正態(tài)性，作出回歸殘差的頻率分布直方圖并得到其核密度估計與按正態(tài)分布擬合的概率密度函數(shù)，核密度估計曲線(Kernel Density Estimation Curve, KDE)與正態(tài)分布概率密度曲線(Norm)的偏差情況反映了殘差e的分布與正態(tài)分布的偏離程度。自變量與因變量的線性相關(guān)性強弱通過回歸決定系數(shù)(R2)來衡量，它是兩變量線性相關(guān)系數(shù)的平方。異方差性的診斷采用等級相關(guān)系數(shù)法，計算等級相關(guān)系數(shù)(rs)首先需要對回歸模型(1)和(6)中的自變量xi和殘差絕對值|ei|分別進行排序，則

(12)

其中di對應(yīng)于xi和|ei|等級的差數(shù)。對rs進行顯著性檢驗得顯著性水平p，一般認為p<0.05時認為回歸模型存在異方差性。表1列出了OLS回歸和BOX-COX回歸的結(jié)果。

表1 OLS回歸和BOX-COX回歸的結(jié)果

由表1結(jié)果可知，當(dāng)以δp、δa和δb作為EDP建立概率地震需求模型時，均存在顯著的異方差性，通過BOX-COX變換后，異方差問題得到消除。同時BOX-COX回歸的決定系數(shù)較OLS回歸平均提升了3.26%，δp的決定系數(shù)增加達到8.2%，說明自變量和因變量的線性相關(guān)程度有所提高。為直觀顯示殘差的非正態(tài)和異方差問題，圖5作出了δp回歸殘差的分布圖和頻率分布直方圖以及KDE曲線和Norm曲線。由圖像可看出OLS回歸的殘差在PGA較小時方差較小，隨著PGA的增大逐漸增加，即OLS回歸將方差作為一個定值過程中在PGA較小時對方差進行了高估，在PGA較大時低估了方差，且殘差的分布與正態(tài)分布具有一定程度的偏離，進行BOX-COX變換后，殘差的正態(tài)性得到改善，對δa和δb分析所得結(jié)果相同。需要說明的是，當(dāng)應(yīng)用云圖法以Dr作為EDP時，檢驗結(jié)果發(fā)現(xiàn)其概率地震需求模型滿足自身的三個基本假設(shè)。

(a) 基于OLS回歸的PSDM

(b) OLS回歸殘差圖

(c) 基于BOX-COX回歸的PSDM

(d) BOX-COX回歸殘差圖

2.4 構(gòu)件易損性曲線比較

Hwang等[29]建議將橋梁的損傷描述為四個極限狀態(tài)，即輕微損傷(LS1)、中等損傷(LS2)、嚴重損傷(LS3)、完全損傷(LS4)。參考前人的研究結(jié)果，對應(yīng)不同損傷極限狀態(tài)，本文采用的能力均值和對數(shù)標(biāo)準(zhǔn)差見表2。

表2 不同損傷極限狀態(tài)的能力均值和標(biāo)準(zhǔn)差

根據(jù)式(11)可求得各構(gòu)件在不同損傷極限狀態(tài)下的易損性曲線，并與云圖法得到的構(gòu)件易損性曲線進行比較，如圖6所示。由圖6(a)～6(c)可知，因為OLS回歸在PGA較小時高估了方差，導(dǎo)致在該PGA范圍內(nèi)云圖法得到的失效概率偏大，最大偏差2.21%；而當(dāng)PGA較大時OLS回歸低估了方差，故云圖法得到的失效概率偏小，最大偏差27.2%。且極限狀態(tài)對應(yīng)的損傷程度越高，采用云圖法得到該極限狀態(tài)下的易損性曲線偏離實際越明顯。當(dāng)PGA=1g時，δp達到LS3狀態(tài)、δa達到LS2狀態(tài)的和δb達到LS2狀態(tài)的概率分別被低估18.1%、524.7%、30.9%，若仍然采用云圖法的三個基本假設(shè)，將導(dǎo)致易損性計算結(jié)果失真。由圖6(d)可知，當(dāng)采用云圖法能夠滿足自身基本假設(shè)時，基于BOX-COX變換的易損性分析方法(采用BOX-COX回歸建立需求模型)所得易損性曲線與云圖法(采用OLS回歸建立需求模型)計算結(jié)果相同，即BOX-COX變換不會導(dǎo)致信息的丟失。

3 結(jié) 論

(1) 提出了結(jié)合BOX-COX變換和蒙特卡洛抽樣的地震易損性分析方法，消除由于云圖法的三個基本假設(shè)帶來的易損性計算誤差，并通過算例驗證了其應(yīng)用于連續(xù)梁橋地震易損性分析的有效性。

(a) 被動土壓力方向橋臺位移

(b) 主動土壓力方向橋臺位移

(c) 滑動支座位移

(d) 固定墩墩頂漂移率

(2) 采用BOX-COX變換建立的概率地震需求模型，無須遵循云圖法所特有的三個基本假設(shè)，且改善了云圖法特有的三個基本假設(shè)帶來的易損性計算誤差，使模型的線性程度提高，正態(tài)性改善，異方差性消除；同時，保持了云圖法非線性時程分析次數(shù)少的優(yōu)點。

(3) 在云圖法適用的情況下，采用本文提出的方法計算的易損性曲線與云圖法結(jié)果一致；當(dāng)概率地震需求模型存在異方差性時，使用云圖法將導(dǎo)致結(jié)構(gòu)易損性曲線計算結(jié)果的失真，且極限狀態(tài)對應(yīng)的損傷程度越高，該極限狀態(tài)下的易損性曲線偏離實際越明顯。