混凝土面板堆石坝地基防渗墙塑性损伤数值分析

温立峰，李炎隆，柴军瑞

（西安理工大学省部共建西北旱区生态水利国家重点实验室，陕西西安 710048）

1 研究背景

随着能源需求的攀升和筑坝技术的进步，大量的土石坝已经、正在或计划修建在具有高压缩性和强透水性的覆盖层地基上。防渗墙是覆盖层地基最有效的防渗措施之一，其可以有效截断地基渗流通道，并具有良好的耐久性［1-3］。经过50余年的发展，我国在防渗墙的设计、施工及复杂地质条件下的建设技术等方面均取得了长足的发展［3-5］。在坝体填筑和蓄水引起的不同外荷载作用下，防渗墙结构通常呈现复杂的力学特性。外荷载可能引起墙体产生显著的塑性应变，进而导致墙体开裂。多个长期运行的土石坝除险加固防渗墙已经观测到墙体的开裂或破坏现象［6］。Rice 等［7］发现，即使防渗墙的裂缝开度小于1 mm，墙体的有效渗透系数也可能增加几个数量级。研究防渗墙的应力变形和损伤开裂特性对墙体的设计优化和安全评价至关重要。

防渗墙工程具有隐蔽性，受力条件及应力变形特性复杂。墙体不仅承受坝体重力及蓄水引起的水压力，同时与相邻土体存在高度非线性的相互作用。Dascal［8］基于实测资料分析了Manic-3 防渗墙变形和应力分布规律，并讨论了防渗墙的受力特点。Brown 等［9］介绍了Arminou 心墙坝塑性混凝土防渗墙的建设信息，分析了影响墙体力学特性及其失效的可能原因。Song 等［10］基于防渗墙力学特性的现场观察资料，对施工过程中接头管的拔管时间和拔管力进行了研究。Rice 等［6］对除险加固防渗墙的力学性状进行统计分析，讨论了防渗墙的失效现象和机制。温立峰等［11］对新建坝基防渗墙的力学性状开展统计分析，总结了防渗墙应力变形规律及其影响因素。此外，相关学者对防渗墙开展了离心模型或振动台试验研究［12-14］。上述实测、统计及试验结果为揭示防渗墙应力变形规律提供了有益认识，但是难以全面描述防渗墙应力变形和损伤开裂特性的演化规律及其形成原因。

部分学者采用数值分析方法进一步揭示和研究防渗墙的力学性状。丁艳辉等［15］、徐建国等［16］、Li等［17］采用有限元方法，分析了不同材料防渗墙的应力变形响应规律。周小溪等［18］和万宇豪等［19］探讨了坝基防渗墙地震反应规律及其抗震安全性。基于数值模拟手段，大量学者研究了墙体材料、河谷地形和爆炸振动等影响因素对墙体力学性状的影响［20-22］。沈振中等［23］、温续余等［24］、介玉新等［25］根据应力变形数值计算结果，研究了防渗墙的受力特点及其与整个防渗体系的连接形式。Yu等［26-28］对沥青心墙坝地基防渗墙开展了线弹性和塑性损伤计算，分析了三维河谷效应对墙体应力变形及损伤开裂特性的影响机制。数值计算结果通常没有对应的实测结果进行验证。此外，目前对心墙坝地基防渗墙力学特性的研究相对较多，而对面板堆石坝防渗墙力学特性，特别是损伤开裂特性的研究不足，有必要结合实测资料和数值分析方法系统揭示混凝土面板堆石坝地基防渗墙的应力变形及损伤开裂特性。

本文结合实测资料和数值分析，研究面板堆石坝防渗墙的力学性状。基于混凝土塑性损伤模型建立防渗墙力学特性三维数值模型。在数值模型中考虑防渗墙与相邻土体的接触效应及坝体和地基的渗流-应力耦合效应。采用实测结果验证数值计算结果的准确性，进而结合实测和数值结果分析防渗墙的受力机制、应力变形及损伤和开裂特性，讨论墙体位置、材料和渗流-应力耦合效应对防渗墙力学性状的影响。

2 工程实例分析

2.1 工程概况本文以苗家坝水电站面板堆石坝为依托开展研究。苗家坝水电站位于甘肃省文县境内的白龙江干流上。枢纽由面板堆石坝、下游地下厂房以及引水隧洞等组成。大坝高度111 m，坝顶长度348.2 m，水库正常蓄水位800 m。大坝上游坝坡1∶1.4，下游综合坡比1∶1.55。大坝修建在覆盖层地基上，主要包括4 个建设阶段，由2009年10月开始施工，2011年7月建设完成。图1为坝体典型剖面图。坝体填筑材料、碾压方法和填筑过程信息也显示在图1（a）中，主堆石和下游堆石的施工参数一致。

图1 苗家坝面板堆石坝坝体剖面以及防渗墙布置

2.2 地基条件大坝所在的V型河谷中分布有44～50 m厚的深厚覆盖层。为了获得覆盖层的物理力学特性开展了大量现场和室内试验。根据材料特性的差异，覆盖层主要分为3层，上层为卵石和砂砾石层（Q4a3），厚度6～20 m；中部为砂砾石层（Q4a2），厚度12～15 m；底部为块石和砾石层（Q4a1），厚度6～20 m。中部砂砾石层为覆盖层的主要组成部分。覆盖层以下基岩为变质凝灰岩，基岩中未观测到大尺寸断层和强风化区。图2为基于现场和室内试验获得的覆盖层物理力学特性。孔隙率随着覆盖层深度的增加逐渐减小，而覆盖层上部干密度波动较大，整体相对较小，底部干密度随深度呈增加趋势，整体相对较大。覆盖层孔隙率为0.2左右，变形模型在35～55 MPa之间，具有较大的剪切模量和变形模量。水平方向渗透系数沿着深度的分布规律与垂直方向基本一致，但是渗透系数总体相对较大。

图2 覆盖层地基物理力学特性沿深度分布

图3为覆盖层的粒径分布曲线。颗粒粒径范围总体较大，最大粒径大约为600 mm 左右。Q4a3、Q4a2、Q4a13 层覆盖层的粒径分布下包线的不均匀系数（Cu）分别为283、200、237，曲率系数（Cc）分别为13.6、4.3和10.6。图3表明，3层覆盖层的粒径分布均不均匀。

图3 覆盖层粒径分布

上述结果表明，苗家坝面板堆石坝坝基覆盖层具有结构松散、粒径分布范围大、高渗透性和物理力学特性不均匀等特点。为了减小覆盖层的压缩变形，筑坝前使用质量为25 000 kg、宽度为2.2 m的振动碾对覆盖层进行了10遍碾压处理。

2.3 防渗墙和监测系统布置为了控制覆盖层地基渗流，在坝基中修建了一道槽孔型防渗墙。防渗墙厚度1.2 m，顶部长度121 m，总面积达2900 m3。防渗墙最大深度为50.5 m，贯入基岩深度0.5 m。图1（b）为沿坝轴线方向防渗墙的纵向断面和详细设计信息。墙体采用C30 混凝土材料，轴向抗压强度为20.1 MPa，抗拉强度为2.01 MPa。该防渗墙采用拔管法进行板对板的分期施工。防渗墙在坝体填筑的前3 个月进行施工。为了观测防渗墙和相邻覆盖层的力学特性，地基和防渗墙中布置了详细的监测系统如图1（b）和表1所示。

表1 地基和防渗墙中监测系统布置

3 防渗墙力学性状实测结果分析

3.1 水平位移和顶部沉降图4（a）为防渗墙顶部测点IN1-1、IN2-1和IN3-1（测点位置见图1（b））水平位移随时间的变化过程，若干防渗墙实例顶部最大水平位移结果也表示在图中。为了使防渗墙实例数据在时间上具有较好的可比性，将不同实例水库蓄水开始的时间对齐，作为参考时间。图4（a）结果表明，蓄水开始前苗家坝防渗墙向上游变形，位移随着时间的增加而增加，但是变形速度逐渐减小。竣工期（EOC）防渗墙向上游总变形为6.8 cm，其中大约70%发生在坝体填筑的前3个月之内。蓄水开始后，在水压力作用下，苗家坝防渗墙逐渐向下游变形。运行期防渗墙向下游总变形为9.9 cm，其中大约90%发生在蓄水阶段，说明水压力是防渗墙向下游变形的主要荷载。图4（b）为蓄水完成时苗家坝及其他若干实例防渗墙顶部轴线实测沉降分布。图4（b）表明，苗家坝防渗墙顶部沉降相对较小，其中中间部位的沉降相对两岸的部位较大，最大沉降值为2.0 cm。

由图4（a）可知，所有面板堆石坝防渗墙的变形规律与苗家坝防渗墙基本一致，此时防渗墙位于上游坝基（简称上游防渗墙）。上游防渗墙在竣工期和蓄水期（EOF）均产生较大水平位移，但是沉降相对较小。心墙坝的防渗墙位于中部坝基（简称中部防渗墙），竣工期呈现相对上游防渗墙较小的水平位移，而产生较大的沉降变形，蓄水期同时产生较大水平位移和沉降变形。上述结果差异由防渗墙的不同受力特点引起。竣工期中部防渗墙相邻土体水平位移较小，墙体上下游面承受基本对称的侧土压力，而上游防渗墙承受来自下游侧土体水平位移引起的较大侧土压力，因此产生较大水平位移。蓄水期在上下游侧水压力差作用下两种位置防渗墙均产生较大水平位移。中部防渗墙承受来自上部坝体的垂直土压力和向下摩阻力，因此沉降变形较大。上述结果表明，防渗墙变形性状与其位置直接相关。

图4 苗家坝及若干防渗墙实例变形分布

3.2 垂直应力图5为苗家坝防渗墙测点SM1-5、SM2-5、SM3-3实测垂直应力及若干防渗墙实例最大垂直应力随时间的变化过程。从图5可见，施工期苗家坝防渗墙上游面（UF）处于受压状态，下游面（DF）垂直应力随着时间逐渐减小，并在测线SM1 和SM3 底部观测到一定的拉应力。竣工期在测点SM1-5 和SM3-3 观测到的最大拉应力分别为1.6 和1.0 MPa。上述拉应力主要由防渗墙受约束的弯曲变形以及防渗墙和相邻土体不一致的轴向变形引起的墙体拖曳效应引起。蓄水开始后防渗墙下游面垂直应力开始增加并于蓄水后逐渐趋于稳定。蓄水阶段，向下游的弯曲变形在上游面引起拉应力。拉应力主要分布在上游面的底部和靠近两岸部位，在测点SM1-5、SM2-5、SM3-3观测到的最大拉应力分别为1.2、0.1和1.5 MPa。

图5 防渗墙典型测点及若干典型防渗墙实例实测垂直应力变化过程

察汗乌苏防渗墙实测垂直应力的演化过程与苗家坝防渗墙结果基本一致，但是中部防渗墙呈现明显不同的演化过程。中部防渗墙垂直应力施工期和蓄水期均随时间逐渐增加。上下游面应力分布和演化规律基本一致。中部防渗墙主要处于受压状态，不同位置防渗墙应力分布的差异由上述分析的不同受力特点引起。

4 数值模型

4.1 混凝土防渗墙塑性损伤模型苗家坝防渗墙实测结果表明墙体一定区域内承受较大拉应力，可能引起损伤或开裂。本文采用数值方法进一步揭示防渗墙的力学特性。数值计算中通常采用线弹性模型模拟混凝土的应力变形关系，但是只有在小荷载作用下混凝土才展现线弹性行为，随着拉伸应变的增加，混凝土将产生损伤和开裂［26］。因此，分析高土石坝混凝土防渗墙的力学性状时有必要考虑材料的损伤行为。本文采用Lee 等［29］提出的塑性损伤模型描述混凝土防渗墙的应力变形关系，该模型可以考虑独立的压缩和拉伸损伤模式以及刚度恢复，已经应用于面板堆石坝面板的有限元损伤分析中［30］。该模型中将应力应变关系表示为：

式中：σ为混凝土的应力；d为一个表示弹性刚度减小的标量退化损伤变量；σˉ为有效应力；Ei为未损伤的弹性刚度；为ε总应变；εp为塑性应变。

模型引入损伤变量kn描述混凝土的损伤状态，退化的损伤变量表示为：

式中：n为不同的状态，n为T表示拉伸状态，n为C表示压缩状态；gn为混凝土的耗散能量密度，可以通过断裂能Gn与单个网格单元的特征长度ln的比值确定；σn为单轴应力；dT、dC分别为拉伸和压缩单轴损伤变量；kT、kC分别为拉伸和压缩为损伤变量；s为权重因子。

采用上述模型模拟已有的混凝土拉伸和压缩加载试验结果如图6所示。从图6可以看出，数值结果与试验结果吻合较好，说明上述塑性损伤模型可以较为合理地模拟混凝土的力学响应。根据相似混凝土材料的已有研究［26］，本文数值计算中，混凝土材料的拉伸强度和断裂能分别取为2.01 MPa 和325 N/m。根据防渗墙单元尺寸，网格单元的特征长度ln取为0.24 m。苗家坝防渗墙实测结果表明防渗墙承受的压应力总体相对较小，因此数值计算中只考虑混凝土的拉伸损伤。为了简化，面板和趾板、连接板均采用线弹性模型模拟，其密度、弹性模量、泊松比及渗透系数分别取为2.45 g/cm3、28 GPa、0.167和1×10-12m/s。

图6 混凝土单轴和双轴荷载拉伸和压缩试验数值模拟与试验结果对比

为了对比不同混凝土模型的结果，同时采用线弹性模型计算防渗墙的力学性状，其密度、弹性模量、泊松比分别取为2.45 g/cm3、26 GPa、0.167，其中弹性模量由单轴压缩试验确定。同时为了比较常规混凝土（OC）与塑性混凝土（PC）防渗墙力学性状的差异，也采用线弹性模型计算塑性混凝土防渗墙力学性状。由于加入了膨润土，塑性混凝土具有相对普通混凝土更强的适应变形的能力。塑性混凝土弹性模量通常在0.5～2.5 GPa之间［2］，根据已有研究［21］，本文塑性混凝土防渗墙的密度、弹性模量和泊松比分别取为2.20 g/cm3、1.5 GPa和0.25。

4.2 堆石和覆盖层材料的弹塑性模型苗家坝堆石和覆盖层材料三轴试验结果表明：（1）试件压缩性随着法向应力的增加而减小；（2）偏应力和轴向应变曲线随围压的增加而增加；（3）随着围压的增加，体积应变从剪胀变为剪缩。本文采用双屈服面弹塑性模型［31］描述堆石和覆盖层材料体积应变行为，该模型采用两个屈服面来判断是否发生塑性应变：

式中：q为主剪应力；p为平均应力；f1、f2为两个屈服面，f1>（f1）max或f2>（f2）max表明材料产生塑性应变，（f1）max和（f2）max分别是两个屈服面函数的历史最大值。该模型的详细信息可以参考文献［31］。基于三轴试验获得的模型计算参数如表2所示。

表2 弹塑性模型计算参数和初始渗透系数

4.3 接触面模型为了真实模拟防渗墙工作状态，采用基于接触力学的无厚度摩擦接触方法［27］模拟防渗墙与相邻土体间的接触效应。该方法基于增强的拉格朗日方法求解接触摩擦问题，通过对罚函数修正项进行反复迭代确定精确的拉格朗日乘子。由于结构与土体之间的滑移不受网格不连续的限制，该方法可以获得土体-结构接触面的接触和不连续性状。法向接触压力P采用下式计算：

式中：ε0为侵入容差；λi为拉格朗日乘子（迭代步i）；Kn为法向接触刚度；μn为两个接触面之间的距离。

接触剪切特性和本构关系采用库仑摩擦模型描述：

式中：τ为等效剪应力；f为摩擦系数；σ为接触压应力；c为黏聚力；τlim为极限抗剪强度。

为了验证接触摩擦方法的准确性，采用该方法模拟了文献［32-33］中紫坪铺和水布垭面板堆石坝面板与垫层的直剪试验。根据上述两个工程直剪试验结果，并类比相关工程的接触摩擦计算参数，紫坪铺直剪试验模拟中摩擦系数取为0.3，接触容差为0.5 mm，水布垭直剪试验模拟中摩擦系数取为0.3，接触容差为1.0 mm。图7为剪切位移-剪切应力关系曲线数值与试验结果的对比。图7表明，试验与数值结果吻合良好，说明接触摩擦方法可以获取土体和结构接触面的主要力学性状。由于缺少试验资料，本文根据经验和已有相关研究确定接触模型的计算参数。数值计算中通常建议覆盖层和防渗墙之间的摩擦系数取为0.2～0.4［21］，因此本文接触面模型的摩擦系数取为0.2，接触容差设置为0.5 mm。

图7 面板和垫层剪切试验数值模拟与试验结果对比

4.4 渗流-应力耦合分析方法土石坝数值计算中通常将水压力直接施加在防渗体上，而不考虑地基和坝体中的渗流作用。土石坝的力学性状依赖于材料的力学和水力特性。渗流-应力耦合过程对坝体的力学性状具有重要影响［33］，特别是修建在覆盖层地基上的土石坝，此时地基和部分坝体均处于饱和渗流状态。苗家坝覆盖层地基为具有强渗透性的典型粗粒土材料，蓄水期覆盖层和部分坝体处于饱和渗流状态。为了考虑渗流-应力耦合效应的影响，采用Chen 等［33］提出的水力耦合分析方法模拟坝体和地基渗流-应力耦合效应。渗流-应力耦合过程受连续介质力学的动量和质量守恒定律控制，其控制方程为：

式中：φ为总水头；v为流速；ρw为水的密度；Sw为水的存储量；u为位移矢量；pw为孔隙水压力；D为切向弹性模量；α为固结系数；fb为体积力矢量；H（φ-z）为Heaviside 罚函数；εv为体积应变；δ为Kronecker三角矢量；t为时间。

渗透系数随变形的变化通过修正的Kozeny-Carman 方程来表征，其中孔隙率演化与体积应变有关：

式中：ni、n0分别为当前和初始孔隙率；k、k0分别为当前和初始渗透系数；β1为相关修正系数。

试验获得的不同材料初始渗透系数如表2所示。本文采用水力耦合分析方法模拟地基和坝体渗流-应力耦合效应及其对防渗墙力学性状的影响，而未考虑墙体本身的水力耦合效应对防渗墙力学性状的影响。因此本文中防渗墙渗透系数的取值主要影响渗流控制效果，而不影响墙体的力学性状。普通混凝土和塑性混凝土均具有良好的防渗性能，只要两种材料防渗墙能取得良好的防渗效果，便不会引起明显的地基和坝体渗流场分布差异，因此也不会通过影响坝体和地基的力学特性影响防渗墙的力学性状。因此，本文常规混凝土和塑性混凝土两种防渗墙的初始渗透系数均取为1×10-9m/s。

采用交叉迭代算法解决渗流-应力的耦合问题。详细的迭代方法在文献［33］进行了介绍。在每个时间中，首先采用变分不等式方法求解渗流过程，然后采用中点增量法求解非线性变形过程。每个时间步后根据变形计算结果更新材料的渗透系数。在某时间步同时满足变形和渗流的收敛标准时，终止迭代完成计算。本文在上述计算程序基础上引入混凝土塑性损伤模型实现考虑渗流-应力耦合作用的防渗墙塑性损伤分析。

4.5 有限元模型数值计算三维有限元网格如图8所示。采用空间八节点等参单元模拟材料的力学特性，模型包含44 870个单元和49 936个节点，其中6425个防渗墙单元和6857个防渗墙节点。为了更加准确描述防渗墙的力学特性，沿防渗墙厚度方向划分5排单元。模型真实考虑坝体不同分区和覆盖层地基地质条件。

图8 三维有限元网格

模型模拟坝体的实际填筑过程和相应的蓄水过程，如图1（a）所示。根据施工过程设置坝体填筑的时间步长。采用26 个时间步模拟坝体的填筑过程，3 个时间步模拟面板施工过程，坝体填筑过程模拟厚度小于5 m。之后，根据水库水位的上升速度，模拟的时间步长设置为10 d。模型的底部和两侧设置为法向约束并且为不透水边界。模型上下游侧设置为不透水边界，大坝上下游水位分别为800和710 m，淹没在水位以下的上下游河谷及坝面设置为总水头边界，其他边界设置为潜在出渗边界。

5 结果分析

5.1 防渗墙受力分析图9为防渗墙最深断面墙体和相邻土体沉降、防渗墙承受的侧土压力和摩阻力沿墙体深度分布规律。施工期，防渗墙下游侧产生压缩变形而上游侧产生向上的隆起变形，计算的顶部最大向上变形为5.0 cm。向上隆起变形是由顶部无约束覆盖层受到向上推力引起。竣工期防渗墙下游侧土体顶部最大计算沉降为6.0 cm，明显大于墙体顶部1.5 cm的沉降变形。蓄水阶段，在水压力作用下防渗墙上下游侧土体均产生向下的压缩变形。上下游侧覆盖层与防渗墙顶部最大沉降的比值（相对沉降）分别为3.7和2.4。计算的防渗墙顶部沉降和相对沉降与实测结果2.0 cm和4.1基本一致。统计结果表明，上游防渗墙顶部沉降和相对沉降的统计范围分别为0.02%H～0.05%H（H为防渗墙深度）和3～6［11］，本文结果在上述范围之内。与上游防渗墙不同，中部防渗墙上下游侧土体竣工期和蓄水期均产生压缩变形，而且压缩变形整体较大［11］。

防渗墙的力学行为主要取决于其所承受的水平向和垂直向荷载。水压力和侧土压力是墙体承受的主要水平荷载。水压力主要出现在蓄水以后，一般呈线性分布，与浸润线位置直接相关。坝体的摊铺和碾压改变土体应力，进而引起覆盖层地基的水平位移。覆盖层较大的水平位移使防渗墙上下游面承受较大的侧土压力差。如图9所示，防渗墙侧土压力随深度方向呈非线性分布。由于相邻土体较大的水平位移，防渗墙中部承受的侧土压力较大。竣工期计算所得上游面和下游面的最大侧土压力分别为0.73 和0.98 MPa。蓄水期在水压力综合作用下，防渗墙上下游面的侧土压力均相对减小。竣工期和蓄水期防渗墙的侧土压力与实测结果均较为接近，说明本文采用的接触摩擦方法较为合理地描述了土体和结构接触面的力学特性。由于覆盖层水平位移的特点，中部防渗墙的竣工期上下游侧土压力较小且基本对称，然而蓄水期下游面相对于上游面承受较大侧土压力［2］。两种位置防渗墙侧土压力差异是引起不同水平位移分布规律的主要原因。

图9 防渗墙与相邻土体沉降、防渗墙最大深度断面上下游侧侧土压力及摩阻力分布

防渗墙承受4类垂直荷载，包括垂直土压力、垂直水压力、墙体自重和摩阻力。由于防渗墙压缩性显著低于相邻土体，在顶部土压力和水压力作用下，覆盖层和防渗墙之间产生较大的沉降差，在墙体两侧引起较大的摩阻力。在摩阻力作用下，墙体中部的垂直应力大于顶部应力［8］。此外，防渗墙相邻土体某一深度处的实测垂直应力小于上覆土压力［9］。上述结果均表明摩阻力显著影响防渗墙的力学性状。如图9所示，受覆盖层和防渗墙相对变形规律的影响，竣工期防渗墙上游面承受向上的摩阻力而下游面承受向下的摩阻力。计算所得竣工期下游面最大摩阻力为0.47 MPa，大于上游面向上的摩阻力。蓄水作用改变防渗墙承受的摩阻力的分布规律，墙体上下游面均承受向下摩阻力，上游面和下游面最大摩阻力分别为0.1 和0.4 MPa，均发生在防渗墙顶部。防渗墙的沉降变形主要由墙体压缩变形和墙体刚体位移组成。在基岩的约束作用下，覆盖层底部的压缩变形较小，因此防渗墙底部的沉降变形可能大于相邻土体的变形。摩阻力在防渗墙和相邻土体沉降变形一致的位置减小为0。该位置通常称为墙体的中性点。该深度以下，由于防渗墙的沉降变形相对较大，墙体承受向上的摩阻力。由图9可知，竣工期防渗墙中性点的位置靠近墙体底部，位于基岩以上大约1.5 m 的位置，蓄水期该位置向上移动大约1.5 m。中部防渗墙竣工期和蓄水期中性点以上均承受向下摩阻力，而且整体摩阻力相对较大，是引起中部防渗墙产生较大沉降变形和承受较大压应力的主要原因。

由图9可知，塑性混凝土防渗墙的沉降明显大于常规混凝土防渗墙的结果。覆盖层与塑性混凝土防渗墙的相对沉降未超过1.5。较小的相对沉降意味着塑性混凝土防渗墙与相邻土体变形一致并共同承担外部荷载，因此塑性混凝土防渗墙承受的摩阻力显著减小。由于刚度差异，竣工期塑性混凝土防渗墙上游面承受的侧土压力大于常规混凝土防渗墙，但是下游面侧土压力相对较小。蓄水期塑性混凝土防渗墙承受的侧土压力相对于常规混凝土的结果恰好相反。在相邻覆盖层水平位移的作用下，塑性和常规混凝土防渗墙上下游面承受的侧土压力没有明显差异。不考虑渗流-应力耦合作用情况下，下游面侧土压力和摩阻力相对较小，而上游面受力差异较小。

5.2 防渗墙变形特性分析图10为计算的防渗墙整体水平位移分布及3 条测线竣工期和蓄水期水平位移分布。结果表明，竣工期墙体弯向上游变形，最大值为7.0 cm，蓄水期墙体弯向下游变形，最大值为10.0 cm，最大值均发生在顶部中间部位。受到基岩约束效应的影响，防渗墙底部变形较小。下游面侧土压力是引起防渗墙竣工期向上游变形的主要原因，而上游面水压力是引起蓄水期向下游变形的主要原因。三条测线计算的变形分布和数值大小与实测结果基本吻合，说明数值模型较好地描述了防渗墙的变形性状。

图10 防渗墙水平位移整体分布及三条测线水平位移分布（EOC指竣工期，EOF指蓄水期）

由图4可知，苗家坝防渗墙的变形规律与其他上游防渗墙实例变形结果基本一致。计算和实测最大水平位移均在竣工期和蓄水期0.07%H～0.20%H和0.05%H～0.35%H的统计范围［11］之内。数值和实测结果表明，苗家坝防渗墙水平位移显著大于顶部沉降。图4和Yu 等［26］数值计算均表明，中部防渗墙，墙体产生较大沉降变形，蓄水过程中进而产生较大水平位移。在不同受力特点作用下，上游防渗墙和中部防渗墙的变形模式存在明显差异。

如图10所示，塑性混凝土防渗墙的水平位移与常规混凝土防渗墙基本一致。这是因为虽然不同材料墙体的刚度存在较大差异，但是均难以抵抗相邻土体的水平位移，而产生与覆盖层水平位移相适应的变形，此时不同材料防渗墙的应力状态将存在较大差异。图10比较了不考虑地基和坝体渗流-应力耦合效应的结果，此时水荷载作为面力施加在防渗体表面。蓄水期，考虑渗流-应力耦合效应情况下防渗墙的水平位移相对不考虑时增大3.0 cm。渗流-应力耦合效应主要通过影响坝体和地基的变形，进而改变防渗墙承受的侧土压力和摩阻力（如图9（b）所示）影响墙体的力学性状。

5.3 防渗墙应力分析图11为3 条测线上游面和下游面的实测和计算垂直应力分布。防渗墙大主应力和小主应力分布如图12所示。竣工期，防渗墙下游面两岸及底部部位存在一定的拉伸区域。计算所得墙体承受的最大压应力和拉应力分别为20.0和1.9 MPa，发生在底部靠近两岸部位。上述结果均未超出材料的允许强度。水库蓄水后，防渗墙下游面逐渐转变为受压状态。蓄水期，防渗墙上游面靠近两岸部位出现拉伸区域。计算的墙体蓄水期最大压应力和拉应力分别为22.0和2.0 MPa。如图11所示，在靠近两岸的断面，防渗墙垂直压应力由顶部到底部逐渐增加。在向下的摩阻力作用下，防渗墙最深断面位置垂直压应力由顶部到中性点的位置逐渐增加，之后向底部逐渐减小。这些结果表明，防渗墙最危险的压缩断面位于靠近两岸的部位，而不是河床中间部位。如图11所示，计算的防渗墙上部垂直应力结果与实测结果较为接近，而靠近防渗墙底部的位置，计算结果相对实测结果偏大。上述差异与数值计算模型的局限性有关，比如几何模型的简化、复杂加载和地质条件的简化等。在水压力和侧土压力水平荷载作用下，面板堆石坝防渗墙主要承受弯曲效应，在竣工期和蓄水期均产生较大的拉应力。上述应力分布与中部防渗墙的应力分布明显不同。Dascal［8］和Yu 等［26］发现，中部防渗墙主要处于受压状态，在摩阻力作用下，最大压应力发生在中性点的位置，只在受基岩约束的墙体上部尖端部位出现一定拉应力。Dasca［8］发现Manic-3 中部防渗墙85%的垂直应力是由摩阻力引起的。

图11 防渗墙上下游面3条测线垂直应力分布（EOC指竣工期，EOF指蓄水期）

防渗墙不同阶段的拉应力主要由受约束的弯曲变形促使的正弯矩和地基与墙体不均匀轴向变形促使的拖曳效应引起。图13为不同阶段墙体典型位置内力的空间分布。在基岩约束作用下，竣工期防渗墙向上游的变形在下游面两岸及底部位置引起正弯矩，而蓄水期向下游的变形在上游面两岸及顶部位置引起正弯矩。正弯矩是引起上下游面拉应力的主要原因。此外，覆盖层和防渗墙不均匀变形对墙体具有拖曳效应，在四周特别是两岸部位引起轴向拉力，进一步在墙体两岸部位引起拉应力。垂直方向上防渗墙承受垂直压缩轴向力，这是墙体底部拉应力比两岸拉应力小的主要原因。

图13 防渗墙典型位置内力分布

如图11所示塑性混凝土防渗墙的垂直应力显著小于常规混凝土防渗墙。蓄水期，塑性混凝土防渗墙下游面最大垂直压应力和拉应力相对常规防渗墙分别减小8.0和1.7 MPa。摩阻力的减小（图9）是拉压应力减小的主要原因。王清友等［2］发现，塑性混凝土防渗墙的垂直应力通常只有常规混凝土防渗墙的1/10～1/15。本文结果与该结果基本吻合，说明塑性混凝土材料可以显著改善防渗墙的应力状态。此外，考虑渗流-应力耦合作用下，由于侧土压力和摩阻力增加，防渗墙的垂直应力整体相对较大。蓄水期下游面压应力和上游面拉应力分别增加4.0和0.5 MPa，说明渗流-应力耦合作用对防渗墙的应力状态存在一定影响。

图12（b）比较了采用混凝土塑性损伤模型和线弹性模型计算的防渗墙小主应力分布结果。采用线弹性模型计算的墙体拉应力明显较大，竣工期和蓄水期最大拉应力均大于3.0 MPa，超过材料的拉伸强度。线弹性模型使用弹性模量和泊松比描述防渗墙的力学性状，无法描述当拉应力或压应力超过材料的允许强度时产生损伤或开裂特性。因此线弹性模型获得的结果只能用于评估墙体潜在的开裂或失效风险。相反，采用塑性损伤模型计算的防渗墙最大拉应力未超过材料的允许强度，而且防渗墙拉伸区域明显较小。塑性损伤模型可以模拟材料刚度退化和应力重分布现象。因此在混凝土防渗墙力学特性的数值计算中，塑性损伤模型可以较为准确地描述墙体的力学特性。

图12 防渗墙大小主应力分布

5.4 防渗墙损伤和开裂分析当防渗墙承受的压应力或拉应力超过材料的允许强度时，墙体将产生开裂或者失效。计算的防渗墙最大压应力大约为混凝土抗压强度的73%左右。线弹性分析结果也表明，防渗墙压应力未超过材料抗压强度。上述结果表明压应力对苗家坝防渗墙安全的影响较小。本节主要分析可能引起防渗墙开裂的拉伸损伤结果。

弹塑性损伤模型中，采用损伤变量kT表示材料的损伤程度。当kT<0.1 时，表明材料几乎不发生损伤；当kT>0.8 时，表明材料发生严重损伤，可能产生开裂失效［29］。图14为防渗墙拉伸损伤分布结果。施工期，防渗墙下游面底部和靠近两岸的部位产生一定的拉伸损伤区。蓄水后，在上游面底部和靠近两岸部位出现拉伸损伤区域，蓄水期防渗墙的拉伸损伤程度比施工期较大。如图12（b）所示，采用两种数值模型计算的最大小主应力也出现在相应的位置。两种结果显示出良好的一致性，相互印证。大部分区域损伤因子小于0.1，说明墙体未损伤或者只是轻微的损伤。计算的防渗墙最大拉伸损伤因子均未超过0.8，说明防渗墙尚不会产生开裂失效。上述防渗墙的损伤分布规律与中部防渗墙存在明显差异。Yu 等［28］对中部防渗墙开展损伤分析发现施工期拉伸损伤区域主要位于顶部尖端部位，蓄水后拉伸损伤区逐渐扩散到下游面靠近基岩的部位。Rice 等［7］发现中部防渗墙上下游面的水压力差足以引起墙体开裂，墙体与基岩的接触部位是最可能开裂的位置。

图14 防渗墙拉伸损伤分布（由损伤变量kT表示）

表3为若干典型工程防渗墙失效或开裂案例。面板堆石坝或斜墙坝防渗墙主要失效模式为拉伸或剪切破坏，失效发生在防渗墙顶部或底部与基岩的接触部位。该开裂失效结果与本文计算的损伤分布结果基本吻合，也与Brown等［9］发现的在弯曲效应作用下上游坝基防渗墙更易于产生水力开裂的结果基本一致。心墙坝中部防渗墙的失效模式主要为压缩失效，并且发生在底部。上覆土压力和来自相邻土体的摩阻力引起的过大压应力是引起上述失效的主要原因。此外，如表3所示，虽然Arminou塑性防渗墙观测到墙体开裂，但其原因主要是材料侵蚀，而非结构应力。由于塑性混凝土防渗墙弹性模量与覆盖层地基相当，因此在塑性混凝土防渗墙中较少观察到结构性失效。

表3 混凝土防渗墙失效或开裂的若干典型实例

6 结论

本文基于实测资料和数值分析研究面板堆石坝地基防渗墙的力学特性。在验证数值计算模型的基础上，讨论了混凝土材料、坝体和地基渗流-应力耦合效应、墙体位置对力学性状的影响。本文主要获得以下几点结论：（1）苗家坝面板堆石坝混凝土防渗墙在水压力和侧土压力作用下呈现弯曲效应，竣工期在下游面底部产生拉应力，蓄水期在上游面底部及两岸部位产生拉应力。防渗墙拉伸损伤区主要发生在底部及靠近两岸部位，尚未引起墙体开裂。苗家坝防渗墙应力变形状态在合理范围之内，墙体运行良好。（2）面板堆石坝上游防渗墙与心墙坝中部防渗墙的力学性状存在明显差异。上游防渗墙主要承受弯曲效应，可能在底部或靠近两岸部位发生拉伸或剪切失效，而中部防渗墙主要承受压缩效应，可能在底部产生压缩失效。防渗墙位置是力学性状的关键影响因素。塑性混凝土可以显著改善墙体的应力状态。坝体和地基的渗流-应力耦合效应对墙体力学性状也有一定的影响。（3）防渗墙应力变形监测数据与数值模拟结果基本一致，表明本文考虑混凝土塑性损伤、坝体地基渗流-应力耦合效应、防渗体与土料接触效应的数值计算模型较为合理地模拟了坝体和防渗体的力学性状，塑性损伤模型可以较为合理地描述防渗墙的损伤力学行为。