基于附加质量法的三维重力坝地震响应研究

王志坤，陈艳江，杨　璐,*

(1.沈阳工业大学 建筑与土木工程学院，沈阳 110870;2.黑龙江大学 水利电力学院，哈尔滨 150080)



基于附加质量法的三维重力坝地震响应研究

王志坤1，陈艳江2，杨璐1,*

(1.沈阳工业大学 建筑与土木工程学院，沈阳 110870;2.黑龙江大学 水利电力学院，哈尔滨 150080)

应用ABAQUS有限元软件二次开发功能，编写三维附加质量单元子程序UEL，实现了三维模型中节点附加质量的快速施加，并模拟Koyna混凝土重力坝地震破坏过程。研究发现：坝体中上部特别是坝面突变处是抗震薄弱部位，坝顶以下40 m左右出现贯穿裂缝；下游折坡点位置易发生应力集中，瞬时竖向拉应力达到12.738 MPa，应尽量采用圆弧等形式避免坝面突变。地震过程中坝体上部震动剧烈，特别是坝顶处最大震动幅度为41.42 mm。计算结果基本符合1967年Koyna大坝震害记录数据。数值结果证明三维附加质量法子程序的开发是正确可靠的，用该二次开发子程序方法解决此类三维水体—结构相互作用问题是可行的。

附加质量法；动水压力；重力坝；时程分析

近年来，我国水利事业蓬勃发展，建设了一批大型水利工程，但我国是一个地震多发国家，不少已建、在建、拟建的高坝大库都位于强震区：双江口堆石坝设计地震加速度0.205g，小湾拱坝设计地震加速度0.31g，大岗山拱坝设计地震加速度0.557g[1-2]，这些高坝大库的发展也引起了专家学者对水坝抗震性能的重视，而其中一个重要的方面便是对地震作用下动水压力的研究。

目前工程界处理地震作用下水坝动水压力的问题采用的是Westergaard提出的不考虑水体可压缩性的附加质量法[3-4]。但目前所有的大型有限元软件中，均没有相应的附加质量单元类型，ABAQUS中仅提供了二维问题的附加质量单元快速施加程序，但二维模型应用范围有限[5]。鉴于此本文基于ABAQUS软件用户子程序接口，编写三维附加质量单元子程序UEL，实现了三维模型中节点附加质量的快速施加，为今后模拟此类三维水体-结构相互作用问题提供了简便快捷的处理方法。

1　计算原理

在地震荷载作用下，加速度激励是随时间不断变化大小和方向的，坝体产生与之对应的往复加速晃动，水体与坝体之间产生大小和方向也不断变化的相对惯性作用力以及相对滑动。水体对坝体的作用可分为两部分压力的叠加：①静水压力；②地震作用下的惯性力，即动水压力[6]。

Westergaard假定坝面最大动水压力沿水深呈抛物线分布，并根据实际动水压力对坝踵力距与近似动水压力图形对坝踵的力距相等的条件导出了韦氏附加质量公式[7]：

(1)

式中ρω为水的密度；hω为设计水深；z为到水面的距离。

由于该方法简单实用、易于计算、偏于安全，美国、日本等国家的建筑物抗震设计规范中至今仍沿用此忽略坝体变形和水体可压缩性的动水压力公式[8-10]，我国《水工建筑物抗震设计规范》亦是采用附加质量模型作为重力坝抗震设计中动水压力模拟方法[11]。

2　模型建立

2.1工程概况

迄今为止，全球范围内的大坝遭受强震震害的实例稀少，遭受Ⅷ度以上强震的百米级重力坝仅有4例，即中国的新丰江大坝、宝珠寺大坝，伊朗的西菲罗(Sefid Rud)大头坝，印度的柯依那(Koyna)重力坝[1]。其中，Koyna大坝作为少数几个在强震中破坏且有比较完整记录的重力坝，一直是重力坝抗震分析中的经典研究对象，已有大量学者[12-15]对该坝的地震破坏过程进行了模拟分析。

1967年12月11日Koyna坝址区域遭受了一次6.5级强烈地震，地震发生时的坝前水位91.75 m[12]。地震作用后，坝顶以下40 m左右产生了多条裂缝，下游面出现了严重的漏水现象。

2.2计算模型

Koyna混凝土重力坝位于印度的Koyna河上，坝长850 m，坝高103 m，坝底宽70 m，坝顶宽度14.8 m，地震发生时坝前水位为91.75 m，具体尺寸见图1。

混凝土重力坝是大体积结构，为防止开裂、满足施工要求往往分缝分块浇注，需沿坝轴线方向按约20 m设置横缝，横缝需在坝体冷却至稳定温度时经灌浆后形成整体，但灌浆的浆体仅能起传递压应力的填充作用，抗拉强度极低以至可以忽略。故而重力坝抗震分析一般取单个坝段进行计算，坝段间横缝的影响可以不计。

依据上述Koyna混凝土重力坝尺寸选取单个坝段进行模拟计算，坝段宽度取为20 m，分别自上、下游坝踵、坝趾位置延伸2倍坝高，同时向下延伸1.5倍坝高，作为坝基，建立有限元模型见图2，采用C3D8R单元，即8节点实体减缩单元，上游坝体与水体接触表面采用附加质量单元，模型中共33 792个节点，29 850个单元，其中29 200个C3D8R单元，650个附加质量单元。

图1　坝体剖面(单位m)Fig.1　Dam section profile(company:m)

图2　计算模型Fig.2　Computational model

2.3计算参数

在模型中坝体材料采用混凝土损伤塑性本构模型，坝基部分采用弹性模型。根据文献[11]的规定，在抗震强度计算中，混凝土动态强度和动态弹性模量的标准值可较其静态标准值提高30%，混凝土动态抗拉强度的标准值可取为动态抗压强度标准值的10%[11]。依据相关文献及Koyna大坝工程简介，坝体混凝土各项力学参数见表1。

表1　坝体力学参数

另外，由于混凝土抗压强度一般较大，故而没有考虑压缩引起的损伤[13]，混凝土拉伸屈服应力、拉伸损伤因子和开裂位移之间的关系见图3。

图3　混凝土塑性指标Fig.3　Concrete plastic properties

2.4地震波

在进行坝体结构动力分析中，同时计入坝体和地基的刚度，但只计入坝体的质量，地基单元只考虑弹性不考虑质量，以消除波的传播效应，避免人为的放大作用，这样的处理方式较为接近实际工程情况。

地震波在空间的传播方向十分复杂，不同方向的地震波对重力坝影响差异较大，模型中在坝基底面边界施加水平和垂直方向地震加速度，选取1967年Koyna大坝实测的柯依那波进行计算，时间间隔0.01 s，地震总时程10 s。地震波时程曲线见图4。

图4　地震加速度时程曲线Fig.4　Seismic acceleration time history curve

3　模型分析

模态分析是各种动力学分析类型中的基础内容，在进行其他动力学分析之前首先要进行模态分析，采用ABAQUS有限元软件中的Lanczos法提取大坝前4阶频率。表2、表3分别为Koyna重力坝空库和满库状态下前4阶频率，通过与文献[14]对比可知误差较小，证明建立的模型及开发的三维附加质量单元用户子程序是正确的。

表2　Koyna重力坝(空库)前4阶频率

表3　Koyna重力坝(满库)前4阶频率

4　时程分析

运用ABAQUS有限元软件对三维Koyna坝体进行动力时程运算并分析计算结果，分别提取坝体上关键位置点，并研究其位移、应力等随时间的变化情况，从而揭示混凝土重力坝在地震作用下的一般反应规律。关键位置点布置见图1，其中A、B分别为为上、下游坝顶，C为坝踵，D为坝趾，E为下游折坡点。

地震持续过程中坝体上各特征点在各个方向上应力和位移的峰值及相应峰值出现的时刻见表4，其中水平位移负值为上游方向，竖直位移负值为垂直向下方向，负应力值代表拉应力。由表4可见，坝顶处应力很小，震后也未出现塑性破坏，只是在地震过程中晃动比较严重，出现较大位移。坝踵和坝趾处在地震中的晃动几乎可以忽略，水平方向的应力也很小，在地震波最剧烈的4～5 s内坝踵处的竖直方向应力较大，产生了塑性应变，这是由于刚性地基没有考虑地基的影响造成的，实际震害中不存在这种状况，故不加分析。

相关震害记录及研究资料表明，重力坝上部尤其是断面突变处，即折坡点位置，是抗震薄弱部位。表4及图5也表明地震对重力坝的破坏作用主要集中在断面突变处，故着重分析折坡点E处应力、位移等在地震中的变化规律。

表4　各特征点峰值应力、位移

4.1应力时程曲线

地震过程中特征点E的应力时程曲线见图5。由图5可见，折坡点处水平和竖直方向应力变化趋势基本相同，只是峰值不同，水平应力峰值为4.522 MPa，竖直应力峰值为12.738 MPa，均超过了混凝土动态抗拉强度，出现了较严重的损伤破坏。两个方向应力峰值均出现在4.8 s左右，这也是坝体裂缝扩展最快的时刻，略滞后于地震波最大加速度出现时间。

在地震持续过程中，前期由于地震加速度小，坝体能够抵御由震动引起的破坏作用，两个方向的应力均低于混凝土的抗拉强度，坝体只发生晃动，并未出现裂缝等不可逆转的破坏现象。随着地震波加强，地震的破坏作用超过了坝体自身抵抗能力，坝面突变处产生应力集中，应力值随着地震的持续急剧增加，尤其是竖向应力值，最大达到抗拉强度的4倍以上。地震作用末期，加速度减弱，拉应力也相应减小，但竖向拉应力依然较大，直到地震结束。

图5　折坡点E应力时程曲线 Fig.5　Stress time history curve of slope point E

4.2位移时程曲线

下游折坡点E的时间位移时程曲线见图6，其中水平方向主要是随着地震过程的持续坝体沿上下游方向晃动，晃动幅度基本相同，均随着地震加速度的加强而增大，但略滞后于加速度，在4.48 s时达到最大位移11.43 mm。随后地震作用减弱，坝体晃动幅度减小，至地震作用结束时坝体基本回到原来位置，仅向上游方向偏移1.45 mm。

在地震作用前期，坝体在竖直方向仅发生上下晃动，晃动幅度在2 mm以内，但在地震加速度峰值过后，下游坝面突变处即E点处出现塑性损伤，4 s后E点的竖直方向出现2 mm左右的位移，并一直持续到地震结束，表明在E处出现了宽度为2 mm左右的裂缝，与1967年的Koyna大坝震害记录吻合。

图6　折坡点E位移时程曲线Fig.6　Displacement time curve of slope point E

4.3损伤分析

Koyna大坝是国外遭受Ⅷ度以上强震的两个百米级重力坝之一，震灾情况为12～18号、24～30号坝段坝顶以下40 m左右位置产生了多条水平裂缝，下游面出现了严重的漏水现象，即裂缝已经贯通了上下游。说明重力坝上部，尤其是断面突变处，是其抗震薄弱部位，在强震作用下，上部坝体易开裂[16]。震后Koyna坝体保持了整体稳定性，坝基面未出现拉裂和剪切损坏[17]。

模拟的地震过程中坝体断面拉伸损伤云纹图见图7。由图7可见，坝体破坏主要集中在下游折坡点处。地震前期地震动加速度较小，未对坝体造成损伤破坏，图7(a)为2 s时的坝体断面拉伸损伤云纹图，损伤值为零。2 s后地震加速度变大，震动剧烈，坝体局部区域主要是下游折坡点位置及坝顶以下40 m左右处，4 s左右坝体开始出现明显的破坏，损伤区域由坝面处向坝体内部扩展，逐渐形成贯穿上下游的水平裂缝。5 s后地震作用减弱，破坏面积没有进一步扩大，主要是已经形成的贯穿裂缝继续发展，严重的贯穿裂缝导致下游面出现漏水现象，符合1967年的震害记录。

坝踵处出现了轻微的损伤破坏，原因是刚性地基没有考虑地基的影响[18]，但震后坝体保持了整体稳定性，亦符合当时的实际情况。说明三维附加质量法的运用是合理可行的，能相对准确地模拟地震作用下水体对坝体的相互作用情况，为今后工程界处理此类问题提供了简便可靠的解决方法。

图7　坝体断面损伤云纹图Fig.7　Damage moire of dam section

5　结　论

依据1967年印度Koyna混凝土重力坝地震相关资料数据，建立坝体三维计算模型，采用Westergaard提出的附加质量法考虑地震作用下动水压力的影响，利用ABAQUS软件用户子程序接口UEL，开发编写三维附加质量单元子程序。进行地震时程计算，并分析坝体损伤、应力、位移等在地震过程中的响应规律，得出以下结论：

1) 计算结果与1967年震害记录基本吻合，说明三维附加质量子程序的开发是合理可行的，为工程界处理类似三维水体-结构动水相互作用问题提供了新的简单快捷的方法。

2) 重力坝在地震过程中的损伤破坏主要集中在坝体中上部，坝面突变处是其抗震薄弱部位，特别是下游折坡点位置易出现贯穿裂缝，应引起足够重视。

3) 下游折坡点处易产生应力集中，应力值远大于混凝土的动态抗拉强度值，最大拉应力达到12.738 MPa，应避免出现坝面突变，折坡点处尽量采用圆弧等形式。

4) 地震中坝顶震动剧烈，主要沿上下游方向晃动，最大振幅达到41.42 mm，坝顶处构造设施应采取相应措施避免发生破坏。

[1]潘家铮. 水电与中国[M]. 北京：中国水利水电出版社,2007.(PAN Jia-Zheng. Hydro-power and China[M]. Beijing: China Water Conservancy and Hydro-power Press,2007.(in Chinese))

[2]陈厚群. 混凝土高坝抗震研究[M]. 北京：高等教育出版社，2011.(CHEN Hou-Qun. Study on seismic behavior of concrete dam[M]. Beijing: Higher Education Press,2011.(in Chinese))

[3]Lee J, Fenves G L. A plastic-damage concrete model for earthquake analysis of dams[J].Earthquake Engineering and Structural Dynamics,1998,27(4)：937-956.

[4]Omidi V S. Seismic cracking of concrete gravity dams by plastic-damage model using different damping mechanisms[J].Finite Elements in Analysis and Desing,2013,63(6)：80-97.

[5]徐磊. 基于ABAQUS的四节点三维附加质量单元开发[J]. 水电能源科学，2011,29(6)：100-102.(XU Lei. Development of four nodes 3-D added mass element in ABAQUS[J]. Water Resources and Power, 2011,29(6)：100-102.(in Chinese))

[6]宋波，刘浩鹏，张国明. 基于附加质量法的桥墩地震动水压力分析与实例研究[J]. 土木工程学报，2010,43(增刊)：102-107.(SONG Bo, LIU Hao-Peng, ZHANG Guo-Ming. Analysis and case study of earthquake-induced hydrodynamic pressure on deep water bridge based on added masses approach[J]. China Civil Engineering Journal,2010,43(s):102-107.(in Chinese))

[7]Westergaard H M. Water pressures on dams during earthquakes[J].Transactions of the American Society of Civil Engineers, 1933,98(3)：418-433.

[8]王铭明，陈健云，范书立. 重力坝地震动水压力试验研究[J]. 水电能源科学，2012,30(5)：51-53.(WANG Ming-Ming, CHEN Jian-Yun, FAN Shu-Li. Experimental study on seismic water pressure of gravity dam[J]. Water Resources and Power, 2012,30(5)：51-53.(in Chinese))

[9]Pailoplee, Santi. Earthquake hazard of dams along the Mekong mainstream[J]. Natural Hazards,2014,74(3)：1813-1827.

[10] Das, R., Cleary, P. W. A mesh-free approach for fracture modeling of gravity dams under earthquake[J].International Journal of Fracture,2013,179(2)：9-33.

[11] 中华人民共和国水利部.水工建筑物抗震设计规范：SL203-97[S].北京：中国水利水电出版社，1997. (Ministy of Water Resources of the People’s Republic of China.Code for seismic design of hydraulic structures：SL203-97[S]. Beijing:China Water Conservancy and Hydropower Press,1997.(in Chinese))

[12] 徐金英，李德玉，郭胜山. 基于ABAQUS的两种库水附加质量模型下重力坝动力分析[J]. 中国水利水电科学研究院学报，2014,12(1)：98-103.(XU Jin-Ying, LI De-Yu, GUO Sheng-Shan. Dynamic analysis of gravity dam under two kinds of reservoir water quality model based on ABAQUS[J]. Journal of China Institute of Water Resources and Hydropower Research, 2014,12(1)：98-103.(in Chinese))

[13] 杨璐，金峰．水工混凝土结构中的数值计算与实例[M]. 北京：科学出版社,2015.(YANG Lu, JIN Feng. Numerical calculation and examples of hydraulic concrete structures[M]. Beijing: Science Press,2015.(in Chinese))

[14] 赵光恒. 结构动力学[M]. 北京：中国水利水电出版社，1996.(ZHAO Guang-Heng. Structural dynamics[M]. Beijing: China Water Conservancy and Hydro-power Press,1996.(in Chinese))

[15] 黄耀英，孙大伟，田斌. 两种库水附加质量模型的重力坝动力响应研究[J]. 人民长江，2009,40(7)：64-66.(HUANG Yao-Ying, SUN Da-Wei, TIAN Bin. Study on dynamic response of gravity dam with two kinds of reservoir water quality model[J]. Yangtze River, 2009,40(7)：64-66.(in Chinese))

[16] 方修君，金峰，王进廷. 基于扩展有限元法的Koyna重力坝地震开裂过程模拟[J]. 清华大学学报:自然科学版. 2008,48(12)：2065-2069.(FANG Xiu-Jun, JIN Feng, WANG Jin-Ting. Simulation of seismic cracking of Koyna gravity dam based on extended finite element method[J]. Journal of Tsinghua University(Science and Technology), 2008,48(12)：2065-2069.(in Chinese))

[17] 邵长江，钱永久. Koyna混凝土重力坝的塑性地震损伤响应分析[J]. 振动与冲击，2006,25(4)：129-131.(SHAO Chang-Jiang, QIAN Yong-Jiu. Analysis of plastic seismic damage response of Koyna concrete gravity dam[J]. Journal of Vibration and Shock, 2006,25(4)：129-131.(in Chinese))

[18] Kartal, Murat Emre,Bayraktar. Non-linear response of the rock-fill in a concrete-faced rock-fill dam under seismic excitation[J]. Mathematical and Computer Modelling of Dynamical Systems,2014,21(1)：77-101.

Research on seismic response of 3-D gravity dam based on added mass method

WANG Zhi-Kun1, CHEN Yan-Jiang2, YANG Lu1,*

(1.SchoolofArchitectureandCivilEngineering,ShenyangUniversityofTechnology,Shenyang,110870,China;2.SchoolofHydrauicandElectricPower,HeilongjiangUniversity,Harbin150080,China)

3-D added mass subroutine UEL based on the redevelop function of ABAQUS is developed，fast application of added mass in 3-D model is realized, and seismic damage process of Koyna concrete gravity dam is simulated. According to the result: middle upper part especially the mutation of dam body is weak part in Seismic, and through cracks appear below the dam crest of 40m around; The position of the downstream slope occurs stress concentration, the instantaneous vertical tensile stress reaches 12.738 MPa, and the circular arc and other forms should be used to avoid the mutation of the dam surface. In earthquake process, the upper shakes violently especially the maximum vibration amplitude of the dam crest reaches 41.42 mm. The calculation results are in accord with the damage record of Koyna dam in 1967, the development of 3-D added mass subroutine is correct and reliable, and the method to deal with the interaction problem of 3-D water-structure is feasible.

added mass method；hydrodynamic pressure；gravity dam；time history analysis

10.13524/j.2095-008x.2016.03.035

2016-05-06

国家自然科学基金资助项目(11102118)

王志坤(1990-)，男，河北昌黎人，硕士研究生，研究方向：水工结构抗震及ABAQUS数值模拟， E-mail：xxzhikun@163.com;*通讯作者：杨璐(1973-),女，山东掖县人，教授，博士，硕士研究生导师，研究方向：混凝土弹塑性损伤本构和ABAQUS的数值模拟，E-mail:yanglu515@163.com。

TV312

A

2095-008X(2016)03-0012-08