徐 強(qiáng) 姚文斌 陳健云

(1.大連理工大學(xué) 海岸與近海工程國家重點(diǎn)實(shí)驗(yàn)室, 遼寧 大連 116024;2.大連理工大學(xué) 工程抗震研究所,遼寧 大連 116024)

我國水能資源豐富,尤其是西南地區(qū),具有較高的開發(fā)利用前景.因此,我國在該區(qū)域修建了數(shù)座百米級高拱壩,在灌溉、防洪、供水、發(fā)電、航運(yùn)及水產(chǎn)養(yǎng)殖等領(lǐng)域發(fā)揮著重要作用.但是,大壩作為主要的擋水建筑物,一旦發(fā)生破壞,將直接影響下游人民的生命和財(cái)產(chǎn)安全.近年來,我國唐山、汶川等地發(fā)生了烈度較高的強(qiáng)地震動,超過了部分建筑物的設(shè)計(jì)地震烈度,使當(dāng)?shù)氐母黝惤ㄖ镌馐芰瞬煌潭鹊钠茐?地震動作為大型結(jié)構(gòu)損傷破壞的重要因素,研究拱壩在地震動作用下的損傷破壞機(jī)理是有必要的.

對于特高拱壩這種大型混凝土結(jié)構(gòu),由于其結(jié)構(gòu)復(fù)雜,探究壩體在不同地震動強(qiáng)度下的非線性響應(yīng)時,需要反復(fù)進(jìn)行加速度時程的調(diào)幅計(jì)算.然而這一過程十分耗時,而且計(jì)算結(jié)果可能無法收斂.為了解決這一問題,Estekanchi[1]提出了一種耐震時程法(endurance time analysis,ETA),詳細(xì)闡述了耐震時程加速度曲線的合成方法,并驗(yàn)證了該方法在保證計(jì)算精度的同時,能有效減少計(jì)算量.白久林、楊樂等[2]基于我國的抗震反應(yīng)譜,合成了ETA 時程曲線,并將其與增量動力分析方法(IDA)進(jìn)行了對比,結(jié)果驗(yàn)證了ETA 方法的可行性.

地震的危險(xiǎn)性分析、易損性分析和災(zāi)害損失評估是預(yù)測地震風(fēng)險(xiǎn)的3個主要方面[3].其中,易損性分析作為評估結(jié)構(gòu)安全性的核心內(nèi)容,近年來在研究大壩的抗震性能方面被廣泛關(guān)注.大壩地震易損性分析是為了預(yù)測大壩受不同強(qiáng)度地震動作用時,能滿足各級性能水準(zhǔn)要求的可能性.地震易損性分析按數(shù)據(jù)資料的來源可分為經(jīng)驗(yàn)法與解析法.經(jīng)驗(yàn)法就是利用現(xiàn)有地震記錄及震害資料分析易損性,解析法主要通過多組數(shù)值計(jì)算模擬地震中結(jié)構(gòu)的反應(yīng).對高拱壩進(jìn)行易損性分析可以發(fā)現(xiàn)結(jié)構(gòu)在設(shè)計(jì)中存在的薄弱環(huán)節(jié),為后期的除險(xiǎn)加固提供依據(jù)[4-6].

國內(nèi)學(xué)者曾對拱壩的地震易損性進(jìn)行了研究.金愛云[7]分別采用了云分析、多條帶分析(MSA)以及增量動力分析方法研究了高拱壩在地震動作用下的易損性曲線,結(jié)果表明多條帶分析方法的可靠度最高;鐘紅等[8]考慮了多種參數(shù)的不確定性,定義了拱壩的5級極限狀態(tài),分別為無損傷、輕微損傷、中度損傷、嚴(yán)重?fù)p傷和破壞.姚霄雯和蔣建群[3]等考慮了拱壩的材料非線性和橫縫接觸非線性,分析了塑性損傷、橫縫數(shù)量對拱壩動力響應(yīng)的影響,以及15種地震強(qiáng)度指標(biāo)與拱壩動力響應(yīng)的相關(guān)性.最終,他們以壩頂位移和橫縫開度為指標(biāo),提出了基于性能的高拱壩抗震安全評價(jià)方法.王進(jìn)廷等[9]以大崗山拱壩為例,考慮材料不確定性,定義了幾種極限狀態(tài),分別為輕微破壞、中等破壞和嚴(yán)重破壞,基于拱壩的損傷程度以及是否發(fā)生貫穿性損傷.章明旭[5]比較了直接擬合法和概率分布法用于擬合地震易損性曲線時的差異.范書立等[10-11]基于向量地震強(qiáng)度參數(shù)和多種工程需求參數(shù),建立了以拱冠位移、橫縫開度和損傷體積比為指標(biāo)的地震易損性曲線和曲面,采用了基于響應(yīng)面的方法進(jìn)行了拱壩的地震易損性分析,研究結(jié)果表明響應(yīng)面分析比傳統(tǒng)的蒙特卡羅模擬法效率更高.這些研究成果為高拱壩抗震設(shè)計(jì)和評估提供了科學(xué)依據(jù).

基于上述理論方法,本文根據(jù)水工結(jié)構(gòu)設(shè)計(jì)規(guī)范[12]中的標(biāo)準(zhǔn)設(shè)計(jì)反應(yīng)譜生成了20 組耐震時程加速度曲線;建立了白鶴灘拱壩-地基系統(tǒng)有限元模型,進(jìn)行非線性時程分析,分析結(jié)構(gòu)的損傷特點(diǎn);選取相對位移、損傷體積比和最大橫縫開度為指標(biāo),提取各指標(biāo)時程,并分別基于MSA 法、回歸曲線法、矩估計(jì)法及截?cái)嘧畲笏迫还烙?jì)法對結(jié)構(gòu)的易損性進(jìn)行了分析,得出易損性曲線.通過對高拱壩進(jìn)行抗震安全評價(jià)及失效概率分析,可為結(jié)構(gòu)的抗震設(shè)計(jì)提供可靠依據(jù).

1 耐震時程法

耐震時程法是一種動力推覆過程,其主要優(yōu)點(diǎn)在于保證預(yù)測結(jié)構(gòu)響應(yīng)精度的同時,可以減少計(jì)算量[13-14].該方法的核心思想是生成ETA 加速度時程曲線,使得不同時間段對應(yīng)的反應(yīng)譜呈現(xiàn)既定的倍數(shù)擬合關(guān)系.

ETA 加速度時程的生成過程是先生成平穩(wěn)加速度時程,然后通過反復(fù)迭代優(yōu)化,使得加速度反應(yīng)譜與目標(biāo)反應(yīng)譜的擬合程度更佳.在地震動持續(xù)時間內(nèi),不同時間段對應(yīng)的目標(biāo)反應(yīng)譜與目標(biāo)反應(yīng)譜的倍數(shù)關(guān)系為:

式中:tTarget為目標(biāo)時間;SaC(T)為標(biāo)準(zhǔn)設(shè)計(jì)反應(yīng)譜;SaT(T,t)為目標(biāo)反應(yīng)譜;t為任意時間.

根據(jù)加速度反應(yīng)譜與位移反應(yīng)譜的對應(yīng)關(guān)系可以推導(dǎo)出位移譜的公式為:

式中:SuT(T,t)為時刻t的目標(biāo)位移反應(yīng)譜.

但是,想讓加速度時程上的任一時刻前對應(yīng)的反應(yīng)譜都能以目標(biāo)譜完美擬合是不可能的,因此,可以采用式(3)對加速度時程進(jìn)行無約束優(yōu)化,使其對應(yīng)的反應(yīng)譜與目標(biāo)譜更加吻合:

式中:ag為初始生成的加速度時程曲線;Tmax是地震動時長;α是位移譜的權(quán)重系數(shù),本文只考慮加速度反應(yīng)譜的影響,因此α值取0;Sa(T,t)和Su(T,t)分別為周期T下0～t時刻的加速度反應(yīng)譜和位移反應(yīng)譜.

基于上述理論,根據(jù)標(biāo)準(zhǔn)設(shè)計(jì)反應(yīng)譜,共生成20組隨機(jī)ETA 地震動時程,每組包含3條地震動時程,分別作為X、Y和Z向的動力輸入.其中,X向和Y向的最大峰值加速度調(diào)幅為0.6g,Z向最大峰值加速度調(diào)幅為0.4g.圖1展示了ETA 時程和IDA 時程的對比圖.從圖中可以觀察到,ETA 時程主要具有峰值加速度不斷增大的特點(diǎn).分別對ETA 時程0～5 s、0～10 s、0～15 s和0～20 s各時段時程進(jìn)行反應(yīng)譜優(yōu)化擬合,擬合效果如圖2所示,可以看出,不同時段的反應(yīng)譜與標(biāo)準(zhǔn)譜呈倍數(shù)關(guān)系.

圖1 ETA 與IDA 加速度時程對比

圖2 ETA 時程加速度反應(yīng)譜

2 結(jié)構(gòu)地震易損性分析方法

結(jié)構(gòu)地震易損性是指結(jié)構(gòu)在某一強(qiáng)度地震動作用下,結(jié)構(gòu)達(dá)到極限狀態(tài)的概率[15],采用下式對其進(jìn)行描述:

式中:F(im)為易損性概率;LS為指定的破壞指標(biāo)極限狀態(tài),本文選取相對位移、最大橫縫開度和損傷體積比作為破壞指標(biāo);IM為表征地震動強(qiáng)度的指標(biāo),本文選取峰值加速度PGA 為地震動強(qiáng)度指標(biāo);P[LS|IM=im]為地震動參數(shù)IM=im時,結(jié)構(gòu)達(dá)到極限狀態(tài)的概率.

2.1 基于回歸曲線的易損性分析方法

假設(shè)地震動需求(Dd)與IM呈對數(shù)相關(guān)性,采用直線擬合lnDd與lnIM之間的關(guān)系:

式中:a、b表示回歸參數(shù);Dd表示地震動需求,也就是本文選取的損傷指標(biāo).

lnDd在每個地震動強(qiáng)度IM水平下,具有恒定的方差,其標(biāo)準(zhǔn)差可表示為:

2.2 基于矩估計(jì)的易損性分析方法

用含兩個參數(shù)的對數(shù)正態(tài)累計(jì)分布函數(shù)來定義易損性函數(shù),其表達(dá)式為[17]:

式中:θ為在超越概率為50%的IM水平下易損性函數(shù)的中位數(shù);β表示lnIM的標(biāo)準(zhǔn)差.通過對任意工況下結(jié)構(gòu)開始發(fā)生破壞的IM值取對數(shù)來對易損性函數(shù)的參數(shù)進(jìn)行估計(jì),其均值和標(biāo)準(zhǔn)差分別為[18]:

式中:n為工況總數(shù)量;IMi為第i個工況下結(jié)構(gòu)開始發(fā)生破壞的PGA 值.

2.3 基于截?cái)嘧畲笏迫还烙?jì)的易損性分析方法

最大似然估計(jì)法的原理是通過使似然函數(shù)達(dá)到最大值來對參數(shù)進(jìn)行估計(jì).通過設(shè)定所考慮的最大地震動強(qiáng)度IMmax,假設(shè)在進(jìn)行的n次地震動分析中,有m次結(jié)構(gòu)在低于IMmax的IM水平下發(fā)生破壞,通過式(12)對θ和β進(jìn)行估計(jì):

式 中:^θ和^β是 使 得 該 函 數(shù) 取 得 最 大 值 的θ和β值;IMi為第i個工況下結(jié)構(gòu)開始發(fā)生破壞的PGA 值.該函數(shù)需要借助MATLAB 中的fmincon函數(shù)進(jìn)行求解.將參數(shù)估計(jì)值代入式(9)中求解易損性曲線.

2.4基于多樣條分析(MSA)的易損性分析方法

MSA 分析方法通過分析每一個IM水平下,損傷指標(biāo)超過設(shè)定閾值的工況所占總工況比例,基于最大似然估計(jì)對參數(shù)進(jìn)行估計(jì)[19]:

式中:m表示不同地震動強(qiáng)度的個數(shù);xj表示逐漸增大的第j個地震動強(qiáng)度;zj和nj分別表示工況總數(shù)和在地震動強(qiáng)度為xj條件下導(dǎo)致結(jié)構(gòu)發(fā)生破壞的工況數(shù).式(13)的求解過程同式(12).

3 數(shù)值模擬

3.1 有限元模型

以白鶴灘拱壩為研究對象,其壩高為289 m,壩頂寬度13 m,最大壩底寬度72 m,壩頂高程為834 m,壩頂中心線弧長為709 m.壩體-地基系統(tǒng)有限元模型如圖3所示;壩體共設(shè)有30條橫縫,橫縫分布如圖4所示.

圖3 白鶴灘拱壩-地基系統(tǒng)有限元模型

圖4 白鶴灘拱壩橫縫分布

分析中考慮的靜態(tài)荷載有壩體自重和靜水壓力.壩體混凝土密度為2 400 kg/m3,地基密度為2 800 kg/m3.正常蓄水情況下,壩體上游靜水壓力水位為825 m,下游靜水壓力水位為604 m.壩體動水壓力按照Westergaard附加質(zhì)量添加,采用附加質(zhì)量的形式模擬不可壓縮水體對結(jié)構(gòu)的動水壓力.地震動采用ETA 生成的耐震加速度時程輸入.

采用接觸邊界模型模擬混凝土拱壩的橫縫,每個橫縫接觸面包含主面和從面,接觸關(guān)系[18]如圖5所示,法向接觸壓力的計(jì)算為:

圖5 主從面接觸模型

式中:初始接觸距離c取100 mm,初始壓力p0取0.3 MPa.

為模擬地基的彈性恢復(fù)能力,采用黏彈性人工邊界[20]來實(shí)現(xiàn)有限地基對無限地基的模擬作用,該人工邊界還能吸收邊界輻射能量,模擬效果良好.示意圖如圖6所示.人工邊界上,節(jié)點(diǎn)上的參數(shù)由式(15)計(jì)算給出:

圖6 黏彈性人工邊界示意圖

式中:Kn和Ks分別表示彈簧法向和切向的彈性剛度;Cn和Cs分別表示法向和切向的阻尼;Cp和Cs表示地震波的傳播速度;ρ為地基密度;r為震源距;λ和G為拉梅常數(shù);a和b為修正系數(shù),a取0.8,b取1.1.

壩體混凝土材料和地基巖石的其他力學(xué)性能參數(shù)分別見表1和表2.

表1 壩體混凝土材料參數(shù)

表2 地基基本參數(shù)

3.2 動力響應(yīng)分析

對白鶴灘拱壩進(jìn)行動力學(xué)分析.圖7展示了壩體結(jié)構(gòu)在峰值加速度水平(PGA)分別為0.2g、0.4g和0.6g時,拱冠梁、下游面和上游面的拉應(yīng)力損傷分布情況.從圖中可以觀察出,隨著PGA 的增大,損傷成明顯的加重趨勢.在PGA 為0.2g時,只有壩踵和下游面下部存在局部輕微拉應(yīng)力損傷;PGA 達(dá)到0.4g時,壩踵處損傷加重,下游面多處出現(xiàn)大面積損傷,但損傷值均為超過0.5,上游面損傷加重不明顯;當(dāng)PGA 達(dá)到0.6g時,壩踵處損傷進(jìn)一步加重,下游面出現(xiàn)大面積嚴(yán)重拉應(yīng)力損傷,上游面也在多處呈現(xiàn)不同程度損傷.拱冠梁的拉應(yīng)力損傷云圖也呈現(xiàn)明顯加重的趨勢.

圖7 不同PGA 水平下拉應(yīng)力損傷云圖

圖8 各指標(biāo)隨PGA 的變化時程

通過對損傷云圖的分析,在遭遇地震動作用時,最先被破壞的是壩踵部位,若PGA 繼續(xù)增大,下游面會出現(xiàn)損傷,上游面需要在更大的PGA 水平地震動作用下才會出現(xiàn)損傷.此時,壩體在中上部出現(xiàn)貫穿損傷,結(jié)構(gòu)存在失效風(fēng)險(xiǎn).

3.3 損傷指標(biāo)分析

采用最大橫縫開度、損傷體積比和相對位移作為損傷評價(jià)指標(biāo).計(jì)算結(jié)果顯示,隨著PGA 值增大,各指標(biāo)呈現(xiàn)明顯的增大趨勢.最大橫縫開度是取各時刻30條橫縫的最大開度作為該時刻的開度值.最大橫縫開度隨PGA 的增大平穩(wěn)增長,PGA 達(dá)到0.6g時,各工況最大橫縫開度最大值為9 cm,最小值約4.7 cm,平均值約為7 cm;損傷體積比[21]可以宏觀地衡量損傷單元所占的比例,損傷體積在PGA 較小時增速較為平緩,但隨著PGA 的增大旗增速逐漸加快,峰值加速度0.6g時,各工況損傷體積比最大值為0.28,最小值約0.09,平均值約為0.17;順河向的相對位移是指拱冠梁底部和頂部之間的順河向相對位移,與損傷體積比相反,在PGA 較小時增速較快,但隨著PGA 的增大旗增速逐漸減緩,峰值加速度0.6g時,各工況順河向相對位移最大值為0.34 m,最小值約0.18 m,平均值約為0.26 m.

為分析各個工況在不同PGA 水平下的離散性,分別計(jì)算了各指標(biāo)在不同峰值加速度水平下的標(biāo)準(zhǔn)差,計(jì)算結(jié)果如圖9所示.從圖中可以觀察到,標(biāo)準(zhǔn)差的整體走勢是隨著PGA 的增大逐漸增大.說明PGA越大,各工況離散性越高.但各指標(biāo)的標(biāo)準(zhǔn)差大小和增減趨勢均有所差異,PGA 在小于0.4g情況下,損傷體積比和最大橫縫開度各工況的標(biāo)準(zhǔn)差較小,但在0.4g后,各工況損傷體積比標(biāo)準(zhǔn)差迅速增大.順河向相對位移各工況標(biāo)準(zhǔn)差始終處于較高水平,說明各工況計(jì)算結(jié)果有較大差異性.

圖9 各指標(biāo)標(biāo)準(zhǔn)差隨PGA 變化曲線

從離散性角度出發(fā),最大橫縫開度始終保持較小的標(biāo)準(zhǔn)差,更適合作為評判拱壩結(jié)構(gòu)損傷破壞程度的指標(biāo).

3.4 易損性分析

本節(jié)以上文介紹的易損性分析方法,基于最大橫縫開度、損傷體積比和順河向相對位移等指標(biāo),根據(jù)以往研究成果[5,10],結(jié)合實(shí)際分析結(jié)果設(shè)定損傷閾值,對拱壩進(jìn)行易損性分析,并評價(jià)了各方法與ETA方法相結(jié)合的分析效果.

3.4.1 最大橫縫開度

以最大橫縫開度為損傷指標(biāo),設(shè)定其損傷閾值為3 cm,分別采用截?cái)嘧畲笏迫还烙?jì)法、MSA、矩估計(jì)法和回歸曲線擬合法分析指標(biāo)在各PGA 水平下達(dá)到該閾值的概率,計(jì)算結(jié)果如圖10所示.

圖10 以最大橫縫開度為指標(biāo)的易損性曲線

分析結(jié)果顯示,矩估計(jì)與MSA 方法具有較高的一致性;通過截?cái)嘧畲笏迫还烙?jì)計(jì)算的易損性曲線在PGA 小于0.3g時偏大,0.35g之后與其他方法計(jì)算結(jié)果較為接近;通過一次函數(shù)擬合方法計(jì)算的易損性曲線在0.2g～0.4g范圍內(nèi)偏小.

3.4.2 損傷體積比

以損傷體積比為損傷指標(biāo),設(shè)定其損傷閾值為0.02,通過各分析方法計(jì)算的易損性曲線如圖11所示.其中,矩估計(jì)法和MSA 法依然保持了較好的一致性,一次函數(shù)擬合和截?cái)嘧畲笏迫还烙?jì)法計(jì)算結(jié)果偏大.

圖11 以損傷體積比為指標(biāo)的易損性曲線

3.4.3 順河向相對位移

以順河向相對位移為損傷指標(biāo),設(shè)定其損傷閾值為14 cm.易損性計(jì)算結(jié)果如圖12所示.同樣,矩估計(jì)和MSA 方法的計(jì)算結(jié)果高度一致,截?cái)嘧畲笏迫还烙?jì)的計(jì)算結(jié)果最大,一次函數(shù)擬合方法的結(jié)果稍大,且曲線局部出現(xiàn)震蕩現(xiàn)象.

圖12 以順河向相對位移為指標(biāo)的易損性曲線

通過對上述3個損傷指標(biāo)進(jìn)行易損性分析,發(fā)現(xiàn)矩估計(jì)和MSA 方法可以較好地與ETA 方法相結(jié)合,對拱壩的損傷和易損性進(jìn)行評估.原因是上述兩種方法均依賴于給定PGA 水平下拱壩結(jié)構(gòu)的響應(yīng),且給定的PGA 水平越多,曲線精度越高.本文采用的耐震時程法可以給定任意峰值加速度水平下結(jié)構(gòu)的響應(yīng),因此,耐震時程法與矩估計(jì)和MSA 法相結(jié)合可以以理想的精度計(jì)算易損性曲線;雖然一次函數(shù)擬合法同樣依賴于給定PGA 水平的數(shù)量,但給定更多PGA 水平下的結(jié)構(gòu)響應(yīng)并不能有效提高易損性曲線的精度,因此,分析結(jié)果往往有一定偏差,且曲線局部易出現(xiàn)震蕩現(xiàn)象;截?cái)嘧畲笏迫还烙?jì)需要預(yù)先給定地震動強(qiáng)度閾值,這一閾值的選取沒有確定方法,且該閾值往往對易損性分析結(jié)果產(chǎn)生影響.

4 結(jié) 論

本文建立了白鶴灘拱壩-地基有限元模型,基于耐震時程法理論,生成了耐震時程加速度曲線,作為動力輸入,作用于白鶴灘拱壩,進(jìn)行動力時程分析.對多個PGA 水平下拱壩結(jié)構(gòu)的拉應(yīng)力損傷分布進(jìn)行分析,以最大橫縫開度、損傷體積比和順河向相對位移為損傷指標(biāo),分別基于截?cái)嘧畲笏迫还烙?jì)法、MSA、矩估計(jì)法和一次函數(shù)擬合法分析了拱壩的易損性,主要得出以下結(jié)論:

1)拱壩在地震動作用下的損傷破壞主要集中在壩踵和下游面中上部,隨著地震動強(qiáng)度的不斷增強(qiáng),上游面出現(xiàn)損傷,并出現(xiàn)貫穿損傷.

2)標(biāo)準(zhǔn)差能反映各工況的離散性,在最大橫縫開度、損傷體積比和順河向相對位移3個損傷指標(biāo)中,隨著PGA 的增大,標(biāo)準(zhǔn)差均呈現(xiàn)增大趨勢,但最大橫縫開度具有相對較小的標(biāo)準(zhǔn)差,更適合作為評價(jià)指標(biāo).

3)由于計(jì)算原理不同,MSA 法和矩估計(jì)法更適合與ETA 法相結(jié)合,計(jì)算結(jié)構(gòu)在不同PGA 水平下的易損性概率.