库水位骤降下非饱和岩质边坡稳定性数值模拟

杨琳

摘要：为研究库水位骤降下岩质边坡渗透稳定特性，推导了基于Hoek-Brown准则下考虑岩质岩体非饱和效应的极限平衡杭剪强度公式，并结合Biship法推导了基于该方法的安全系数表达式，基于Geostudio与Python平台实现了库水位骤降下某岩质边坡渗透稳定性的分析，以4种不同Hoek-Brown参数表征不同类型的岩质边坡并对其进行了敏感性分析。结果表明：库水位下降速率较小时，边坡内部浸润线呈现“先密后疏”的规律，库水位下降速率较大时，边坡内部浸润线呈现“先疏后密”的规律，库水位下降速率越大，前期浸润线的“弯折”程度越大;上部监测点孔压在深度2m以下随深度呈线性分布，而中部监测点与下部监测点孔压随时间呈线性分布。库水位下降速率越快，监测点孔压下降越快;不同工况安全系数随时间呈现先减小后保持不变的趋势，单轴杭压强度σ_c、岩体性状m_i、地质强度指标GSI与安全系数成正相关，而岩体损伤因子D与安全系数成负相关;不同 Hoek-Brown参数对安全系数的敏感程度大小排序为mi≥GSI≥D≥σ_c。

关键词：库水位骤降;渗透稳定;Hoek-Brown准则;Biship法;岩质边坡

中图分类号：TV139.14;TU432 文献标志码：A

库水位骤降是引发边坡失稳的重要因素，大型水库在库水放空的情况下，极易诱发滑坡灾害，如1963年的意大利Vajoint水库，在库水位骤降的情况下发生了大规模的滑坡[1];再如我国三峡水库，2003年建成以来，库水位变动诱发多处滑坡，造成巨大的经济损失。库水位骤降对边坡的影响主要体现在2个方面：①库岸浸润线以下的部分受到水的浮托力作用，有效应力减小;②库水位骤降导致距离库岸较远处浸润线下降出现“滞后”现象[2-3]，使得渗流力指向库岸，从而加大下滑力。

针对库水位骤降下的边坡稳定分析，较多集中对土质边坡的研究，然而大型水库库区存在不少的岩质边坡[4]，国内外学者对岩质边坡的研究较少，且存在以下三方面问题：①仅仅将岩质岩体参数等效为土质岩体参数，以考虑岩体边坡非饱和渗流的方法对岩质边坡进行有限元分析[5-6]，忽略了岩质岩体的本身特有属性;②通过岩体渗流力增量离散的方法来考虑岩质边坡的非饱和效应，忽略了岩质岩体存在损伤效应的事实[7-8];③考虑了岩质岩体的损伤效应，也考虑了岩质岩体的非饱和特性，但忽略了岩质岩体的非饱和抗剪强度[9-10]。

Hoek E.和Brown E.T.[11-i2]在1980年结合岩石性状方面的理论研究和实践经验，基于Griffith的脆性断裂理论，创造性地提出了现今应用极为广泛的岩石强度准则——Hoek-Brown强度准则，综合考虑了岩体各方面的特性[13-17]。基于此，笔者推导了考虑岩质岩体损伤效应和非饱和效应的极限平衡安全系数公式，并结合算例分析了某岩质边坡在库水位骤降下的渗流特性与稳定性变化规律，为相应工况下岩质边坡稳定性的计算提供一定参考。

1 计算理论

1.1 基于Hoek-Brown准则的强度参数推导

2002年，Hoek E.等[16]提出了最新的Hoek-Brown准则，即广义Hoek-Brown经验强度准则：式中：σ₁和σ₃分別为岩体破坏时的最大和最小有效主应力，kPa;σ_c为岩体单轴抗压强度，kPa：m_b与s为与岩体特性有关的参数;m;为岩体性状;GSI为地质强度指标，取值范围为0～100;D为岩体损伤因子，取值范围为0～1;a为曲线参数。

为方便参数转换，绘制Hoek-Brown准则的包络线，不同应力水平下的包络线切线方程等同于相同应力下的Mohr-Coulomb准则[17]，见图1。

Mohr-Coulomb准则可以表示为式中：φ为内摩擦角;c为黏聚力。

在应力水平为σ_a处的切线方程为式中;k为图1所示Mohr-Coulomb切线的斜率;b为Mohr-Coulomb切线与σ₁主应力轴的交点纵坐标值。

本文中GSI≥45，为推导方便，可以认为式（4）中a近似为0.5，则可得到Mohr-Coulomb切线的k与b：

将式（7）、式（8）对照式（5）可得对应于Mohr-Coulomb准则下的c、φ值：

式（9）、式（10）即为基于Hoek-Brown准则下内摩擦角与等效黏聚力的表达式。

1.2 考虑非饱和效应的Hoek-Brown抗剪强度推导

考虑非饱和效应的抗剪强度，多采用Drnevich[18]提出的非饱和抗剪强度理论公式计算：式中：τ_f为岩体实时抗剪强度;σ为法向总应力与孔隙气压力的差值;u_a为孔隙空气压力，这里u_a=0;u_w为孔隙水压力;φ_b为由负孔隙水压力而提高的强度参数。

将式（9）、式（10）代入式（11），便可以得到基于Hoek-Brown准则下的非饱和强度公式：式中：τ_总为基于Hoek-Brown准则下的非饱和抗剪强度;c_总为基于Hoek-Brown准则下的非饱和黏聚力。

1.3 基于Hoek-Brown准则非饱和安全系数公式

根据边坡稳定性理论，采用Biship公式[19]，结合式（13），可以得到基于Hoek-Brown准则考虑非饱和效应的极限平衡安全系数公式：其中式中：FS为安全系数;β为滑块的底边长度;W为土条的重力;α为土条底面的倾斜角;φ为岩体的内摩擦角。

1.4 基于Hoek -Brown准则非饱和安全系数程序实现

基于Hoek-Brown准则非饱和孔压以及安全系数的计算是通过Geostudio和Python平台来实现的。首先在Geostudio里进行边坡饱和一非饱和渗流计算，通过导出各个节点的孔压值，然后在Python平台基于式（14）编制相应的程序实现不同工况安全系数的计算。

2 计算模型及计算参数

2.1 计算模型及边界条件

计算模型选择文献[20]中的计算模型，如图2（a）所示。模型初始条件为左边界abcd以及右边界gh均为13m水头的稳定渗流场：①abcd为库水位变动边界，Seep/W中规定[20]cb边界上每一节点的水头值随时间不断变化，软件每运行一步便判断节点水头是否小于节点高程，若条件成立，则该节点的流量Q便被赋予0;②ah、dedg糖为不透水边界;③gh为13m定水头边界。模型网格如图2（b）所示，计算单位全局网格尺寸约为1m，一共剖分为553个节点668个单元。

2.2 材料参数及计算工况

材料基础参数根据文献[20]取值，饱和渗透系数取k_sat=0.02m/d、重度取γ=26kN/m³，边坡土的土水特征曲线见图3。岩体的基礎Hoek-Brown参数σ_c=150kPa、m_i=14、GSI=65、D=0.6，相应工况为改变Hoek-Brown参数来反映不同类型岩质边坡，从而计算边坡的渗透稳定性。库水位下降速率为0.5、1、2m/d，从正常蓄水位13m下降到死水位8m，计算时间为100d，相应的工况见表1。

3 计算结果分析

3.1 渗流特性分析

3.1.1 岩体内部浸润线变化

不同库水位骤降速率下岩体内部浸润线变化见图4。库水位下降速率较小工况下，边坡内部的浸润线呈现“先密后疏”的规律，即在库水位下降时，岩体内部浸润线下降较为缓慢，表现为不同时间的浸润线较为密集，而在库水位下降结束后，岩体内部浸润线开始逐渐下降;库水位下降速率较大工况下，边坡内部的浸润线变化规律则刚好相反，呈现“先疏后密”的规律，库水位下降速率越大，浸润线向下“弯折”程度越大。

3.1.2 岩体内部孔压变化规律

为反映库水位骤降情况下边坡内部的渗流特性，选取图2（a）中的3个监测面（上部监测面、中部监测面、下部监测面）来实时监测边坡内部的孔压变化，不同工况的孔压变化见图5。

上部监测面孔压变化情况见图5（a）。在库水位下降阶段，监测面上部孔压变化幅度小于下部孔压的变化幅度，当监测面深度大于2m时，孔压随深度成线性变化。相同库水位骤降速率工况下，孔压随时间的变化曲线呈现“先疏后密”，即在库水位骤降的过程中孔压变化幅值较大，孔压线较为稀疏，而在库水位骤降后期孔压变幅较小，即孔压线较为密集。对于不同库水位骤降速率下的孔压变化，库水位骤降速率越大，孔压整体下降越快;对于库水位下降到死水位后的孔压变化，不同库水位骤降速率下的差异较小。

中部监测面孔压变化情况见图5（b）。与上部监测面不同，中部监测面的孔压随深度几乎呈现线性变化，最大孔压变幅约为23.8kPa，大于上部监测面的平均最大孔压变幅。库水位下降速率越大，整体孔压下降越快，但在库水位下降结束后不同库水位下降速率下的孔压差别并不大。

下部监测面的孔压变化情况见图5（c）。孔压随深度的变化呈现线性分布，但是不同库水位骤降情况下的孔压变化幅度差异较大，库水位下降速率 0.5m/d情况下的孔压变化明显缓于库水位骤降速率为1m/d与2m/d的情况，在库水位下降至死水位后，不同库水位下降速率下的孔压变化较小。

3.2 不同工况下的安全系数变化

不同工况下的安全系数变化规律见图6。

图6（a）为不同单轴抗压强度σ_c（即工况A）下的安全系数随时间的变化规律，库水位骤降下的边坡安全系数随时间呈现先减小后趋于稳定的趋势。相同库水位下降速率的情况下，岩体单轴抗压强度越大，整体安全系数也越大;单轴抗压强度σ_c相同、不同库水位下降速率情况下，库水位下降速率越大，安全系数的前期降幅越大，但是最终（即第10d安全系数稳定后）不同库水位下降速率下的安全系数趋于一致。

图6（b）为不同岩体性状m_i（即工况B）下安全系数随时间的变化规律，库水位骤降下的边坡安全系数随时间呈现先减小后趋于稳定的趋势。相同库水位下降速率情况下，m_i越大，整体安全系数也越大;m_i值相同、不同库水位下降速率下，最终安全系数趋于一致。

图6（c）为不同地质强度指标GSI值（即工况C）下的安全系数随时间的变化规律，库水位骤降下的边坡安全系数随时间呈现先减小后趋于稳定的趋势。相同库水位下降速率情况下，GSI值越大，整体安全系数越大;GSI值相同、不同库水位下降速率下，最终安全系数趋于一致。

图6（d）为不同岩体损伤因子D值（即工况D）下的安全系数随时间的变化规律，库水位骤降下的边坡安全系数随时间呈现先减小后趋于稳定的趋势。相同库水位下降速率情况下，D值越大，整体安全系数越小;D值相同、不同库水位下降速率下，最终安全系数趋于一致。

以最终稳定的最小安全系数作为评判不同Hoek-Brown参数的敏感性指标，对于不同单轴抗压强度σ_c来说，σ_c为140、145、150、155kPa时的最小安全系数分别为1.595、1.611、1.627、1.645，随着σ_c的增大，最小安全系数的平均增幅为1.033%;对于不同岩体性状m_i来说，m_i为2、8、14、20时的最小安全系数分别为0.869、1.379、1.628、1.747，平均增幅为28.01%;GSI为45、55、65、75时的最小安全系数分别为1.133、1.341、1.588、1.940，平均增幅为19.65%;对于不同岩体损伤因子D来说，D为0.2、0.4、0.6、0.8时的最小安全系数分别为1.847、1.701、1.586、1.480，平均降幅为7.66%，可见最小安全系数对不同Hoek-Brown参数的敏感程度大小排序为m_i≥GSI≥D≥σ_c。

4 结语

（1）本文推导了一种基于Hoek-Brown准则同时考虑岩质岩体非饱和效应的极限平衡安全系数计算公式，并利用该公式基于Python平台针对不同Hoek-Brown参数进行了渗透稳定敏感性分析，对岩质边坡在库水位骤降下的稳定性评价提供参考。

（2）库水位下降速率较小时，边坡内部浸润线呈现“先密后疏”的规律;库水位下降速率较大时，边坡内部浸润线呈现“先疏后密”的规律，同时库水位下降速率越大，前期浸润线的“弯折”程度越大。

（3）上部监测面孔压在深度2m以下随深度呈线性分布，而中部监测面与下部监测面孔压一直随时间呈线性分布。库水位下降速率越快，监测面孔压下降越快。

（4）库水位骤降情况下不同工况安全系数随时间呈现先减小后保持不变的趋势，单轴抗压强度σ_c、岩体性状m_i、地质强度指标GSI与安全系数成正相关，而岩体损伤因子D与安全系数成负相关;最小安全系数对不同Hoek-Brown参数的敏感程度大小排序为m_i≥GSI≥D≥σ_c。

（5）本文仅仅对表征岩体的4个参数进行了敏感性分析，对于岩体非饱和参数的敏感性有待进一步研究。

参考文献：

[1]徐永強，祁小博，张楠.基于降雨与库水位联合的三舟溪滑坡渗流模拟及稳定性分析[J].水文地质工程地质，2016，43（5）：111-118.

[2]卢博，郭永成，赵二平，等.库水位变化和降雨条件下边坡渗流特性及稳定性分析[J].三峡大学学报（自然科学版），2017，39（2）：54-59.

[3]宋亮.库水位变动条件下滑坡体变形特性的离散元研究[J].交通科技，2017（2）：15-17.

[4]曾铃，邱祥，付宏渊，等.考虑损伤及非饱和效应的炭质泥岩路堤稳定性分析[J].中南大学学报（自然科学版），2016，47（10）：3546-3554.

[5]OH S，LU N.Slope Stability Analysis Under UnsaturatedConditions：Case Studies of Rainfall-Induced Failure of CutSlopes[J].En矛neefing Gedogy，2015，184：96-103.

[6]曾铃，付宏渊，何忠明，等.饱和-非饱和渗流条件下降雨对粗粒土路堤边坡稳定性的影响[J].中南大学学报（自然科学版），2014，45（10）：3614-3620.

[7]荣冠，王思敬，王恩志，等.强降雨下元磨公路典型工程边坡稳定性研究[J].岩石力学与工程学报，2008，27（4）：704-711.

[8]吴李泉，张锋，凌贤长，等.强降雨条件下浙江武义平头村山体高边坡稳定性分析[J].岩石力学与工程学报，2009，28（6）：1193-1199.

[9]傅鹤林，刘运思，李凯，等.裂隙损伤岩体在渗流作用下的边坡稳定性分析[J].中国公路学报，2013，26（4）：29-35.

[10]常来山，李绍臣，颜廷宇.基于岩体损伤模拟的露井联采边坡稳定性[J].煤炭学报，2014，39（增刊2）：359-365.

[11]EBERHARDT E.The Hoek-Brown Failure Criterion[J].Rock Mechanics & Rock Engineering，2012，45（6）：981-988.

[12] WM R Judd.Underground Excavations in Rock[J].Engi-neering Geology，1983，19（3）：244-246.

[13]HOEK E，MARINOS P，BENISSI M.Applicability of theGeological Strength Index （GSI） Classification for VeryWeak and Sheared Rock Masses： the Case of the AthensSchist Formation[J].Bulletin of Engineering Geology &the Environment，1998，57（2）：151-160.

[14]HOEK E，BROWN E T.The Hoek-Brown FailureCriterion-A 1988 Update[J].Journal of Heuristics，1988，16（2）：167-188.

[15] SHAH S.A Modified Hoek-Brown Failure Criterion forJointed Rock Masses[J].Thomas Telford，1992，14（3）：134-138.

[16] HOEK E，CARRANZA-TORRES C.Hoek-Brown FailureCriterion： 2002 Edition[J].Proceedings of the Fifth North A-merican Rock Mechanics Symposium，2002，22（2）：18-22.

[17]贾雪慧，任旭华，张继勋.Hoek-Brown准则在深埋隧洞围岩稳定性分析中的应用[J].三峡大学学报（自然科学版），2016，38（3）：28-31.

[18] DRNEVICH V P，FREDLUND D G，WONG D.Calibration ofThermal Conductivity Sensors for Measuring Soil Suction[J].Geotechnical Testing Journal，1989，12（3）：7.

[19]GEO-SLOPE International Ltd.Seepage Modeling withSEEP/W 2007[R].Calgary，Alberta，Canada： Geo-SlopeInternational Ltd.，1983：244-246.

[20]刘博，李江腾，王泽伟，等.非饱和土渗流特性对库岸边坡稳定性的影响[J].中南大学学报（自然科学版），2014（2）：515-520.