舟曲江顶崖滑坡的早期判识及风险评估研究

韩旭东，付 杰，李严严，王高峰，曹 琛

（1.中国地质调查局水文地质环境地质调查中心/自然资源部地质环境监测工程技术创新中心，河北 保定 071051；2.中国地质调查局水文地质环境地质调查中心/自然资源部三峡库区地质灾害—重庆巫山野外科学观测研究基地，重庆 404700；3.北京工业大学城市建设学部，北京 100124；4.吉林大学建设工程学院，吉林 长春 130026）

随着国家对地质灾害防治工作更高的要求及期许，如何提前判识滑坡变形，并对该滑坡进行早期风险评估已成为地质灾害防治领域最关注的焦点及难点。在滑坡变形早期判识工作方面，国外发展比较迅猛，欧洲尤其是意大利已经基本实现了基于InSAR 技术的全国范围滑坡变形早期普查[1]。近年来，国内研究者积极探索用于滑坡变形早期判识的InSAR 技术，在四川茂县新磨村滑坡，金沙江白格滑坡及上下游滑坡，丹巴五里牌滑坡等众多滑坡案例中得到了成功应用[2−4]。在滑坡变形早期风险评估方面，传统地质分析可以有效地预测滑坡变形失稳。在此基础上，利用斜坡稳定性力学分析计算方法，如极限平衡计算方法，有限元数值计算方法及离散元数值计算方法[5−6]，可获取滑坡变形失稳时的滑体几何信息。滑坡变形失稳后，对滑体的动力过程进行评估时，国内外研究者采用了雪橇经验模型、等效流体模型和离散元模型等进行分析，获取滑体的动力特征[7−9]。

2018年7月12日上午8 时左右，白龙江甘肃舟曲南峪段江顶崖区域发生滑坡，滑坡体向白龙江河道错动，推挤白龙江河道，造成白龙江上游河道水位上涨，国道345 线中断，南峪乡段众多房屋被淹。本文以该滑坡为分析案例，综合研究滑坡变形的早期判识及风险评估。首先通过小基线集技术（SABS-InSAR）对该滑坡进行早期形变及范围判识，然后基于地质-力学联合分析，判断该滑坡的变形失稳机制及滑体几何形态，最后通过基于深度积分的连续介质等效流体数值模型计算该滑坡的动力过程，评估其堵江风险。文章拟通过对江顶崖滑坡进行综合全面的分析，建立一套可行的滑坡变形早期判识及风险评估全流程分析模式，为舟曲白龙江沿岸类似滑坡变形的早期判识及风险评估提供借鉴。

1 舟曲江顶崖滑坡

甘肃舟曲县南峪乡段江顶崖滑坡位于白龙江北岸。滑坡体呈较典型“簸箕”状，滑坡体物质主要由碎石和黄土混合而成，碎石含量为65%～75%，碎石主要为灰色的炭质板岩和千枚岩岩块，局部夹有较大的灰岩岩块。此外，滑坡体中下部物质掺杂了大量的灰黑色黏土与炭质板岩以及千枚岩碎屑，遇水极易软化，抗剪强度低[10]。受连续降雨影响，2018年7月12日上午8 时左右，江顶崖区域发生滑坡，滑坡体长度约为680 m，宽度约为210 m，平均厚度约为35 m，滑坡体总量达到了5×106m3，为大型滑坡，见图1（a）。

图1 舟曲江顶崖滑坡和江顶崖滑坡区域遥感影像Fig.1 Jiangdingya ladslide and remote sensing image in the Jiangdingya landslide area in Zhouqu

结合区域遥感影像可以发现，江顶崖滑坡为典型的白龙江江岸堆积层滑坡，其发育于巨型老滑坡堆积体前端，老滑坡堆积体长度约1 600 m，宽度约1 400 m，见图1（b）。

2 InSAR 技术解译及分析

2.1 研究区域数据

本文使用16 期欧洲航空局发布的哨兵-1A SAR影像数据及POD 精密定轨星历数据。数据拍摄时间为江顶崖滑坡失稳前的2018年1月9日—7月8日，数据间隔为12 d。区域地形数据采用日本宇宙航空研究所2014—2015年发布的较高精度的ALOS PALSAR观测数据，见表1。

表1 哨兵-1 影像数据和ALOS PALSAR 地形数据信息Table 1 Information of Sentinel-1A images and ALOS PALSAR topography data

2.2 区域地表解译

SBAS-InSAR 是Berardino 等[11]提出的一种时序InSAR 技术解译方法。该方法通过设置时间和空间两个维度的基线阈值构建InSAR 影像数据间的两两组合，较PS-InSAR 技术解译方法可获取更多高相干性的SAR 影像数据组合，提高了解释过程中对InSAR 影像数据在时间和空间2 个维度上的采样频率，有利于减轻解译过程中InSAR 影像数据间的时空失相关和大气相位影响[12]。SBAS-InSAR 作为一种毫米精度的InSAR 后处理算法已经被广泛认可。本文采用此技术，首先获取研究区域的地表形变速率及地表位移[13]。解译流程如下：

（1）生成16 期哨兵-1A SAR 影像数据间的空间和时间基线关联表，空间和时间基线阈值分别设置为80 m和40 d。42 对哨兵-1A SAR 影像数据组合情况，见图2。

图2 16 期哨兵-1A SAR 影像数据的空间与时间基线关联表Fig.2 Association of spatial and temporal baselines of 16 periods of sentinel-1A SAR data

（2）生成42 对哨兵-1A SAR 影像数据组合的差分干涉图，采用16 期哨兵-1A POD 精密定轨星历数据和ALOS PALSAR 地形数据去除差分干涉时的平地相位，采用经典的Goldstein 滤波方法滤除干涉相位噪声，提高信噪比，最终获取了研究区域42 幅质量较好的SAR 影像数据的干涉图，其中江顶崖滑坡区域2018年6月14日影像与2018年7月8日的影像差分干涉图结果示例见图3。

图3 江顶崖滑坡区域SAR 影像数据干涉相位分布Fig.3 Interference phase distribution of SAR images in the Jiangdingya landslide area

（3）获取研究区域地表形变速率和时序地表形变，选取40 个稳定的地表控制点用于轨道精炼及去除地形残差相位的影响，进行时空过滤计算，移除大气延迟相位，采用经典的最小费用流方法进行干涉相位解缠，最后采用奇异矩阵变换，获取研究区域雷达视线方向的地表变形速率及时序变形。江顶崖滑坡2018年1月9日—7月8日期间平均变形速率及地表累积变形曲线，见图4 和图5。

图5 江顶崖滑坡体雷达视线方向累积变形量曲线Fig.5 Cumulative deformation curves of the Jiangdingya landslide in radar line of sight

根据SBAS-InSAR 技术解译结果，基于差分干涉图和地表变形速率图（图3 和图4），可以清晰辨别江顶崖滑坡的范围。从江顶崖滑坡时序累积变形曲线（图5）看，江顶崖滑坡时序累积变形曲线是一条典型的滑坡变形斋藤曲线[14]，具有典型的3 个阶段演化特征：（1）滑坡体初始变形阶段。江顶崖滑坡体开始变形，变形曲线变现为减速变形特征，时间大约发生于2018年3月10日之前；（2）滑坡体匀速变形阶段。江顶崖滑坡体近匀速变形，变形曲线表现为等速特征，发生时间段大约从2018年3月10日—6月2日；（3）滑坡体加速变形阶段。江顶崖坡体变形速率不断增加，发生时间为2018年6月2日以后。

3 地质-力学分析

SBAS-InSAR 解译结果表明，舟曲江顶崖滑坡失稳前已存在形变加速趋势，结合现场地质调查资料可以明显看出（图6），该滑坡发育于老堆积体前缘，为典型的老堆积层滑坡前缘局部变形，坡脚长期受白龙江江水冲刷，降雨可触发其牵引式滑移破坏。从滑坡体内部物质的成分看，中下部物质中掺杂了大量的灰黑色黏土和炭质板岩碎屑，遇水极易软化，抗剪强度低，滑坡体很可能会在该深度范围内滑动。为进一步获取江顶崖滑坡失稳破坏时的滑体几何数据，本文采用极限平衡分析方法对该滑坡进行力学分析计算，估算滑体的几何信息。根据文献[10]及该区域地质调查的实验数据，滑坡体物质的物理力学参数取值见表2。本次滑坡体的滑面搜索过程，采用GeoStuido 软件内嵌的常规剪入—剪出滑面搜索方法。为使滑面的搜索结果科学客观，滑面搜索时，不固定剪入点和剪出点，而是将剪入口设置为SBAS-InSAR 技术解译分析过程中圈定的滑坡范围后部附近区域，将剪出口设置为江顶崖滑坡前缘坡脚处的小陡坎区域，稳定性系数计算方法采用常规的Bishop 和Janbu 法[15]。经计算，江顶崖滑坡的稳定性参数在1.02～1.04，整体处于欠稳定状态，存在失稳破坏的可能，江顶崖滑坡滑体平均厚度约为35 m，滑床整体坡度较缓，见图7。

图7 江顶崖滑坡稳定性计算模型Fig.7 Stability calculation model for the Jiangdingya landslide

表2 滑坡体物质物理力学性质Table 2 Physical and mechanical properties of the landslide body

图6 江顶崖滑坡地质剖面Fig.6 Geological profile of the Jiangdingya landslide

4 堵江风险评估

江顶崖滑坡为典型的白龙江江岸滑坡，通过SBAS-InSAR 技术解译及地质-力学联合分析，江顶崖滑坡失稳前已存在变形失稳趋势，具有堵江风险，本文采用基于深度积分的连续介质模型对其动力过程进行计算分析，为堵江风险评估提供分析数据，计算模型如下：

式中：t—时间；

x—笛卡尔坐标系下的x方向坐标；

q—守恒变量；

fx— 方向通量；

s—源项；

h—滑体厚度；

ux— 方向的滑移速度；

g—重力加速度；

Sgxx— 方向滑床坡度项；

Sfxx— 方向滑坡体基底摩擦项。

数值计算前，采用2 m 间距对江顶崖滑坡的滑床及前部的白龙江河道地形进行剖分，白龙江河道的地形数据参考实地调查数据及文献[16]中的数据。数值计算时采用Godunov 型中心格式有限体积方法离散计算模型，采用HLLC 格式黎曼求解器进行通量求解，采用MUSCL 重构法和Hancock 预测-校正格式提高Godunov 型中心格式的时空求解精度，采用隐式格式处理基底摩擦值，数值求解过程见文献[8,17]。此外，舟曲江顶崖滑坡的滑床整体坡度较小，滑体不具有高速流动性特征，基底摩擦模型采用较贴合这一特征的库伦摩擦公式[18]：

式中：φ—滑体对滑床的摩擦角，本次取16°，数值计算结果见图8 和图9。

图8 江顶崖滑体滑后几何形态Fig.8 Post-sliding geometry of the Jiangdingya landslide body

图9 江顶崖滑体滑移速度监测曲线Fig.9 Velocity monitoring curve of the Jiangdingya landslide body

数值计算结果表明，江顶崖滑坡体运动方式大体表现为整体向前错动，前缘推挤白龙江河道，计算时间20 s 后，江顶崖滑坡体前缘错动完成，从速度监测点1，2，3 位置处的滑坡体速度监测曲线看：当滑坡体前缘错动完成后，滑坡体滑移速度从前缘到后缘快速降至为0，表现为牵引式运动特征，并且滑坡体在向白龙江河道错动过程中，滑移速度并不是很大，最大值约为2.2 m/s，堵江风险很低。

5 结论

针对滑坡变形的早期判识及风险评估研究，本文提出了小基线集雷达干涉（SBAS-InSAR）技术解译分析，地质-力学联合分析，动力过程数值模拟分析相结合的滑坡变形早期判识与风险评估全流程分析模式，并选取甘肃舟曲白龙江流域江顶崖堆积层滑坡案例进行了综合研究分析，得出以下结论：

（1）SBAS-InSAR 技术解译及分析结果能准确地判识江顶崖滑坡的范围及早期形变特征，滑坡的变形破坏模式为牵引式，滑坡体长度约680 m，宽度约210 m，与江顶崖滑坡实际数据吻合。

（2）基于早期识别信息，江顶崖滑坡的地质-力学联合分析结果表明，江顶崖滑坡为典型的老堆积层滑坡前缘局部变形，滑坡体平均厚度约35 m，滑床整体坡度较缓，失稳后运移速度不大，分析结果与现场地质调查匹配。

（3）江顶崖滑坡的动力过程计算分析结果表明，滑坡体滑移速度不大，最大值约为2.2 m/s，运动方式表现为推挤白龙江河道，堵江可能性较小，当江顶崖滑坡体前缘错动完成后，滑坡体的滑移速度从前缘到后缘快速降至为0，表现为牵引式运动特征。动力过程数值计算时采用的库伦摩擦模型符合江顶崖滑坡的滑体滑移摩擦特征，模型选择合理。

针对江顶崖滑坡案例，本文的分析结果与江顶崖滑坡实际相符，吻合度较高，采取的综合分析方法及研究模式可用于舟曲白龙江沿岸类似滑坡的早期判识及风险评估。