劉依松 陳燈紅

（1.三峽大學(xué) 水利與環(huán)境學(xué)院，湖北 宜昌 443002；2.三峽大學(xué) 土木與建筑學(xué)院，湖北 宜昌 443002）

系統(tǒng)評(píng)價(jià)重力壩抗震安全性的關(guān)鍵之一是解決大壩與地基的動(dòng)力相互作用的問(wèn)題.Clough對(duì)此提出了一種無(wú)質(zhì)量地基模型，認(rèn)為大壩附近一定范圍的地基是線彈性、無(wú)質(zhì)量的、地震激勵(lì)均勻作用于截?cái)噙吔缟?這種模型雖在工程中應(yīng)用較廣，但存在局限性：①地基實(shí)際上是有質(zhì)量的半無(wú)限介質(zhì)，地震波動(dòng)能量向無(wú)限遠(yuǎn)處散逸，即無(wú)限地基的幅射對(duì)散射地震波動(dòng)能量將起到一種吸能作用；②深山峽谷中的高壩，其基底延伸較長(zhǎng)，沿壩基交界面的地震動(dòng)振幅和相位差異明顯，高壩空間整體作用效應(yīng)明顯，忽略這種不均勻性將難以反映壩體的實(shí)際地震反應(yīng).為了克服這些不足，近年來(lái)，Lysmer等［1］提出了黏性邊界，Deeks等［2］提出了時(shí)域解耦的黏彈性邊界，廖振鵬等［3］提出了透射邊界等，比例邊界有限元法［4－5］也在發(fā)展之中.但黏性邊界僅考慮邊界吸能阻尼，未計(jì)入介質(zhì)彈性恢復(fù)影響，易致低頻漂移.透射邊界常會(huì)高頻失穩(wěn)，需增設(shè)許多邊界節(jié)點(diǎn)和單元.黏彈性人工邊界條件基于衰減散射波的影響，在人工邊界上設(shè)置并聯(lián)彈簧和阻尼器，使得低頻精度好，高頻穩(wěn)定性好.采用基于衰減散射波的黏彈性思想，建立了重力壩－地基系統(tǒng)［6］二維平面黏彈性模型，采用ANSYS中的APDL參數(shù)化設(shè)計(jì)語(yǔ)言進(jìn)行重力壩動(dòng)力響應(yīng)分析.

1 研究方法

1.1 黏彈性人工邊界

將波陣面方程引入幾何擴(kuò)散因子，Deeks［2］、劉晶波［7－9］推衍出反映無(wú)限介質(zhì)彈性恢復(fù)性能的黏彈性人工邊界條件.將平面波和遠(yuǎn)場(chǎng)散射波混合經(jīng)驗(yàn)疊加且考慮多角度透射，杜修力［10－11］改進(jìn)了黏彈性邊界公式，對(duì)于波源問(wèn)題、衰減散射波問(wèn)題均有較好的計(jì)算精度.對(duì)于二維黏彈性邊界，剛度系數(shù)和阻尼系數(shù)為

式中，KBN、KBT分別表示法向和切向彈簧系數(shù)，rb為散射波源至人工邊界的距離，λ為拉梅常數(shù)，G為介質(zhì)剪切模量，ρ為介質(zhì)質(zhì)量密度，并具有G＝和λ＋2G＝的關(guān)系；CBN、CBT分別為表示法向和切向阻尼系數(shù)，cp、cs分別為P波和S波波速，α、β為無(wú)量綱參數(shù)，分別取0.8、1.1.

黏彈性人工邊界作為一種應(yīng)力連續(xù)分布的人工邊界條件，可采用與普通單元類似的形函數(shù)進(jìn)行離散，結(jié)點(diǎn)i處的形函數(shù)為

位移形函數(shù)矩陣為

黏彈性人工邊界上分布的彈性矩陣為

單元的剛度矩陣公式為

經(jīng)式（5）計(jì)算得到，黏彈性人工邊界單元的剛度矩陣如式（6），阻尼矩陣具有相同的形式，只需將剛度系數(shù)K用相應(yīng)的阻尼系數(shù)C替換即可.

其中，l為黏彈性人工邊界單元的長(zhǎng)度.

1.2 地震荷載輸入方法

黏彈性邊界能夠吸收由近場(chǎng)向遠(yuǎn)域散射的外行波，將外源地震波引入計(jì)算區(qū)域時(shí)，只需在人工邊界處施加自由場(chǎng)荷載即可，可以把近域地基作為半無(wú)限空間自由場(chǎng)地基的子結(jié)構(gòu)，由其在入射地震波作用下的邊界相互作用力給出［10－11］

圖1 壩－地基系統(tǒng)黏彈性邊界及自由場(chǎng)示意圖

自由場(chǎng)應(yīng)力可由彈性理論求得，即

當(dāng)在底邊界垂直入射P波時(shí)，底部的邊界條件為u＝0，v＝v0（t），則自由場(chǎng)在沿高度h處，自由場(chǎng)位移、速度可表示為入射波和反射波的疊加，即

式中，H為底邊界到自由表面的距離.

并且對(duì)于底部邊界，有h＝0，n＝［0－1］T，則將式（10）代入式（9）中，可得自由場(chǎng)應(yīng)力σfb為

將式（10）、（11）代入式（8）中，等效荷載FB可求得，即

當(dāng)在底邊界垂直入射SV波時(shí)，底部的邊界條件為u＝u0（t），v＝0，按照同樣的方法，等效荷載FB為

對(duì)于側(cè)邊界的自由場(chǎng)應(yīng)力也可按類似的方法推導(dǎo)出.

2 數(shù)值實(shí)現(xiàn)與算例驗(yàn)證

黏彈性人工邊界具有較多優(yōu)點(diǎn)，最主要的是和有限元的結(jié)合實(shí)現(xiàn)起來(lái)比較方便，只需在有限元模型人工邊界上結(jié)點(diǎn)的法向和切向分別設(shè)置并聯(lián)的彈簧單元和阻尼器單元.在ANSYS程序中，彈簧－阻尼器元件可以用COMBIN14單元來(lái)模擬，該單元具有彈簧和阻尼選項(xiàng).值得一提的是，在實(shí)際的有限元實(shí)現(xiàn)過(guò)程中，各向邊界物理元件參數(shù)的實(shí)際取值等于式（1）給出的量值與對(duì)應(yīng)有限元結(jié)點(diǎn)所代表網(wǎng)格面積的乘積.利用ANSYS中的APDL參數(shù)化設(shè)計(jì)語(yǔ)言，編寫了實(shí)現(xiàn)黏彈性人工邊界的宏文件VSB.mac，程序流程圖如圖3所示.利用Fortran語(yǔ)言編寫了式（12）、（13）所示的加載子程序Load.for，可以方便地施加等效地震荷載.

圖2 彈性半空間網(wǎng)格圖

圖3 黏彈性人工邊界程序流程圖

為驗(yàn)證上述黏彈性邊界單元及波動(dòng)輸入方法程序的正確性，現(xiàn)采用均勻彈性半空間［12］受底部垂直入射SV波進(jìn)行驗(yàn)證.波動(dòng)輸入采用的力學(xué)計(jì)算模型如圖2所示.介質(zhì)的力學(xué)參數(shù)為：彈性模量E＝13.23 GPa，泊松比ν＝0.25，質(zhì)量密度ρ＝2 700kg／m3，波速為cs＝1 400m／s，cp＝2 425m／s，計(jì)算區(qū)的范圍為762m×381m，有限元單元尺寸為19.05m×19.05 m，時(shí)間步長(zhǎng)取為0.01s，入射波的方程見(jiàn)式（14）.

經(jīng)過(guò)計(jì)算得到 A（－381，381），B（0，381），C（381，381），D（381，0），E（0，0）等5個(gè)點(diǎn)的x向位移時(shí)程，如圖4所示.

圖4 頂部和底部位移時(shí)程

從A、B、C 3個(gè)點(diǎn)的x向位移時(shí)程圖可以看出：A、B、C 3點(diǎn)的位移值相同，說(shuō)明是由入射波和反射波共同引起的位移時(shí)程，且頂部的最大值為接近SV入射波幅值的2倍.從D、E兩個(gè)點(diǎn)的位移時(shí)程圖可以看出，D、E兩點(diǎn)位移值相同，位移曲線的前半段為入射波引起的位移時(shí)程，后半段為反射波引起的位移時(shí)程.5個(gè)點(diǎn)的解均與理論解相吻合.因此，在ANSYS編寫的黏彈性邊界單元及實(shí)現(xiàn)波動(dòng)輸入子程序是正確的.

3 工程應(yīng)用

某水電站具有發(fā)電、防洪、航運(yùn)等綜合效益，屬Ⅰ等大（1）型工程，大壩為碾壓混凝土重力壩，與泄水建筑物同為1級(jí)建筑物，防洪標(biāo)準(zhǔn)按500年一遇設(shè)計(jì)，10 000年一遇校核.7孔溢流壩布置于河床中央，溢流壩段右側(cè)為擋水壩段，壩段寬22m，壩頂寬18m，下游壩坡1∶0.73.上游壩坡在高程270.00m以上為直立坡、以下為1∶0.25的斜坡，壩體橫剖面見(jiàn)圖5所示.壩址的地震基本烈度為7度，大壩按8度地震控制，設(shè)計(jì)烈度為8度，設(shè)計(jì)地震加速度為0.2g.按照本文分析方法，結(jié)合波動(dòng)輸入子程序，在正常蓄水位條件下，對(duì)大壩－地基系統(tǒng)進(jìn)行考慮地基輻射阻尼效應(yīng)的抗震分析，有限元網(wǎng)格如圖6所示.

圖5 壩體位移和應(yīng)力的關(guān)鍵部位

圖6 壩體網(wǎng)格圖

在計(jì)算過(guò)程中，地震荷載不再按照慣性力輸入，而是將合成的人工地震波及反應(yīng)譜（見(jiàn)圖7）加速度轉(zhuǎn)化為位移、速度，然后求得自由場(chǎng)應(yīng)力，進(jìn)行等效地震荷載輸入.輸入的各向位移、速度如圖8所示.動(dòng)力時(shí)程分析各特征點(diǎn)的位移、應(yīng)力結(jié)果見(jiàn)表1，不同地震輸入模型下壩體位移時(shí)程比較如圖9所示，壩體關(guān)鍵部位的應(yīng)力比較如圖10所示.

圖7 合成的人工地震波及反應(yīng)譜

圖8 輸入的各向速度和位移

圖9 不同地震輸入模型位移值比較

圖10 不同地震輸入模型應(yīng)力值比較

由表1可知，考慮了黏彈性人工邊界后的無(wú)限域地基輻射阻尼影響的地震響應(yīng)分析所得到的大壩－地基動(dòng)力響應(yīng)峰值比固定邊界無(wú)質(zhì)量地基模型的所得結(jié)果減小了4%～38%.

表1 特征點(diǎn)位移、應(yīng)力峰值比較

4 結(jié) 論

1）在有限元模型人工邊界上結(jié)點(diǎn)的法向和切向分別設(shè)置一系列并聯(lián)的彈簧和阻尼器單元，以考慮基于衰減散射波的地基輻射阻尼的影響，利用ANSYS中的APDL參數(shù)化設(shè)計(jì)語(yǔ)言，編寫實(shí)現(xiàn)黏彈性人工邊界的宏文件，對(duì)重力壩－地基系統(tǒng)進(jìn)行了動(dòng)力響應(yīng)分析.

2）與無(wú)質(zhì)量地基模型相比，考慮了黏彈性人工邊界后的無(wú)限域地基輻射阻尼影響的地震響應(yīng)分析的大壩－地基動(dòng)力響應(yīng)峰值降低了4%～38%.

3）無(wú)質(zhì)量地基模型并非是有限的、線彈性的、無(wú)質(zhì)量的均勻作用于截?cái)噙吔缟?，而是夸大了結(jié)構(gòu)的動(dòng)力響應(yīng)，所得成果偏于保守，需要改進(jìn)和優(yōu)化.在進(jìn)行動(dòng)力響應(yīng)分析時(shí)，考慮無(wú)限域地基的輻射阻尼影響更符合工程實(shí)際.

［1］ Lysmer J，Kuhlemeyer R L.Finite Dynamic Model for Infinite Media［J］.Journal of Engineering Mechanics，ASCE，1969，95（1）：759－877.

［2］ Deeks A J，Randolph M F.Axisymmetric Time－domain Transmitting Boundaries［J］.Journal of Engineering Mechanics Division，ASCE，1994，120（1）：25－42.

［3］ Liao Z P，Wong H L，Yang B.A Transmitting Boundary for Transient Wave Analyses［J］.Scientia Sincia Methematica，1984，27：1063－1076.

［4］ Wolf J P，Song Ch.The Scaled Boundary Finite－element Method－aprimer： Derivations［J］.Computers and Structures，2000，78：191－210.

［5］ Song Ch，Wolf J P.The Scaled Boundary Finite－element Method－aprimer：Solution Procedures［J］.Computers and Structures，2000，78：211－225.

［6］ 郝明輝，張艷紅，陳厚群.基于ABAQUS的黏彈性人工邊界在重力壩分析中的應(yīng)用［J］.中國(guó)水利水電科學(xué)研究院學(xué)報(bào)，2012，10（2）：120－126.

［7］ 劉晶波，谷 音，杜義欣.一致黏彈性人工邊界及黏彈性邊界單元［J］.巖土工程學(xué)報(bào)，2006，28（9）：1070－1075.

［8］ Liu J B，Du Y X，Du X L.3DViscous－spring Artificial Boundary in Time Domain［J］.Earthquake Engineering and Engineering Vibration，2006，5（1）：93－102.

［9］ 谷 音，劉晶波，杜義欣.三維一致黏彈性人工邊界及等效黏彈性邊界單元［J］.工程力學(xué)，2007，24（12）：31－37.

［10］杜修力，趙密，王進(jìn)廷.近場(chǎng)波動(dòng)模擬的人工應(yīng)力邊界條件［J］.力學(xué)學(xué)報(bào)，2006，38（1）：49－56.

［11］趙建鋒，杜修力，韓 強(qiáng)，等.外源波動(dòng)問(wèn)題數(shù)值模擬的一種實(shí)現(xiàn)方式［J］.工程力學(xué)，2007，24（4）：52－58.

［12］賀向麗.高混凝土壩抗震分析中遠(yuǎn)域能量逸散時(shí)域模擬方法研究［D］.南京：河海大學(xué)，2006.