基于FLAC 3D的高边坡安全稳定静动力分析

高佳东，刘逸辉，张 巍

(1.广东省水利电力勘测设计研究院有限公司，广东 广州 510635；2.华南农业大学水利与土木工程学院，广东 广州 510640)

我国重大水利工程大坝建设处于较为领先的水平，而水库大坝作为重要基础设施之一，广泛应用于防洪、供水、泥沙控制、改善水环境和解决水资源分布不均等问题[1-2]。但由于工程质量缺陷及洪水、地震、冻害现象、气候变化等极端事件的不利影响，以及人为因素影响、管理不当等引发的非常规运行工况，均会影响到水库大坝的运行安全[3-7]、经济损失和环境影响[8-10]。

甘肃玉门抽水蓄能电站属于一等大(1)型工程，工程区地震基本烈度为Ⅷ度。其上水库布置于照壁山山顶台地，东南外侧陡崖高差大，坡面陡峭，卸荷强烈，裂隙发育，岩体破碎，其稳定性直接影响到上水库大坝的运行安全。

随着现代计算土力学的发展，为更好地研究边坡安全稳定性，数值计算已成为主要计算手段。任东兴等11]、董杰华等[12]以边坡工程实例为背景，利用ABAQUS软件模拟边坡破坏机理，研究其应力变形特征及其稳定性；何彬[13]、吴旭敏等[14]针对边坡失稳的问题，利用Geo-studio软件进行数值模拟分析边坡稳定性的影响因素；徐宏等[15]运用MIDAS GTS NX软件计算不同力学参数下的边坡安全系数，研究边坡安全系数影响因素的敏感性；匡义等[16]、毛正君等[17]通过对不同地区的边坡采用FLAC 3D软件进行数值模拟计算，研究边坡破坏模式和稳定性影响因素。目前，大部分三维有限元法中采用隐式求解方案，求解非线性复杂问题所花费时间较长，且存在累计误差；在模拟物理上极限状态的不稳定过程时采用静态方程求解，存在数值上的障碍。而FLAC 3D软件基于显式有限差分法计算原理[22]，在求解非线性问题时，不需通过迭代来满足本构关系，避免了累计误差以及解决了花费时间较长的问题，应力通过应力-应变关系，跟随应变实时变化，进而跟踪系统的演化过程，使得计算时间几乎与线性本构关系相同；同时采用动态运动方程描述地质材料在达到强度极限或屈服极限时发生破坏的力学行为，可用于模拟岩体开挖卸荷效应问题。

因此，本文利用FLAC 3D软件对上水库外侧高边坡进行计算分析，评价上水库在不同工况下边坡的稳定性和抗震安全性，并结合分析结果，提出边坡防护工程措施及上水库工程布置合理性建议。

1 工程概况

上水库库区为北、东及南侧三面开挖，西侧筑坝成库。库坝轴线环库全长2863m，西北侧窄，东南侧稍宽。主坝布置在西南侧和南侧，长1228m，最大坝高65.8m，副坝位于水库东南侧，长92m，最大坝高12m，坝型均采用沥青混凝土面板堆石坝，库坝线距离东南端侧山体的最近水平距离约200m。上水库外侧高边坡示意如图1所示。

图1 上水库外侧高边坡示意图

照壁山四周现状地形均为陡崖，东南外侧陡崖边缘与沟底最大高差约900m，平面距离1350～1900m。陡崖段平面长度350～600m，高度400m左右，山体坡度以50°～70°为主，局部近直立或倒坡；在高程2350～2400m形成宽度100～300m的缓坡平台，2350m平台向沟底形成一连续的斜坡段，平面长度约1km，地形坡度整体以小于30°为主。受断裂构造、岩体卸荷及冻融作用等影响，陡崖段近山顶处岩体形成不稳定块体，崩塌、坐落现象强烈，在坡脚及缓坡段堆积有崩积物。

2 计算模型及参数

2.1 模型建立

2.1.1三维整体模型

选取上水库外侧照壁山区域建立三维整体计算分析模型。模型的x坐标轴与大地坐标系的东西方向重合，总长2740m；模型的y坐标轴与大地坐标系的南北方向重合，总长1130m；模型的z坐标轴与铅直向重合，向下覆盖的最低高程为1260m，向上则建模至地表，覆盖范围内最低的地表高程为2020m，最高的地表高程为2900m，故模型在高程方向上的高度为760～1640m。

计算分析模型共考虑覆盖范围内的强风化、弱微风化奥陶系微泥晶灰岩等岩层和F1断层，以及上水库临边坡一侧的堆石坝和坝后堆渣体。其中，与堆石坝平行布置的坝后堆渣体以实体单元形式建出(如图2所示)，其他坝后堆渣体则通过在计算时施加覆盖范围内地表等效荷载的方式来考虑。

图2 三维整体模型网格图

2.1.2准三维模型

在主坝最大坝高位置沿垂直坝轴线建立准三维模型。该模型长1930m，高690～1320m，覆盖范围的地表高程从2320m增加至2950m，网格划分考虑强风化、弱风化和微风化地层，以及F1断层和L1结构面。进一步给出了准三维模型的边坡坡表监测点布置，如图3所示。

图3 准三维模型的模型图

2.2 本构模型及计算参数

边坡岩体的本构模型采用FLAC 3D的考虑拉伸修正的摩尔库伦模型和遍布节理模型。根据地勘结果可知，上水库外侧陡崖地层岩性主要为微泥晶灰岩，其边坡岩体物理力学性质和透水性指标见表1、构面力学指标见表2。采用动力时程计算地震作用时，岩体的动态物理力学参数与静态值一致。考虑工程主要关注边坡岩体的稳定性，初始地应力通过自重应力场方法计算。

表1 边坡岩体物理力学性质和透水性指标

表2 边坡岩体构面力学指标

2.3 计算分析工况

计算工况选取天然状态工况、库水渗漏工况和设计地震工况。其中设计地震工况是在坝后堆渣工况的基础上，计算设计地震荷载采用100年超越概率2%设防水平、基岩水平峰值加速度为514.2g工况。

以上工况均采用整体三维模型进行计算分析，其中设计地震工况采用拟静力法计算，并进一步采用准三维模型，通过动力时程法进行计算分析。根据地震安全性评价建议的基岩设计规准反应谱，采用人工合成的方法，生成加速度时程曲线，如图4所示。并根据规范考虑计算分析模型的底部到地表的距离超过50m，因此直接将加速度时程折半后，作为设计地震荷载输入模型底部。

图4 加速度时程曲线

3 边坡稳定分析

3.1 天然状态工况

通过探究上水库外侧陡崖边坡的初始地应力场与边坡应力分布的对应规律，可为后续边坡稳定性分析提供可靠的初始应力场条件。采用三维强度折减方法，对处于天然状态下的边坡求取强度储备安全系数。如图5所示，当强度折减系数为1.50时，边坡关键点变形的变化曲线开始出现拐点，显示边坡开始进入临界失稳状态。因此，天然状态下的边坡安全系数可取为1.50。

图5 边坡关键点变形随强度参数折减系数的变化曲线

对于处于临界状态的边坡，其失稳模式为剪切滑移+后缘拉裂破坏。其中，前缘剪出口位于2520m高程左右的边坡坡表，后缘为拉裂破坏，拉裂面在2800m高程附近的边坡坡表出露。边坡变形和边坡剪切滑动带分布如图6—7所示。

图6 处于临界状态的边坡变形图

图7 处于临界状态的边坡剪切滑动带分布

3.2 库水渗漏工况

假定正常蓄水工况条件下，库盆防渗措施失效，引发上水库库水持续渗入边坡，且上水库水源可得到稳定补给，使上水库库水水位保持不变，使计算分析模型最终形成稳定渗流场。

在上水库库水持续入渗边坡且维持库水水位不变的条件下，地下水在边坡坡表溢出，最终形成稳定渗流场时，边坡坡表溢出区域的高程在2600m附近，其孔隙水压力如图8所示。

图8 库水渗漏入渗边坡后孔隙水压力在三维模型内的分布(单位：MPa)

根据边坡渗流场计算渗流荷载，并分析由渗流荷载引起的边坡岩体变形。由库水入渗的渗透荷载引起的边坡变形在18.1mm以内，变形矢量总体为沿着边坡坡表指向下，且边坡岩体的变形分布受到F1断层的一定影响，在断层上下盘出现一定的不连续分布特征。边坡变形和边坡剪切滑动带分布如图9—10所示。

图9 处于临界状态的边坡变形图

图10 处于临界状态的边坡剪切滑动带分布

如图11所示，库水入渗工况条件下的边坡安全系数可取为1.48，对应的边坡失稳模式与天然状态工况下的基本一致。

图11 边坡关键点变形随强度参数折减系数的变化曲线

3.3 设计地震工况

3.3.1拟静力法

首先采用拟静力法计算设计地震荷载，如图12所示，设计地震工况条件下的边坡安全系数可取为1.26，对应的边坡失稳模式与天然状态工况下的基本一致。

图12 边坡关键点变形随强度参数折减系数的变化曲线

3.3.2动力时程法

通过准三维计算分析模型，采用动力时程分析方法考虑设计地震荷载效应，计算分析边坡的稳定性。采用图4给出的加速度记录输入模型。通过观察边坡岩体变形，分析边坡岩体动位移分布；通过A、B、C、D四个部位监测点反应边坡的动位移时程曲线特性，以及地震动作用完成时刻的残余动位移特征；在完成动力时程计算分析的基础上，进一步采用强度折减法求取边坡安全系数。

(1)边坡岩体的动位移分布和边坡监测点的动位移特性在地震动作用过程中，边坡岩体总体上呈现位移一致响应的特征，边坡的整体性较好。边坡水平向动位移分布在2.19～24.50cm。边坡岩体的动位移在绝大部分时刻均指向坡外侧，且在断层部位出现动位移不连续现象，如图13所示。

图13 边坡岩体水平向动位移分布

如图14所示，边坡坡表不同高程的监测点，在地震动荷载作用下的动位移时程曲线形态相同，变化规律也较为一致，仅是量值上存在一定区别。监测点A～D在地震动作用完成时的残余动位移分布在6.74～12.62cm，见表3，其中量值相对较大的部位为结构面附近的测点，其他部位的残余动位移相对较小。

表3 边坡监测点的动位移特征参数统计 单位：cm

图14 边坡监测点的动位移时程曲线

(2)安全系数。如图15所示，在加速度记录输入的设计地震工况条件下的边坡安全系数可取为1.15，对应的边坡失稳模式与天然状态工况下的基本一致。

图15 边坡关键点变形随强度参数折减系数的变化曲线

对比拟静力法计算结果发现，在整体三维模型计算时，潜在失稳边坡受到两侧未失稳区域的约束，即存在边坡变形的三维效应，使得计算得到的安全系数大于准三维模型计算得到的安全系数。

4 结语

研究采用FLAC 3D软件，对玉门抽水蓄能电站上水库外侧高边坡的稳定性开展了静动力数值模拟分析，主要结论为：

(1)针对上水库外侧高边坡进行稳定分析，各计算工况的边坡潜在失稳区域和失稳模式一致，安全系数均大于容许值，上水库布置方案基本合理。考虑因断层构造、岩溶溶隙等不良地质现象影响，产生的局部不稳定块体，需进行局部的锚喷加固处理。

(2)地震动作用对边坡岩体的结构面影响相对显著，结构面附近的浅表层边坡发生一定量值的错动且在临空面出露区域更为明显。但总体位移绝对量值较小，对边坡岩体的影响有限。

(3)研究通过动力时程法分析地震工况时输入了单一的时程曲线，为提高计算结果准确性，需进一步考虑采用不同的相角对应的时程曲线进行计算。