沟道松散物质起动模型试验及冲出量预测——四川省以北川青林沟为例

林 斌，张 友 谊，罗 珂，黄 成，袁 亚 东

(1.西南科技大学 土木工程与建筑学院，四川 绵阳 621010； 2.武汉市市政建设集团有限公司，湖北 武汉 430023)

泥石流冲出固体物质是泥石流灾害中直接作用于受灾体的物质，同时该冲出量的大小直接反映出泥石流灾害的危害性，也是防灾减灾工程设计中的重要参数。对于泥石流的冲出规模，国内外相关学者已经做了大量的研究工作，并取得很多的研究成果，目前大多数学者采用统计学经验公式和现场调查的方式开展研究。

国外学者Hampel就将泥石流的流域面积和堆积扇的坡度作为推导泥石流活动规模预测模型的影响因子[1]；Inners选取单因子集水区面积为参量建立了泥石流冲出体积的经验公式；Ikeya,MizuYama选取泥石流的沟道长度和宽度为参数建立了预测泥石流最大冲出量模型[2-3]；池谷浩于1979年建立了流域面积与泥石流冲出量的经验统计公式[4]；Takei建立了泥石流冲出量与流域面积的经验统计公式[5]。Bovis,Jakob选取与泥石流活动规模密切相关的流域内的物源、流域高差作为参数[6]，建立了泥石流活动规模的预测模型；美国学者Joseph E统计了科罗拉州大量的泥石流沟事件[7]，结合流域面积、流域内滑坡面积、累积降雨量3个因子，对泥石流的内外因做了充分的考虑，并得出了冲出量与三者的关系。

中国学者对此也有相关研究：王继康经过研究得出泥石流的冲出规模为暴发时间和流量乘积的19.72倍[8]；Ma Chao通过对汉川震区泥石流分析得出泥石流冲出体积与松散物源量的数学关系[9]，并且冲出量与流域面积的相关性大于台湾、美国的火烧云地区、意大利的阿尔卑斯山的东部地区；我国台湾地区游繁结等人研究陈有兰溪的泥石流发现冲出规模与集水区面积成正比；台湾学者Chang C W选取59个小流域的地形参数[10]、地质指数和降雨因素来建立研究台湾地区极端降雨条件下暴雨泥石流的冲出量暴的预警模型；田运涛等考虑5种因素来评价各泥石流沟的冲出量及其影响[11]；翟淑花等分析了松散物源堆积体在不同含水率下的强度衰减特征,获得了含水率与物源堆积体强度参数之间的相关关系[12]，建立了松散物源堆积体失稳时临界含水率模型，为目标区域泥石流启动辅助预警提供技术支持；宋兵等运用RAMMS预测泥石流在不同洪水频率下的冲出量[13]；方群生考虑总物源量、动储量总量、流域面积、沟道纵向长度、相对高差等自然因素和人为因素建立适用于汶川震区暴雨沟谷型泥石流一次性冲出量的预测的模型[14-15]；王纳纳采用线性回归分析27条泥石流沟得到泥石流冲出量与影响因子的数学计算模型[16]；蒋志林通过物源的滑坡面积、流域面积、流域高差3个因子建立模型，预测研究区泥石流一次最大冲出量[17]。

综上研究成果，不难看出对于泥石流冲出量的研究都是针对泥石流暴发时的某一方面进行研究，并没有对相关因素各方面综合考虑，它虽为泥石流活动规模预测指明一个新的方向，但并不全面。同时，研究成果并没能够在泥石流爆发之前给予现场居民一定的预见性，并且部分统计成果与现场的调查结果相差较大，难以体现其准确性。因此本文以北川青林沟为例，考虑不同沟床粗糙程度，沟道坡度和降雨强度下的沟道松散物质起动，并针对起动完成之后形成3处具有代表性的沟道宽度与沟道坡度、降雨和沟床糙率进行函数拟合，最终得到可以利用当地降雨预报即可预测青林沟泥石流一次冲出固体物总量的模型。

1 室内模型实验设计

1.1 相似比选取

(1)

(2)

根据上关系式可以得到本次模拟实验中降雨强度相似比为1∶10，但是由于现实条件的限制，无法严格按照相似比将实际的降雨强度缩小，但是可以在计算汇流时，根据公式(1)计算值进行缩小。

1.2 试验模型参数选取

土体的宽度：根据现有的四川省地质工程勘查院提供的《北川县陈家坝乡场镇地质灾害勘查》中工程地质平面图及其他相关勘察资料，同时考虑室内模型缩尺后的尺寸，选取青林沟高程为714.5 m处的沟道横断面作为室内模型参考原始横断面，根据该平面图得到该处横断面(D横原型)为50 m，即：

(3)

计算得到模型中横断面尺寸为D横模型=0.5 m。

土层厚度：按照游勇提出的泥石流最佳排导矩形断面为土体厚度与断面宽度比为1∶2，因此试样的厚度选择为25 cm[20]；土样堆积及孔隙水压力计布置见图1。

图1 土样堆积及孔隙水压力计布置(单位：cm)Fig.1 Soil sample accumulation and pore water pressure gauge layout

土体材料：试验土体采用参考实验原型青林沟准泥石流体进行室内配制，先对试样进行颗粒筛分待用。为了减小尺寸效应的影响，筛选粒径小于2.5 cm的颗粒作为实验材料。根据(SL237-1999)《土工试验规程》提供的剔除法，计算并得到剔除后颗粒级配曲线，如图2所示。

图2 试验土样级配曲线 Fig.2 Test sample grading curve

1.3 试验参数的选取

结合前人的研究成果及青林沟调查结果[21]，选取主要的影响因子为沟道坡度、降雨强度和沟道糙率。根据实际情况，主要实验参数取值如下。

(1) 沟道坡度。现有的勘察资料显示，青林沟松散堆积体处坡度为8.2°，因此综合勘察资料及庄建琦、崔鹏的研究成果[22]，确定沟道内准泥石流堆积体的坡度为8°，12°，16°，20°。

(2) 降雨强度。由于理论的实验受场地、降雨喷头规格、大风天气等不利因素的影响，无法严格按照相似比来将雨强进行缩小，但是可以根据青林沟现有的降雨资料以及前人研究的激发泥石流的降雨强度，设计4种降雨强度：30，50，70，90 mm/h。

径流设计：由于现实过程中因降雨产生的地表径流对泥石流的起动的影响较大，因此，根据缩尺模型的比例，以及现有的青林沟的尺寸对地表的径流进行相似原理计算，获得上述4种降雨强度下地表径流量。

其中径流量与雨强的关系一般通过经验公式计算得出。根据铁道部第二勘测设计院的研究成果(TBJ107-1992)《铁路小桥涵设计》中关于西南地区小流域暴雨计算公式，可以计算：

Ωm=0.278×F×C×am×ym

(4)

式中，Ωm为暴雨洪峰流量，m3/s；C为产流系数系数，依据《铁路小桥涵设计》中表1查得为0.5；am为设计频率的最大小时暴雨强度,mm/h，分别为30，50，70，90 mm/h；F为流域面积,km2，其值为23.7；ym以设计频率的最大小时暴雨强度,mm/h，分别为30，50，70，90 mm/h。

计算得到，降雨强度分别为30，50，70，90 mm/h对应全流域径流量为57.2，92.7，127.2，160.8 m3/s。同时根据缩尺比得到即可算的分别为206,33,458,579 L/h。

根据：

γ=0.36×am0.4τ

(5)

(6)

式中，τ为流域最远点到达的径流运动时间，h；L为流域分水岭沿流程只计算横断面处之沟道长度,km，取值为9.94；A1为阻力系数，从《铁路小桥涵设计》中表2取得，该值为1.5；I为流域平均坡度(小数)，其值为0.105。计算得到降雨强度分别为30，50，70，90mm/h对应时间分别为4 080，3 609，3 320，3 134 s，根据公式(2)得到模型试验中，对应时间分别为408，361,332,313 s。

(3) 沟底糙率。根据已有的西南科技大学结构试验中心的橡胶防滑垫，且已知4种规格的防滑垫的糙率分别为0.3,0.5,0.6,0.7。

1.4 实验设计方法

本试验涉及三因素四水平，因此采用18组试验。具体试验参数设计见表1～2。

2 室内模型试验

本次实验系统包括高速摄影系统、孔隙水压力测定系统、人工模拟降雨、沟槽系统等5部分组成，而沟槽系统包括地表径流系统、门式升降机和冲出物收集箱。实验装置见图3。

表1 模型试验参数对照

Tab.1 Test parameters

沟道坡度/(°)降雨强度/(mm·h-1)糙率汇流时间汇流流量 (L·h-1)8300.36'48″20612500.5 6'1″33316700.65'23″45820900.75'1″579

表2 实验设计方案

Tab.2 Orthogonal experimental design schemes

第1组 沟道坡度/(°) 降雨强度/(mm·h-1)沟底糙率 第2组 沟道坡度/(°) 降雨强度/(mm·h-1)沟底糙率8300.312500.58500.512700.68700.612900.712300.516500.612500.616700.712700.316900.516300.620500.716500.320700.516700.520900.6

图3 实验装置示意Fig.3 Schematic diagram of the experimental device

高速摄影系统。型号为日本产MEMRECAMH- X-4高速摄像机。

孔隙水压力测定系统。采用由南京丹陌电子科技有限公司设计、生产的应变式孔隙水压力计和由四川德翔科创仪器有限公司生产的DEX3811型号应变箱。

人工降雨系统。流量表、喷头、水管、和支架等部分组成。

模型槽。长2.5 m，宽0.5 m，深0.5 m，两侧采用钢化玻璃制作，骨架采用钢材焊接而成，并在玻璃上利用颜料画纵向、横向间距为5 cm的线，以便观察土样变化。

地表径流系统。钢质水槽，鉴于水管出水具有初始能量，出水口距离水槽底部具有一定高度，且进水口在底部附近，水槽进水不会对出水口造成影响，因此能够减少原始水管出水的初始能量对泥石流沟道堆积体起动的影响。

门式升降机。于门式钢架上安装1.2 t起重葫芦。

3 准泥石流堆积体起动冲出物预测

根据表2进行试验，每次试验后对沟道松散物质堆积体的前缘、中部、后缘宽度进行测量，测量部位见图4。各部位测量数据见表3。

图4 测量部位示意 Fig.4 Measuring point

Tab.3 Measuring result of each part

项目沟道坡度/(°)降雨强度/(mm·h-1)沟底糙率前缘宽度/cm中部宽度/cm后缘宽度/cm第1组8300.340.623.416.38500.533.015.79.58700.641.219.717.712300.522.516.015.212500.629.418.013.012700.319.012.68.016300.631.522.613.116500.320.013.07.016700.527.016.08.7第2组12500.524.013.59.512700.628.616.812.912900.743.331.627.516500.630.017.410.916700.745.036.825.516900.522.013.47.220500.736.029.623.920700.521.013.45.920900.628.517.49.4

3.1 模型前缘、后缘、中部形成沟道宽度确定

现认为前缘及后缘沟道宽度(DQ模型、DH模型)与沟道坡度(S)、降雨强度(Q)及糙率(n)有关，即宽度与之存在某函数关系D=g(S，Q，n)；假设中部某部位处存在沟道宽度为DM模型，见图5。

假定：沟道宽度为DM模型的部位沟道前、后宽度变化率较大，即：(DM模型-DH模型)/L2与(DQ模型-DM模型)/L1相比较而言，(DQ模型-DM模型)/L1相对变化较大，且DM模型为堆积体坡度转折点处。

图5 冲沟宽度示意 Fig.5 Schematic diagram of the gully width

图中,DQ,DM,DH分别为堆积体形成冲沟的前缘宽度，中部宽度，后部宽度；H为堆积体厚度；L1，L2分别是堆积体前缘纵向长度和后缘纵向长度。

根据该地区现有的气象水文资料，结合沟道坡度范围、降雨强度范围定义以下参数：相对沟道坡度SR=S/14°；相对降雨强度QR=Q/60 mm/h。

现拟公式D=g(S，Q，n)=A×(SR)B+C×(QR)D+E×(n)F+G,其中：A，B，C，D，E，F，G为不确定常数，A，B为与沟道坡度有关常数，C，D为与降雨强度有关的常数，E，F为与糙率有关的常数。先通过Origin Pro软件对该表格中前缘宽度、后缘宽度进行拟合，得到如下关系式。

前缘宽度DQ模型：

(7)

通过回归分析结果显示，复相关系数R2=0.904，方程回归效果较好。

后缘宽度DH模型：

(8)

通过回归分析结果显示，复相关系数R2=0.972，方程回归效果较好。

现利函数关系式(7)～(8)，认为中间某部位部沟道宽度DM模型=f(DQ，DH模型)=K×DQ模型+I×DH模型+J，其中，K,I,J为常数。经过Origin Pro软件对该表格中前缘宽度、后缘宽度进行拟合，得到如下关系式：

DM模型=0.222×DQ模型+0.762×DH模型+2.380

(9)

式中，DQ模型,DH模型,DM模型单位均为cm。

3.2 野外原型冲沟体积的确定

由于根据模型是缩尺模型，有：

(10)

现将函数关系式(7)～(9)代入公式(10)中可以得到：

前缘宽度DQ模型

(11)

(12)

DM原型=22.2×DQ+76.2×DH+238

(13)

式中，DQ原型，DH原型，DM原型单位均为cm。

现利用图5冲沟示意图可以计算得到：

(14)

式中，V前，V后，V原型单位均为cm3。现将公式(11)～(13)中单位转化为m，即可得到：

(15)

(16)

DM模型=0.222×DQ模型+0.762×DH模型+2.380

(17)

式中，DQ原型，DH原型，DQM原型单位均为m。

利用图5即可得到：

(18)

式中，V原型单位为m3。

3.3 验证模型

现以青林沟“9·24”泥石流作为验证标准，分别采用以下3种方法验证。

(1) 现场调查。据现场调查，该泥石流沟堆积区与物源及流通区没有严格的分界线，“9·24”泥石流冲出物堆积分布范围高程约在660～750 m之间，沿沟长约1.85 km，上游呈现沟道条状，下游至沟口成扇区，最宽处690 m，堆积面积15万m2，堆积体厚3～30 m不等，由上游到下游逐渐变薄，平均厚4.7 m，堆积体积70.5万m3。

(2) 根据现有(DZT 0220-2006)《泥石流灾害防治工程勘测规范》推荐计算以下参数。

一次冲出泥石流流体总量：通过断面的一次泥石流总量可根据泥石流历时T和断面峰值流量Qc按下式进行计算：

Q=KTQc

(19)

式中，Q为一次泥石流流体总量，m3；T为泥石流历时，s，按“9·24”爆发的泥石流所经历的时间计，2.5 h，即为9 000 s；Qc为通过计算断面的最大流量，m3/s；K取值0.083 5 (取自《泥石流灾害防治工程勘查规范》1.1.3条，用内插法取得)。

一次冲出泥石流固体物质总量：

QH=Q(γc-γw)/(γH-γw)

(20)

式中，γc为泥石流流体重度，t/m3，取值为1.71；γw为清水重度，t/m3，取值为1.0；γH为固体颗粒重度，t/m3，取值为2.57。计算结果为：K为0.083 5，T为30 600 s，Qc为697.09 m3/s，Q为52.39万m3，QH为58.40万m3。

(3) 利用章节3.2中野外原型冲沟体积计算公式计算。根据现场调查，列出如下数据：“9.24”泥石流降雨强度为20 a一遇，降雨强度为78 mm/h；现场泥石流流体容重测试，泥石流平均容重值为1.71 t/m3，属于黏性泥石流，同时查阅(DZT0220-2006)《泥石流灾害防治工程勘测规范》表Ⅰ-5，可以得到青林沟堆积区内糙率为0.125；青林沟松散堆积体区域由堰塞湖至下游沟口处长度为1 850 m，其中前缘长度L1为150 m，后缘长度L2为1 700m(自上游堰塞湖至沟口处长度)；沟道平均坡度为8.2°，沟床糙率为0.125；将上述数据代入公式(17～20)得到V原型=54万m3。

现将3种方法的预见性、准确性、便捷性作对比，见表4。由表4可以看出，青林沟一次冲出物的预测数值比现场实测及规范推荐方法偏小，相对误差为-23.4%，该误差范围与方群生(-43.05%～-11.12%)[14-15]、王纳纳(-49.21%～400.17%)[16]、蒋志林(-99%～28%)等人的研究模型相比[17]，其误差范围可以接受。该方法利用地区气象局对降雨预报，再进行简单的现场调查，包括沟道松散堆积体长和沟床坡降，再根据沟床状况查询泥石流灾害防治工程勘测规范中沟床糙率n，即可计算出泥石流冲出量，提前作好防灾、减灾准备。

同时，根据表格中不难发现，现场调查与规范法和本文中所提模型计算结果也有一定的差异。由于2008年“5·12”地震造成青林村滑坡，该滑坡体前缘将青林沟中部堵塞，大部分堆积于沟道的松散物质在短时间内随着地表径流被带走，少部分堆积体在沟道内形成堰塞体。当“9·24”暴雨到来时，洪水对堰塞体冲刷致使堰塞体溃决时，形成溃决型泥石流。根据成都理工大学地质灾害防治与地质环境保护国家重点实验室研究成果，溃决型泥石流冲出量较一般型泥石流要大出几倍甚至几十倍[15]。

4 结 语

(1)通过对沟道松散物质起动室内模型试验完成后形成沟道的前缘、中部及后缘宽度进行测量，并通过OriginPro软件对其数据进行分析、拟合，得到所形成的沟道前缘、中部、后缘的宽度的函数关系式，并通过继续简化可以计算沟道松散物质起动过程固体物质冲出量。计算结果可以提前预测泥石流起动过程冲出规模及其危害程度，并提前对沟口周边居民进行疏散处置，并可以根据起动过程冲出量预测对后续相关防治工程进行设计指导。

(2)本研究提出的泥石流冲出物预测模型较其他模型准确度较高，并且具有预见性、便捷性。