徐金英，李德玉，郭勝山

（中國水利水電科學(xué)研究院工程抗震研究中心，北京100048）

基于ABAQUS的兩種庫水附加質(zhì)量模型下重力壩動力分析

徐金英，李德玉，郭勝山

（中國水利水電科學(xué)研究院工程抗震研究中心，北京100048）

庫水附加質(zhì)量反映了壩體與庫水間的相互作用。本文以Koyna重力壩為例，基于ABAQUS并輔以部分MATLAB程序計算了兩種不可壓縮性庫水附加質(zhì)量模型下壩體的自振特性。分析Koyna各計算條件下的自振頻率，結(jié)果表明：與庫水有限元模型相比，W estergaard模型夸大了壩體與庫水間的相互作用，應(yīng)當(dāng)進(jìn)行折減，折減系數(shù)取0.8較為適宜。考慮折減的Westergaard模型下壩體地震響應(yīng)計算結(jié)果與實際震害更符合，從而驗證了折減系數(shù)的準(zhǔn)確性。

ABAQUS；重力壩；動力分析；庫水有限元；Westergaard；折減系數(shù)

1 研究背景

重力壩作為擋水建筑物，地震時水庫將產(chǎn)生動水壓力，討論動水壓力的問題實際是一個庫水耦合問題［1］，因此，水體與壩體間的動力相互作用是其地震動力響應(yīng)分析中的一個主要特點［2］。1933年，美國學(xué)者Westergaard通過研究剛性直立壩面動水壓力，提出了不考慮庫水可壓縮性［3］的附加質(zhì)量模型。許多研究者針對Westergaard研究中的假設(shè)條件開展了多方面的研究，但都假定壩面是剛性的。毫無疑問，這一假定與大壩是可變形彈性體的實際情況相差較大。隨著有限元技術(shù)的發(fā)展，人們開始對壩體-庫水間的動力相互作用展開研究。陳厚群等［4］通過白山拱壩模型實測以及三向電模擬試驗求得振型動水壓力，發(fā)現(xiàn)有限單元數(shù)學(xué)模型給出的結(jié)果與試驗結(jié)果很符合，而Westergaard模型給出的結(jié)果偏大，建議對其進(jìn)行折半修正。李德玉等［5］認(rèn)為庫水有限元模型能更好模擬庫水對壩體動力特性和地震動力反應(yīng)的影響，而Westergaard模型夸大了庫水的動力影響，建議在適當(dāng)折減的基礎(chǔ)上采用，初步建議取0.7。雖然，庫水有限元模型更符合實際的壩體-庫水相互作用機理，但基于概念簡單、易于采用以及偏于安全等諸多因素的考慮，我國現(xiàn)行的抗震規(guī)范［6］在計算重力壩附加質(zhì)量時仍推薦采用Westergaard模型。所以，對簡單的Westergaard模型進(jìn)行折減，求解合理的拆減系數(shù)是一件非常重要的工作。由于陳厚群、李德玉計算的折減系數(shù)均采用振動臺試驗數(shù)據(jù)得到的，本文作者嘗試采用數(shù)值編程來實現(xiàn)，并與李德玉結(jié)果進(jìn)行比較，驗證計算程序的準(zhǔn)確性。

本文利用ABAQUS軟件和MATLAB程序，通過計算兩種附加質(zhì)量模型下壩體的自振特性，得到Westergaard模型相比于庫水有限元模型下的折減系數(shù)，并利用ABAQUS軟件分析考慮Westergaard模型折減后的重力壩動力響應(yīng)。

2 兩種附加質(zhì)量模型的基本原理

2.1 W estergaard模型假定壩面是直立剛性的，庫水上游方向無窮遠(yuǎn)，庫底為剛性水平面，庫水作無旋小變形運動并忽略表面波影響，則動水壓力表達(dá)式為：

式中：Pi為上游壩面的動水壓力，以壓為正；αi為節(jié)點i的附加質(zhì)量；ρ為水體密度；Hi為包含節(jié)點i的鉛直截面上的水深；Zi為包含節(jié)點i的鉛直截面上從壩基算起到節(jié)點i的高度；y?ni為在節(jié)點i的壩面法向加速度。

2.2 庫水有限元模型假設(shè)庫水為不可壓縮性流體，以動水壓力p表示的庫水運動方程為：

式中：Δ2為拉普拉斯（Laplace）算子。

由于本文只考慮水平向的附加質(zhì)量，所以壩體-庫水系統(tǒng)數(shù)值模擬的邊界條件如下：

不考慮庫水可壓縮性，經(jīng)有限元離散的庫水運動方程為

式中：[H]為庫水的剛度矩陣；{B}為壩面運動引起的“荷載”向量。

采用高斯法作靜力凝聚，使庫水運動方程只有壩體迎水面上的節(jié)點，即：

式中：[HS]=[HSS]-[HSW][HWW]-1[HWS]，其中，[HSS]為壩體自身節(jié)點間的剛度矩陣，[HWW]為水體自身節(jié)點間的剛度矩陣，[HSW]、[HWS]分別為壩體與水體節(jié)點的剛度矩陣。

令：

式中：[T]為由壩面動水壓力到壩面節(jié)點力的轉(zhuǎn)換矩陣；X?G（t）、X?S（t）分別為地面加速度和壩面節(jié)點上的相對加速度。

由此得到庫水有限元附加質(zhì)量表達(dá)式

3 MATLAB編程計算重力壩的自振特性

由于庫水有限元模型計算出來的壩面附加質(zhì)量矩陣是滿陣，不能直接應(yīng)用到ABAQUS中去，所以自行編寫了MATLAB程序來計算壩體的自振特性，從而實現(xiàn)兩種模型下壩體自振特性的比較。

3.1 M ATLAB計算庫水有限元模型下的附加質(zhì)量用MATLAB編程計算庫水有限元模型下的附加質(zhì)量矩陣時取3倍壩高的水庫域進(jìn)行有限元離散［7］，過程如下：（1）形成庫水單元剛度陣，從而形成整體剛度陣，并通過靜力凝聚得到僅與壩-庫接觸面上節(jié)點有關(guān)的等效剛度矩陣，考慮壩-庫邊界條件，進(jìn)一步形成簡化剛度陣；（2）計算壩面動水壓力與等效節(jié)點力向量之間的轉(zhuǎn)換矩陣[T]；（3）利用式（12）計算出庫水有限元附加質(zhì)量矩陣[Mp]1；（4）由（3）計算得到的附加質(zhì)量矩陣是滿陣，為了方便計算同時根據(jù)附加質(zhì)量矩陣元素的物理意義，將[Mp]1中的每列元素相加，作為相應(yīng)節(jié)點上的等效節(jié)點質(zhì)量，從而將[Mp]1“對角化”為[Mp]2。

3.2 M ATLAB編程計算兩種模型下重力壩的自振特性用MATLAB計算重力壩的自振特性，主要是利用子空間迭代法來實現(xiàn)［8］。過程如下：（1）形成單元剛度陣和單元集中質(zhì)量陣，從而形成整體剛度陣和整體質(zhì)量陣；（2）考慮壩體-庫水系統(tǒng)的相互作用，需要在（1）中得到的整體質(zhì)量陣的基礎(chǔ)上加上附加質(zhì)量陣，從而形成最終的等效整體質(zhì)量陣；（3）利用子空間迭代法，對由（1）得到的整體剛度陣和由（2）得到的等效整體質(zhì)量陣進(jìn)行迭代計算，完成壩體自振特性的求解。

4 兩種附加質(zhì)量模型下重力壩動力分析

本文是以Koyna重力壩為例進(jìn)行計算。Koyna重力壩位于印度西南部，1967年發(fā)生了里氏6.5級地震。震害情況為12～18號、24～30號壩段壩頂以下40 m左右產(chǎn)生了多條水平裂縫，下游面出現(xiàn)了嚴(yán)重的漏水現(xiàn)象，但是水庫水位無明顯下降。壩高103m，寬70m。地震發(fā)生時，水庫水位為91.75 m。壩體材料參數(shù)：彈性模量為31 027 MPa；泊松比為0.15；密度為2 643 kg/m3；剪脹角為36.31°；初始壓縮屈服預(yù)應(yīng)力為13.0 MPa；混凝土抗壓強度為24.1 MPa；混凝土抗拉強度為2.9 MPa。

圖1給出的是混凝土拉伸屈服應(yīng)力、拉伸損傷因子和開裂位移之間的關(guān)系。另外由于混凝土抗壓強度一般較大，故沒有考慮壓縮引起的損傷。

4.1 計算兩種附加質(zhì)量模型下Koyna重力壩的自振頻率圖2是ABAQUS和MATLAB計算Koyna重力壩自振頻率的有限元計算網(wǎng)格圖。

圖3（a）為利用MATLAB程序計算“對角化”后的庫水有限元附加質(zhì)量與Westergaard模型計算的結(jié)果圖；圖3（b）為圖3（a）中節(jié)點上兩種附加質(zhì)量的比例及其均值圖。

表1和表2分別為Koyna大壩空庫及滿庫條件下的前4階自振頻率結(jié)果。

由上述計算結(jié)果可知：（1）用軟件ABAQUS和MATLAB程序計算Koyna壩體的自振頻率在空庫和同種模型滿庫時的結(jié)果都十分接近，說明MATLAB程序是正確的；（2）以基頻為例，Westergaard模型下的結(jié)果相比于庫水有限元模型降低了5.59％，相比于庫水有限元模型“對角化”降低了2.66％。同樣的，對于第四階模態(tài)而言，則分別降低了12.09％和1.96％。由此可見，Westergaard模型和“對角化”后的庫水有限元模型都夸大了庫水相互作用，對Westergaard模型進(jìn)行折減是十分必要的；（3）由圖3（b）可知，Westergaard模型取折減系數(shù)為0.8較宜。

4.2 重力壩動力響應(yīng)分析圖4為Koyna實測地震波加速度時程曲線。

計算得到的折減系數(shù)，利用ABAQUS中的混凝土塑性損傷模型分析壩體動力響應(yīng)，Westergaard模型折減前后的損傷情況分布見圖5。

由圖5可以得知，考慮Westergaard模型折減后的計算結(jié)果顯示壩體損傷情況降低了，而且在下游折坡處得到了與實際震害較符合的結(jié)果。

圖6為按照文獻(xiàn)［5］中提出的0.7折減后的損傷分布情況，顯然與實際震害有出入，原因是文獻(xiàn)［5］中計算的重力壩的上游壩面不是直立的。兩圖的壩踵處均出現(xiàn)了橫向損傷，但在實際震害中是沒有的，原因是剛性地基沒有考慮地基的影響［11］。

4 結(jié)論

本文基于ABAQUS軟件并輔以部分MATLAB程序進(jìn)行了兩種庫水附加質(zhì)量模型下重力壩的動力分析，分析結(jié)果表明，相比于庫水有限元模型，Westergaard模型夸大了庫水相互作用，應(yīng)當(dāng)進(jìn)行折減，折減系數(shù)取0.8較為適宜，同時將折減后的Westergaard模型應(yīng)用到壩體動力響應(yīng)分析中，得到與實際震害相近的結(jié)果，從而驗證了折減系數(shù)取值的正確性。

［1］涂勁.混凝土大壩抗震數(shù)值分析理論與工程應(yīng)用［M］.北京：中國水利水電出版社，2007.

［2］王海波，李德玉.拱壩抗震設(shè)計理論與實踐［M］.北京：中國水利水電出版社，2005.

［3］陳厚群，吳勝興，黨發(fā)寧，等.高拱壩抗震安全［M］.北京：中國電力出版社，2012.

［4］陳厚群，侯順載，楊大偉.地震條件下拱壩庫水相互作用的試驗研究［J］.水利學(xué)報，1989（7）：29-39.

［5］李德玉，張伯艷，王海波，等.重力壩壩體—庫水相互作用的振動臺試驗研究［J］.中國水利水電科學(xué)研究院學(xué)報，2003，1（3）：216-220.

［6］DL5073-2000，水工建筑物抗震設(shè)計規(guī)范［S］.

［7］李瓚，陳興華，鄧建波，等.混凝土拱壩設(shè)計［M］.北京：中國電力出版社，2000.

［8］徐斌，高躍飛，余龍.MATLBA有限元結(jié)構(gòu)動力學(xué)分析與工程應(yīng)用［M］.北京：清華大學(xué)出版社，2009.

［9］趙光恒.結(jié)構(gòu)動力學(xué)［M］.北京：中國水利水電出版社，1996.

［10］黃耀英，孫大偉，田斌.兩種庫水附加質(zhì)量模型的重力壩動力響應(yīng)研究［J］.人民長江，2009，40（7）：64-66.

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

The dynam ic analysisof gravity dam under tw o added w ater m assm odels based on ABAQUS

XU Jin-ying，LI De-yu，GUO Sheng-shan

（Earthquake Engineering Research Center，IWHR，Beijing 100048，China）

Added water mass reflects the interaction between dam and water.Taking the gravity dam Koyna as an example in this paper，the natural vibration characteristics of the dam are calculated by using ABAQUS and MATLAB，under two different added water mass models which assume that the water is in?compressible.The result shows that the vibration frequency from the Westergaard model exaggerates the in?teraction between the dam and water compared with that of the finite element model and should be re?duced with a reasonable reduction factor of 0.8.The result of Westergaard model by considering the reduc?tion factor is consistent with the real damage of earthquake，thus verifying the accuracy of reduction factor values.

ABAQUS；gravity dam；dynamic analysis；finite element for reservoir；Westergaard；Reduc?tion factor

TV642

A

10.13244/j.cnki.jiwhr.2014.01.017

1672-3031（2014）01-0098-06

（責(zé)任編輯：王冰偉）

2013-03-04

徐金英（1989-），女，山東濟(jì)南人，碩士生，主要從事水工結(jié)構(gòu)抗震研究。E-mail：908464507@qq.com