基于分形理论和模型试验的沟道物源动储量评价模型

2022-10-25张友谊王云骏袁亚东

中国地质灾害与防治学报 2022年5期

张友谊，王云骏，袁亚东

（1.西南科技大学土木工程与建筑学院，四川绵阳 621010；2.四川正天创意建筑设计有限公司，四川成都 610041）

0 引言

2008年汶川地震以来，各类物源参与泥石流活动的贡献度随着时间的推移而改变，崩滑物源的供给不断衰减，预计将在2025年达到震前水平[1]。但由于沟道中仍有足够的物源储量，泥石流活动将持续更长的时间——沟道物源逐渐成为震后泥石流暴发的主要物源之一[2]。

传统的泥石流物源统计方法是根据实地调查来估算河道网络中储存的物源量[3]。近年来，通过数学统计方法，有学者[2,4]选取震区多条泥石流的物源及各类物源样本数据建立回归模型，得到总的物源量与动储量具有线性关系。面对极端降雨情况，台湾地区的学者[5]从地貌参数、地质指标和降雨因素出发，构建多元回归模型，以此估算泥石流流域的总物源量。针对沟道物源，亦有学者[6-7]通过人工降雨模型试验来研究其起动模式，并以降雨强度为主要影响因子预测其一次冲出量。

在实际的工程中，工程设计所需的泥石流动储量仍采用定性分析、比例统计、形态判别、稳定分析等统计方法。这些方法没有统一的标准衡量，往往不同的方法计算结果差异较大。目前，对泥石流动储量的研究主要是从地形地貌等方面出发，室内物模型试验方法的研究相对缺乏。

沟道内的松散物质是流域侵蚀演化的产物[8-9]，其动储量不仅与地形地貌和降雨因素有关，粒径组成亦是影响其起动量大小的关键因素之一，但土体复杂的粒度组成导致其缺乏相适的定量参数来进行分析，而分形理论为此提供了新的理论方法[10]。

因此，通过室内物理模型试验，研究沟道内不同分维的堆积体在各频率降雨条件下的侵蚀规律，构建雨强和分维值为双影响因素的动储量评价模型，可为震后泥石流治理工程计算沟道动储量提供新思路。

1 分形理论及分维值

近年来，有国内学者提出应加强针对固体物质岩土特性随时间的动态演变规律的研究[11]。因此，利用分形理论中具有自相似性的线性分形，局部按照整体成比例缩小来定量描述沟道物源的粒度演化，这种自相似性的参数即是分维值。

在土体工程分类实例计算中[12]，分维值作为一个综合定量指标的特征参数，可反映各种粒径大小在土体中的组成特点。

土体工程分类中，设土的粒径为z，N(z)为大于粒径z的数量，有：

式中：D——分维；

z′——粒径自变量；

P(z′)——粒径z的分布密度函数。

粒径组成的分维亦可用对应的质量关系来讨论，设m(z)为小于粒径z的质量，M为样本土体的总质量，有：

则：

dM∝za-1dz

对式（1）求导得：

dN∝z-D-1dz

由于粒径的增加与质量的增加是相对应的，则有：

dM∝z3dN

也即：

za-1dz∝z3·z-D-1dz

即分维值：

2 试验方案

2.1 模型相似比

泥石流室内模型试验的动力特征、几何形态以及边界条件等相似条件都与泥石流的物质组成、流变特性息息相关，根据其运动特性，主要引用以下相似理论。

（1）相似比例

将七盘沟主沟下游作为试验模型原型，结合试验场地条件和经济性，确定缩尺比例为：

（2）初始条件

人工逐层夯实堆土并静置1d，使其初步固结，并达到野外现场采样所测5%含水率，以确保试验土体与野外相符。

（3）降雨条件

由相似原理的第二定理推导得到降雨时长和降雨量的缩尺比例为：

2.2 试验模型几何参数

研究段沟道平均坡降为176 ‰，将堆积体坡度定为10°，按照相似比，试验模型宽度为1 m，土体堆积长度为2.3 m。将堆积体分为前缘、中部、后缘三部分以便对其侵蚀规律分析。按照泥石流排导最佳矩形断面尺寸为土体厚度与断面宽度比为1∶4～1∶2[13]，试验选取土体厚度为0.3 m。对于双层堆积体，一般表层较薄，因此将土层分为0.1 m 的表层和0.2 m 的底层，如图1所示。

图1 堆积体示意图Fig.1 Schematic diagram of deposit

2.3 试验影响因子

（1）降雨

降雨强度按照《四川省中小流域暴雨洪水计算手册》求得研究区各雨频下的小时降雨量。再结合铁道部第二勘测设计院的《铁路小桥涵设计》中关于西南地区小流域暴雨计算公式求得径流流量，并按模型相似比对径流流量进行1∶106折减，得到相应试验径流量，降雨用时按式（5）折减[14]，最终得表1。

表1 模型试验降雨条件Table 1 Rainfall conditions of the model test

（2）粒径组成及堆叠

汶川地震后七盘沟多次爆发泥石流，其中2013年“7·10”和2019年“8·20”两次规模较大，经资料收集及现场筛分试验（图2），得到“7·10”（表2）及“8·20”（表3）泥石流后沟道物源颗粒累积百分含量。

图2 七盘沟沟道物源野外筛分试验Fig.2 Field screening test of materials source in Qipan gulley

表2 “7·10”后沟道物源颗粒累积百分含量Table 2 Cumulative percentage of particles in the channel after “7·10”

表3 “8·20”后沟道物源颗粒累积百分含量Table 3 Cumulative percentage of particles in the channel after “8·20”

试验土体使用以上两种土体为原型土，按其分维大小相应分为粗粒土和细粒土。鉴于试验水槽及仪器设备的尺寸限制，试验土体粒径最大界限为20 mm[15]，再按土工试验等量替代法处理得到试验土级配曲线（图3—4）。

图3 “7·10”试验土级配曲线Fig.3 “7·10” test soil grading curve

由于当前沟道物源呈明显的分层堆叠形式，试验中堆积体按粗粒土、细粒土、上细下粗型及上粗下细型四种组成形式考虑。通过分维理论对堆积体粒度成分特征进行指标界限值的分类，得表4。

图4 “8·20”试验土级配曲线Fig.4 “8·20” test soil grading curve

表4 试验土分维值Table 4 Fractal dimension of test soil

2.4 试验设计

试验设计为两因素四水平单因素对照试验，试验装置包括高速摄影系统（型号为日本产MEMRECAMHX-4 高速摄像机）、人工降雨系统（流量表、喷头、水管和支架等部分组成，流量表型号：SJLWGYB10L）、模型槽系统（长3 m，宽1 m，深1 m，两侧采用钢化玻璃制作，骨架采用钢材焊接而成，并在玻璃用有色胶带纵向划分堆积体的分层界限，横向划分堆积体的范围）等3 部分组成，而沟槽系统包括地表径流系统和冲出物收集箱，如图5所示。共进行16 次试验，具体方案见表5。

图5 试验布置Fig.5 Test layout

表5 单因素对照试验设计方案Table 5 Single factor-controlled trial design scheme

3 侵蚀分析

径流形成之前，前期降雨使堆积体土体趋近于饱和，黏聚力与内摩擦角的减小使得土体抗剪强度减弱，为沟道物源的的起动提供了条件。

降雨越大，沟道径流的流量、流深以及速度等参数越大，相应的流体的侵蚀能力、侵蚀范围以及挟沙能力均呈正相关变化[16-18]，而堆积体组成粒径大小及组成也不同程度地影响着侵蚀能力。

3.1 单层堆积体

（1）细粒土

由于颗粒粒径较小易起动，雨强较小时（33.2～38.1 mm/h）前缘陡坡段因流速增大发生下切，随后演变成溯源侵蚀；雨强较大时（44.4～49.1 mm/h）中部缓坡形成多条径流展开侵蚀，并逐渐汇聚，前缘陡坡段以下切侵蚀为主，随着水力梯度不断增大，逐渐发生溯源侵蚀，形成下切→陡坎→溯源的循环侵蚀模式（图6）。

图6 细粒土49.1 mm/h 试验过程现象Fig.6 Test process phenomenon in 49.1 mm/h (fine-grained soil)

侵蚀深度受该模式的影响其峰值出现在距后缘0.5 m 的缓坡段，深度约为0.25 m。边坡随着侧缘侵蚀（以及掏蚀）的加剧，坡脚失稳垮塌，大量物源并入泥石流的发生，侵蚀宽度增大并在一定范围内波动（图7）。

图7 细粒土侵蚀宽度与深度Fig.7 Erosion of width and depth (fine-grained soil)

（2）粗粒土

粗粒土由于内摩擦角较大难以启动，但由于渗透系数较大，来流加速渗流作用并发生潜蚀现象，内部细小颗粒的迁移导致堆积体中部沉降，水流汇聚。而前缘陡坡段水力梯度较大，随着渗流作用的加剧，孔隙水压力增大，前缘发生牵引式破坏（图8）。

图8 粗粒土49.1 mm/h 试验过程现象Fig.8 Test process phenomenon in 49.1 mm/h (coarse-grained soil)

雨强较大时（44.4～49.1 mm/h）这种滑塌使得溯源侵蚀显著且剧烈，物源起动的规模成倍放大；但较小的雨强（33.2～38.1 mm/h）只能使前缘陡坡粗粒土滑动堆积或者发生有限的溯源侵蚀，难以发生大规模侵蚀（图9）。

图9 粗粒土侵蚀宽度与深度Fig.9 Erosion of width and depth (coarse-grained soil)

（3）单位侵蚀量变化

由图10和图11可知，粗粒土雨强较小时不易发生侵蚀或侵蚀规模较小，而在极端降雨下（49.1 mm/h）一旦有充足的水动力条件使其起动，由于粒径大，黏粒少，土体骨架结构一旦被破坏，侵蚀效应成倍放大，发生泥石流的规模增大。

图10 细粒土单位侵蚀量Fig.10 Erosion amount per unit (fine-grained soil)

图11 粗粒土单位侵蚀量Fig.11 Erosion amount per unit (coarse-grained soil)

3.2 多层堆积体

（1）上细下粗土

在地形坡度突增，表层细粒土在该位置率先起动，溯源侵蚀逐渐向后缘发展，而前缘陡坡因流速增大，粗粒土受下切作用不断被挟带滚落。

在中部缓坡段，由于底层颗粒较粗且受表层土体自重应力影响，不易向两侧发生侧缘侵蚀，主要是表现为下切侵蚀和溯源侵蚀，边坡两侧出现张拉裂缝（图12）。

图12 上细下粗土49.1 mm/h 试验过程现象Fig.12 Test process phenomenon in 49.1 mm/h (fine-grained on coarse-grained soil)

对于图13中各降雨条件下，上细下粗土在前缘段侵蚀宽度的突增，是由于颗粒分选效应导致大量粗颗粒在前缘形成堆积，其较大的抗剪强度，使得泥石流下深度有限，为满足泥石流的泄流，侧缘侵蚀加剧，侵蚀宽度增大。

图13 上细下粗土侵蚀宽度与深度Fig.13 Erosion of width and depth (fine-grained on coarse-grained soil)

（2）上粗下细土

随着径流到达前缘，逐渐出现管涌现象，细小颗粒的大量流失使土体颗粒骨架结构被破坏，堆积体从前缘开始发生塌陷。

塌陷后粗颗粒发生沉降，管涌现象向四周扩散，堆积体不断滑塌加入泥石流活动，溯源侵蚀急剧且强烈，发生大规模的物源起动（图14）。

图14 上粗下细土49.1 mm/h 试验过程现象Fig.14 Test process phenomenon in 49.1 mm/h (coarse-grained on fine-grained soil)

在该堆叠模式中，亦会出现颗粒分选效应使前缘粗颗粒堆积导致的侵蚀宽度突增这一现象。除此之外，一旦强降雨（44.4～49.1 mm/h）使得该次泥石流的侵蚀深度超过土体分层界限（0.1 m），其发生的侵蚀规模将显著增加，如图15所示。

图15 上粗下细土侵蚀宽度与深度Fig.15 Erosion of width and depth (coarse-grained on fine-grained soil)

（3）单位侵蚀量变化

图16所示上细下粗的多层堆叠土在四个类别的堆积体中，其侵蚀规模和变化幅度也最小的。即使是百年一遇的降雨条件（49.1 mm/h），其单位侵蚀量最大值也不足50×10-3m3，这种堆叠模式有益于沟道堆积体的稳定。

图16 上细下粗土单位侵蚀量Fig.16 Erosion amount of per unit (fine-grained on coarse-grained soil)

图17所示上粗下细土的多层堆叠土与图11所示粗粒土类似，在雨强较小时，侵蚀作用相对较弱，仅在前缘陡坡段发生较小规模的侵蚀，但强降雨作用下（44.4～49.1 mm/h），上层粗粒土对下层细粒土的铲刮以及细粒土土体失稳造成的溯源侵蚀扩大使得侵蚀规模增大。这种堆叠模式与粗粒土侵蚀现象相似的同时，大规模泥石流发生的降雨阈值下降。

图17 上粗下细土单位侵蚀量Fig.17 Erosion amount of per unit (coarse-grained on fine-grained soil)

可见，除了降雨因素外，颗粒分选效应及粒序分布形式对沟道物源的起动有着显著影响。由试验分析可知，对沟道侵蚀作用的效应排序为：溯源侵蚀＞下切侵蚀＞侧缘侵蚀＞潜蚀。溯源侵蚀是沟道物源起动加入泥石流的重要侵蚀过程。

将试验单次泥石流过后堆积体每隔0.1 m 所测的侵蚀宽度及深度按照式(7)求和得到总的侵蚀量（表6）。

式中：V侵——总侵蚀量/m3；

b——侵蚀宽度/m；

h——侵蚀深度/m。

测量计算数据汇总得表6。

表6 沟道侵蚀数据统计Table 6 Data statistics of the channel erosion

4 动储量评价模型

采用Origin Pro 2018 对表6的数据进行拟合，得到四种不同分维值的堆积体以雨强为影响因子的侵蚀量曲线如图18所示。

图18 试验模型侵蚀量与雨强的曲线关系Fig.18 Curve relationship between erosion amount and rainfall intensity in the test model

将上述四种模型分别无量纲化处理，通过SPSS 软件数据处理得到关于降雨强度和分维值的双变量拟合公式，并模型试验缩尺比例对应放大得到（R2=0.84）：

式中：V拟——侵蚀量/m3；

L——流通堆积区沟道长度/km；

D——分维值；

Q——降雨系数，Q=q/q0；

q——降雨强度/（mm·h-1）；

q0——单位降雨强度/（mm·h-1）。

上述公式不仅考虑了不同的设计雨强，亦考虑了堆积物的物质组成结构影响，更方便于工程应用。

5 验证

选取宽缓型沟道七盘沟—“7·10”和锄头沟—“8·20”泥石流分别验证其动储量，并以方群生的统计模型[2]、实际调查统计得到的动储量与文章拟合的公式相互验证、分析讨论。

锄头沟主沟长约8.9 km，平均纵坡降约184 ‰，计算各雨频下的降雨强度得表7，并对流通堆积区进行现场颗分试验得表8。

表7 雨强计算Table 7 Calculation of rain intensity

表8 锄头沟“8·20”后沟道物源颗粒累积百分含量Table 8 Cumulative percentage of particles in Chutougou after “8·20”

（1）文章的拟合公式

根据上述条件，由式（8）得到模型试验各降雨频率下的侵蚀量（表9）。

表9 本文拟合计算结果Table 9 The results of the fitting calculations in this paper

由于实际泥石流过程中，一次泥石流导致的物源起动量和最终冲出量并不能完全代表其动储量。因此对于上述式（8）的计算结果，并不能直接代表沟道物源的动储量，并且在实际工程应用中会考虑不同的设防标准来计算拦排固体物质总量，常采用“5年一遇3 次+10年一遇2 次+20年一遇1 次”、“3 次20年一遇”或者按“50年一遇”以此对文章拟合公式进行改进，得到动储量计算公式：

式中：V——沟道动储量/104m3；

μ——常数，μ=1.152。

计算得表10 如下。

表10 动储量计算Table 10 Dynamic reserve calculation

该优化公式结合了现有的泥石流防治工程设计规范，考虑了设防标准，以便实际工程中运用定量计算沟道物源的动储量。

（2）文献[2]的统计模型

统计拟合出沟道物源动储量与总物源量的相关性公式，即：

式中：VO(g)——沟道物源动储量/104m3；

Vg——沟道总物源量/104m3。

汇总得到表11。

表11 计算结果比较Table 11 Comparison of calculation results

文献[2]其拟合结果与实际调查得到的动储量存在一定差异（误差＞20%）。文章拟合公式从沟道堆积体的粒度成分和降雨条件出发，得到与实际调查更加贴合的沟道物源动储量评价模型（误差＜10%）。

文章拟合结果误差讨论：与实际调查结果相比偏小，但实际调查结果会由于调查人员的工作经验存在不同程度的测量误差。除此之外，在试验本身上，缩尺效应的影响难以避免，并且文章的动储量计算公式是在尺寸固定的沟道模型上建立，而在实际的沟道地形中，其每段沟道的宽度和堆积体厚度都不是定值，存在着固有的试验误差。对于-4.59～9.38%的误差值，该公式在允许范围之内，仍能满足于宽缓型沟道泥石流的动储量计算。