天花板水电站田坝村堆积体地质-力学机理-形变耦合分析预警

刘昕昕， 汤明高， 王李娜， 向育才， 吴辉隆

(成都理工大学，地质灾害防治与地质环境保护国家重点实验室， 成都 610059)

滑坡的孕育和诱发由于具有复杂性、多样性及随机性的特点，其预警预报一直是相关领域研究的难点和热点[1]。滑坡的地质-力学机理-形变(geology-mechanism-deformation，GMD)预警预报模型是以斜坡的地质结构为基础，结合多源监测数据揭示斜坡内部的地质—力学作用模式，再通过数值和物理模拟将其反映至宏观表现(变形)的综合GMD耦合模型[2]。GMD数值预警模型从深层次上将斜坡变形过程和内部力学作用机理及滑动面的行程演化过程相统一，突破了以往数学定量模型对中长期预报和事件预报的限制[3]，目前在中国已得到广泛的应用。李秀珍[4]建立了泄滩滑坡及金龙山II区蠕滑体和石榴树包滑坡的GMD数值预警模型，使用FLAC3D模拟了泄滩滑坡和金龙山II区蠕滑体的变形过程并使用时间延拓的方法预测了滑坡未来的变形趋势，对石榴树包滑坡不同水库蓄水位的变形稳定性进行模拟计算，界定其所处变形阶段并对其未来的演化趋势进行预测。丁秀美[5]使用GMD数值模型预警方法对泄滩滑坡蓄水条件下的变形及稳定性进行预测。马旭[6]建立了西山村滑坡的GMD预警模型，运用时间延拓结合滑坡变形阶段对西山村滑坡进入加速阶段的时间进行了预测。冉佳鑫[7]建立了垮梁子滑坡GMD预警模型，针对后缘拉陷槽不同冲水高度下滑坡的变形情况进行了UDEC模拟，得到了不同拉陷槽水位的滑坡预警判据。李焕焕[8]使用强度折减法对渭玉高速公路 K32+690～K32+840 段滑坡的破坏过程进行模拟，并对其使用抗滑桩治理后的发展趋势进行了预测。

借助FLAC3D软件平台建立田坝村堆积体GMD数值预报模型，使用自编命令强度折减法仅折减堆积层的强度参数，并采用位移收敛率作为计算收敛标准，相比软件自带强度折减命令以体系不平衡力比率收敛标准，在堆积层滑坡中更具准确性，且可以更直观反映滑坡的变形破坏过程，从而找出滑坡破坏的临界状态下各监测点的变形值作为预警依据。随后以典型监测点的临界变形值为依据提出堆积体预警阈值，使用动态流固耦合同步时间步模拟的方法，得到了不同库水位下降速率下堆积体变形达到预警阈值所需的天数，可作为天花板堆积体“库水位下降速率”的事件预报判据，也为在FLAC3D软件中对水位变动这一事件的模拟提供了一种新的思路。

1 堆积体概况

田坝村堆积体位于牛栏江左岸云南省昭通市巧家县老甸乡田坝村，距天花板水电站坝址1.0～2.3 km，属近坝库段(图1)。前缘顺河岸最大宽度约1 km，分布高程1 020～1 662 m，总体积约4 295×104m3。堆积体属低纬度山地季风气候，受高山峡谷影响，立体气候显著，故有寒、温、热兼有的立体气候特征，年平均降水量约881.8 mm，降雨多集中在5—9月。水库正常蓄水位1 071.0 m，死水位1 050.0 m。

图1 研究区域示意图

堆积体物质组成主要为碎块石土夹孤石。堆积体地形复杂，由高处到低处整体上表现为顶部平缓、中部较陡、下部较缓的阶梯状地形。区内横向地形起伏大，有两条NNE向冲沟及同向山梁在岸坡分布。构成了由平缓台地、陡缓斜坡、冲沟及凸起山梁组成的复杂地形。堆积体内地下水类型为裂隙地下水和孔隙地下水两类。岸坡处地下水水位高于牛栏江水面，因此地下水补给主要来自大气降水，并以渗入、地表径流等形式排向牛栏江。基岩裂隙水主要赋存于岩体中的节理裂缝及断层带中，受断裂构造形成的裂缝网络所控制[9]。

2010年12月18日，天花板水电站正式开始下闸蓄水，堆积体前缘随即出现塌岸，蓄水后一个月左右，堆积体表面开始出现裂缝，并不断增加至数十条，裂缝最长可达500 m。2011年起堆积体开始建立GPS地表变形监测系统，先后共布置了6个断面，共计26个监测点(图2)，此后经历的几次库水位较大幅度变动，均监测到堆积体发生一定程度的变形。

图2 堆积体分区及监测点布置

2 田坝村堆积体变形特征及机理分析

2.1 水库蓄水后库水位变化特征

天花板水库为不完全季调节水库，2010年12月18日开始下闸蓄水，库水位正常调动区间为1 050～1 071 m。图3为天花板水电站库水位监测曲线，最早的库水位监测数据为2011年3月1日，库水高度为1 049.1 m，此后至2011年9月1日期间库水位一直维持在1 050～1 060 m小幅度涨落，到2011年10月17日，水库蓄水突然至正常蓄水位(1 071 m)附近。之后天花板水电站经历有四次库水位较大幅波动(单次连续水位变动幅度超过10 m)，具体变动情况如表1所示，其中第4次库水位变动为配合8·3鲁甸地震抗震救灾工作，期间水库最低水位1 035 m，超过水库正常调动区间。

表1 天花板水电站库水位变动分析

①～④为波次

2.2 蓄水后堆积体变形特征

2.2.1 堆积体裂缝发展

图4为田坝村堆积体变形及裂缝分布，大致经历：出现—贯通—后延—稳定4个阶段。

图4 堆积体裂缝统计

(1)初期出现阶段(2010年12月—2011年3月)：2010年12月水库蓄水后，堆积体前缘随即出现塌岸裂缝，次年1月，I区湾子社房屋开始产生裂缝(L4、L5)，3月I-2区及II-2区后缘陡坡产生拉裂缝(L28—L32)。

(2)裂缝贯通阶段(2011年3月—2012年8月)：前一阶段产生的裂缝进一步发展扩大并逐渐贯通。

(3)裂缝后延阶段(2012年8月—2013年11月)：2012年8月—10月，II区前部裂缝(L31、L32)急剧增大，黄家坪子房屋开裂，随即裂缝L35～L38出现。

(4)稳定发展阶段(2013年11月至今)：此后，田坝村堆积体裂缝经长期发展，堆积体中上部位裂缝大多已被填塞，仅新出现三条规模较小的裂缝，分别为I-2的 L18、L19及II-4区后缘的L39。

2.2.2 堆积体变形发展

图5为田坝村堆积体经历四次大幅度水位变动后的监测位移等值线图，堆积体主要变形区域为I区及II区前部，自建立监测以来，该区域变形量已超过4 000 mm， 2013年11月前，II区后部，即II-3及II-4区变形不明显，此后至2014年11月期间，II-3及II-4区变形量急剧增加。此外，自水库蓄水以来，III区变形发展情况不明显。

图5 经历4次大幅泄水后堆积体监测位移变化

结合田坝村堆积体裂缝发展特征及变形监测来看，前堆积体缘产生较大塌岸导致后部土体失去承载力，从而导致斜坡前部(I区及II区前部)先产生变形，并随着变形的不断发展，II-3区前缘产生临空，随即II-3区和II-4区产生变形，具有渐进后退式滑坡[10]的典型特征。

2.3 堆积体变形与库水位变化的关系

根据堆积体典型监测点监测曲线显示(图6)，堆积体监测位移的增加与库水位大幅下降关系密切，此前历史上4次库水位大幅变动时期的库水位调动均引起堆积体位移产生明显的响应。分析4次天花板水电站库水位大幅波动时期各代表监测点的合位移变化量与库水位变化量之间的关系(图7)可以更明显地反映其相关性。以第三次库水位骤降为例[图7(c)]，2013年11月5日起库水位开始连续下降，持续8 d，最大日降幅超过4 m/d，11月10日，各监测点变形量开始相应并不断增加，11月13日，库水位降到最低值，此时监测点变动速率达到峰值附近，此后库水位开始抬升，监测点变形速率开始逐渐降低，此次库水位变动充分体现了监测点变形速率与库水位之间的负相关及滞后性关系，且滞后期约在3～5 d，其他几次水库位变动及堆积体位移响应稍复杂，但同样大致体现了负相关和滞后性的关系，因此田坝村堆积体是典型的受库水位下降控制的渗透压力型滑坡[11]。

图6 典型监测点位移监测曲线

图7 典型监测点位移变化量与库水位变化量关系曲线

2.4 堆积体变形模式

综前所述，田坝村堆积体主要变形模式如下。

(1)水电站水库蓄水改变了坡体前缘的水文地质条件，降低了堆积体的力学参数，引起坡体应力重分布。在库水骤降时又增加了库岸土体内的静水压力和动水压力，致使岸坡下滑力增加，超过岸坡土体的抗剪强度，引起库边岸坡变形、坍塌，进而牵引至上方坡体出现拉裂变形。

(2)根据堆积体分区位置，I区、II-1区、II-2区变形模式为岸坡下部受库水升降影响首先变形或产生圆弧形破坏，进而牵引上部产生较大变形；II-3区、II -4区变形模式是在II -1区、II -2区出现较大变形时，由于失去顶托作用，沿基覆界面(基岩层面)产生向上游偏坡外的变形滑动破坏。

(3)III区崩塌堆积体底面平缓，破坏模式为圆弧形，但不受库水影响，未发生较大变形。

因此田坝村堆积体是受库水位骤降控制的渐进后退式斜坡，即库水位下降导致斜坡前部产生塌岸和变形，随后在斜坡中部产生裂缝、临空，后部土体失去支撑从而产生变形(图8)。

图8 田坝村堆积体变形模式

3 堆积体GMD预警模型研究

3.1 预警模型

岩土体失稳(尤其是滑动失稳)是一个变形从量变的累积到质变的过程，可以看作其最薄弱滑动面不断贯通的过程，此过程中随着滑动面塑性区的增加，滑动面整体承受剪切力的能力不断降低，当累积的强度降低维持不住岩土体的稳定，岩土体就会失稳，该过程存在3个因子，其相互关系如图9所示，因此，稳定性度量K可以看作是位移U的函数。对于一个特定的斜坡，不论斜坡结构形式如何，其失稳(K<1)必然存在一个所能承受的最大变形量Umax，这是GMD(geology-mechanism-deformation)数值模型预报的基本出发点，而其分析计算的关键就是最大变形量Umax的确定[3]。强度折减法可以通过折减岩土体参数的方法近似模拟岩土体从稳定状态向不稳定状态演化的过程，找出岩土体失稳的临界状态，就可以建立其临界失稳状态的最大变形量预警判据[12]。随后，只要考虑将斜坡失稳的事件因素(库水位下降速率)，计算某一事件下堆积体的变形量是否达到最大变形量预警阈值，就可以得到事件预报的预警判据(图10)。

图9 斜坡失稳过程各因素相互关系

图10 田坝村堆积体GMD预警模型

3.2 变形量预警阈值

图11为田坝村堆积体FLAC3D三维计算模型，1 071 m水库正常高水位浸润面基于钻孔水位监测数据建立；坐标轴选取为：x轴指向长江下游，y轴垂直指向坡外，z轴垂直向上；边界条件固定模型底面及前、后、左、右5个边界作为数值计算的位移边界，I区及II区前缘水面以下的临空部分设置为透水边界，其余所有边界为不透水；变形量监测点与前文各分区典型监测点相对应；岩土体计算参数基于室内土工实验、相关文献及类似地质条件滑坡的岩土体参数类比确定(表2)[9]。

表2 数值计算岩土体参数

图11 数值计算模型

3.2.1 最大变形量预警判据

在FLAC3D中通过自编命令流强度折减的方法，对田坝村堆积体在高库水位运行状态下(1 071 m库水位)的失稳过程，之所以选择高库水位作为强度折减时的潜水面，是因为考虑到高库水位状态下斜坡中处于饱水状态下的土体范围最大，其强度参数静水压力均为最不利状态，因此计算的结果偏保守，更符合预警的原则，该状态下堆积体孔隙水压力分布情况如图12所示。传统强度折减法采用的是力不平衡比率小于某一临界值作为收敛标准，其临界值难以确定在学术界尚存争议[13]，强度折减计算收敛标准取相邻10个时步位移监测点的变形收敛率小于0.01%，折减过程如图13所示。

图12 初始孔隙水压力云图

紫色符号表示自定义变量；蓝色符号为用户指定参数

图14为堆积体强度折减计算得到的时步-变形曲线，得到田坝村堆积体1 071 m库水位状态运行下的稳定性系数为1.046。临界破坏状态的堆积体变形情况如图15所示，堆积体前缘在水压力的作用下产生大量塌岸，变形量值大于10 m，除前缘塌岸以外，堆积体主要形区域为I区，变形量在4～10 m，II区变形量均在2 m以上，数值模拟变形区域与前文基于监测数据的堆积体变形区域基本一致，因此可以认为结果可信。各分区典型监测点在临界状态下的变形量如表3所示。

表3 临界破坏状态各监测点变形值

图14 强度折减计算监测点位移曲线

图15 临界状态堆积体变形云图

3.3 库水位调动预警判据

在FLAC3D中调用渗流模块，对比模拟收敛速度与实际位移监测曲线得到本文模型计算时间步与真实时间之比为20步∶1 d。使用FLAC3D内嵌fish语言分时步调整堆积体前缘水面以下的临空部分的透水情况，并同步开启力学模块，以实现动态流固耦合的效果，库水位下降模拟实现流程如图16所示。考虑到天花板水电站低水位运行高程为1 050 m，因此模拟最大降幅为21 m。预警判据以典型监测点C103作为预警点，达到其最大变形量的85%(5.466 m)为红色预警阈值，60%(3.858 m)为橙色预警阈值。

图16 库水位下降模拟过程

图17为监测点C103在1～5 m/d库水位下降作用下的变形值-时间步曲线，参考天花板水电站历史变形监测数据，2014年8月4日—2014年8月15日，天花板水电站库水位在12 d连续下降了34.03 m，平均日降幅为2.84 m，监测站C103变形值增大了1.59 m，结果显示监测点C103在库水位下降33 m时的变形值为1.43 m，误差为0.16 m，可以认为模拟结果可信。

图17 不同库水位下降速率C103时间(步)-位移曲线

堆积体在库水位以1 m/d及2 m/d的下降速率下，监测点C103变形值趋于稳定，在库水位下降到1 050 m时变形值分别为0.46、0.79 m，未达到预警阈值；堆积体在3 m/d的库水位下降速率下，在库水位下降到1 050 m时变形值为0.73 m，未达到预警阈值，但是库水位继续下降，监测点变形速率会逐步增加；4、5 m/d的库水位下降速率下预警点位移呈指数型上升，其中库水位以4 m/d下降速率下C103会在第6天达到橙色预警阈值，库水位以5 m/d下降速率下降C103会在第2天、第3天分别达到橙色预警及红色预警阈值，相对应的库水位下降幅度分别为8.7、11.8 m。

4 结论

(1)通过分析田坝村堆积体的变形及裂缝发展历史，总结得出堆积体为前缘产生较大塌岸导致后部土体失去承载力发生形变，并不断向后发展的渐进后退式滑坡。

(2)田坝村堆积体为受库水位骤降因素控制的渗透压力型滑坡，以往4次单次连续降幅超过10 m较大幅度库水位下降均引起了堆积体变形的响应，并体现了一定的滞后性。

(3)通过强度折减法计算得到堆积体在1 071 m高库水位运行状态下的安全系数为1.046，建立了各典型监测点的最大变形量预警判据，可为将来的滑坡相关预警工作提供一定的参考。

(4)以C103监测点为依据，通过流固耦合数值模拟得到田坝村库水位下降速率预警判据：在库水位变动范围内，监测点C103的变形量在1、2、3 m/d的库水位下降速率不会达到预警阈值；在4 m/d的库水位下降速率下在第6天达到橙色预警阈值；在5 m/d库水位下降速率下分别在第2天和第3天达到红色预警阈值。