陳燈紅 曾怡陽 趙藝園 劉云龍 岳 夢

(1.防災(zāi)減災(zāi)湖北省重點(diǎn)實(shí)驗(yàn)室(三峽大學(xué)), 湖北 宜昌 443002;2.三峽大學(xué) 土木與建筑學(xué)院, 湖北 宜昌443002)

在大壩-地基動力相互作用的研究中,工程上最常用的是無質(zhì)量地基模型[1-2].但地基本身是有質(zhì)量的半無限體,并沒有固定邊界,地震波的能量會被地基輻射吸收,若使用無質(zhì)量地基則無法考慮遠(yuǎn)域地基的輻射阻尼效應(yīng)[3].在地震波輸入方面,現(xiàn)階段用于模擬結(jié)構(gòu)地震響應(yīng)時一般采用垂直入射的方式[4],這種方式對于模擬結(jié)構(gòu)的遠(yuǎn)場地震響應(yīng)較為合適.但是,地震波傳播到壩體結(jié)構(gòu)時的入射波類型及角度存在一定的不確定性,對于距離結(jié)構(gòu)較近的淺源地震來說一般不是垂直入射的.

為了更加真實(shí)地反映壩體在地震作用下的響應(yīng),國內(nèi)外學(xué)者針對地震波斜輸入問題進(jìn)行了研究.在地下結(jié)構(gòu)的抗震設(shè)計與計算分析中,杜修力等[5-6]研究了不同地震波斜入射下地下結(jié)構(gòu)的時域地震動反應(yīng),建立了半空間場地的平面波傾斜輸入方法.趙密等[7-8]提出并采用了一種結(jié)構(gòu)-地基動力相互作用的地震波輸入方法,以此模擬P波或SV 波以不同角度在成層介質(zhì)中的輸入.周曉潔等[9]實(shí)現(xiàn)了地震波斜入射下層狀場地地下綜合管廊地震反應(yīng)分析,建立了不同場地條件下地下綜合管廊分析模型,分析了SV 波斜入射下的綜合管廊結(jié)構(gòu)地震響應(yīng).陳燈紅等[10]利用Fortran語言對Abaqus軟件進(jìn)行二次開發(fā),將黏彈性邊界有效地嵌入到Abaqus中,通過數(shù)值算例驗(yàn)證了程序的正確性.孫靖云等[11]研究了地震荷載作用下隧道內(nèi)油氣管道應(yīng)力響應(yīng),說明了軸向地震作用下管道應(yīng)力、位移較大.劉述虹等[12]進(jìn)行了典型綜合管廊體系的抗震性能數(shù)值模擬研究,分析了綜合管廊這一淺埋地下構(gòu)筑物的特性及其震害特點(diǎn).Heymsfield[13]采用直接邊界積分方程法,研究了不同角度下SH 波斜入射時傾斜基巖對地表位移的影響.周鵬等[14]建立了海水-沉管隧道-海床的整體分析模型,基于該模型分析了地震P波斜入射作用下,不同入射角度對沉管隧道結(jié)構(gòu)的動力響應(yīng)的影響.

在水工結(jié)構(gòu)抗震領(lǐng)域中,苑舉衛(wèi)等[15]將地表地震動時程分量分解為斜入射的平面SV 波和P波,證明了地震波斜入射對重力壩結(jié)構(gòu)影響顯著,尤其在壩-基交界面上.朱代富等[16]研究了不同角度地震波入射下的壩體動力響應(yīng),證明了動力分析時考慮地震波同時斜入射的必要性.岑威鈞等[17]分別研究了SV波和P波斜入射下高面板堆石壩的地震反應(yīng)特性,證明了常規(guī)垂直入射結(jié)果相較斜入射偏于保守.杜修力等[18]結(jié)合小灣拱壩研究了SV 波斜入射下拱壩的地震響應(yīng),表明入射角度對低頻區(qū)幅值放大系數(shù)影響顯著.孫奔博等[19]研究了平面SV 波斜入射下重力壩的地震響應(yīng),證明入射角度的存在對壩體位移和應(yīng)力有顯著影響.Maeso等[20]基于動力分析的邊界有限元模型,研究了不同的簡諧波斜入射下拱壩的動力響應(yīng).趙密等[21]提出了一種用于近場波動有限元分析的高階精度人工邊界條件,表明了提出的高階精度人工邊界條件具有精確性、穩(wěn)定性和高效性.王剛等[22]將時間顯式積分動力有限元方法與黏彈性動力人工邊界理論相結(jié)合,采用動力學(xué)方法求解靜力作用效應(yīng),求解了含接觸非線性問題的重力壩靜、動力荷載組合作用的響應(yīng).值得注意的是,上述學(xué)者大多針對地震波斜入射情況下的結(jié)構(gòu)動力響應(yīng)進(jìn)行了重點(diǎn)研究,但在重力壩靜動力響應(yīng)影響的規(guī)律方面還有待深入研究.

基于上述學(xué)者對地震波斜入射的研究,本文對該問題進(jìn)行深入研究,探究P波和SV 波斜入射對重力壩靜動力響應(yīng)影響的規(guī)律.在對結(jié)構(gòu)進(jìn)行地震響應(yīng)分析時考慮初始應(yīng)力場[23-24](不僅考慮結(jié)構(gòu)的動力響應(yīng)情況,還綜合考慮了靜-動力分析問題),分別推導(dǎo)了P波和SV 波傾斜入射時的等效地震力,構(gòu)建了考慮地震波斜入射的黏彈性靜-動力統(tǒng)一型人工邊界.通過用戶子程序接口編寫二維地震波的波動輸入程序,基于有限元軟件Abaqus完成地震動輸入,以黃登混凝土重力壩為研究對象,研究P波和SV 波在不同斜入射角度下混凝土壩的動力響應(yīng).

1 研究方法

1.1 黏彈性人工邊界

黏彈性人工邊界是一種局部人工邊界,它有效地解決了黏性邊界的低頻失穩(wěn)問題,能夠較好地模擬無限地基的彈性恢復(fù)能力和輻射阻尼效應(yīng)[25].基于Abaqus軟件平臺在計算模型的人工邊界節(jié)點(diǎn)的每個方向施加一端固定的單向并聯(lián)的彈簧和阻尼器元件,以此完成人工邊界的施加,對應(yīng)的公式如下所示:

二維黏彈性人工邊界等效物理系統(tǒng)的彈簧系數(shù)和阻尼系數(shù)分別為[4]:

切向邊界

法向邊界

式中:R為散射波源至人工邊界點(diǎn)的距離;cs、cp為SV 波和P 波波速:G為介質(zhì)剪切模量;ρ為介質(zhì)密度;αT、αN分別為黏彈性人工邊界切向、法向參數(shù),其中αT取0.5,αN取1.施加彈簧-阻尼元件時只需將式(1)和式(2)乘以邊界結(jié)點(diǎn)影響面積并施加在該結(jié)點(diǎn)上即可(如圖1所示).

圖1 黏彈性邊界示意圖

對于黏彈性人工邊界的外源波動輸入問題,可將總波場uT分解為散射波場uS和自由波場uF:

由位移平衡條件與力學(xué)平衡條件,人工邊界上任意一點(diǎn)n的運(yùn)動方程為[26]:

式中:A為結(jié)點(diǎn)影響面積;(t)為自由場應(yīng)力;K nj、C nj分別為邊界結(jié)點(diǎn)黏彈性組件的彈簧系數(shù)和阻尼系數(shù).

1.2 基于黏彈性邊界的二維斜入射波動輸入方法

1)P波斜入射等效荷載

假設(shè)P波以任意角α入射,在入射區(qū)域內(nèi)的任意一點(diǎn)(x,y)的各向位移時程可記為up(t)(計算區(qū)域如圖2所示),圖中反射P 波與反射SV 波的反射角分別為α、β,對應(yīng)的波幅系數(shù)分別為A1、A2,左邊人工邊界高度為H,下邊界和地表面長度為L.

圖2 P波斜入射示意圖

由Snell定律可得如下關(guān)系式[27]:

式中:cs、cp分別為P波和SV 波波速,可由介質(zhì)特性參數(shù)求得.左邊界處入射和反射的P 波、反射的SV波相對于入射波up(t)零時刻波陣面的延遲時間分別為

式中:H為左邊界的高度;y為豎向坐標(biāo)值.

左側(cè)邊界的內(nèi)行場位移由直接入射的P波、反射P波和反射SV 波構(gòu)成,邊界節(jié)點(diǎn)的位移公式如下:

式中:u lx(t)、u ly(t)分別為水平位移、豎直位移;Δt1、Δt2、Δt3分 別 為 入 射P 波、反 射P 波、反 射SV 波 相對于波陣面的時間延遲.由式(8)求導(dǎo)或差分計算可得對應(yīng)的速度場.

計算自由波場傳播應(yīng)力需引入局部坐標(biāo)系(ξ,η),其中ξ為波的傳播方向;η為ξ的法線方向,在均勻彈性介質(zhì)中平面P波ξ方向的傳播應(yīng)力可表示為[28]:

對于平面SV 波在均勻彈性介質(zhì)中η方向的傳播應(yīng)力可以表示為:

由式(9)、(10)局部坐標(biāo)系將波場的傳播應(yīng)力轉(zhuǎn)化為全局坐標(biāo)系,求得左側(cè)人工邊界上的應(yīng)力.結(jié)合式(7)～(10)得出左邊界的等效節(jié)點(diǎn)力如式(11):

地基截斷邊界的右邊界、底邊界處節(jié)點(diǎn)的等效荷載計算方法與左邊界類似.

2)SV 波斜入射等效荷載

假設(shè)SV 波以任意角α1入射,在入射區(qū)域內(nèi)的任意一點(diǎn)(x,y)的各向位移時程可記為usv(t)(計算區(qū)域如圖3所示),圖中反射SV 波與反射P 波的反射角分別為α1、β1,對應(yīng)的波幅系數(shù)分別為B1、B2,左邊人工邊界高度為H,下底邊界和地表面長度為L.

圖3 SV 波斜入射示意圖

左邊界處的入射和反射SV 波、反射P波相對于入射波usv(t)零時刻波陣面的延遲時間分別為:

式中:H為左邊界高度;y為左邊界豎向坐標(biāo)值.

左邊界內(nèi)行場位移由入射SV 波、反射P波和反射SV 波構(gòu)成,邊界結(jié)點(diǎn)上的位移公式如下:

由式(9)、(10)將波場的傳播應(yīng)力從局部坐標(biāo)系轉(zhuǎn)化為全局坐標(biāo)系,求得在左人工邊界上的應(yīng)力.結(jié)合式(12)、(13)得出左邊界的等效節(jié)點(diǎn)力如式(14):

地基截斷邊界的右邊界、底邊界處節(jié)點(diǎn)的等效荷載計算方法與左邊界類似,只需注意應(yīng)力方向即可.

1.3 考慮靜動力邊界的統(tǒng)一分析方法

在結(jié)構(gòu)動力響應(yīng)分析中,通常要考慮結(jié)構(gòu)動力分析的初始條件,稱之為靜動力綜合分析問題.該分析方法綜合考慮了靜動荷載對結(jié)構(gòu)地震響應(yīng)的影響,解決了因靜動邊界不能共存而只能考慮其中一種邊界的問題.靜動邊界不能共存的原因在于二者邊界條件不同,在靜力邊界中通常側(cè)邊界采用滾軸邊界、底邊界采用固定邊界(如圖4所示);本文動力邊界采用黏彈性人工邊界(圖1所示).

圖4 靜力分析邊界

因靜、動力邊界條件不同,若同時考慮會導(dǎo)致結(jié)果失真,對結(jié)構(gòu)進(jìn)行地震響應(yīng)分析時不能忽略該問題.文章采用一種黏彈性靜-動力統(tǒng)一人工邊界,該方法可實(shí)現(xiàn)靜動力綜合分析問題中的邊界轉(zhuǎn)換,其求解步驟如下:

1)地應(yīng)力平衡.進(jìn)行動力分析時,為使重力荷載作用時地基仍保持初始狀態(tài),要進(jìn)行地應(yīng)力平衡.因此先對結(jié)構(gòu)進(jìn)行靜力分析,即添加靜力邊界條件和重力荷載,導(dǎo)出地基的應(yīng)力場.

2)求解支座反力.將上一步導(dǎo)出的應(yīng)力場作為初始邊界條件,施加靜力邊界條件、重力、水荷載等靜力荷載進(jìn)行求解.邊界結(jié)點(diǎn)的受力可分解為靜力約束對地基的作用和地基對靜力邊界作用的疊加,如圖5所示,最后導(dǎo)出支座水平向反力f1和豎直向反力f2.

圖5 求解支座反力

3)動力分析.將得到的支座反力作為靜力約束的等效替代施加于邊界,如圖6所示.施加動力邊界條件(黏彈性邊界)以及靜力(包含重力荷載、水壓力等)、地震荷載等進(jìn)行動力分析.

圖6 動力人工邊界

2 算例驗(yàn)證

為了驗(yàn)證文中波場分解方式下地震動輸入方法的計算精度,考慮二維彈性地基模型,將該模型的左、右、底邊界全部設(shè)置黏彈性邊界,上方為半空間自由表面,計算時分別輸入入射波角度為0°和30°的P或SV 波.選用動力隱式算法,計算區(qū)域大小為700 m×350 m,網(wǎng)格大小取為5 m,時間增量步為0.01 s,計算時長2 s.計算中采用的土體介質(zhì)參數(shù)E=15 GPa,μ=0.24,ρ=2 700 kg/m3,cs=1 496.71 m/s,cp=2 558.93 m/s.入射波位移函數(shù)如下:

2.1 P波程序?qū)崿F(xiàn)與驗(yàn)證

P波在0°、30°斜入射下地表中心點(diǎn)水平或豎直位移時程曲與解析解對比圖如圖7～8所示.圖9～10分別為P波在0°、30°斜入射下不同時刻的波場位移云圖.

圖7 P波0°斜入射地表中心點(diǎn)位移時程曲線

圖8 P波30°斜入射地表中心點(diǎn)位移時程曲線

圖9 P波0°斜入射位移云圖(t=0.3 s)

圖10 P波30°斜入射位移云圖(t=0.16 s)

2.2 SV波程序?qū)崿F(xiàn)與驗(yàn)證

SV 波在0°、30°斜入射下地表中心點(diǎn)水平位移時程曲與解析解對比圖如圖11～12所示.圖13～14分別為SV 波在0°、30°斜入射下不同時刻的波場位移云圖.

圖11 SV 波0°斜入射地表中心點(diǎn)位移時程曲線

圖12 SV 波30°斜入射地表中心點(diǎn)位移時程曲線

圖13 SV 波0°斜入射位移云圖(t=0.27 s)

圖14 SV 波30°斜入射位移云圖(t=0.23 s)

對上述P或SV 波不同斜入射角度位移時程曲線與解析解的對比圖與位移云圖分析得出:P 或SV波斜入射二維算例計算結(jié)果的位移時程曲線與解析解的位移時程曲線吻合良好.所以整體看來本節(jié)基于Fortran語言自編的二維P或SV 波斜入射子程序計算結(jié)果與解析解基本吻合良好.

3 地震波斜入射作用下混凝土壩地震響應(yīng)分析

本節(jié)基于FORTRAN 語言編寫的二維P波斜入射子程序,并利用黃登混凝土重力壩-庫水-地基模型進(jìn)行工程算例分析,因SV 波入射時存在臨界角,在此由該工程模型的基本參數(shù)及公式(16)可計算得出臨界角α1為35.8°.

其中,SV 波臨界角的公式為[29]

在此,對于P波斜入射分以0°、5°、10°、15°、20°、25°、30°、45°、60°;對于SV 波斜入射分以0°、5°、10°、15°、20°、25°、30°進(jìn)行工程算例分析.并對壩體關(guān)鍵點(diǎn)(壩踵、上游折坡處、壩趾)的有限元計算結(jié)果進(jìn)行提取與分析.

3.1 模型參數(shù)及作用荷載

3.1.1 模型參數(shù)

壩體不模擬分區(qū),只采用一種混凝土材料,并假定為線彈性,動彈模取為靜彈模的1.5倍,阻尼比取5%;地基假定為均質(zhì)線彈性介質(zhì),且忽略材料阻尼,材料參數(shù)見表1,計算模型如圖15所示.

表1 模型材料力學(xué)特性參數(shù)

圖15 重力壩-庫水-地基系統(tǒng)計算模型

此外,在實(shí)際的動態(tài)分析中還考慮了壩體結(jié)構(gòu)的阻尼,本節(jié)使用瑞利阻尼的標(biāo)準(zhǔn)形式來計算,即兩兩參數(shù)模型見式[11]:

式中:ω1、ω2分別為壩體-地基結(jié)構(gòu)自振狀況下的前兩階模態(tài)頻率;ξ為阻尼系數(shù),對于重力壩ξ=5%.

3.1.2 作用荷載

本工況選取黃登混凝土重力壩12號擋水壩段,考慮上游正常蓄水位情況進(jìn)行數(shù)值模擬計算,作用荷載如下:

1)自重荷載;

2)靜水壓力:上游正常蓄水位1 619.00 m,下游尾水位1 422.00 m;

3)動水壓力:采用Westergaard公式計算的動水附加質(zhì)量[7];

4)地震荷載:采用Koyna地震時程,總歷時為10 s,選擇固定時間步長Δt=0.01 s進(jìn)行計算.加速度時程如圖16所示.

圖16 地震動加速度時程曲線

3.2 P波斜入射下壩體動力響應(yīng)分析

3.2.1 壩體關(guān)鍵點(diǎn)位移響應(yīng)

如圖17為不同P波入射角度的壩頂壩踵點(diǎn)相對位移時程曲線,為了更清晰直觀地看出壩頂壩踵點(diǎn)相對位移隨著角度增加的變化規(guī)律,如圖18所示分別提取了不同入射角度的相對位移峰值.從圖中可以得出:隨著入射角度的增加壩頂壩踵點(diǎn)水平、豎直相對位移時程曲線整體上逐漸增大,其中60°較0°的峰值相對位移,水平向增加了110.11%,豎直向增加了48.52%,這說明P波斜入射在0°～60°內(nèi)入射角度為60°時斜入射時對地震響應(yīng)影響較大.

圖17 P波斜入射下壩頂壩踵點(diǎn)豎向相對位移曲線

圖18 P波斜入射下壩頂壩踵點(diǎn)相對位移峰值曲線

3.2.2 壩體關(guān)鍵點(diǎn)應(yīng)力響應(yīng)

圖19～20分別為壩踵點(diǎn)、壩趾點(diǎn)兩個壩體易損害點(diǎn)的不同入射角度的最大或最小主應(yīng)力時程曲線,對應(yīng)分別提取了不同入射角度的關(guān)鍵點(diǎn)主應(yīng)力峰值如圖21～22所示.從圖中可以得出:隨著入射角度的增加壩踵點(diǎn)、壩趾點(diǎn)的應(yīng)力時程曲線整體上隨著入射角度的增加最大、最小主應(yīng)力均呈現(xiàn)出逐漸增加的趨勢.對比兩個關(guān)鍵點(diǎn)不同入射角度的主應(yīng)力峰值曲線發(fā)現(xiàn)隨著入射角度的增加峰值應(yīng)力均逐漸增加,其中壩踵處的峰值應(yīng)力60°為0°的11.15倍;壩趾處最小主應(yīng)力增加了10.42%.這說明P 波斜入射在0°～60°內(nèi),入射角度為60°時對壩體動力響應(yīng)的影響程度較大.

圖19 P波斜入射下壩踵最大主應(yīng)力時程曲線

圖20 P波斜入射下壩趾最小主應(yīng)力時程曲線

圖21 P波斜入射下壩踵點(diǎn)最大主應(yīng)力峰值

圖22 P波斜入射下壩趾最小主應(yīng)力峰值

3.3 SV波斜入射下壩體動力響應(yīng)分析

3.3.1 壩體關(guān)鍵點(diǎn)位移響應(yīng)

圖23為SV 波斜入射時,壩體3 個不同入射角度的壩頂壩踵點(diǎn)相對位移時程曲線,因相對位移隨著角度的增加變化較為復(fù)雜,同時為了更清晰直觀地表達(dá)壩頂壩踵點(diǎn)相對位移隨著SV 波入射角度增加的變化規(guī)律,如圖24所示提取了不同入射角度的相對位移峰值.對比不同角度的相對位移峰值發(fā)現(xiàn)在15°后位移峰值的變化規(guī)律出現(xiàn)了不同,整體來看入射角度為0°(即垂直入射)與其他入射角度相比無論是橫向還是縱向相對位移是最大的,其中較30°的峰值相對位移,水平向降低2.31%,豎直向降低1.18%,這說明SV 波0°(即垂直入射)入射對地震響應(yīng)影響變化不大.可以看出:隨著入射角度的增加壩頂壩踵點(diǎn)水平、豎直相對位移時程曲線整體上變化不大.

圖23 SV 波斜入射下壩頂壩踵點(diǎn)相對位移曲線

圖24 SV 波斜入射下壩頂壩踵點(diǎn)相對位移峰值曲線

3.3.2 壩體關(guān)鍵點(diǎn)應(yīng)力響應(yīng)

圖25～26分別為上游折坡處、壩趾點(diǎn)的不同入射角度的最大或最小主應(yīng)力時程曲線,為了更清晰直觀地看出隨著角度增加關(guān)鍵點(diǎn)的主應(yīng)力變化規(guī)律,如圖27～28所示,分別提取了不同入射角度的關(guān)鍵點(diǎn)主應(yīng)力峰值.

圖25 SV 波斜入射下上游折坡處最大主應(yīng)力時程曲線

圖26 SV 波斜入射下壩趾最小主應(yīng)力時程曲線

圖27 SV 波斜入射下上游折坡處最大主應(yīng)力峰值

圖28 SV 波斜入射下壩趾最小主應(yīng)力峰值

從圖中可以得出:入射角度為0°與其他入射角度相比無論是上游折坡處還是壩趾的主應(yīng)力都是最大的.其中30°相較0°,上游折坡處最大主應(yīng)力降低了4.85%,壩趾處最小主應(yīng)力降低了5.94%.總的來說,對于SV 波斜入射情況下,隨著入射角度的增加上游折坡處、壩趾點(diǎn)兩個關(guān)鍵點(diǎn)的應(yīng)力時程曲線整體上影響不大.

4 結(jié) 論

以黃登混凝土重力壩為研究對象,基于黏彈性靜-動統(tǒng)一人工邊界,研究P 波和SV 波在不同斜入射角度下混凝土壩-地基-庫水系統(tǒng)的地震響應(yīng),得出以下結(jié)論:

1)不同地震波斜入射下的動力響應(yīng)規(guī)律與垂直入射相比變化差異較大.P波斜入射時的地震響應(yīng)規(guī)律與垂直入射相比有明顯的不同,P波與SV 波斜入射時的地震響應(yīng)規(guī)律也存在差異;在重力壩動力響應(yīng)分析中很有必要考慮地震波的不同入射角度.

2)P波斜入射時,隨著入射角度的增加壩體關(guān)鍵點(diǎn)的水平、豎直相對位移與主應(yīng)力值逐漸增加,其中30°與0°相比,峰值相對位移水平向增加了110.11%,豎直向增加了48.52%;壩踵處最大拉應(yīng)力為0°時的11.15倍;壩趾處最小主應(yīng)力增加了10.42%.這說明P波斜入射在0°～60°,內(nèi)入射角度為60°時斜入射時壩體的地震響應(yīng)較大.

3)SV 波斜入射時,隨著入射角度的增加壩體關(guān)鍵點(diǎn)的水平、豎直相對位移與主應(yīng)力值變化不大,雖在15°后位移峰值的變化規(guī)律出現(xiàn)了不同,但整體來看入射角度為0°(即垂直入射)與其他入射角度相比無論是水平還是豎直相對位移均是最大的,這說明相較其他角度,SV 波0°入射時壩體的地震響應(yīng)較大.