不同材料非线性组合对重力坝抗震能力的影响

刘 智，赵兰浩，刘勋楠，胡国平，周清勇

（1.江西省水利科学研究院，南昌 330029；2.河海大学 水利水电学院，南京 210098）

重力坝工程库水-坝体-地基作为一个整体系统在地震过程中的存在强烈的相互作用，各类非线性行为的耦合机制异常复杂，给我国西南强震地区的重力坝抗震设计提出了严峻挑战。因此，围绕无限地基辐射阻尼效应、库水-坝体动力相互作用、坝体混凝土劣化破坏以及抗震钢筋-混凝土耦合作用等非线性行为的动力响应分析是近几十年来我国重力坝抗震安全研究的主要内容，学者们开展了大量的研究工作并取得了较大的研究进展[1]。

当前重力坝结构抗震性能研究中，针对坝体混凝土损伤开裂的材料非线性和局部配筋后钢筋混凝土相互作用的研究成果已经十分丰富[2-9]，但仍存在一些不足。循环荷载下混凝土材料表现出明显的滞回现象，传统损伤模型通常将软化段滞回效应简化为线性表达，这种通过线性形式描述非线性现象的假定无法完整反映材料损伤积累过程，滞后的损伤积累无法及时反馈劣化部位的应力转移，势必将影响后续的仿真结果。地震荷载下的钢筋混凝土相互作用同样存在滞回现象，循环往复的荷载特性使得钢筋对混凝土的刚度贡献并不是一成不变的，传统抗震分析中采用的分布式模型和嵌入式模型假设钢筋与混凝土之间不产生相对滑移，因此无法精准描述两者之间黏结滑移性能的变化过程，从而无法反映黏结性能变化对整体结构或局部构件承载能力和裂缝发展的影响。

相较而言，坝基的劣化破坏过程更是重力坝抗震研究中的薄弱环节，传统抗震分析中通常假定地基为线弹性材料或弹塑性材料，无法准确描述坝基岩体开裂后的软化现象，导致仿真计算中往往会在坝踵部位出现应力集中现象，使得计算结果与实际监测资料无法匹配[10]。规范中建议通过有限元分析求得坝底应力结果，随后推求相应的截面内力，再按平截面假定原理计算坝踵“等效应力”作为强度安全评价指标，但这种做法并不适用于坝体水平截面不垂直于中和轴的坝工结构[11]。

汶川地震之后，学者逐渐关注坝基开裂后的破坏演化。郭胜山等[12]采用D-P模型与损伤模型对Koyna重力坝地震过程中地基系统破坏模式进行对比分析，研究结果认为，地基岩体与混凝土同属非均质准脆性材料，材料中的微裂缝在荷载作用下不断演化最终导致宏观破坏，相较于传统弹塑性模型，损伤模型用以模拟地基岩体材料破坏更为切合实际震情，损伤由坝踵处向地基深处延伸，坝基交界面并未出现破坏现象；李峰等[13]选用Lee-Fenvas模型中的屈服准则，基于自由能函数建立混凝土损伤本构模型，并将模型应用于鲁迪拉重力坝坝基损伤分析，研究结果认为考虑地基为损伤材料时，仿真结果更接近实际震灾。

本文将基于库水-坝体动力相互作用、坝体-地基材料非线性以及钢筋-混凝土耦合作用等非线性行为，组合设立三类仿真分析条件，建立循环荷载下考虑滞回效应的混凝土损伤模型和考虑动态黏结滑移性能变化过程的钢筋混凝土分离式模型，以金安桥重力坝某非溢流坝段为工程依托，采用超载法将输入的地震动荷载按比例放大，分析不同材料非线性组合对重力坝动力响应以及极限抗震能力的影响，以期为重力坝抗震设计提供理论依据和技术支撑。

1 坝体-地基损伤模型

根据已有研究成果，本文中坝体混凝土与地基材料均采用四参数损伤模型进行研究。

1.1 四参数损伤模型

根据韦未等[14]提出的基于应变空间的四参数破坏准则，以及在此基础上的四参数等效应变计算方法[15]，将复杂的多轴问题转化为简单的单轴问题。损伤模型骨架线选取过镇海提出的应力应变全曲线[16]，该曲线已得到国内外科研工作者的认可。

混凝土单轴受拉应力应变曲线表达式为

式中:σ和ε分别为混凝土的应力与应变；ft，r为混凝土单轴抗拉强度；εt，r为ft，r对应的应变；dt为单轴受拉损伤变量；Ec为混凝土弹性模量；at为混凝土受拉应力应变曲线软化段参数。

混凝土单轴受压应力应变曲线表达式为

式中:fc，r为混凝土单轴抗压强度；εc，r为fc，r对应的应变；dc为单轴受压损伤变量；ac为混凝土受压应力应变曲线软化段参数。

当结构处于多轴状态时，需用等效应力应变σ*和ε*替代真实应力应变σ和ε，通过应变张量第一不变量I′1判断单元拉压状态。

1.2 软化段卸载残余应变量值

残余应变是混凝土达到软化段后客观存在不可恢复的塑性变形。本文选用Vecchio公式[17]

式中:εp为塑性残余应变；εr为抗拉或抗压峰值强度对应应变；εun为卸载点处应变。

1.3 循环荷载下的滞回效应

采用Otter等[18]通过试验观测数据建立的数学经验公式推求软化段滞回效应的加卸载特征点，公式对素混凝土及钢筋混凝土均有良好的适用性，且不依赖于应力应变骨架线的形状。

1.3.1 完全加卸载循环

滞回规则中完全加卸载循环，如图1所示。滞回循环由卸载路径AC和重新加载路径CB组成，滞回期间的损伤积累程度在模型中体现为骨架线上应变由卸载点A发展至重新加载点B。

图1 损伤模型完全加卸载示意图Fig.1 Diagram of complete unloading and reloading of damage model

根据卸载点A处卸载应变εun计算重新加载点B处应变εre

式中，kr为重新加载系数，一般取建议值0.1。

滞回规则中卸载曲线的曲率反映了刚度变化过程，割线模量由大到小不断改变。重新加载曲线可简化为线性，期间割线模量保持恒定，期间损伤值不发生变化。滞回规则完全卸载曲线采用José等[19]通过试验拟合的经验表达式，式中计入损伤变量，能够反映卸载过程中损伤积累与刚度变化过程

重新加载阶段割线模量将保持恒定，残余应变点处的损伤值与重新加载点处的损伤值相同，因此卸载产生的损伤积累变化量为dre-dun，卸载时损伤计算公式为

当卸载至残余应变点时，d=dre。

1.3.2 局部重新加载循环

滞回规则中局部重新加载循环，如图2所示，滞回循环由卸载路径AD和重新加载路径DE组成。此时荷载未完全卸载至残余应变点C，相应重新加载点E与完全加卸载循环中的重新加载点B不同，其应变值将是点A与点B之间的插值

图2 损伤模型局部重新加载示意图Fig.2 Diagram of local reloading of damage model

式中:εrx为重新加载点E处应变值；σu为局部卸载最低点D处应力值；npu为插值参数，经试验拟合与敏感性分析，取建议值8。

卸载后的局部重新加载公式与完全卸载时相同

式中:εu为局部卸载最低点D处应变值；σrx为重新加载点E处应力值。

1.3.3 局部卸载循环

滞回规则中局部卸载循环，如图3所示，滞回循环由卸载路径FG和重新加载路径GH组成。图中局部加载最高点F对应的骨架线卸载点应变值将是A点与B点之间的插值

图3 损伤模型局部卸载示意图Fig.3 Diagram of local unloading of damage model

式中:εux为F点对应骨架线卸载点应变值；σx为F点应力；npr为插值参数，取建议值8。

局部卸载循环中卸载路径与完全加卸载循环类似

2 基于单弹簧联结单元法的钢筋混凝土分离式模型

建立基于单弹簧联结单元法的钢筋混凝土分离式模型模拟钢筋与混凝土之间的相互作用。模型通过设置切向弹簧单元模拟两者之间的切向作用，约束法向强制混凝土单元与钢筋单元间的法向自由度协同，其中切向刚度则通过黏结-滑移本构关系确定。如图4所示，单弹簧联结单元法采用实体单元模拟混凝土，使用杆单元或梁单元模拟钢筋，在混凝土与钢筋接触部位设置重复结点分别作为混凝土与钢筋插值点，并建立弹簧单元作为两者间的联结单元[20]。

图4 单弹簧联结单元法求解示意图Fig.4 Diagram of single-spring joint element method

在单调轴向加载作用下，钢筋混凝土黏结-滑移本构关系骨架线与其相互作用机理关系密切，可分为胶结阶段、强度上升阶段、强度峰值阶段、强度退化阶段和摩擦稳定阶段。依据该理论建立的徐有邻五段线模型被我国GB 50010—2010《混凝土结构设计规范》[21]所采用，本文模型中骨架线选用规范模型。

循环荷载下钢筋混凝土动态黏结-滑移本构关系的建立依托于动态试验中测得的滞回加卸载路径，即卸载曲线，重新加载曲线及反向加卸载曲线。本文将循环荷载下钢筋混凝土黏结-滑移滞回关系简化为多线性形式进行表达，如图5所示，其中卸载特征参数选用滕智明等[22]试验拟合参数，卸载刚度取295 N/mm3，卸载摩擦段的黏结强度按试验经验公式计算，τup=-afτun，其中黏结摩擦因数af取0.12。

图5 动态黏结-滑移本构关系示意图Fig.5 Diagram of dynamic bond slip constitutive relationship

3 金安桥重力坝非线性动力响应分析

3.1 工程概况与仿真参数

金安桥水电站位于云南省丽江市境内，工程枢纽主要挡水建筑物为碾压混凝土重力坝，坝顶长度640 m，最大坝高160 m，坝底最大宽度150 m，坝顶高程1 424 m，水库正常蓄水位高程1 418 m。某非溢流坝段断面二维模型示意图如图6所示，地基范围各方向均选取2倍坝高以便考虑坝基损伤，坝体断面体型突变区与近场地基等潜在开裂区域网格做局部细化处理，建立薄层单元模拟地基建基面。设立概化的防渗帷幕，帷幕中心线距坝踵处25 m，帷幕深度取1/2坝前水位高度。

图6 金安桥重力坝模型示意图Fig.6 Diagram of Jin’anqiao gravity dam

选用黏弹性人工边界作为地基边界条件考虑无限地基辐射阻尼效应，采用时程法在静力分析的基础上进行动力分析，通过广义Newmark法确定每一时刻坝体与地基的应力分布及变形情况。建立水体单元，通过流固耦合算法考虑地震荷载下库水-坝体的动力相互作用。设置上游面坝顶点A与坝踵处点B为特征点以分析仿真计算结果。

金安桥坝址区域属于强震区，场地烈度为Ⅷ度，场地100年基准期超越概率2%与1%的基岩水平峰值加速度为0.399 5g与0.475g。地震荷载采用规范标准反应谱生成的人工地震波，水平向峰值加速度基准值为0.475g，竖直向峰值加速度基准值为0.317g。

表1 算例材料参数Tab.1 Material parameters of an example

材料动态参数在静态情况下根据GB 51247—2018《水工建筑物抗震设计标准》[23]将弹性模量提升50%，抗压强度提升10%。

图7 标准设计反应谱Fig.7 Standard design response spectrum

3.2 仿真条件与仿真流程

使用超载法将输入的地震动荷载按比例放大，分析不同仿真条件下金安桥重力坝的极限抗震能力，以坝体出现损伤贯通作为结构失效破坏依据，对其最终破坏模式做出相应评价。

基于不同材料非线性组合，将仿真条件设为三类。第Ⅰ类仅考虑坝体材料非线性，地基选用线弹性模型；第Ⅱ类在第Ⅰ类仿真条件的基础上，考虑坝体局部配筋后钢筋与混凝土的相互作用；第Ⅲ类在第Ⅱ类的基础上，计入地基材料非线性。具体计算条件与计算流程如图9所示。

图9 计算流程示意图Fig.9 Diagram of calculation flow

第Ⅱ类与第Ⅲ类仿真分析时需对金安桥重力坝地震荷载作用下易损区域进行局部配筋，由于本文研究重点在于分析不同仿真条件下重力坝抗震能力差异，因此配筋量及钢筋布置方式未按照工程实际布置情况并作详细分析，仅遵循钢筋混凝土设计中的基本原则。局部配筋简图如图10所示，上下游面各布置2排直径32 mm的钢筋，其保护层厚度为200 mm，钢筋间距为300 mm，钢筋横河向间距为500 mm。钢筋材料计算参数采用抗震钢筋HRB400E的力学性能特征值，即屈服强度 σs=400 MPa，弹性模量E=200 GPa，抗拉强度ft=540 MPa。

图8 人工地震波速度与位移时程曲线Fig.8 Time history curve of velocity and displacement of artificial seismic wave

图10 金安桥重力坝局部配筋示意图Fig.10 Model reinforcement diagram of Jin’anqiao gravity dam

3.3 仿真结果分析

本文仿真计算时采用黏弹性人工边界，未对地基边界进行固定约束，因此时程图中的位移数值均为坝顶特征点A处与坝踵特征点B处的差值。同时结构失效将致使脱离部分位移出现发散现象，因此图中仅截取位移发散前的时程曲线进行对比分析。

图11～图13给出了不同仿真条件下大坝特征点A处的位移时程曲线，图14则给出了结构失效破坏后的损伤分布情况。如图所示，在所给定的相同标准谱人工地震波荷载作用下，当仅考虑坝体材料非线性时，大坝极限承载能力约为基准值地震荷载作用的1.54倍；当坝体局部配筋后，大坝极限承载能力得到小幅提升，约为基准值地震荷载作用的1.66倍，此时坝体与坝基处的损伤破坏面积增大，钢筋在上游坝坡折角处、坝头部位以及下游坝坡面均表现出较为明显的限裂作用；当计入地基材料非线性后，大坝极限承载能力进一步提升，约为基准值地震荷载作用的1.68倍，坝踵处损伤演化将向坝基底部发展，与传统计算中沿坝基向下游面发展并致使帷幕贯穿有明显差异，同时坝体损伤分布面积相比前两类仿真条件显著降低。

图11 仿真条件Ⅰ时特征点A位移时程图Fig.11 Time history diagram of characteristicpoint A displacement under simulation conditionⅠ

图12 仿真条件Ⅱ时特征点A位移时程图Fig.12 Time history diagram of characteristic point A displacement under simulation conditionⅡ

图13 仿真条件Ⅲ时特征点A位移时程图Fig.13 Time history diagram of characteristic point A displacement under simulation conditionⅢ

图14 不同仿真条件下大坝失效破坏模式示意图Fig.14 Diagram of dam failure modes under different simulation conditions

图15与图16对不同仿真下的位移时程曲线进行对比。当大坝处于基准值地震荷载作用时，三种仿真条件下的位移变化规律相似，其中仿真条件Ⅱ局部配筋后结构整体质量增大，位移响应较仿真条件Ⅰ有所增加；仿真条件Ⅲ地基进入损伤后，大坝-地基系统整体刚度降低，水平向位移进一步放大，而竖直向位移受此影响较小。

图15 基准值地震荷载作用下不同仿真条件特征点A位移时程对比图Fig.15 Time history comparison of characteristic point A displacement by different simulation conditions under the action of reference value earthquake load

图16 大坝失效破坏荷载作用下不同仿真条件特征点A位移时程对比图Fig.16 Time history comparison of characteristic point A displacement by different simulation conditions under failure load

当大坝处于失效破坏荷载作用时，考虑钢筋混凝土相互作用对位移响应的影响明显。仿真条件Ⅱ结构出现损伤贯通后，在10 s左右才出现明显的位移发散现象，仿真条件Ⅰ和仿真条件Ⅲ分别在6 s和8.5 s左右即出现位移发散。而计入地基损伤将使结构更快进入失效破坏状态，约在4.5 s左右结构失效破坏，出现位移发散趋势，而其他两种仿真条件的结构失效破坏均出现在5 s后。

4 结 论

在所给定的相同标准谱人工地震波荷载作用下，考虑不同材料非线性组合将对重力坝抗震能力以及失效破坏模式造成一定影响。本文设立三类仿真条件，结合循环荷载作用下考虑滞回效应的混凝土损伤模型和考虑动态黏结滑移性能变化过程的钢筋混凝土分离式模型，以金安桥重力坝为工程依托，研究不同仿真条件对重力坝抗震能力及破坏模式的影响。研究结果表明，当考虑坝体局部配筋后，重力坝极限承载能力有所提升，钢筋限裂作用显著，同时钢筋混凝土的相互作用将有效抑制结构失效破坏后的灾变发展速度；在计入地基材料损伤演化后，大坝极限承载能力进一步提升，坝体损伤演化路径与损伤分布面积有明显改善；坝体-地基系统整体质量和刚度的变化对大坝位移响应幅值的影响十分明显，其中质量与位移响应正相关，刚度与位移响应负相关。

鉴于影响重力坝抗震能力的非线性行为众多，且耦合作用机制复杂，本文比选的三类仿真条件中仍有部分影响因素未被顾及，特别是计入地基材料损伤后未考虑坝基渗流场对损伤演化路径及损伤分布面积的影响。后续将重点研究地基材料损伤效应与应力场和渗流场的耦合影响机制，进一步探索重力坝的真实抗震能力。