基于FLO-2D数值模拟的结底岗村公路泥石流危险性评价

刘福臻， 张浩韦， 肖东升

(西南石油大学土木工程与测绘学院, 成都 610500)

西藏自治区工布江达县因为其特殊的地理位置，新构造运动强烈，基岩风化侵蚀严重，公路沿线山体沟谷纵横，崩塌滑坡现象明显，物源储量较高，遇降雨、冻融等极端工况条件下，中小型泥石流频发。对于该地区公路沿线泥石流而言，普遍具有活动频率高、多点暴发、暴发规模以中小型居多的特性，暴发时冲击力强、破坏力巨大，沟道内的物源在水流的挟裹下短时间内大量冲出，对公路造成巨大威胁。

对公路泥石流危险性评价被认为是减少泥石流危害的有效途径[1]。目前中外学者对公路泥石流灾害已进行了一定程度的研究：邓恩松等[2]以中巴公路奥布段为研究对象，运用主成分分析法建立了公路沿线泥石流灾害危险性评价模型；陈维等[3]运用模糊综合评判法对公路泥石流的危险性进行半定量评价；尚彦军等[4]利用层次分析法结合地理信息系统(geographic information system，GIS)叠加分析评价奥布公路沿线泥石流危险性；王晨晖等[5]为预测泥石流危险性，提出了基于粗糙集理论的粒子群算法优化支持向量机模型。专家学者已采用不同的方法对公路泥石流的危险性进行了评价研究，这些评价方法多是基于某种数学方法采用半定量的分析方法建立评价模型，可以定性公路的危险性，但对公路的影响缺乏实际定量的描述，例如泥深、冲击力、堆积范围等方面，有进一步研究的必要。

根据研究范围，泥石流危险性评价可以分为点评价和面评价[6]。目前公路泥石流多从流域、区域等宏观角度进行评估灾害危险性，能够反映区域灾害的整体特征以及灾害的空间异质性，但受评价单元的限制，评价精度有限[6]。因此可从评价单沟泥石流危险性出发，展开到整体灾害危险性评估，以提高评价精度。

FLO-2D模型主要模拟洪水与泥石流在复杂的地形下的运动与堆积情况，对泥石流有较好的重现性，能明确地展示泥深、流速和反映流体流变时空特性[7]。目前中外已有学者运用FLO-2D模型对泥石流进行研究与评价：方群生等[8]以瓦窑沟为例，采用FLO-2D模拟典型急陡型泥石流在非工程和工程措施下的冲出量；Lin等[9]运用FLO-2D数值程序模拟了中国台湾松河地区泥石流的流动条件，并按不同的危险程度生成了风险分布图，均取得了较好的结果。

现以野外现场调查及光学遥感影像为数据基础，借鉴基于数值模拟评估泥石流危险性的方法和理念，模拟结底岗村公路泥石流在暴雨频率为1%、2%、5%、10% 4种不同降水工况下的堆积特征。选取50年一遇的模拟结果与泥石流沟实际冲淤特征比较，确定FLO-2D对于该地区的适用性；根据泥深、流速确定泥石流强度进而评价典型公路泥石流的危险性，以期为工布江达县公路沿线泥石流的调查评价提供参考。

1 公路泥石流概况

结底岗村泥石流位于工布江达县城北侧的结底岗村，沟口坐标：东经93°14′37.56″，北纬29°53′40.45″，属尼洋河左岸一级支沟。该沟位于强烈挤压、碰撞的冈底斯陆块、雅鲁藏布江结合带内，构造活动强烈；地处嘉黎-然乌、多其木-东久强震带附近，地震活动相对频繁，研究区地震基本烈度为Ⅶ度。

受地形控制径流由北西向南东，流域平面形态为桃叶状，两侧支流(沟)呈树枝状发育并向尼洋河排泄，整个沟谷为深切沟谷形态，中上游为“V”形槽谷，沟口为典型扇状堆积，该沟流域面积为7.64 km2，主沟长约5.97 km，相对高差1 950 m，上陡下缓，沟谷平均纵坡降326‰。沟域内两侧山高坡陡，平均坡度在49°以上，沟谷纵坡较大，特别是主沟上游段及两侧冲沟纵坡多在400‰以上，为泥石流启动提供了良好的沟道条件，泥石流概况图如图1所示。

泥石流堆积区沿沟口呈扇形展布，扇地完整性80%，扇顶至扇缘主轴纵坡降115‰。扇面发展趋势以堆积为主，扇体内冲沟深度较浅，主沟槽宽1.5～3.0 m，切割深1.0～2.5 m，沟道较顺直，堵塞程度轻微。一旦形成较大的洪水或泥石流，将沿扇体低洼处改道形成漫流状。沪聂线(G318)工布江达县城延线段从堆积扇前缘经过，处于泥石流的威胁范围内，为典型的公路泥石流。结底岗村位于泥石流堆积扇上方，工布江达县城分布于泥石流堆积扇前缘，承灾体较为复杂。

根据遥感解译以及实地调查，结底岗村泥石流的物源区主要位于沟谷上游的寒冻风化崩塌区、坡面侵蚀区和下游沟底的洪积物堆积区，物源类型主要分为沟岸8处崩塌物源、2处滑坡堆积体、6处坡面侵蚀物源和2处沟床堆积物源。通过对物源量进行分析估计，结底岗泥石流沟沟床堆积固体物源总量为1.95×104m3；崩滑堆积物物源总量为3.08×104m3；坡面侵蚀物源总量为1.25×104m3，松散固体物源量6.28×104m3，如图2所示。

结底岗泥石流为暴雨型泥石流，形成主要取决降雨的动力条件。根据资料显示，结底岗泥石流所在的工布江达县城多年平均降水量640.1 mm，日最大降水量45.2 mm，小时最大降水量16.5 mm。沟谷常年流水。枯水期流量为0.4 m3/s，地表水流量受降雨影响变化大，5—9月为工布江达的雨季，占全年降水量的85%[10]。

图1 结底岗村泥石流概况图Fig.1 Overview map of debris flow in Jedigang Village

图2 不同物源的储量与动储量Fig.2 Reserves and momentum reserves of different sources

据实地调查，该沟曾于1959年左右发生大型泥石流，大规模淤埋农田，冲毁房屋三十多间，最高泥位达4 m，沟口结底岗村受灾严重；近年于2013年8月暴发过较大规模的泥石流，冲出泥砂近6 000 m3，国道318(工布江达县城延线段)受损严重，造成几十多间民房不同程度地受损。

2 研究方法

2.1 FLO-2D理论模型

FLO-2D软件是1988年由O’Brien所提出的，将数字高程模型(digital elevation model，DEM)划分为相同大小和规则的地形格网，利用非牛顿流体与中央有限差分法求解泥石流运动的控制方程，以数值定量的方法来模拟泥石流冲出的流动过程、堆积范围以及评价危险区域。受理论模型限制，过程中需要满足以下限制和假设条件：①假定流体为静水压力分布；②假定格网内各参数(高程值、粗糙系数)保持一致；③假定为浅水波模式；④假定差分时间间隔内为固定恒流；⑤不考虑沟道侵蚀现象；⑥不考虑流动过程中的跳跃和震荡现象；⑦不考虑泥石流对于工程结构的损毁现象[11-13]。

通过运动方程以及连续方程，可以计算出每个网格中x方向上和y方向上的流体流速、堆积深度，进而得知流体的运动范围。通过连续方程控制泥石流的质量守恒，利用运动方程计算出相邻网格间流体的速度变化。模型方程如下。

(1)运动方程[11]。

(1)

式(1)中：i为降雨强度，mm/h；h为流体流深，m；t为泥石流流域降雨历时，h；λ、θ分别为x轴、y轴方向上的平均流速。

(2)连续方程[11]。

(2)

(3)

式中：sox、soy分别为x、y方向的泥石流沟床坡降；sfx、sfy分别为泥石流沟x、y方向上的摩擦坡降。

若流体内部固体物质浓度较大，则流体运动过程中，固体物质的颗粒间碰撞增多，流体的扩散应力变大，需要考虑泥石流运动时颗粒之间的碰撞对泥石流流动阻力的影响。

(3)流变方程[11]。

(4)

式(4)中：sf为摩擦坡降；sy为屈服坡降；sv为黏性坡降；std为紊流分散坡降；τy为屈服应力；ym为流体比重；k为层流阻力系数；h为泥深；η为流体黏滞系数；n为曼宁系数；v为流体流速[7]。

3 公路泥石流数值模拟

3.1 数据来源

根据所划定泥石流流域高分辨率的DEM作为基础数据，在ArcGIS中转化为FLO-2D可以识别的ASCⅡ码。依据泥石流沟的特征建立10 m×10 m的计算网格，并对网格进行高程插值。结底岗泥石流所需要的主要数据如表1所示。

表1 数据来源及数据类型Table 1 Data sources and data types

3.2 基础参数设定

据野外地质调查，结底岗村以往暴发的泥石流主要为水石流，现场配制泥石流浆体难度较大，本次实验泥石流流体重度采用查表法所得泥石流重度为1.641 t/m3，为稀性泥石流。

泥石流是一种固液两相体，通常泥石流体夹杂着大量的固体，且沿着深度方向和前进方向呈不均匀分布，固体物质在泥石流中所占比例很大程度上决定着泥石流的流变特征，同时也对泥石流堆积体的前缘形态有一定程度的影响，通过式(5)计算泥石流体积浓度[14-15]。

(5)

式(5)中：γc为泥石流重度，取1.641 t/m3；γw为水的重度，取1.0 t/m3；γh为泥石流固体物质的重度，取2.65 t/m3,可得泥石流浓度Cv为0.388。

其他模拟参数如屈服应力与黏滞系数(α1、β1、α2、β2)、曼宁系数(n) 和层流阻力系数(K) 依据FLO-2D手册[11]的经验数值结合结底岗泥石流冲出模拟实验进行修正，最终参数取值见表2。

工布江达县中小流域暴雨洪水的特征与四川西南地区相似，故可以参照《四川省中小流域暴雨洪水计算手册》，以工布江达县1970—2020年的降雨数据作为基础，按手册给出的建议公式计算出不同降雨频率下的暴雨洪峰流量(表3)。一般的泥石流流量过程线为单峰型涨落曲线，所以通过五边形法确定流域清水流量过程曲线。泥石流运动过程中会产生一定的放大效应，故最终输入的泥石流流量过程数值为泥石流清水流量乘以体积膨胀系数BF=1/(1-Cv)(图3)[16-17]。

表2 泥石流模拟参数Table 2 Debris flow simulation parameters

表3 暴雨洪峰流量Table 3 Rainstorm peak flow

图3 泥石流流量过程曲线Fig.3 Debris flow flow process curve

在运用FLO-2D进行模拟泥石流时，集水点的选取至关重要，大部分的泥石流起动点在松散物源聚集处、在沟道物源堆积处和水动力条件充足的位置[18]。综合考虑结底岗泥石流水流条件和物源条件，将泥石流启动点大致选在形成区和流通区交界点附近。再根据遥感解译的大规模崩塌滑坡体，结合实地踏勘泥石流沟的灾害发育特点，最后确定主支沟交汇处下游150 m处作为本次实验的集水点。

3.3 模拟结果

为进一步了解结底岗泥石流对于典型公路的冲淤特征，通过对不同降雨频率下的泥石流冲出过程的模拟，模拟泥石流在沟道的运动和流动过程，获取冲出后堆积扇的流深、泥深、淤积分布，模拟结果见图4和图5[19]。同时总结泥深、流速、冲出范围、冲出量以及与典型公路的拓扑关系，如表4所示，分析公路泥石流对于公路影响。

根据模拟结果显示，结底岗泥石流冲出过程中，受沟道影响，最大堆积厚度主要集中在沟道内转角以及沟口处，最大流速位于沟道内地形急剧变化处。堆积区地形东高西低，受到地形影响，数值模拟所形成的堆积区并不是狭义上的扇形。

结底岗泥石流在降雨频率P=10%时，降雨强度较小，水动力条件较差，地表径流难以汇聚，绝大部分物源仍停留在沟道内，未能形成有效的冲出，对典型公路未形成影响；在降雨频率为5%、2%、1%的情况下，均对公路造成了不同程度的影响，如表4所示；泥石流堆积区的泥深、流速、堆积面积以及冲出量与降雨频率呈正相关。

3.4 模拟情况验证

为了检验FLO-2D模型对于高海拔高寒地区泥石流沟的适用性，本文以结底岗村泥石流50年一遇降雨工况下的模拟计算结果为例，将模拟泥石流冲出的横向最大距离、纵向最大距离和堆积扇面积与1959年暴发的大型泥石流实际堆积情况进行分析比较，结果如表5所示。

根据王高峰等[20]提出的精度系数Ia来验证模拟的合理性：以泥石流模拟堆积扇面积Am、泥石流实际堆积扇面积Ar与重叠面积A0作为数据基础，利用式(6)计算模拟结底岗村泥石流精度，精度系数越高，越能体现实验的合理性。通过公式结合表5数据计算可得FLO-2D模拟结底岗泥石流的精度为86%，满足模拟精度要求，模拟结果与实际情况基本相符。

(6)

图4 不同降雨频率下公路泥石流泥深Fig.4 Depth of mud-rock flow on highway under different rainfall frequency

表5 模拟结果与实际结果对照表Table 5 Comparison table of simulated results and actual results

3.5 模拟情况统计

根据4个不同重现期的泥石流模拟结果，统计泥石流在威胁公路距离、公路路面平均堆积泥深和平均流速3个方面的情况，确定不同重现期泥石流对公路的影响程度，如图6所示。

泥石流对公路的影响与降雨工况(暴发频率)呈正相关，随着降雨的增加，泥石流威胁公路的距离、公路路面平均泥深、平均流速均有不同程度上升。

4 危险性评价

目前泥石流危险性分级标准的影响因子主要有泥深、流速以及暴雨暴发的频率和强度等。通过泥深(h)以及流速与泥深乘积结果(vh)划分泥石流强度(表6)，结合不同暴雨重现周期按泥石流冲出的区域分为高危险区、中危险区、低危险区3种类型，运用ArcGIS对各危险区进行赋值，构建泥石流危险性逻辑关系(图7)，并绘制结底岗泥石流不同分级的危险区划图(图8)[21]。

根据结底岗泥石流危险区划图结果显示：受泥石流特征以及古堆积区地形影响，古堆积扇东南区域以及堆积扇边缘为低危险区，由堆积边缘依次向内侧沟口处危险性逐渐增加；国道318(工布江达县城延线)与危险性分区为相交拓扑关系，共威胁公路335 m，高危险性路段占公路受威胁距离的6%，为20 m，中危险性占比9%，为30 m，低危险性占比85%。

表6 结底岗村公路泥石流强度划分表Table 6 Debris flow intensity division table of Jiedigang Village highway

图6 不同重现期下泥石流对公路影响Fig.6 The impact of debris flow on highways in different return periods

图7 公路泥石流危险性评定标准Fig.7 Village highway debris flow risk assessment standard

图8 结底岗村公路泥石流危险性分区图Fig.8 The dangerous zoning map of the debris flow on the highway in Jiedigang Village

5 结论

(1)模拟了结底岗村公路泥石流在4种不同降雨工况下(10%、5%、2%、1%)的堆积特征。结果显示：泥石流泥深、流速随着降雨量的增大而增大；在暴雨频率为10%时，公路泥石流堆积面积为3.84×104m2，未对公路造成影响；暴雨频率为1%的堆积面积为16.89×104m2，是暴雨频率为2%的1.1倍，暴雨频率为5%的2.1倍，均对公路造成不同程度的影响。

(2)对FLO-2D模型模拟高寒、高海拔以及地质构造强烈地区公路泥石流可行性进行了检验。检验结果表明：典型公路泥石流泥深、流速及堆积面积冲淤特征与泥石流实际冲出相符，模拟精度为0.86，证明了该模拟方法对该地区的适用性。

(3)通过数值模拟得到了10、20、50、100年4个重现期的泥石流强度，并根据灾害等级分类绘制结底岗沟泥石流公路危险性分布图。结果表明：随着重现期的增加，泥石流强度逐步上升；结底岗公路泥石流典型公路范围内属于高危险性区域占6%，中危险性区域占9%，低危险区域占85%；公路泥石流高危险性区域主要分布受影响路段中部；根据该危险分布图可以确定威胁公路具体位置，能够为公路泥石流防灾预案以及紧急避难点提供有效的参考。