采用Midas軟件的土石壩動力非線性有限元分析

孫海瑞

（新疆云灃水利設計咨詢有限公司，烏魯木齊 830002）

水庫大壩在強震作用下易遭受破壞， 尤其是松散土石壩，在地震作用下易出現陷落、坍塌等現象，結構穩定和安全運行受到威脅， 因此開展土石壩地震動力響應分析具有重要的現實意義［1］。本文基于土石料等價線性黏彈性模型和動力平衡計算公式［2］，建立有限元模型， 并以烏斯通溝水庫瀝青混凝土心墻砂礫石壩為例，進行動力分析計算，旨在為土石壩的設計、施工及運行安全提供依據。

1 工程概況

烏斯通溝水庫是一座具有綜合效益的水庫樞紐工程，主要任務是灌溉和工業供水，遏制灌區地下水超采，水庫總庫容1440 萬m2，工程等別為Ⅲ等，工程規模為中型。

大壩為瀝青混凝土心墻砂礫石壩， 最大壩高為73.0 m ，為2 級建筑物。壩頂長度為261.90 m，大壩壩頂寬9.0 m，上游壩坡1∶2.25、1∶2.2，上游壩坡采用250 mm 厚的現澆混凝土板護坡； 下游壩坡1∶2.0、1∶1.8，下游壩坡采用300 mm 厚的干砌石+混凝土網格梁護坡。壩體防滲采用瀝青混凝土心墻，左、右岸壩下基巖采用帷幕灌漿防滲；河床段采用混凝土防滲墻+帷幕灌漿防滲；左壩頭防滲線延長段長度47.0 m，右壩頭防滲線延長段長度49.1 m，大壩防滲線總長363.0 m。

工程區地震動峰值加速度為0.22 g，地震動反應譜特征周期為0.40 s， 對應的地震基本烈度為Ⅷ度，所以，本工程地震設防烈度按Ⅷ度設防，大壩采用乙類抗震設防。

本文按照規范的標準和方法， 利用Midas 有限元分析軟件，對該大壩進行動力計算，計算地震工況下的壩體動力反應及地震永久變形， 驗證結構設計的合理性和工程的安全可靠性。

2 計算參數及方法

2.1 計算參數

2.1.1 靜力參數

動力計算是在靜力計算分析的基礎上進行，初始靜應力對土石壩地震荷載作用下的動力反應存在較大影響，靜力分析時，壩體壩基材料采用鄧肯等人提出的非線性彈性的雙曲線E～V 模型。

該土石壩材料靜力計算參數如表1［3］。

2.1.2 動力參數

El Infiernillo、Takimi、碧口、紫坪鋪等已建大壩遭受的實際地震災害表明： 現代重型振動碾壓土石壩在地震過程中主要表現為振縮特性， 這一方面使大壩震后更為密實，抵御地震破壞的能力得到加強，同時密實的壩體地震反應更為劇烈， 大壩特別是壩頂的加速度反應更大，鞭梢效應更強，大壩局部抗震設計措施必須得到加強。

壩體動力反應計算采用國內外廣泛使用的等價線性黏彈性模型， 地震永久變形采用等價節點力法進行計算。

根據地勘試驗資料， 該土石壩動力分析計算的各項參數如表2。

表2 動力計算參數

2.2 計算公式

2.2.1 動力平衡計算公式

由相關研究成果，動水壓力對建于基巖上的100 m級土石壩影響較?。?］，本文忽略動水壓力的影響，動力平衡計算公式如下：

其中：

故

式中 [M ]為整體質量矩陣； [K ]為整體勁度矩陣；[C ]為整體阻尼矩陣； {δ¨g}為地震加速度列陣； {δ¨}為相對加速度列陣； {δ˙}為相對速度列陣；{δ }為相對位移列陣；[G ]為地震加速度3 個分量到n 個自由度體系的n 維空間的轉換矩陣。

動力計算采用Wilson-θ 法， 將整個地震運動過程按每時段1～2 s 劃分成多個時段，進行時程逐步數值積分， 求解動力平衡方程式。迭代精度控制條件要求為：

式中 Gi和λi分別為單元第i 次迭代計算時的動剪模量和阻尼比；ε 值在10%左右。

2.2.2 壩體震后永久變形計算公式

粗粒料的動力性質研究成果比較豐富， 一般采用振動硬化模型， 該模型可反映粗粒料振動硬化對大壩抗震的不利影響［5］。

通過動力非線性有限元分析計算壩體各單元在地震過程中的殘余應變，計算公式如下：

式中Δεγp為殘余剪切應變；Δενp為殘余體積應變；Pa為大氣壓力，為平均有效主應力；k1，k2，k3，nGM，n1為試驗參數。

依據Prandtl-Reuss 流動法則，根據式（8）將殘余應變換算成直角坐標系下的應變：

式中p 為平均主應力，q 為廣義剪應力。

計算出直角坐標下的殘余應變增量后， 可得等效結點力為：

式中 [B ]為應變轉換矩陣；[D ]為彈性矩陣； {Δ εp}為直角坐標系下的殘余應變增量。

將所求的等效結點力加載于大壩結構， 再按照靜力非線性有限元計算， 得到的壩體變形結果就是動力計算的地震永久變形。

3 有限元模型

3.1 壩體有限元計算網絡剖分

本次計算采用Midas 軟件， 以河床典型剖面為基準， 沿壩軸線方向將大壩分為43 個計算剖面，計算坐標系X 軸為河道下游指向上游，Y 軸為垂直地面方向以向上為正，向下為負；壩軸線方向設為Z 軸方向，以向右岸為正，向左岸為負。

整個計算模型共劃分得到總結點數10792 個，總單元數10367個， 其中包括75個防滲墻單元和728 個心墻單元，網格劃分圖如圖1［6］。

圖1 大壩剖分網格立視圖

3.2 壩體填筑及蓄水加載過程模擬

根據擬定施工順序和大壩結構單元劃分情況，將整個大壩施工（包含地基覆蓋層）和蓄水過程共分為23 級荷載進行加載，由此動態模擬大壩的填筑過程及蓄水過程，為了較好地模擬加載，每一級荷載級均采用一次性加載方式， 動力計算工況為設計規范波地震，計算水位為正常蓄水位905.00 m 高程。

4 壩體動力非線性有限元計算成果分析

4.1 壩址基巖地震動輸入

土石壩抗震計算分析中， 基巖地震動輸入規范反應譜地震波（以下簡稱規范波）地震動加速度時程曲線進行大壩動力計算［4］。

查閱《中國地震動參數區劃圖》，可確定烏斯通溝工程場地特征周期Tg=0.40 s，根據《水電工程水工建筑物抗震設計規范》，標準設計反應譜平臺最大值βmax=1.60，最小值βmin=0.478，大于βmax的20%，滿足規范要求。查閱《中國地震動參數區劃圖》,可確定烏斯通溝地震峰值加速度為0.20 g，因此本文計算過程中以0.20 g 的地震峰值加速度作為設計規范波地震。

動力分析需確定該地區的地震加速度時程曲線,采用SIMQKE_GR 程序生成規范譜地震波，水平向（順河向和壩軸向）地震動加速度時程曲線的生成結果如圖2。 其中設計地震規范波水平向地震動峰值加速度為2.0 m/s2，豎向地震動加速度取水平向地震動加速度的2/3，峰值加速度為1.33 m/s2。

圖2 規范譜地震波合成（水平向）

4.2 計算成果與分析

4.2.1 壩體絕對加速度反應

大壩在設計地震規范波作用下順河向的絕對加速度反應極值分布如圖3（a）。由圖可知，在河床0+137 壩段附近的壩頂部位出現順河向地震加速度反應極值，放大倍數為2.61 倍，最大反應為522 gal。

圖3 設計地震規范波作用下壩體加速度反應極值分布 單位：gal

大壩在設計地震規范波作用下豎直向的絕對加速度反應極值分布如圖3（b），豎直向地震加速度反應極值出現在河床0+159 壩段附近的壩頂部位，放大倍數為1.47 倍，最大反應為293 gal。

大壩在設計地震規范波作用下壩軸向的絕對加速度反應極值分布如圖3（c），同順河向類似，在河床0+137 壩段附近的壩頂部位出現壩軸向地震加速度反應極值，放大倍數為1.94 倍，最大反應為387 gal。

由圖3 可知， 大壩3 個方向的加速度反應極值最大值均發生在壩頂部位，且幾乎保持對稱分布，符合實際規律。

4.2.2 壩體最大振幅

大壩在設計地震規范波作用下壩體的最大振幅分布如圖4，圖中數據為3 個方向的合位移。由圖4 可知，壩體最大振幅出現在壩頂，最大達5.4 cm。

圖4 設計地震規范波作用下壩體最大振幅 單位：cm

4.2.3 動剪應力極值

大壩在設計地震規范波作用下壩體最大橫剖面的動剪應力極值分布如圖5。Txy方向的動剪應力在壩體心墻與反濾料相接區約850.00 m 高程出現較大值， 最大動剪應力為144 kPa； 在心墻與反濾料相接區域，Txz方向的動剪應力較大， 最大動剪應力為122 kPa；Tyz方向的動剪應力極值主要存在于心墻與基座相接部位，最大動剪應力為75 kPa。

圖5 設計地震規范波作用下壩體動剪應力極值分布 單位：kPa

4.2.4 地震永久變形

大壩在設計地震規范波作用下壩體的順河向、豎直向和壩軸向的永久變形分布如圖6。 由圖6 可知，壩體順河向的永久變形基本指向下游，最大值為4.3 cm，出現在0+137 剖面；豎直向永久變形最大值為12.4 cm， 位于0+137 剖面的壩頂；壩軸向永久變形為兩岸向河谷內發展，基本成對稱趨勢分布，最大值左岸指向右岸為3.5 cm，右岸指向左岸為3.5 cm。

圖6 設計地震規范波作用下壩體震后永久變形 單位：cm

5 結語

（1）大壩在設計地震規范波作用下3 個方向的絕對加速度反應極值基本出現在壩頂部位，最大加速度反應的放大倍數約2.6 倍，存在“鞭梢效應”。

（2）壩體動剪應力峰值主要分布在心墻與反濾料相接的兩側區域，最大值為144 kPa；壩體最大振幅出現在壩頂，最大達5.4 cm。

（3）計算豎直方向地震永久變形最大值為12.4cm，位于0+137 剖面的壩頂偏上游區域。

（4）綜上所述，最大加速度反應、最大永久變形都位于壩頂，其分布規律合理，烏斯通溝瀝青混凝土心墻壩設計方案合理，采取一定抗震措施后，能滿足抗震要求。