王 敏，李 毅，尚 钦，郑华康，周佳庆

(1.长江勘测规划设计研究院，湖北 武汉 430010；2.长沙理工大学水利工程学院，湖南 长沙 410076)

0 引 言

在水库运行期，库水位随着水库的调节不断地发生改变，引起水库周围地下水变化，使得库岸岩土体的性质随着水环境的变化而发生调整，进而影响到库岸边坡的稳定性[1，2]。近几年，随着三峡水库的运行，国内外对库水作用下岸坡稳定性分析开展了广泛的研究工作[3，4]，在非稳定渗流计算过程中，地下水浸润面的确定与计算方法的选择有着密切关系，在缺乏实测资料的条件下，准确地模拟库水涨落条件下坡体内浸润面的变动过程，才能对库岸边坡稳定性的演化进行有效评价。

作为国际通用的三维岩土力学软件之一，FLAC3D在求解流-固耦合问题有着明显优势[5]，主要体现在两个方面：一是可采用饱和-非饱和的方法，对非稳定渗流自由面进行模拟，二是可以利用fish语言，对FLAC3D进行二次开发，实现流-固耦合的间接模拟。基于上述优势，本文以西沱变形体为工程背景，采用FLAC3D模拟库水变动条件下岸坡变形体的非稳定渗流场，应用强度系数折减法研究边坡变形体的稳定性演化过程。

1 变形体稳定性分析

1.1 水库运行对岸坡的影响

水库在运行期水位的涨落对库区岸坡的作用是一种典型的水-岩相互作用过程，如图1所示，主要体现在岩土体力学参数以及所受荷载的变化上。

图1 水-岩作用示意图

以三峡水库为例，其运行期水位变动的范围在145.00～175.00 m之间，如图2所示，在水库蓄水和泄水过程中，库岸地下水的抬升和降落一般是滞后库水的变动，滞后效应一般由水库水位消涨的速率、岩土体的渗透参数及山体地下水的补给等因素综合决定，它们将直接影响到库岸边坡的渗透水压力的变化，对岸坡的稳定性可能造成不利的影响，导致边坡失稳。

1.2 强度折减法

在强度折减法计算过程中，把库岸岩体的强度参数(抗剪强度值c和φ)除以一个折减系数Ft，得到新的强度参数c1和φ1，进行重新计算，并不断增大折减系数，直到边坡失稳，此时的折减系数Ft便为岸坡对稳定安全系数Fs[6]。于是有

(1)

本文岩土材料考虑理想弹塑性模型，采用摩尔-库伦屈服准则，其形式如下：

-cosφ=0

(2)

式中：I1为应力张量的第1不变量；J2为应力偏张量的第2不变量；θ为罗德角。

图2 水库运行条件下岸坡水位变化示意图

2 数值模拟

2.1 西沱变形体地质条件

重庆市石柱县城北约47.00 km的西沱镇位于三峡水库右岸，在库区145.00～175.00 m库水消涨区间存在一松软滑坡变形区，其工程地质平面图如图3所示，其中地质剖面3-3’如图4所示。西沱变形体前缘高程为112.00～116.00 m，后缘高程180.00～219.83 m，直至坎崖，平均高差107.80 m，主滑方向318°，长215.00 m，宽305.00 m，一般厚度在2.40～12.00 m，具后缘厚、前缘薄之特点，面积约6.56×104m2，体积约39.3×104m3。变形体地表广泛分布第四系堆积物，为松散土层，多为粘性土，通过现场压水实验测得该松散粘性土的透水性在10～50 Lu之间，约为1.32～4.12×10-4cm/s，不利于降水或库水的快速入渗。

从地质平面图和剖面图可以看出：变形体前缘部分潜在滑移面较缓，后缘较陡，在库水的动、静压力作用下，容易造成变形体的荷载加重，受推移下滑力的影响，后缘的粉质粘土容易出现土体拉裂，可能导致变形体整体失稳，造成岸坡发生大规模滑移，危害库区安全，需要对其进行稳定性评价。

图3 西沱变形体工程地质平面图

图4 3-3’地质剖面图

2.2 有限元模型和计算条件

选择3-3’地质剖面作为典型剖面，建立如图5所示的数值模型，计算模型共划分单元1 492个，节点3 176个。采用三峡水库运行计划，分四个阶段模拟水库水位的变化过程：(1)145.00 m稳定水位；(2)145.00 m升至175.00 m水位，水位上升速度为1 m/d；(3)175.00 m稳定水位；(4)175.00 m下降至145.00 m水位，水位下降速度为1 m/d。根据水文地质调查得到山体侧的水位在库水上升期和下降期分别为160.00 m和175.00 m，水库蓄水前，长观孔ZK401常年水位约为145.33 m。岩土体的物理力学参数取值结合现场勘察和室内实验统计[7]综合确定，非饱和渗透参数根据工程类比[8]得到，如表1所示。

图5 计算模型图

岩土体饱和渗透系数/(cm/s)给水度μ单位贮存量/(10-6m-1)重度/(kN/m3)天然强度参数饱和强度参数天然饱和c/kPaφ/(°)c′/kPaφ′/(°)变形体(Qcol4+ Qdl4)2×10-40.1803.020.720.915.512.511.49.3基岩(J2s)1×10-60.0201.026.527.0870.035.0870.035.0

2.3 非稳定渗流场模拟

采用饱和-非饱和的方法计算岸坡的非稳定渗流场。对FLAC3D进行二次开发，运用修改渗透系数法精确求解渗流的自由面和溢出位置，其具体方法和实现过程见李毅等[5]。首先计算145.00 m稳定水位条件下的边坡渗流场，得到长观孔ZH401的水位为145.36 m，与实测值(145.33 m)相差较小，说明渗流参数取值及初始水文边界合理。初始渗流场如图6所示，在低水位条件下，山体侧地下水位较高，自由面均匀的向坡面下降。在145.00 m水位稳定渗流场模拟的基础上进行蓄水过程模拟，通过编辑时间函数，改变边界水位，并将自由面以上材料的渗透系数重新定义为饱和渗透系数，进行下一水位的计算，得到不同时间(对应不同水位)岸坡的渗流场。蓄水刚完成，水位上升到正常蓄水位175.00 m时，岸坡渗流场等水头云图如图7所示，由图可以看出，由于坡面水位的提升，岸坡周围地下水来不及调整，形成坡面水向后延山体补给的趋势。不同时刻岸坡自由面位置如图8所示，当坡面水位高于山体侧水位时，随着库水位不断的抬升，由于滑坡变形体较基岩渗透性较大，变形体材料区域自由面下降较慢，水位较高，而基岩内自由面下降较快，滑坡体底面水位抬升较慢，渗透压力方向指向坡内，有利于滑坡变形体的稳定。

图6 145.00 m初始水位时等水压云图

图7 水位上升至175.00 m时等水头云图

图8 水位上升过程中不同时期的边坡的自由面

与蓄水期工况类似，编辑时间函数，计算水库水位消落工况，得到不同时间岸坡的非稳定渗流场。水位下降过程中各水位状态下的岸坡浸润面如图9所示，随着库水位的下降，地下水不断由山体侧向河道补给，并随着降落时间的增加，滑坡变形体坡底和坡面水位差不断增大，从而形成较大动水压力，将不利于滑坡变形体的稳定性。

图9 水位下降过程中不同时期的边坡自由面

2.4 稳定性演化

通过计算得到非稳定条件下水库蓄水和降水的渗流场，考虑孔隙水压力变化对变形体的影响，对FLAC3D进行二次开发，采用强度折减法求解滑坡变形体的安全系数。表2和表3分别给出了库水上涨和消落过程中滑坡变形体的安全系数计算结果。值得说明的是，在初始稳定水位状态下，该安全系数代表滑坡变形体初始条件下的稳定性，为稳定渗流场条件下的计算结果，其余水位对应的安全系数则为非稳定渗流条件下的结果。

水位上涨过程中，初始低水位(145.00 m)条件下滑坡变形体处于基本稳定状态，安全系数为1.269，随着库水不断上涨，滑坡变形体的安全系数逐渐增大，当达到正常蓄水位175.00 m时，安全系数也增至最大值为1.313，表明滑坡变形体一直处于较为稳定状态。从非稳定渗流场分析可以得到，在水位上升过程中，滑坡底面水位抬升较慢，变形体内的动水压力主要指向坡内，增加了滑坡变形体的抗滑力，从而增加了其稳定性。但总体而言，安全系数增加不大。

水位消落过程中，初始高水位(175.00 m)条件下，滑坡变形体的安全系数较大为1.297，表明变形体比较稳定，随着库水的消落，滑坡变形体的安全系数也逐渐下降，当库水位降至150.00 m时，滑坡变形体开始出现了局部位移快速增大，出现失稳现象，失稳时滑坡变形的剪应变率云图如图10所示。众多的实验研究和工程实践表明[9]，当边坡失稳时，会产生明显的局部化剪切变形，这种局部化现象一旦发生，表明边坡已经发生了明显的失稳破坏。由非稳定渗流分析可以得到，水位下降过程中，滑坡变形体内动水压力方向指向坡面，不利变形体的稳定安全，且水位下降越大，动水压力越大，滑坡变形体稳定性越差。

表2 库水上涨过程中变形体安全系数

表3 库水消落过程中变形体安全系数

图10 滑坡变形体失稳时剪应变率云图

3 结 论

本文以西沱滑坡变形体为分析对象，采用数值模拟方法分析变形体的稳定性，得到的主要结论如下：

(1)数值计算结果表明，变形体的稳定性随库水上升而增强，随库水下降而降低；

(2)水库水位的周期性涨落带来水环境的改变，影响库岸岩土体的水岩相互作用，可能导致岩土体因不能继续保持与周围环境的平衡而发生失稳等灾变；

(3)建议加强对变形体潜在滑裂区的位移监测，并在滑裂面出口新建挡土墙或抗滑桩；

(4)本次计算未考虑库水对岩土参数的弱化作用，对于遇水软化的岩土材料，需要进一步考虑这一特性。