譙 雯

(華南農(nóng)業(yè)大學(xué)水利與土木工程學(xué)院，廣州 510642)

0 引言

在進(jìn)行壩體—地基系統(tǒng)地震動(dòng)力反應(yīng)分析時(shí)，一個(gè)值得關(guān)注的問題就是地基輻射阻尼的影響。解決這個(gè)問題最有效的方法就是引入人工邊界。人工邊界的概念最早是由Alterman等[1]于1968年提出的。隨后眾多人工邊界方法被提出，其中由Deeks[2]提出的黏彈性人工邊界，因既能夠模擬地基輻射阻尼的作用，又能夠考慮地基的彈性恢復(fù)力，同時(shí)還具有良好的穩(wěn)定性而得到了廣泛的研究和應(yīng)用。劉晶波等[3]在此基礎(chǔ)上相繼提出了靜—?jiǎng)恿y(tǒng)一人工邊界、一致黏彈性人工邊界及黏彈性人工邊界單元的概念[4]。苑舉衛(wèi)等[5]、李浩然等[6]在大型有限元商用軟件的平臺(tái)上通過二次開發(fā)技術(shù)，實(shí)現(xiàn)了黏彈性人工邊界的添加，取得了良好的模擬效果。黏彈性人工邊界單元因便于編程實(shí)現(xiàn)，被廣泛應(yīng)用于結(jié)構(gòu)—地基動(dòng)力有限元問題分析中。

在壩體—地基相互作用問題中，地震動(dòng)的輸入是另一個(gè)關(guān)鍵問題[7]。對(duì)于大型水工結(jié)構(gòu)而言，在考慮地基輻射阻尼的情況下，地震動(dòng)輸入方式通常按波動(dòng)輸入問題來考慮。劉晶波等[8]通過將地震波轉(zhuǎn)化為人工邊界上的等效節(jié)點(diǎn)荷載來實(shí)現(xiàn)地震動(dòng)的輸入。這種方法既適用于波的豎向輸入，又適用于多點(diǎn)輸入。

綜上原因，本文基于黏彈性人工邊界及其相配套的地震動(dòng)輸入方法，引入無厚度Goodman單元模擬黏彈性人工邊界，推導(dǎo)了地震動(dòng)輸入方式在有限元中的實(shí)現(xiàn)過程。該方法在保證計(jì)算精度的前提下，不增加計(jì)算規(guī)模，保證了計(jì)算效率，并且在有限元中易于實(shí)現(xiàn)。通過經(jīng)典數(shù)值算例驗(yàn)證了所編程序的正確性，最后將無厚度Goodman黏彈性人工邊界單元應(yīng)用于混凝土重力壩—地基—庫水系統(tǒng)的地震動(dòng)力響應(yīng)分析中，并與傳統(tǒng)的無質(zhì)量地基模型進(jìn)行比較，進(jìn)一步驗(yàn)證了無厚度Goodman黏彈性人工邊界單元及其與之相適應(yīng)的地震動(dòng)輸入在實(shí)際工程應(yīng)用中的強(qiáng)實(shí)用性、施加的方便性和計(jì)算結(jié)果的合理性。

1 黏彈性人工邊界及地震動(dòng)輸入

基于球面波動(dòng)理論推導(dǎo)出的三維黏彈性人工邊界，通過在模型邊界節(jié)點(diǎn)設(shè)置連續(xù)分布的并列彈簧—阻尼器來模擬無限介質(zhì)。彈簧剛度系數(shù)和阻尼系數(shù)根據(jù)式(1)和式(2)取得：

(1)

(2)

式中：KT、KN分別為彈簧的切向和法向剛度系數(shù)；C為阻尼系數(shù)；G、ρ分別為介質(zhì)的剪切模量、密度；cs、cp分別為剪切波速與壓縮波速；αT為彈簧切向修正系數(shù)(取值范圍0.5～1.0，推薦取0.67)[9]；αN為彈簧法向修正系數(shù)(取值范圍1.0～2.0，推薦取1.33)[9]；R為人工邊界至波源的距離。

黏彈性人工邊界是將輸入地震動(dòng)轉(zhuǎn)化為等效節(jié)點(diǎn)荷載：

(3)

(4)

2 黏彈性人工邊界及其地震動(dòng)輸入在有限元中的處理

2.1 黏彈性人工邊界單元在有限元中的處理

為了方便數(shù)值模擬，基于矩陣等效原理，本文采用文獻(xiàn)[4]推導(dǎo)的利用等效實(shí)體邊界單元替代彈簧—阻尼器原件，并且這種單元對(duì)于不規(guī)則邊界具有更強(qiáng)的適應(yīng)性。通過比較證明，單元厚度對(duì)計(jì)算結(jié)果影響甚微，因此本文利用無厚度Goodman單元模擬黏彈性人工邊界。三維模型單元類型如圖1所示。具體三維黏彈性邊界模型的剛度矩陣、阻尼矩陣的推導(dǎo)過程可參見文獻(xiàn)[10]。

(a)無厚度六面體單元

(b)無厚度五面體單元圖1 Goodman單元的單元類型

2.2 地震動(dòng)輸入在有限元中的處理

假定黏彈性人工邊界單元形函數(shù)為：

(5)

(6)

將式(6)代入式(3)可得到黏彈性人工邊界側(cè)邊面上的等效結(jié)點(diǎn)荷載列陣為：

(7)

同理，將式(6)代入式(4)可得黏彈性人工邊界底邊上的等效節(jié)點(diǎn)荷載為：

(8)

根據(jù)已知輸入地震波的位移和速度時(shí)程曲線，按式(7)～(8)即可求得人工邊界上結(jié)點(diǎn)的等效結(jié)點(diǎn)荷載，從而實(shí)現(xiàn)黏彈性人工邊界上地震動(dòng)輸入。

3 Lamb表面源問題的算例驗(yàn)證

為驗(yàn)證自編有限元程序中三維黏彈性人工邊界實(shí)現(xiàn)的有效性，采用經(jīng)典Lamb問題對(duì)其進(jìn)行驗(yàn)證。模型長L=1，寬B=1，深度H=0.5，根據(jù)對(duì)稱性原則，在水平方向上取原模型的1/4作為計(jì)算模型，如圖2所示。詳細(xì)計(jì)算參數(shù)見表1。通過在模型的四周和底面設(shè)置厚度為0的三維Goodman單元作為等效一致黏彈性人工邊界單元，并將人工邊界單元的最外層結(jié)點(diǎn)固定。在半空間表面處施加垂直集中荷載，荷載激勵(lì)如圖3所示。

圖2 Lamb表面源問題計(jì)算模型示意圖

表1 計(jì)算參數(shù)

圖3 輸入荷載時(shí)程

取距加載位置0點(diǎn)為0.2和0.4的兩點(diǎn)作為觀測(cè)點(diǎn)，計(jì)算結(jié)果與對(duì)比結(jié)果如圖4所示。

(a)r=0.20 m處應(yīng)力時(shí)程曲線

(b)r=0.40 m處應(yīng)力時(shí)程曲線圖4 計(jì)算結(jié)果與對(duì)比結(jié)果圖

由計(jì)算結(jié)果圖4顯示，采用固定邊界模擬時(shí)，脈沖波往復(fù)反射，觀測(cè)點(diǎn)的位移產(chǎn)生了大幅度的波動(dòng)，模擬結(jié)果與解析解相差甚大。而在黏彈性人工邊界模型中，當(dāng)脈沖波到達(dá)黏彈性人工邊界后，能量幾乎完全被邊界的黏滯阻尼吸收，模擬結(jié)果與解析解非常一致。證明了本文基于無厚度Goodman單元的黏彈性人工邊界能夠很好地模擬無限地基輻射阻尼，同時(shí)也驗(yàn)證了所編程序的可靠性和有效性。

4 工程應(yīng)用

4.1 有限元計(jì)算模型

由于實(shí)際工程中，重力壩在沿壩軸線方向設(shè)置了橫縫，強(qiáng)震作用下橫縫傳遞的剪力遠(yuǎn)小于慣性力。加之相鄰壩段會(huì)由于壩高、自振周期的不同而形成單個(gè)壩段的振動(dòng)。因此，對(duì)于重力壩的動(dòng)力分析，一般不考慮沿著壩軸線方向的振動(dòng)，常見做法是選取單個(gè)壩段進(jìn)行動(dòng)力分析。因此，本文選取福建省水口水電站19th壩段作為分析對(duì)象，對(duì)重力壩—地基—庫水系統(tǒng)進(jìn)行地震動(dòng)力反應(yīng)分析，分別對(duì)地基被施加黏彈性人工邊界(模型Ⅰ)和固定邊界(模型Ⅱ)兩種情況建立了三維有限元模型，其中，施加固定邊界的模型的地基范圍較大。分析時(shí)壩體和地基均取為線彈性模型。假定庫水是不可壓縮的無阻尼理想流體，忽略振動(dòng)中庫水表面的影響，不計(jì)庫底豎向位移對(duì)水體的影響。采用Westergard附加質(zhì)量法實(shí)現(xiàn)庫水動(dòng)水壓力對(duì)大壩地震響應(yīng)的影響，采用直接濾頻法求解特征方程。壩體地基網(wǎng)格劃分如圖5所示，庫水位網(wǎng)格圖如圖6所示，大部分壩體、地基均采用三維八結(jié)點(diǎn)六面體等參單元離散，少部分壩身用六結(jié)點(diǎn)五面體單元過渡。在地基左右和底側(cè)截?cái)噙吔缣幐餮由煲粚雍穸葹?的Goodman單元作為等效邊界單元來模擬黏彈性人工邊界。對(duì)于施加了黏彈性人工邊界的模型Ⅰ，共劃分為5 852個(gè)結(jié)點(diǎn)，2 890個(gè)單元，而施加固定邊界的模型II則劃分為了6 972個(gè)結(jié)點(diǎn)，3 364個(gè)單元。

4.2 計(jì)算參數(shù)

重力壩的壩體以及地基假設(shè)為均勻各向同性體。在進(jìn)行動(dòng)力分析之前，需對(duì)重力壩進(jìn)行靜力計(jì)算，以確定壩體和壩基的初始應(yīng)力狀態(tài)。本文結(jié)合壩體和地基材料參數(shù)以及文獻(xiàn)[11]的建議，參數(shù)取值見表2。

(a)模型Ⅰ：壩—地基系統(tǒng)有限元網(wǎng)格

(b)模型Ⅱ：壩—地基系統(tǒng)有限元網(wǎng)格圖5 地基系列有限元網(wǎng)絡(luò)

圖6 庫水位61.00 m的水體有限元網(wǎng)格圖

4.3 地震波的選擇與輸入

這里選用位于Ⅲ類場(chǎng)地的EI Centro波實(shí)際地震記錄，根據(jù)本工程相關(guān)資料，福州地區(qū)場(chǎng)地類別為Ⅲ類，處于抗震設(shè)防烈度7度區(qū)，現(xiàn)將EI Centro波峰值加速度調(diào)整到0.1g，時(shí)間步距0.02 s，調(diào)整后的EI Centro波加速度時(shí)程曲線如圖7所示。

圖7 調(diào)整后的EI Centro波加速度時(shí)程

4.4 動(dòng)力計(jì)算成果分析

由于篇幅限制，本文只給出了壩頂65號(hào)結(jié)點(diǎn)和壩底3 109號(hào)結(jié)點(diǎn)的動(dòng)位移時(shí)程，壩頂4號(hào)單元和壩底543號(hào)單元第一主應(yīng)力和第三主應(yīng)力時(shí)程，并同時(shí)給出了傳統(tǒng)無質(zhì)量地基固定邊界的結(jié)果作為比較。

從圖8～13對(duì)比可以看出：采用本文方法的計(jì)算模型進(jìn)行計(jì)算，重力壩壩頂和壩基處的幾個(gè)主要部位的順河向動(dòng)位移響應(yīng)和動(dòng)主應(yīng)力響應(yīng)與采用傳統(tǒng)無質(zhì)量地基固定邊界條件下的計(jì)算模型對(duì)比，其計(jì)算的結(jié)果有明顯的降低，其降低幅度大致在19%～55%之間。主要原因是由于無質(zhì)量地基未考慮地基輻射阻尼的作用，導(dǎo)致反射波的存在，從而使得動(dòng)力響應(yīng)值偏大。從圖中還可以看出：壩頂結(jié)點(diǎn)動(dòng)位移響應(yīng)較大，沿建基面位置附近的動(dòng)力響應(yīng)相對(duì)較小，其主要原因是地基及壩體對(duì)地震波的放大作用造成的。

表2 計(jì)算參數(shù)

圖8 壩頂65號(hào)結(jié)點(diǎn)順河向位移時(shí)程曲線

圖9 壩底3 109號(hào)結(jié)點(diǎn)順河向位移時(shí)程曲線

圖10 壩頂4號(hào)單元第一主應(yīng)力時(shí)程曲線

圖11 壩底543號(hào)單元第一主應(yīng)力時(shí)程曲線

圖12 壩頂4號(hào)單元第三主應(yīng)力時(shí)程曲線

圖13 壩底543號(hào)單元第三主應(yīng)力時(shí)程曲線

從以上結(jié)果比較可以看出，本文方法地基邊界處對(duì)地震波動(dòng)具有很好的吸收作用，波動(dòng)在邊界處的反射作用得到了明顯的減弱，能夠考慮地基輻射阻尼的影響，使得計(jì)算模型的計(jì)算結(jié)果更加符合實(shí)際，再次驗(yàn)證了本文動(dòng)力有限元計(jì)算程序的可靠性和黏彈性人工邊界施加的正確性。

5 結(jié)語

黏彈性人工邊界不僅能夠吸收外行散射地震波能量，而且能夠模擬人工邊界半無限地基介質(zhì)的彈性恢復(fù)性能，被廣泛地應(yīng)用于實(shí)際工程。本文探討了無厚度Goodman單元模擬黏彈性人工邊界的可行性，推導(dǎo)了地震動(dòng)輸入轉(zhuǎn)化為等效節(jié)點(diǎn)荷載在有限元中的實(shí)現(xiàn)過程。得出以下結(jié)論：

1)經(jīng)典Lamb問題算例表明：使用一致黏彈性人工邊界時(shí)，在接近加載點(diǎn)附近的結(jié)點(diǎn)和靠近邊界的結(jié)點(diǎn)上，其位移計(jì)算結(jié)果與解析解很接近，計(jì)算的精度高。同時(shí)也可以看出，在對(duì)黏彈性人工邊界的模擬上，采用本文的無厚度Goodman單元代替彈簧—阻尼器的方法，將使黏彈性人工邊界的施加方式更為簡單，保證計(jì)算效率。

2)重力壩—地基—庫水系統(tǒng)進(jìn)行地震動(dòng)力反應(yīng)分析實(shí)例表明：合理選取和模擬人工邊界條件將直接影響動(dòng)力分析結(jié)果的準(zhǔn)確性，黏彈性人工邊界條件原理簡單、易于運(yùn)用，具有良好的模擬效果，穩(wěn)定性較好，而且可以在邊界不規(guī)則的情況下使用，從而可為各類工程的抗震分析提供參考。