考虑横缝接触耦合坝体材料非线性下的拱坝损伤开裂分析

张建伟，刘鹏飞，王 涛，焦延涛，李兆恒，赵建军

(1.华北水利水电大学水利学院，河南 郑州 450046；2.水资源高效利用与保障工程河南省协同创新中心，河南 郑州 450046；3.河南省水工结构安全工程技术研究中心，河南 郑州 450046；4.广东省水利水电科学研究院，广东 广州 510635)

对大坝等大型水工结构进行抗震分析计算时，相关规范中规定了水工结构相应的拉应力与压应力的控制标准，其思想是以线弹性理论为基本准则进行控制，黄熠辉等[1]校核检验了坝体结构在线弹性影响因素下的杨房沟拱坝的抗震能力。考虑到坝体混凝土材料在极限荷载下所呈现的非线性应力-应变软化特征，特别是在地震加速度瞬时增大的情况下，坝体混凝土结构会出现损伤，结构刚度下降进而出现结构失效的现象时，线性控制准则不再适用。故应选择合理的塑性本构模型研究地震荷载作用下坝工结构的薄弱区域分布及结构的失稳问题，以达到还原真实结构地震应力场的目的[2]，同时杜小凯等[3]指出对于高坝结合混凝土材料非线性的动力试验依旧是未来研究方向。

拱坝需要用大体积混凝土进行浇注，作为混凝土粘合剂的水泥在硬化的过程中会不断产生水化热，如果这些热量不及时排出坝体，会在拱坝的内部结构产生细小的微裂隙，微裂隙在坝体内部扩展、汇集，最终形成宏观裂纹损害坝体结构。为及时将拱坝在浇筑过程中产生的热量散射出去，通常在浇筑时对坝体进行分块施工，部分沿径向铅垂布置的施工缝被称为横缝[4]。同时，拱坝在结构特性上被定义为高次超静定结构，无论温升或温降的情况下都会沿拱坝轴线产生不利于结构的温度应力，横缝的存在一定程度上缓解了温度应力带来的沿径向变形损害。但是在地震荷载作用下，横缝的不断张开闭合破坏了坝体结构的完整性，削弱了拱的作用，横缝是拱坝结构在地震分析中不可忽视的因素。

涂劲等[5]认为拱坝的安全问题都将通过坝体的位移突变和不断增长以致丧失承载能力反映出来，所以可据此提出拱坝的破坏准则。但该准则需要预设有初始强度的接触缝(如：坝基交界面接触缝、坝肩滑块接触缝等)，这就存在人为的主观性和经验性，因此对于横缝的非线性接触行为理应慎重选择更为合理的接触模型。胡昱等[6]提出了缝面蓄能概念并研究了施工期拱坝横缝增开现象及成因；赵兰浩等[7]提出了摩擦接触问题中有初始间隙的有限元混合法，并分析了不同库水模型对横缝开度的影响；郭胜山[8]建立了较为完整的拱坝-地基-库水损伤模型，着重研究重力坝在地震灾害下的坝体-地基损伤状况。

本文基于Lee和Fenve[9]提出的弹塑性损伤模型耦合Chaudhary和Bathe[10]横缝的非线性接触的动接触力模型，根据实际资料构造地基-大坝-库水三维有限元模型，着重考虑了拱坝在地震之初由坝体自身重力层层堆叠平衡所形成的初始应力场，并分别进行了不同峰值加速度(PGA=0.2g、0.4g、0.6g)的地震动响应分析，以便获取不同地震工况下的拱坝损伤结果。

1 基本原理

1.1 塑性损伤原理

由于准脆性材料表现出抗压不抗拉的特点，特别是以混凝土材料为代表承受荷载时呈现出的塑性损伤的特性，Lee和Fenves[9]根据其特点演化出一种各向同性，兼顾拉伸与压缩的塑性损伤模型——CDP本构模型，能够充分演示材料弹塑性损伤的不可逆过程。

弹塑性应变理论根据应变的增量特性可分解为弹性应变部分和塑性应变部分：

(1)

(2)

式中，E0为初始弹性模量。

为了将混凝土出现损伤后所受应力与有效应力结合起来，特引入损伤因子的概念代为表述：

(3)

式中，d为损伤因子，可分解为拉伸部分损伤因子dt和压缩部分损伤因子dc。对应同一结构在地震荷载作用下瞬时受拉或者受压的特点。同时从数值上表示为当d=0时，代表结构完好，未出现损伤；当d接近1时，表示结构失效，严重部位退出结构工作。于是d的假定表示为

1-d=(1-stdc)(1-scdt)

(4)

式中，st、sc分别为拉、压应力荷载作用下，结构刚度复原应力状态函数，可表示为

(5)

式中，wt、wc分别为拉、压应力荷载作用下结构刚度的恢复系数。

综上可得混凝土的刚度削弱结果，即

Dc=(1-dc)E0

(6)

Dt=(1-dt)E0

(7)

式中，Dc为结构的单轴受压损伤演化参数；Dt为结构的单轴受拉损伤演化参数。

1.2 横缝非线性接触模型原理

(8)

图1 主从面接触关系

(9)

(10)

(11)

2 工程概况及模型简介

2.1 工程主体建模

本文以拉西瓦双曲拱坝为例建立拱坝有限元模型。工程建基面高程2 210 m，顶部高程2 460 m，坝顶基本宽度为10 m，坝底基本宽度为49 m，最大坝高250 m，坝前正常蓄水位下的工作水深为242 m[11]。限于坝体材料的非线性与横缝接触的非线性的复杂性，徐轶慷[12]曾提出沿拱圈径向设置3条横缝就可比较合理模拟出实际情况下的各种物理指标，故沿径向均匀设置5条横缝，如图2所示，图3为建立的整体有限元模型及坝体中部切片。

图2 拉西瓦双曲拱坝坝体横缝布置示意

图3 坝体有限元模型及坝体中部切片

2.2 工程坝址处地震波选取原则

结合振型分解反应谱法，查询得到罕遇地震荷载情况下的百年基准期超越概率为2%的坝址处的地基基岩峰值加速度为0.23g[13]，同时水工结构抗震规范的规定拱坝结构的最大抗震放大系数为2.5，对顺河向水平地震波进行处理，设置峰值加速度分别为0.2g、0.4g、0.6g的地震工况，其他2个方向的地震波各自乘以相应的换算后的放大系数，处理地震波后输入有限元模型进行地震动力计算。图4为所选地震波三方向反应谱与对应拱坝的规范反应谱。

图4 所选地震波三方向反应谱与对应拱坝的规范反应谱

该条记录的天然地震波持续时间为101.11 s，并按相关标准选取前10 s地震波数据进行拉西瓦拱坝地震动力时程计算，计算步长设置为0.02 s。

3 试验设计及结果分析

3.1 坝体地震动仿真试验设计

为能模拟坝体从初建到正常蓄水位下的工作状态再到遭遇地震的全过程，必须先要进行坝体地震前的初始应力场的模拟，故采用有限元软件的生死单元技术还原拱坝完建时期到正常蓄水位工作状况下的初始应力场的静力场问题，具体为：

(1)将坝体沿拱圈径向均匀设置5条横缝，并暂时关闭非线性接触的相关属性，再将坝体均匀分成若干拱圈层，逐层激活施加自重，计算得到考虑施工过程中形成的层层拱圈应力平衡状态，最后在激活顶拱单元集后，激活横缝非线性接触的相关属性。

(2)施加静水压力荷载，得到正常蓄水位工况下地震发生前的坝体初始应力分布场。

(3)进行地震动荷载分析步的计算，获取坝体横缝开合效应的结果及坝体损伤开裂结果。

图5 横缝开合度动力时程

3.2 坝体地震动结果分析

3.2.1横缝张合结果分析

图5为模型沿拱圈径向均匀布置的5条横缝在地震动荷载下所呈现的张开闭合动力时程。其中，在0～10 s为模型的静力分析，0～9 s内采用有限元生死单元技术层层叠加水平向的坝体拱圈，最终在激活顶层拱圈后施加坝体结构的自重。在9～10 s，激活前期因采用生死单元技术层层叠加水平向拱圈而暂时关闭的横缝接触关系，同时施加静水压力等静力荷载；在10～20 s，进行坝体的地震动力时程分析。从结构力学的角度来看，因拱形承受压力并由顶点扩散至坝端的力学特性，故选取拱坝中部横缝C作为地震荷载张合效应的主要研究对象。

由图5可知，在9～10 s，因施加静水压力等静力荷载，中缝C略微张开了2 mm，但在激活设置横缝间的接触关系后，中缝C的开度又回归至0，说明未设置横缝接触关系前，正常蓄水位下的静水压力使坝体微微变形，横缝呈现略微张开状态，设置接触关系后又回归闭合状态，与现实中正常蓄水位工作下拱坝横缝呈现闭合状态的特征一致，说明了试验设置横缝非线性接触因素的正确性。对比邻近缝B、D在同时段内略微呈现出相似的特征，如缝B、D分别张开0.7、1 mm，而缝A、E未见明显特征，也符合拱形承受压力是由拱中部传递至拱端的一般规律，并且可看出缝A、E由于靠近拱端，即使在未设置横缝间接触关系的情况下，在静水荷载的压力下，拱坝肩部附近的坝段与坝段间所受角度以及应力大小原因紧密压实，所以未呈现明显的开合效应。

进入10～20 s的地震动力时程后，3种工况下的地震加速度在未超越一定峰值的前期，横缝的开合状况基本一致，未出现明显分化。因拱坝的坝段间所受应力不相同，但大体可以确定在10～12 s时段内，3种工况下因横缝接触因素作用的影响，5条横缝开合效应基本一致。随着后期峰值加速度的不断增大，3种工况下的横缝开合效应出现明显的差异性，其中0.6g地震工况下的横缝开合效应最为突出，0.4g、0.2g工况下的开合度依次递减。横缝A、B、C、D、E的在0.6g峰值加速度下的最大开度依次达到了8、25、38、31、5 mm。同时因耦合了坝体混凝土材料的非线性软化因素，0.2g和0.4g工况下的坝体横缝在地震持时结束后趋于闭合状态，而0.6g工况下坝体横缝则出现了一定的残余开度，是由于坝体出现严重损伤开裂造成的，可与后面提取的坝体损伤开裂结果综合进行比较分析。

图6针对横缝A、B、C、D、E在不同地震工况下的开合所体现出的差异性，提取的每条横缝极值开合包络图。包络图对横缝在整个地震持时中横缝张开的最大值进行了合理的归纳统计，是检验了横缝设置的合理性的一种手段。因横缝沿拱圈径向均匀布置，限于展示角度的因素，图6中可对比缝A、B、C在3种地震工况下呈现的差异性。每一段坝段在地震荷载作用下除去坝段间的拱向作用可看成是一根根伫立在地基上的悬臂梁，梁端的位移最大，同时也是横缝最大张开值所在的区域，且缝A、B、C的张开值依次增大。在3种地震工况下包络图显示横缝整体张开极值分别为14、27、38 mm。

图6 横缝开合极值包络图示意

3.2.2坝体损伤开裂结果提取

坝体损伤变量是因混凝土的材料刚度不可恢复的退化所表现的非线性损伤行为，在整个地震动试验时程中是一个不断累积的物理量，故坝体损伤的结果只需提取地震时程结束后最后一个分析步所呈现的坝体损伤结果。

图7为在0.2g峰值加速度的地震荷载下坝体损伤情况。坝体出现了局部损伤开裂区域，大体上都沿横缝边缘出现，除去坝踵与坝基交接面处的因角缘效应引起的应力集中现象，每一坝段的损伤行为也符合混凝土构件的损伤开裂行为，即裂缝先从边缘处出现，随后不断的由边缘向内汇聚扩展，最后形成宏观的裂缝。对比于0.2g工况下下游面未出现明显的损伤，可判定地震中形成的裂缝并未贯穿坝体。同时由于地震中横缝的开合行为，削弱了拱的作用，从而导致梁向作用显著增加，坝体出现的中下部出现的损伤行为也说明此种情况。

图7 坝体损伤示意(0.2g工况)

图8为在0.4g峰值加速度的地震荷载下坝体损伤情况。此时，坝体的损伤区域进一步扩展，特别是沿横缝垂直方向的分布有汇聚成片的趋势，且不同横缝间的损伤区域有汇聚在一起的趋势。同时在0.4g地震工况下坝体下游面也出现了损伤，而相对应的上游面相应区域却未见明显损伤，说明坝体上下游损伤区域还未贯通。但下游面的损伤主要集中在坝顶高程的中上部，且对比于0.2g地震工况下的下游面损伤，损伤区域有明显增加，同时说明了坝体上下游面的损伤是各自开展的，最终是通过上下游的裂缝扩展汇聚，最后裂缝贯通坝体，造成坝体的完全破坏。

图8 坝体损伤示意(0.4g工况)

图9为在0.6g峰值加速度的地震荷载下坝体损伤情况，坝体损伤区域进一步加剧，坝体损伤区域基本连接成一片，且对比于下游损伤开裂的区域，出现了上下游贯通的现象，表示坝体结构已经被损坏，破损较严重的区域已经退出工作。特别是下游面的损伤更为集中，主要出现在拱坝的中上部，贯穿裂缝最先出现的区域也可判定为拱坝结构的薄弱区域，在以后的设计中应着重予以考虑并对薄弱区域用特殊材料进行补强。

图9 坝体损伤示意(0.6g工况)

4 结 论

本文利用ABAQUS软件建立坝体-地基耦联体系有限元模型，取正常蓄水位的静水压力并兼顾地震荷载作用产生的坝前水体的动水压力；采用混凝土塑性损伤本构模型模拟坝体混凝土材料非线性力学特性并耦联拱坝横缝非线性接触的相关因素，同时根据坝体结构的场地信息选择符合规范反应谱具有天然频谱特性的San Fernando地震波处理后进行0.2g、0.4g、0.6g3种地震峰值加速度工况下的坝体损伤分析，重点研究坝体横缝的非线性接触行为以及复杂因素耦合场下拱坝的损伤开裂行为。可知，0.6g工况下的横缝开合效应最为突出，横缝A、B、C、D、E在0.6g工况下的最大开度依次为8、25、38、31、5 mm。同时因坝体混凝土材料的非线性刚度降低的影响，0.2g和0.4g工况下的坝体横缝在地震持时结束后趋于闭合状态，而0.6g工况下坝体横缝出现了一定的残余开度。同时统计了3种地震工况下包络图显示横缝整体张开极值分别为14、27、38 mm。

根据坝体损伤结果，在3种地震峰值加速度工况下，坝体沿横缝布置部位出现损伤区域并逐渐增大，符合随着地震峰值的不断增大，坝体出现损伤的程度也在不断加剧的一般规律。特别是在0.6g工况下，坝体结构出现了贯穿性的裂缝，且主要集中在拱坝中部附近，说明在超强地震荷载作用下，拱坝中部为坝体结构的薄弱区域，应在以后的相似坝型结构对应区域采取加筋加固措施或对应区域采用改性混凝土材料进行补强，研究结论对坝型的优化设计及薄弱部位的判定与加固具有指导意义。