“9·5”泸定震后落井沟泥石流灾变过程预测分析

李林泽, 常 鸣, 李宏杰, 胡锦凤, 徐恒志, 刘沛源

(1. 成都理工大学 地质灾害防治与地质环境保护国家重点实验室, 四川 成都 610059;2. 山东正元数字城市建设有限公司, 山东 烟台 264670)

0 引言

我国是世界上震后地质灾害活动最频繁、损失最严重的国家之一,其中泥石流灾害是震后山区地质灾害中危险性最大的灾害之一。强烈的地震作用在震后短期内具有诱发大量山体滑坡并产生松散堆积物分布在沟道内及沟道两岸的能量,对山体造成的内部损伤影响亦可达数百年之久[1-2]。与此同时,地震也大大降低了诱发泥石流灾害的降雨阈值,导致灾难性的泥石流事件频繁发生,规模加大,人民的生命财产安全受到极大地危害和威胁[3-4]。如2008年5月12日汶川特大地震后在地震的影响范围内接连暴发了大规模泥石流,该地震10年后虽然仍有大规模泥石流暴发,但在震后短期内泥石流启动降雨阈值普遍降低的情况下雨季暴发泥石流的频次更高,规模更大,受灾更严重。因此,对震后短期内可能因强降雨引发泥石流的运动过程和冲出规模需要更为深入的研究。

近二十年来,国内外大量学者对泥石流运动过程及冲出规模开展的研究可分为三类:基于体积方量的预测模型、动态预测模型和基于地形参数的预测模型[5-9]。随着科技水平的提高,AI、数值计算等方法不断出现并更新迭代,其中泥石流数值计算是理论基础和计算数学的结合,将现实中的泥石流运动全过程简化为离散方程后进行求解,加之数值计算平台搭建完善后使用成本低且操作简捷,近年来被大量运用到泥石流研究中[10-13]。但如FLO-2D、FLO-3D、VolcFlow、PFC等数值计算模型均是以获得泥石流流深及流速为目的的单一冲出模拟,不能很好地耦合降雨、启动、侵蚀、运移、堆积整个动态过程,并且忽略了泥石流运动过程中尤为重要的侵蚀冲刷。

OpenLISEM数值计算模型是一个集降雨输入、径流、启动、侵蚀、运移、冲出、堆积于一体的综合水文模型,其原理是耦合分布式流域水文模型与二相动量守恒定律,特点是把流域划分成若干网格区域,对每个网格区域分别输入不同植被、土壤、岩土参数和高程等属性,输入实时降雨条件后分别计算产流量;通过比较相邻网格的高程确定各网格的流向,演算降雨过程中泥石流的启动、运移、堆积整个过程[14]。该模型不但考虑了地形、植被、物源和水文过程,同时也避免了以往常规模拟软件需要的泥石流流量过程线计算与输入、启动点位置选取等问题的不准确性,能更加真实地重现泥石流暴发时的场景。

受泸定6.8级地震影响,四川省泸定县落井沟流域新增大量物源,未来雨季,在泥石流激发降雨阈值降低的条件下研究区极有可能暴发泥石流,威胁沟口联合村隧道、在建的泸定—石棉高速与沟口居民生命财产安全。本文以落井沟为研究区,通过分析其物源分布特征,选取可以较好地模拟泥石流运动侵蚀全过程的OpenLISEM数值计算模型对不同降雨强度下的落井沟泥石流运动全过程进行模拟预测,为该区域的泥石流危险评价、预警预报及工程治理提供具有一定适用性和参考性的方法。

1 研究区概况

1.1 区域环境地质条件

落井沟位于四川省泸定县得妥镇中部,沟口坐标102°10′42″E,29°33′33″N,其形成区距2022年9月5日泸定6.8级地震震中约1.2 km,位处Ⅸ级烈度区(地理位置如图1)。流域较小,总面积约3.4 km2,地势复杂,东高西低,地形起伏范围为1 089～2 572 m,主沟长约2.7 km,流向自东南向西北,平均纵坡降约为390‰,沟谷狭长,两岸坡度较陡,为明显的“V”字型,属典型的高山峡谷区。流域内发育有两条小型支沟,于主沟中部交汇并最终流入与沟口正交且流向为东北至西南的大渡河。

图1 研究区地理位置Fig.1 Geographical location of the study area

图2 泸定县2022-09-05—2023-04-15地震概况Fig.2 General situation of earthquakes in Luding County during 2022-09-05—2023-04-15

图3 落井沟坡度分级与统计Fig.3 Slope grading and statistics in Luojing Gully

流域年平均气温15.7 ℃,年最低气温-4.5 ℃,年平均降雨量664.4 mm。从山巅至谷底可观察到明显的气候垂直分带现象,海拔低处气候较温暖湿润,高处较严寒干燥,属典型的亚热带季风气候。研究区属华南地层,形成区主要发育花岗岩、砂岩;流通区及堆积区主要发育二叠、三叠系变质岩,以板岩为主,同时整个流域发育较多沉积岩,第四系堆积物杂乱且分布不均。得妥乡位于川滇南北向构造带北段,为多组构造的交叉部位,构造环境复杂,新构造运动强烈,主要断裂为鲜水河断裂东南段(也称磨西断裂)及大渡河断裂,且大渡河断裂横穿落井沟堆积区[15]。

1.2 地震概述

2022年9月5日12时52分,四川省甘孜州泸定县发生6.8级地震,震源深度16 km,震中位于29.59°N,102.08°E,震源机制解以走滑型为主。对四川省地震局提供的震讯进行统计,发现从2022年9月5日6.8级地震截止到2023年4月15日,泸定县共发生32次地震,其中包括3次5级及以上地震,即2022年9月5日6.8级、2022年10月22日5.0级及2023年1月26日5.6级地震,为典型的“主震-余震型”地震。

经现场调查,“9·5”泸定地震直接造成落井沟内产生大量滑坡,沟道内堆积较多的松散堆积物,后续接连暴发的31次地震更使得沟道两侧山体内部损伤程度及沟道堆积物松散程度逐步增加。泥石流物源在震后得到大量补给,一旦暴发将对沟口联合村隧道及在建的泸定—石棉高速的使用和建造造成巨大影响,泥石流若发生堵江也可能对沟口东南部约11 km处正在运营的大岗山水电站造成难以估量的损失。

2 落井沟泥石流特征

2.1 地形特征

陡峻的沟道地形条件是形成泥石流至关重要的一环[16]。落井沟流域发育有两条小型支沟,主沟顶部至沟口平均纵坡降为390‰,地势陡峭,为明显的“V”字型,属高山峡谷深度切割区域。其沟谷狭长,两岸坡度较陡,统计发现近70%的流域面积分布在坡度为30°～50°范围内,滑坡松散物质运移到沟道的地形条件极好。这样的地形条件对沟道与坡面物源稳定程度及降水形成径流的水动力强度均具有正向的交互效应,利于泥石流发生。

2.2 物源特征

近年来对泥石流物源的识别技术愈加丰富,如通过遥感高精度影像进行滑坡的人工解译[17]与通过合成孔径雷达干涉测量(Interferometric Synthetic Aperture Radar,ISAR)技术进行潜在滑坡的隐患识别[18],这些识别技术无一例外都需通过现场复核以检验其适用性及准确性,而一些极为困难且危险的野外调查及测绘工作难以正常进行,因此无人机低空数字摄影野外测量技术开始广泛应用[19]。

泸定县6.8级地震后,通过遥感解译技术获取落井沟流域新增滑坡37处,面积约0.125 km2。震后新增滑坡位置分布见图4(b)。

图4 滑坡震前震后影像对比Fig.4 Comparison between landslide images before and after the earthquake

震后不久对其进行现场调查并使用搭载长测程激光雷达的无人机进行仿地飞行采集数据,通过大量点云数据处理等方式获得了分辨率为0.2 m的高精度机载激光雷达(Light Detection and Ranging,LiDAR)数据,包括数字高程模型(Digital Elevation Model,DEL)及数字正射影像(Digital Orthophoto,DO)数据,用于辅助识别落井沟流域内因植被覆盖未能识别的崩滑体,并作为基础数据进行后期泥石流运动全过程数值计算。通过图5中b1与b2对比,发现去植被后的LiDAR数字高程模型能识别出影像不易察觉的滑坡堆积物,同时,由于较多滑坡碎屑流仅是浅层物质剥落,故通过LiDAR航测处理得到的山体阴影数据进行精准识别,剔除未产生与产生极少物源的滑坡,得到实际滑坡堆积面积约7.798 万m2(初始滑坡分布见图1,处理后的滑坡分布见图5)。

图5 LiDAR数字正射影像及去植被的数字高程模型(图a1与b1分别是落井沟两处去植被后 的沟道与坡面滑坡LiDAR数据;图a2与b2分别是落井沟两处未去植被的LiDAR数字 正射影像数据;c为滑坡取样点)Fig.5 LiDAR digital orthophoto and digital elevation model for vegetation removal

结合滑坡体积-面积经验公式[20][式(1)]得到震后落井沟沟道及坡面新增的实际滑坡松散堆积物方量为15.52 万m3,且主要集中在主沟中下段(如图5)。沟道堆积物易被强降雨作用下形成的地表径流冲刷侵蚀,沟岸坡脚堆积物被冲刷带走后形成临空面又会不断失稳直至运移到沟道,增加沟道物源堆积,且随海拔下降物源侵蚀会不断增加,进而暴发泥石流灾害。

v=1.315*A1.201 8

(1)

式中:A为同震滑坡面积(m2);v为同震滑坡体积(m3)。

2.3 水文地质特征

震后大量崩滑体在强降雨作用下极易转化为泥石流[21]。泸定县年平均降雨量664.4 mm,雨季降雨充沛,最高月降雨量可达167 mm,受“V”型峡谷地貌影响,该流域汇水形成径流的能力较强,加之大渡河从沟口流经,地下水丰富,流域内主要为松散堆积层孔隙水和基岩裂隙水,补给方式通常为大气降水。

泥石流灾害的发生是地形、物源及降雨等多因素共同作用的结果[22]。通过对落井沟地形、物源及降雨三个影响因素分析发现,近70%的流域面积分布在坡度为30°～50°范围内,“9·5”泸定震后沟道及坡面的实际新增滑坡松散堆积物方量为15.52 万m3,且主要集中在主沟中下段,加之研究区雨季降雨充沛,地下水丰富,在未来的强降雨条件下具备暴发泥石流的初始条件,故对落井沟泥石流启动、物质运移及堆积规模借助数值计算工具进行深入分析。

3 泥石流运动演化分析

3.1 数值计算基本原理

OpenLISEM数值计算模型通过建立栅格单元,融合SWATRE渗透模型、理查兹方程、圣维南二维非恒定流方程、二相动量守恒定律及无限边坡理论,实现灾害启动、侵蚀堆积和控制流体的微分方程求解,以此模拟地表水沙平衡,其具有很明显的时空差异,是一个离散的二维数值模型。该模型可用于模拟流域在降雨期间或降雨后短时间内的物源变化、土地利用类型变化、防治工程等制约因素对泥石流的冲出趋势、物源侵蚀、洪水等动态流动趋势的影响。

3.2 数据获取

模拟所需数据主要包括物源数据、地形数据和降雨数据,其中地形数据利用QGIS平台的转换功能将0.2 m的机载LiDAR数字高程模型(DEM)转换成模拟所需map格式文件,其他地形数据利用PCraster 软件和Arcgis平台制作成map格式文件,包括流向网、地势梯度及表面粗糙度等。物源的物理力学特性指标包括渗透系数、初始含水率、内摩擦角、孔隙率、内部黏聚力、平均粒径及特征粒径D50和D90等,需通过野外调查及现场取样进行室内试验。降雨数据仅需输入泥石流暴发期间的降雨过程,但由于震后落井沟尚未进入雨季且未暴发泥石流,无法获取降雨启动泥石流后运动过程中的实际降雨历时数据,故结合5%、2%、1%降雨强度得到落井沟洪水流量过程线,输入到OpenLISEM数值计算模型中进行模拟。

(1) 物源数据

通过0.2 m分辨率的LiDAR航拍数据处理得到的落井沟沟道与坡面滑坡松散堆积物体积方量约15.52 万m3。此外,还需获取物源的土体级配特征D50及D90进行建模,通过筛除大于20 mm的现场样品并进行室内颗粒分析试验最终得到土体级配曲线,获取了D50及D90分别为5 mm、11 mm(图6)。取样点在流域中的位置分布见图5中c点。

图6 落井沟滑坡堆积物级配曲线(<20 mm)Fig.6 Grading curve of landslide deposit in Luojing Gully (<20 mm)

通过不同土体的初始含水率、孔隙度及渗透系数等参考值[23]加上现场调查与室内试验综合得到如表1所列的物源指标。

表1 物源特征指标及取值

(2) 洪峰流量

OpenLISEM进行泥石流数值计算的降雨数据不需要通过计算获取泥石流流量过程,仅需提供实际降雨过程便能结合物源及地形参数进行自动化模拟,得出最符合实际的泥石流流量过程。但由于落井沟尚未暴发泥石流,故采用四川省中小流域暴雨洪水计算手册进行计算,获得落井沟泥石流20年、50年、100年一遇暴雨重现期的洪水流量水文图。采用以下公式(2)～(9)进行计算:

(2)

(3)

(4)

(5)

m=0.221θ0.204

(6)

(7)

(8)

(9)

式中:QB为暴雨流量(m3/s);ψ为洪峰径流系数;i为最大平均暴雨强度(mm/h);F为流域面积(km2);s为暴雨雨力(mm/h);n为暴雨公式指数;τ为流域汇水时间(h)。K6、K24为不同重现期的模比系数,可从皮尔逊Ⅲ型曲线模比系数表中获取;θ为流域特征参数;m为汇流参数;τ0为当ψ=1时的流域汇流时间(h);μ为产流参数(mm/h)。相关取值列于表2,洪水流量过程线见图7。

表2 清水洪峰流量相关参数

图7 落井沟洪水流量过程线Fig.7 Flood flow hydrograph in Luojing Gully

3.3 模拟步骤、结果及分析

模拟过程需提前设置四个基础全英文命名的文件,包含地形和物源的map文件、模拟运行的命令流run文件、降雨数据rain文件与存储输出结果的res文件,设置完毕后对落井沟泥石流基于OpenLISEM模型进行运动全过程模拟。

结果显示:5%、2%、1%降雨强度下,落井沟泥石流运动过程的最大泥深分别为6、8、11 m,最大流速分别为5、7、10 m/s,整个运动过程的泥沙堆积显示沟道内淤积较多泥石流固体物质,推测受落井沟流域地形影响,于拐弯处发生弯道超高,降低了泥石流流速,导致泥石流淤积沟道的固体物质较多(图8)。泥沙随空间的变化关系呈现为逐步从主沟上端运移至沟口(图9)。

通过分析落井沟泥石流模拟的最终冲出规模发现:5%、2%降雨频率下落井沟泥石流虽冲出沟口但对大渡河的堵江程度较小,而1%降雨频率下落井沟泥石流造成大渡河半堵江状态。

4 结语

通过遥感解译、现场调查及机载LiDAR航拍,对2022年9月5日泸定6.8级地震后落井沟泥石流新增滑坡的发育特征进行分析,结合OpenLISEM数值计算模型对落井沟泥石流在5%、2%、1%降雨强度下的运动全过程进行模拟,得到运动过程中的最大堆积深度、最大流动速度及最终冲出范围,并给出1%降雨强度下泥石流的泥沙演化全过程,综合反映出落井沟泥石流的运动侵蚀堆积情况,并对大渡河的堵塞程度给出建议。主要结论或建议如下:

(1) 受泸定6.8级地震直接影响,流域新增滑坡37处,滑坡面积约0.125 km2。通过现场调查及0.2 m高精度分辨率机载LiDAR数据对滑坡进行优化识别,并结合滑坡体积-面积经验公式得到落井沟沟道及坡面的实际滑坡松散堆积物方量为15.52万m3,且主要集中在主沟中下段。

(2) 基于OpenLISEM数值模拟平台,结合精细化地形、精准物源及5%、2%、1%降雨强度对落井沟泥石流进行运动全过程模拟,获得运动过程最大泥深分别为6、8、11 m,最大流速分别为5、7、10 m/s。通过对沟口最终的堆积状况分析,5%、2%降雨频率下落井沟泥石流虽冲出沟口但堵江程度极小,1%降雨频率下落井沟泥石流造成大渡河半堵江状态。

(3) 落井沟泥石流沟域面积较小,即便在雨季暴发对沟口联合村隧道和在建的泸定—石棉高速的使用和建造在大渡河下游库岸相连的大岗山水电站的正常运营及沿线居民的生命财产安全造成的威胁仍然较小。但值得注意的是,自2022年9月5日泸定地震以来共发生大小地震32起,包括3次5级以上地震,这种小而多的地震对沟道及坡面滑坡物源的各方面影响(尤其是潜在滑坡临界失稳状态)有着肉眼不易察觉的特征,雨季丰沛的降雨及形成的径流对该流域的潜在物源也有着较大影响,对泸定震区此类物源丰富、地势陡峻的沟道仍需密切关注。