库水位变化对库岸滑坡稳定性影响的预测研究

黄 栋，乔建平，张小刚

（1.长江科学院a.工程安全与灾害防治研究所；b.水利部水工程安全与病害防治工程技术研究中心；c.国家大坝安全工程技术研究中心，武汉 430010；2.中国科学院地表过程与山地灾害重点实验室，

成都 610041；3.中国科学院水利部成都山地灾害与环境研究所，成都 610041）

库水位变化对库岸滑坡稳定性影响的预测研究

黄 栋1a，1b，1c，2，3，乔建平2，3，张小刚2，3

（1.长江科学院a.工程安全与灾害防治研究所；b.水利部水工程安全与病害防治工程技术研究中心；c.国家大坝安全工程技术研究中心，武汉 430010；2.中国科学院地表过程与山地灾害重点实验室，

成都 610041；3.中国科学院水利部成都山地灾害与环境研究所，成都 610041）

以三峡库区晒网坝滑坡为研究对象，在充分考虑双场耦合的条件下，结合极限平衡分析法，得到了库水位变化过程中考虑位移增量和孔隙水压力增量的库岸滑坡稳定性系数的计算方法，并对库水位变化作用下晒网坝滑坡的长期稳定性进行分析。结果显示：随着时间的变化，稳定性系数逐步趋于一定值1.172，且坡体内浸润线也随着库水位的调度逐步保持稳定，趋于不变。

水库滑坡；渗流－应力场耦合；极限平衡法；稳定性分析

据不完全统计，自1982年以来，三峡库区发生滑坡、崩塌、泥石流多达70多次，其中规模较大的就有40余次。大坝修建后所形成的狭长型水库于2003年6月正式蓄水，随着蓄水水位的不断提升，水库库容不断增大。水库蓄水后，因水库季节性或年际性水位变化，对库岸进行反复冲刷、浸泡。古、老滑坡将会受到河流地质作用和库水动力作用，使得原已稳定的滑坡再度复活。2010年10月26日，三峡工程正式达到175 m正常蓄水位，随后库水位将长期在145～175 m之间运行，水位变幅高达30 m。因此，对库水位长期变化下库岸滑坡稳定性的研究显得十分重要［1－2］。

目前，水库运行条件下的滑坡稳定性分析已有大量研究成果，理论和应用都逐渐趋于成熟。但是，对库水位长期变化作用下库岸滑坡稳定性的预测研究鲜有报导。因此，本文以三峡库区晒网坝滑坡为研究背景，通过利用渗流－应力场耦合的分析方法结合极限平衡法，对三峡库区库水位长期周期性调节作用下晒网坝滑坡稳定性进行预测研究。

1 理论分析

1.1 应力－渗流场耦合方程

1.1.1 渗流场有限元方程

对渗流的偏微分方程采用伽辽金的加权残余方法，得到渗流方程有限元形式：

式中：［B］为梯度矩阵；［C］为单元水力传导系数矩阵；｛H｝为节点水头矢量；〈N〉为插值函数矢量；q为穿越单元边的单位流量；τ为单元厚度；t为时间；δ为渗流方程瞬态分析时的存贮项；A为区域单元的综合；L为单元边的综合。

渗流方程的有限元形式可简写为

式中：［K］为单元本征矩阵；［M］为单元质量矩阵；｛Q｝为单元水流量。方程是进行瞬态渗流分析时的通用有限元方程，对于稳定状态分析，水头不是时间的函数，项｛H｝，t消失，方程简化为

此式是基本渗流方程的简短有限元形式。

1.1.2 应力场有限元方程

采用虚功原理获得体系平衡方程［3－5］，其表达式为

式中：｛ε*｝为虚位移；｛δ*｝为虚应变；｛Δσ｝为系统内应力；｛F｝为外部点荷载。

建立有限元方程：体模量；［B］为应变矩阵；［D］为排水条件下有效的应力－应变本构矩阵；［K］为刚度矩阵；［Ld］为耦合矩阵；｛Δδ｝为位移增量；｛Δuw｝为孔隙水压力增量向量。

其方程的简化形式可写为

1.1.3 耦合方程

从以上推导可以得出用位移增量和孔隙水压力增量为方程域变量的饱和－非饱和岩土体耦合固结控制方程，即

由于渗流－应力场存在空间上的差异性和在时间上的非稳态性，进行耦合场有限元分析就比一般弹塑性有限元分析要困难得多。在进行各矩阵计算的过程中，要区别饱和区和非饱和区，即要考虑单元与自由面的关系。如单元位于水下，则要考虑地下水的渗流和静水压力等问题；如单元位于水上，处于非饱和状态，则可认为单元内各点渗透流速为零，不存在渗流力和静水压力；另外，在岩土体体积改变时，其体积改变量并不等于单元体中水量的改变量。

1.2 稳定性理论

Morgenstern-Price法既满足力平衡又满足力矩平衡条件，该法由Morgenstern和Price（1965）［6］以微分方程组形式提出。

图1 条块受力分析Fig.1 Forces on the stick

将边坡划分为n个条块，取第j条块为隔离体（图1）。该条块所受到的作用力有：重力Wj；水压力合力Uj＝γuWj，γu为孔隙水压力比；外力Qj，与竖直方向成ωj角；滑动面上有效法向力Nj′；调用的抗剪强度Tj／Fs，Fs为安全系数；条块间法向作用力Ej和Ej＋1，与底面的垂直距离分别是zj和zj＋1；条块间的剪切力λfjEj和λfj＋1Ej＋1。于是有

其中：hj和bj为条块高度和宽度；αj为底面倾角。

求解上述方程组应满足边界条件：

式中：Ea和Eb为端部条间力；Ma和Mb为力矩。

这样，式（8）和式（10）组成的方程组可简化为

其中En＋1和Mn＋1分别称为剩余力和剩余力矩，可递推求得。

另外，式（8）和式（9）条间力函数f（s）为

式中：A为表征条间力函数形状的系数；B（s）为表征峰值点位置的函数。

其中：s为滑体归一化后的水平坐标；BBF为峰值点的位置。于是，通过Newton-Raphson法迭代求解方程组（11），便可求得λ和Fs。

2 数值计算与分析

2.1 模型的建立

为了分析库水位的变动对晒网坝滑坡体稳定性的影响，结合三峡库区晒网坝滑坡勘察报告，采取如图2所示的计算剖面（晒网坝滑坡剖面Ⅲ）和计算模型。由于计算渗流场的瞬态变化计算点上不能一次达到平衡，因此，如果计算中某点的计算水头大于该点坐标，则下一次迭代中定义该点水头为竖直方向坐标值，直到水头值与坐标值相符为止。有限元网格划分为四边形法，划分为7个区域，共剖分为1 544个结点，1 346个单元。计算模型的应力场边界条件为：水平方向固定，竖直向自由；底边界竖直向和水平向均固定。

图2 稳定性计算模型Fig.2 Stability calculation model

2.2 参数的选取及工况

仍采用晒网坝滑坡稳定性分析模型，岩土体参数也与之相同，模型及边界条件设置见图2。模型仅考虑库水位对晒网坝滑坡的影响，不考虑岩土体物理力学参数的变化。岩土体物理力学参数取值如表1。

表1 物理力学参数Table1 Physical and m echanical param eters

根据三峡库区调度（图3），确定模型的计算工况：模拟三峡库区5个水文年的调度过程，其中水位缓慢升、降速度为0.5，0.25 m／d，快速升、降速度为1，1.5 m／d，时间步长2 d一步，共609步，具体见图4。

图3 三峡库区库水位调度图（长江水利委员会）Fig.3 Water level fluctuation in the Three-Gorges reservoir area（provided by Changjiang Water Resources Comm ission）

2.3 结果分析

根据文献［7－9］最小稳定性系数一般发生在库水位下降时。因此，我们只需对5个水文年的库水位下降时的稳定性系数进行分析就可以得出最小稳定性系数的变化规律。5个水文年库水位下降时的稳定性系数如图5。

图5 库水位下降时稳定性系数（5个水文年）Fig.5 Landslide stability factor during water level decrease in the reservoir（in five hydrological years）

图4 库水位调度图（5个水文年）Fig.4 W ater level fluctuation（in five hydrological years）

对图5分析可知，最小稳定性系数变化趋势符合一阶指数衰减趋势y＝y0＋A1·e－（x－x0）／t1。因此，对最小稳定性系数进行非线性拟合得到晒网坝滑坡最小稳定性系数非线性拟合曲线如图6。

图6 最小稳定性系数拟合图（5个水文年）Fig.6 Fitted curve ofm inimum landslide stability factor（in five hydrological years）

晒网坝滑坡最小稳定性系数一阶指数曲线为

式中：x为时间（d）；y为最小稳定性系数。

对式（14）进行分析可得，当x趋向于无穷大（＋∞）时，最小稳定性系数为定值1.172 45。因此，根据拟合曲线可知，随着时间的变化稳定性系数逐步趋于一定值1.172 45。

5个水文年库水位缓慢下降时的浸润线研究表明，坡体内浸润线随着库水位的调度逐步保持稳定，趋于不变。

3 结 论

（1）以多孔介质渗流理论、Biot固结方程、质量守恒定律等方法为基础，从流固耦合理论出发，推导了基于虚功原理的渗流场和应力场的有限元空间域及时间域的离散方程，得到了用位移增量和孔隙水压力增量为方程域变量的饱和－非饱和土体耦合固结控制方程。

（2）在充分考虑双场耦合的条件下，结合极限平衡分析法，对晒网坝滑坡稳定性进行了数值计算，得到了库水位变化过程中考虑位移增量和孔隙水压力增量的库岸滑坡稳定性系数。

（3）对晒网坝滑坡5个水文年作用下的最小稳定性系数进行拟合得到晒网坝滑坡最小稳定性系数非线性拟合曲线公式，公式显示晒网坝滑坡在蓄水初期2～3 a将会产生变形破坏，随后将保持稳定。

［1］ 乔建平.滑坡减灾理论与实践［M］.成都：科学技术出版社，1997.（QIAO Jian-ping.Theory and Practice of Landslide Disaster Reduction［M］.Chengdu：Science and Technology Press，1997.（in Chinese））

［2］ 中村浩之，王恭先.论水库滑坡［J］.水土保持通报，1990，10（1）：53－64.（NAKAMURA H，WANG Gongxian.On Reservoir Landslide［J］.Bulletin of Soil and Water Conservation，1990，10（1）：53－64.（in Chinese））

［3］ 毛昶熙.渗流计算分析与控制［M］.北京：中国水利水电出版社，2003.（MAOChang-xi.Seepage Calculation Analysis and Control［M］.Beijing：China Water Power Press，2003.（in Chinese））

［4］ MATSUI T，SAN K C.Finite Element Slope Stability A-nalysis by Shear Strength Reduction Technique［J］.Soils and Foundations，1992，32（1）：59－70.

［5］ KULHAWY F H.Finite Element Analysis of the Behavior of Embankments［D］.California：University of California，1969.

［6］ MORGENSTERN N R，PRICE V E.The Analysis of the Stability of General Slip Surfaces［J］.Geotechnique，1965，15（1）：79－93.

［7］ BERILGENM M.Investigation of Stability of Slopes under Drawdown Conditions［J］.Computers and Geotechnics，2007，34（2）：81－91.

［8］ 刘新喜，夏元友，张显书，等.库水位下降对滑坡稳定性的影响［J］.岩石力学与工程学报，2005，24（8）：1399－1444.（LIU Xin-xi，XIA Yuan-you，ZHANG Xian-shu，et al.Effects of Drawdown of ReservoirWater Level on Landslide Stability［J］.Chinese Journal of Rock Mechanics and Engineering，2005，24（8）：1399－1444.（in Chinese））

［9］ 廖红建，盛 谦，高石夯，等.库水位下降对滑坡体稳定性的影响［J］.岩石力学与工程学报，2005，24（19）：3454－3458.（LIAO Hong-jian，SHENG Qian，GAO Shi-ben，et al.Influence of Drawdown of Reservoir Water Level on Landslide Stability［J］.Chinese Journal of Rock Mechanics and Engineering，2005，24（19）：3454－3458.（in Chinese） ）

（编辑：曾小汉）

2.2 锚杆应力监测

地下厂房共布置6个锚杆应力监测断面，每个断面内设27～34支1点式或3点式锚杆应力计对支护锚杆应力和岩体应力进行监测，见图3。

图2 地下厂房典型变形监测断面图Fig.2 Typical section of deformation m onitoring for the underground powerhouse

图3 地下厂房典型锚杆应力监测断面图Fig.3 Typical section of anchor bar stressmonitoring for the underground powerhouse

2.3 锚索锚固力监测

地下厂房（F8断层范围内）上游侧岩锚吊车梁上方设置2台锚索测力计，下游侧岩锚吊车梁上方设置1台锚索测力计，对预应力锚索锚固力进行监测。

图4 地下厂房典型钢筋应力监测断面图Fig.4 Typical section of steel-bar stressmonitoring for the lattice beam of the underground powerhouse

2.4 格构梁的钢筋应力监测

地下厂房（F8断层范围内）顶拱设置3个格构梁钢筋应力监测断面，每个断面内布置10支钢筋计对格构梁钢筋应力进行监测，见图4。

3 监测成果与围岩稳定性分析

3.1 围岩变形

地下厂房围岩表面变形最大值为15.79 mm，发生在F8断层上游侧边墙47.63 m高程处。围岩表面变形平均值为4.37 mm。岩体质量较好部位的围岩变形均＜2 mm。F8断层较为破碎的部位，围岩表面变形基本在5 mm以上。

地下厂房F8断层典型部位的围岩变形随时间变化过程如图5和图6所示。地下厂房围岩变形主要发生在开挖施工期间。开挖完成后，围岩进入蠕动变形期，变形呈缓慢增加的趋势，2009年12月以后地下厂房各部位的围岩变形趋于稳定，围岩变形已收敛。

图5 地下厂房围岩顶拱变形过程线图Fig.5 History of rock deformation in the arch roof of the underground powerhouse

图6 地下厂房围岩边墙变形过程线图Fig.6 H istory of rock deformation in the side-wall of the underground powerhouse

围岩表面及内部变形空间分布如图7所示。桩号厂0＋120和厂房西端等部位的岩体较为完整，变形规律符合一般岩石开挖变形规律，变形主要发生在距围岩表面2～3 m的松动圈范围内，表面变形最大，沿孔深方向变形越来越小。F8断层穿越的部位，岩体质量稍差，该区段内顶拱部位的围岩变形发生在距围岩表面约5 m的范围内，基本呈现从孔口至孔底逐渐减弱的趋势。

图7 地下厂房围岩变形空间分布图（2011年6月30日）Fig.7 Spatial distribution of surrounding rockmass deformation of the underground powerhouse（Jun.30，2011）

地下厂房围岩边墙表面变形大于顶拱表面变形，这与地下洞室围岩变形的基本规律相符。围岩顶拱变形主要发生在开挖初期（I—IV层开挖），第IV层开挖完成后，围岩顶拱变形趋于平稳，第V—VIII层开挖对顶拱变形影响较小，见图8。

图8 开挖初期顶拱围岩变形过程线图Fig.8 H istory of rock deformation in the arch roof during initial stage of the excavation

3.2 支护锚杆应力

锚杆最大拉应力为150.5 MPa，发生在厂0＋120断面上游侧42.0 m高程处。地下厂房平均支护锚杆应力为22.84 MPa。经统计，地下厂房79.3%的观测锚杆其应力测值小于50.0 MPa，约5.4%的锚杆应力位于50.0～100.0 MPa区间，约10.0%的锚杆应力大于100.0 MPa。锚杆应力增量主要发生在地下厂房开挖期间，之后逐渐趋于平稳，见图9。2009年12月后地下厂房内的锚杆应力总体上稳定，仅随围岩温度变化而变化，与温度变化趋势相反。锚杆应力的变动范围较小，平均变幅为5 MPa左右。

图9 地下厂房锚杆应力过程线图Fig.9 History of anchor bar stress in the underground powerhouse

3.3 锚索锚固力

锚索锚固力最大值为1 143.29 kN。地下厂房开挖期间，锚索张力随围岩变形有少量增加，增加幅度约为50 kN。3台锚索测力计均有少许松驰现象，松弛应力约为15 kN。目前各台锚索测力计的测值均已趋于稳定，见图10。

与围岩变形和锚杆应力相比，锚索锚固力变化较小，这主要是因为锚索位于地下厂房拱脚部位，其安装时间为2006年3月，此时地下厂房顶拱的围岩变形已基本趋于稳定（见图5）。

图10 地下厂房锚索锚固力过程线图Fig.10 History of anchor force in the underground powerhouse

3.4 格构梁钢筋应力

地下厂房格构梁最大钢筋拉应力63.76 MPa，发生在厂0＋95桩号顶拱下层钢筋上；最大钢筋压应力75.56 MPa，发生在厂0＋104桩号顶拱与拱脚之间的下层钢筋上。平均钢筋拉应力为16.02 MPa，平均压应力为25.08 MPa。

格构梁钢筋应力变化过程线如图11所示。由图中可见，钢筋应力变化主要发生在地下厂房开挖期间，当顶拱岩体发生变形时，岩体对格构梁拱圈产生压应力，从而导致钢筋应力改变。开挖完成后，围岩变形趋于稳定，钢筋应力也趋于稳定，但随厂房内温度变化而呈正弦规律变化，在一定范围内波动。

图11 地下厂房格构梁钢筋应力过程线图Fig.11 H istory of stress on the steel-bar of lattice beam of the underground powerhouse

从钢筋应力过程线来看，目前格构梁钢筋应力主要受洞室气温的影响。当气温升高时，格构梁混凝土温度升高快于岩体温度，混凝土膨胀，导致钢筋压应力增加，拉应力减少；当气温降低时，格构梁混凝土温度降低快于岩体温度，混凝土收缩，导致钢筋压应力减小，拉应力增加。

3.5 围岩稳定性评价

从围岩变形、支护锚杆应力、锚索锚固力及格构梁钢筋应力监测的成果来看，地下厂房围岩各项监测效应量正常，均小于设计指标，且进入电站试运行期后一直处于比较稳定的状态。

变形、锚杆应力及钢筋应力的变化过程也反映出了三者之间的相互影响关系：当围岩变形增加时，

Prediction of Landslide Stability of Reservoir Bank under the Influence of W ater Level Fluctuation

HUANG Dong1，2，3，QIAO Jian-ping2，3，ZHANG Xiao-gang2，3
（1.Yangtze River Scientific Research Institute，Wuhan 430010，China；2.Key Laboratory of Mountain Hazards and Surface Process，Chengdu 610041，China；3.Institute of Mountain Hazards and Environment，Chinese Academy of Sciences，Chengdu 610041，China）

Landslide stability affected by long-term reservoir water level fluctuation is predicted by finite equilibrium method considering the coupling between seepage field and stress field.Shaiwangba landslide in the Three Gorges Reservoir area is taken as a case study.Through numerical calculation，the stability factor considering the displacement incremental variable and pore water pressure incremental variable is acquired.By analyzing the stability factor during the decrease of reservoirwater level in 5 hydrological years，the nonlinear fitted curve ofminimum stability factor is obtained.Results show that the stability factor tends to be steady at1.172，and the landslide phreatic line tends to be constant with the variation of reservoir water level.It’s concluded that deformation will take place at the Shaiwangba landslide during the first2-3 years after the impoundment，and will remain stable in subsequence.

reservoir landslide；coupling between seepage field and stress filed；limit equilibrium method；stability analysis形进行监测，见图2。

TU457

A

1001－5485（2012）09－0103－04

10.3969／j.issn.1001－5485.2012.09.024

2012－07－05

中国科学院知识创新工程重要方向项目（KZCX2－YW－Q03－5）；国家自然科学基金（41001007）；国家科技支撑计划课题（2012BAK10B04）

黄 栋（1980－），男，湖北黄石人，助理研究员，博士，从事滑坡机理研究，（电话）027－82829878（电子信箱）hd1017＠sohu.com。