0.   引言

西藏自治区地质环境脆弱、地质结构复杂。区内超过85%的面积在海拔4000 m以上，高山峡谷纵横、高陡斜坡发育。受地震活动、降雨和气候变暖等因素影响，近年来特大高位远程地质灾害频发，造成了重大的人员伤亡和财产损失^[1-3]。西藏自治区是我国冰川面积最大的省区，冰川总面积28664 km²，占全国冰川总面积的48%^[4]。气候变化导致冰川消融加速，冰崩灾害频发。冰崩启动时携带或在运动时铲刮裹入冰屑而形成特殊的含冰碎屑流，动力过程极其复杂，破坏能力较大，易引发碎屑物堵塞河道、冰湖溃决等次生灾害^[5-8]。2018年10月17日西藏自治区雅鲁藏布江西岸色东普沟(N29°45′1.47″，E94°56′14.37″)发生冰崩岩崩，冲刷侵蚀下部冰积堆积体和沟床松散堆积物，形成碎屑流冲入雅鲁藏布江，堆积后形成堰塞坝堵塞河道。上游回水严重威胁加拉村、直白村和赤白村公路水利设施、电力通讯设施和耕地等。堰塞坝溃决严重破坏了下游墨脱县、墨脱亚让水电站，16600人受灾，0.34 km²耕地受灾，7 km公路受损，初步估算经济损失3×10⁸元以上。

目前国内外对于色东普沟崩滑-碎屑流的研究主要集中在现场调查和滑坡机理研究等^[9-11]，也有学者用色东普沟崩滑灾害引发的地震波信号分析灾害诱发历史及强度特征^[12]。上述相关研究成果为本文提供了理论基础，但是针对色东普沟崩滑-碎屑流运动全过程的动力学研究较少。因此，本文在现场勘查基础上，利用DAN^3D软件反演了色东普沟崩滑-碎屑流运动全过程，分析其动力学特性，确定合理的数值分析流变模型和参数；通过物探解译结果验证数值模拟的可靠性，为数值分析流变模型及参数的校验提供了新途径。

1.   色东普崩滑-碎屑流基本概况

1.1   研究区环境条件

色东普沟崩滑-碎屑流位于米林县派镇加拉村，雅鲁藏布江左岸色东普沟。该沟流域面积约66.8 km²，坡度大于30°的流域面积占总面积的61.6%；地形上呈陡-缓的地势，区域整体为南低北高，最高点位于流域北东侧，海拔7283 m，最低点位于沟道沟口处，海拔2750 m，沟域相对高差4533 m；冰川极为发育，终年冰川面积23.6 km²，约占流域面积的35%。根据地形和沟域情况，共可划分15座冰川山峰，雪线海拔在4200 m左右；沟域形态呈“树叶”状，由主沟和多条支沟组成，支沟中有4条较大的支沟，整个流域卫星图见图1。米林县气候属高原温带半湿润季风气候区，崩滑-碎屑流暴发期间，日最高温度在15 ℃左右，最低温度在−5～0 ℃(图2a)，该温度范围有利于泥石流的暴发。区域内降雨多集中在6—9月份，色东普流域气象站显示，在2018年10月17日崩滑-碎屑流发生前2—4天出现过集中降水(图2b)。

图  1  色东普沟崩滑-碎屑流地理位置

a-滑体运动路径边界；b-崩滑区；c-流通区；d-堆积区红色箭头表示滑体的滑动方向

Figure  1.  Location of the ice avalanche-debris flow in Sedongpu Basin

下载: 全尺寸图片 幻灯片

图  2  色东普流域气温、降水变化曲线

Figure  2.  Temperature and precipitation variation in Sedongpu Basion

下载: 全尺寸图片 幻灯片

2018年10月17日该沟发生大规模崩滑-碎屑流灾害，堵断雅鲁藏布江干流并形成堰塞湖，泥石流冲出堆积量约2.05×10⁷ m³。10月19日13时30分左右堰塞体局地出现过流，20日9时50分堰塞湖坝前水位基本恢复至灾前水位。

1.2   崩滑-碎屑流分区

根据色东普沟崩滑-碎屑流运动和堆积特征，将其分为崩滑区、流通区和堆积区（图3）。

图  3  色东普沟主沟剖面图

Figure  3.  Longitudinal profile of the main gully in Sedongpu Basin

下载: 全尺寸图片 幻灯片

(1)崩滑区：崩滑区面积约52.42 km²，海拔多在4000 m以上，分布于冰川顶至下部冰斗、冰舌沟段，距离主沟约1～3 km。上部冰川纵坡较陡，坡度多在40°～50°，坡顶可达60°，该区长年冰雪覆盖，基岩出露，风化冻融冻胀作用强烈，易产生岩质崩塌，堆积在坡脚处；下部冰斗、冰舌处坡度在15°～25°，局部可达8°～10°，长期的冰水堆积作用形成厚度超过10 m的冰积和冰碛物源。在上部产生大规模的冰崩、雪崩或岩崩时，这类堆积物质被大量刨铲带走。该区主要为崩滑-碎屑流的形成汇集水源和提供水动力条件、物源条件。

(2)流通区：该区位于沟道形成区以下至沟口堆积区以上，面积约14.18 km²，沟道长度6.79 km，平均坡降182‰。根据地形和沟道可分为上、下两段。

上段：为冰舌末端至冰蚀凹地前缘，长度3.48 km，平均坡降144‰，以冰水堆积为主。早期冰川退缩后，在冰蚀凹地处堆积大量的冰碛物，厚度超过30 m，并形成小规模的冰湖。1～3支沟及主沟上部的冰积物源，被冰雪及降雨携带，堆积在流通区上段冰蚀凹地内，形成上部宽度约2.0 km，下部宽度约1.0 km的沟道堆积物源。近期，随着全球变暖冰川融化和降雨冲刷作用，原来堆积的固体物质大量被冲刷带走，亦形成泥石流的沟道固体物源。

下段：该段长度3.31 km，地形较陡，平均坡降222‰，为深切V型沟谷，沟道已可见基岩出露。在暴雨和泥石流的冲刷、侧蚀作用下，两岸坡脚堆积体不断被掏空，斜坡稳定性降低，最终发生垮塌，参与碎屑流活动。在经历过2017年和2018年两次特大规模的泥石流后，沟道两侧150～250 m范围内的绝大多数的固体物质，已被冲刷带走，两侧坡体已可见明显的基岩出露，沟床后期以冲刷底蚀作用为主。

(3)堆积区：堆积体主要分布沟口附近，长度2.0 km，宽度0.6 km，堆积厚度20～40 m，堆积体总体积2.40×10⁷m³。目前，大部分仍堆积在沟口，未被水流冲刷带走。

2.   崩滑-碎屑流动力学分析

2.1   DAN^3D简介

DAN(Dynamic Analysis)是加拿大学者HUNGR^[13]开发的等效流体动力分析软件，是滑坡动力分析的重要工具之一，在高速远程滑坡模拟中已有较多的应用^[14-16]，其基本思想是将滑体等效为连续介质流体模型，通过不同流变关系，设定滑坡的运动路径，模拟运动速度、时间、路程以及堆积体等特征。多项研究表明，Frictional模型结合Voellmy模型能较好的模拟滑坡动力特征^[17-19]，因此本文选用Frictional-Voellmy复合模型反演色东普沟崩滑-碎屑流的运动过程。

Frictional方程假设抗剪强度(τ)与正应力(σ)有如下关系：

式中：r_u−孔隙水压力系数(即孔隙水压力与基底正应力之比)；

$\phi$−动摩擦角/(°)。

其中，孔隙水压力系数和动摩擦角之间有如下关系：

孔隙水压力系数(r_u)和动摩擦角($\phi$)可用综合摩擦角(${\phi }_b$)表示，即：

Voellmy模型假设滑体受到的阻力来源于摩擦力和湍流造成的阻力之和，表达式为：

式中：f−摩擦系数，等同于$\tan {\phi }_b$；

ρ−滑体密度；

g−重力加速度；

v−滑体平均速度；

$\xi$−紊流系数，与碎屑流的运动速度和密度有关。

2.2   模型及参数选取

根据1∶10000的比例尺生成的色东普沟崩滑-碎屑流数字高程模型如图4所示。通过试错反演，在崩滑区采用Frictional模型，流通区和堆积区采用Voellmy模型，反演得到一组最优计算参数如表1所示。

图  4  色东普沟崩滑-碎屑流数字高程模型

Figure  4.  Digital elevation model of the ice avalanche-debris flow in Sedongpu Basin

下载: 全尺寸图片 幻灯片

表  1  流变模型和计算参数

Table  1.  Rheological model and simulation parameters

	模型	容重/ (kN·m−3)	摩擦角/ (°)	摩擦 系数	湍流系数/ m	内摩擦角/ (°)	最大铲刮 深度/m
崩滑区	Frictional	25	28	−	−	36	0
流通区	Voellmy	21	−	0.04	1200	24	35
堆积区	Voellmy	21	−	0.05	1000	24	17

下载: 导出CSV 

| 显示表格

2.3   DAN^3D运动过程模拟结果

利用DAN^3D计算的色东普沟崩滑-碎屑流运动全过程如图5。

图  5  崩滑-碎屑流过程不同时刻滑坡堆积形态

Figure  5.  Time-lapse image of deposit distribution during the ice avalanche-debris flow in Sedongpu Basin

下载: 全尺寸图片 幻灯片

模拟结果表明：色东普沟崩滑-碎屑流运动时长约480 s，最远运动距离10000 m，平均运动速度20.83 m/s。

0～30 s阶段，上部冰川发生大规模的冰崩，并混杂岩崩。崩塌体启动后，沿NE-SW方向运动，铲蚀下部陡坡，坡面可见明显的冰川擦痕和沟槽。30 s时滑体离开崩滑区，进入流通区。60～150 s阶段，滑体运动于流通区上段，此段堆积有大量的冰碛物质，伴有小型冰水湖泊，滑体运动到该处时，将其全部冲刷带走。150～250 s阶段，滑体运动于流通区下段，冲刷冰碛沟道物质和流通区两侧的固体物质，于300 s时运动至雅鲁藏布江河道。350 s时，滑体前缘运动至河道对岸。至480 s，滑体后缘停止运动，堆积形态不再变化，滑体停止运动。最终堆积形态如图5（i）所示，堆积体大致呈扇形分布，最大堆积厚度30 m，平均厚度20 m，右岸堆积体最宽1200 m。DAN^3D模拟得到的滑坡堆积范围与实际较为吻合。

崩滑-碎屑流模拟铲刮分布如图6所示，高程3200～4200 m位置处，铲刮深度较大，平均铲刮深度20 m，最大铲刮深度35 m。现场调查发现，色东普沟高程4000～4200 m处，灾后坡面呈沟域形态，有明显切割，可见冰蚀擦痕，测量得实际铲刮深度为20～40 m。在高程3100～4000 m之间，沟内无植被发育，堆积有大量泥、泥沙夹碎块石的沟道堆积物，灾后实际铲刮深度为30～40 m。与图6DAN^3D模拟结果对比发现，模拟得到的铲刮分布与实际较为符合。

图  6  崩滑-碎屑流铲刮分布

Figure  6.  Scraping distribution of the ice avalanche-debris flow

下载: 全尺寸图片 幻灯片

最大速度分布模拟结果如图7所示。崩滑区高程较大，崩滑-碎屑流启动后，巨大的势能转化为动能，滑体获得较快的速度，在地形转折处最大速度增加至24 m/s (运动距离2000 m)。滑体进入流通区上段，地形坡度变缓，同时冲刷切割大量冰碛物质，运动速度降低，最大速度降至8 m/s(运动距离4500 m)。滑体进入流通区下段，坡度增加，运动速度再次增大，运动距离6500 m处，最大速度达24 m/s。进入堆积区后，沟道变宽，坡降变缓，泥石流流速逐步降低，直至冲入雅鲁藏布江停积。将DAN^3D模拟得到的速度分布与利用色东普沟崩滑-碎屑流地震波信号反演得到的运动速度^[12]对比，发现二者较为吻合，证明了本次模拟的可靠性。

图  7  崩滑-碎屑流最大速度分布图

Figure  7.  Maximum velocity contour of theice avalanche-debris flow

下载: 全尺寸图片 幻灯片

3.   崩滑-碎屑流堆积物特征和内部构造的高密度电阻率法探测

物探是滑坡勘察的重要手段之一，为详细了解滑坡堆积物深度和内部结构，本文利用高密度电法(Electrical Resistivity Tomography)测量了3组纵向测线(ERT1-ERT3)和4组横向测线(ERT4-ERT7)，根据电阻率的大小、形态、变化趋势，较为准确地确定了堆积物与基岩的边界。

3.1   测线布置

如图8所示，物探测线位于色东普沟崩滑-碎屑流的流通区与堆积区，共3条纵向剖面，4条横剖面。物探测线采用Wenner电极阵列进行电阻率测量，利用RES2Dinv软件对电阻率数据进行反演分析^[20]。

图  8  崩滑-碎屑流物探测线分布

Figure  8.  ERT survey lines distribution of the ice avalanche-debris flow

下载: 全尺寸图片 幻灯片

3.2   物探解译结果

横向测线(ERT4-ERT7)物探成果见图9。结果显示：不同剖面的滑坡沉积物厚度不同，同一剖面不同位置处也不相同。ERT7剖面表层电阻率较低，在150 Ω·m以内，推测为地表覆盖层、堆积层。地表低阻层以下岩体电阻率逐渐升高，整体在2000 Ω·m以上，推测为较完整基岩。整条剖面堆积物深度变化不大，为5～18 m，最大堆积物厚度位于距测线起点400 m处(图9a)；ERT6剖面测量结果表明，该剖面表层电阻率在300 Ω·m以内，地表低阻层以下电阻率大于2000 Ω·m，堆积物的深度范围为10～20 m，最大堆积物厚度位于距测线起点180 m处(图9b)；ERT5剖面表层电阻率在200 Ω·m以内，地表低阻层下岩体电阻率整体在2000 Ω·m以上，堆积物的深度范围为5～15 m，最大堆积物厚度位于距测线起点60 m处(图9c)；ERT4表层电阻率在300 Ω·m以内，地表低阻层下岩体电阻率整体在1500 Ω·m以上，堆积物的深度范围为10～12 m，最大堆积物厚度位于距测线起点50 m处(图9d)。

图  9  横向测线(ERT4-ERT7)物探成果图

Figure  9.  Transverse ERT profiles along the lines ERT4 to ERT7

下载: 全尺寸图片 幻灯片

图10为纵向测线(ERT1-ERT3)的物探成果，ERT1剖面位于河道堆积区内，堆积物深度较大，深度范围为10～29 m，最大堆积物厚度位于距测线起点200 m处(图10a)；ERT2剖面位于流通区下段和堆积区，高程范围2750～3000 m，堆积物深度范围为5～40 m，最大堆积物厚度位于距测线起点1420 m处(图10b)；ERT3剖面堆积物深度范围为10～19 m，最大堆积物厚度位于距测线起点420 m处(图10c)。

下载: 全尺寸图片 幻灯片

图  10  纵向测线(ERT1-ERT3)物探成果图

Figure  10.  Longitudinal ERT profiles along the lines ERT1 to ERT3

下载: 全尺寸图片 幻灯片

3.3   物探解译与 DAN^3D模拟结果对比

DAN^3D模拟得到的堆积物最终分布如图8所示，堆积物的平均深度20 m，最大深度30 m。物探解译与DAN^3D模拟结果对比见表2和图11。如表2所示，DAN^3D模拟与物探解译得到的最大碎屑流堆积厚度基本一致，在ERT2测线处有一定差异，相差15 m，但在其他测线处较为吻合。如图11所示，DAN^3D可以较好地模拟此次崩滑-碎屑流堆积形态分布，尽管由于DAN^3D模型将滑坡体视为了等效流体^[21]，导致ERT4测线处，DAN^3D模拟值偏大2～4 m(图11b)，堆积物深度被高估，但模拟的准确性仍然较高。

表  2  物探解译与DAN^3D模拟堆积物最大深度对比（单位: m）

Table  2.  Comparison of the maximum deposition depth between the ERT method and DAN^3D simulation (Unit: m)

	ERT1	ERT2	ERT3	ERT4	ERT5	ERT6	ERT7
物探解译	29	40	19	12	15	20	18
DAN3D模拟	30	25	20	15	15	10	10

下载: 导出CSV 

| 显示表格

图  11  物探解译和DAN^3D模拟堆积物厚度对比

Figure  11.  Comparison of the landslide deposits depth from the ERT interpretation and DAN^3D

下载: 全尺寸图片 幻灯片

4.   结论

本文针对2018年10月17日色东普崩滑-碎屑流，基于真实遥感影像建立三维数值模型，利用DAN^3D分析其动力学特征，在野外探测基础上，结合高密度电阻率法分析了滑体的堆积特征，得到如下结论：

(1)利用动力学模型DAN^3D反演得到碎屑流堆积特征与滑体运动速度、铲刮深度等动力特征参数，发现Frictional-Voellmy复合模型具有很好的模拟效果，并通过反演试算给出最佳流变参数。结果表明，色东普沟崩滑-碎屑流运动时长约480 s，最远运动距离10000 m，平均运动速度20.83 m/s。通过对比，发现DAN^3D模拟得到的运动持续时间和速度变化与地震波反演结果基本一致，最终堆积分布与现场勘察结果基本一致，证明模拟结果具有较好的准确度。

(2)在野外调查的基础上，结合高密度电阻率成像(ERT)技术，分析了滑体的堆积特征，对比DAN^3D模拟和物探解译得到的碎屑流堆积厚度，发现两种方法得到的滑坡碎屑流堆积厚度基本一致，模拟的准确性较好，研究结论及方法对于此类灾害的机理研究及未来防治工程提供了新的研究思路。