过路涵洞设计中的泥石流模拟计算（吴保生）
吴保生
（清华大学）
关键词：泥石流模拟; 过路涵洞; BRI-STARS; 流管模型
作者简介：吴保生(1959-)，男，清华大学博士。
收稿日期：2000-11-28
1 引言
泥石流是公路交通中经常遇到的重大自然灾害之一。在泥石流多发地区，泥石流携带的大量泥沙颗粒常常会在路面突然堆积，造成道路交通的中断和堵塞；特别是一些较大的泥石流，还可以凭借其巨大的能量将道路和桥梁冲毁，造成严重的交通事故和巨大的经济损失。位于美国科罗拉多州斯奴马西峡谷段（Snowmass Canyon）的第SH-82号高速公路，便经常受到来自位于其靠山坡一侧流域的大小不等的泥石流的威胁，图1是现有SH-82号高速公路穿过典型泥石流堆积扇时的情形。
由于公路交通发展的需要，科罗拉多州高速公路局计划在现有SH-82号高速公路的内侧再增加一条道路，使之成为双向封闭的单行高速公路。为了避免路面的正常运行受到泥石流的频繁威胁，拟对来自靠山坡一侧大小不等流域的泥石流采取必要的工程措施加以控制。初步设计采用了在路面下设置过路涵洞的方法，将来自流域的泥石流通过涵洞输送到公路下侧的谷地[1]。
在传统的过路涵洞设计中，对于泥石流的水力计算并没有现成的方法可寻。本文采用一维泥石流数学模型对过路涵洞设计中的泥石流进行了模拟计算，提出了可靠的过路涵洞与渠道水力设计参数。
2 流域基本特征
第SH-82号高速公路的斯奴马西峡谷段位于美国科罗拉多州的Basalt和Aspen之间。该峡谷段两侧的大小流域可以分为两种类型[2]：高阶流域和低阶流域。高阶流域的河网通常大于3阶，流域面积在21.8ha到36.0ha之间，流域平均坡降在36% 和57%之间。底阶流域的河网通常只有1～2阶，流域面积在4.4ha到40.5ha之间，流域平均坡降在53%和78%之间。 位于斯奴马西峡谷段两侧的流域，在河谷和坡脚地带大多有大量的泥石流堆积物。这些堆积物呈喇叭状由从流域出口向河谷延伸，形成典型的冲积扇外形，见图1。堆积物的表面一般是比较均匀的；沿堆积扇垂向的取样分析表明，泥石流在这一地区是十分活跃的，呈周期性发生。 历史上在该地区所发生的泥石流，多是由6至9月份之间强度较大的暴雨和春季5月份迅速融化的雪水所触发的。年平均降水量约为490.2mm，其中53.4% 为7～4月份的融雪、46.6%为6～10月份的降雨所产生。 |
图1 现有高速公路穿越典型泥石流堆积扇鸟撖图 |
3 泥石流模拟
泥石流是由泥沙颗粒和水组成的泥浆，其中水的体积比一般小于50%（NRC，1982）[3]，通常是由某一时段的暴雨或迅速融雪所产生。泥石流作为非牛顿体，它的运动基本上是以间歇的高速层流运动为主。若不考虑泥石流内部的不均匀性等的影响，则描述泥石流运动的一维基本方程组可以表示为
| (1) |
| (2) |
其中 式（1）为连续方程，式（2）为能量守恒方程。式中A为水流过水断面；g为重力加速度；Q为流量;Sf为由紊动剪切力所产生的摩阻坡降；Se为由河道断面局部扩张或收缩引起的坡降；Si为由内部粘性耗散所产生的摩阻坡降；t为时间；x为沿河道的距离；Y为泥石流的水（表）面高度；α为非均匀流速分布校正系数。
从数值模拟的角度出发，河道水流可以分为洪水（颗粒的体积比含量小于20%）、泥石流洪水（颗粒的体积比含量在20%和40%之间）和泥石流（颗粒的体积比含量在40%和50%之间）[2]。根据其粘性颗粒的含量，泥石流又可以分为粘性泥石流和颗粒泥石流。粘性和颗粒泥石流的模拟方法一般包括（Jin and Fread，1999）[4]：（a）适合模拟粘性泥石流的粘塑技术；（b）适合模拟颗粒泥石流的颗粒滑动技术；和（c）适合模拟粘性和颗粒泥石流的综合阻力系数技术。
（a）粘塑技术 对于粘性泥石流，通用的指数流变方程可以表示为
τi(y)=γmSi(h-y)=τy+μ( |  | )η | (3) |
式中 h为水深；y为距河床的距离；u为沿水流方向的点流速；γm为泥石流水沙混合体的容重；μ为动力粘滞系数；η为指数方程中的粘性剪切力指数(对宾汉流体，η=1)；τi和τy分别为内部剪切力和屈服应力。
如果将Si表示为屈服应力、粘性和有关水力要素的函数，可以得到[4]
Si=(τy/γmD)[1+(m+1)(m+2)Q/0.74+0.656m(τy/μ)mDA]1/m+0.15] | (4) |
式中 m=1/η；D为水力半径。有关参数一般表示为体积含沙量的函数。一般来讲，粘性泥石流要求较小的流速和较大的水深来克服水沙混合体的屈服应力τy。
（b）颗粒滑动技术
颗粒滑动技术是建立在床面颗粒的稳定性分析基础之上的，考虑了河床坡降、孔隙水压力、泥沙颗粒在床面的摩擦角及涨水和落水过程的影响。相应的内部粘性耗散所产生的摩阻坡降可以表示为[4]
Si=hcosθ(x-x1/x2-x1)e[tanφbed/h+Ka/p-1](e/x-x1+1/h |  | ) | (5) |
式中e为系数；x1和x2分别为孔隙水压力恢复区和泥石流前坡的位置；Ka/p为常数系数；θ为以度数表示的河床坡降；φbed为内摩擦角。
（c）综合阻力系数技术 由内摩擦、粘性及紊动产生的剪切力可以用一个综合剪切力来代替，并表示为容重、阻力坡降和水力半径的函数，即
τt=γmRSt | (6) |
式中 St为综合能坡；τt为泥石流床面总剪切力。
由此，所有的阻力项便可以用一个与曼宁糙率系数相似的综合阻力系数nt来表示
St=n2tf(Q,R,A) | (7) |
根据Jin和Fread（1999）[4]的研究，粘塑技术一般适合于模拟河床比降小于φbed（20°～30°）的泥石流，而颗粒滑动技术则适合模拟河床比降大于φbed的泥石流。对于本研究来说，尽管有些河段的河床比降略大于20°，但多数河段、特别是所关心问题河段的河床比降均小于17°，因而本研究选用粘塑技术。另外一个在选择模拟方法时应该考虑的因素是泥沙颗粒的大小及粘性颗粒泥沙含量的多少。在该地区泥石流冲积扇及沉积物取样表明，其组成部分平均有55%是细颗粒泥沙（Mock，1999）[2]。较多的粘性细颗粒泥沙含量也要求选用粘塑技术比较合适。对于综合阻力系数技术来说，只要虚拟的曼宁糙率系数的大小选择合适，其使用范围似乎并没有什么具体的限制。
采用何种泥石流模拟技术和方法，还应根据所研究流域和河网的具体特点来确定。对于本研究来说，有关流域和河网的主要特点可以概括如下：
(1)在所研究的区域内，特别是在靠近冲积扇的上游地带具有良好发育的河网。对设计过路涵洞来说，其来流河段要么已具有一定的冲刷河槽，要么将开挖一定的导流河槽。由收缩和扩散引起的局部水头损失是可以忽略不计的，并且不存在建筑物对水流的阻塞现象。在一般情况下，水流运动可以用一维或准二维模型来描述。
(2)就洪水波的传播和变形而言，河床的坡降是比较陡的。在保证有足够精度的条件下，水流的运动完全可以用运动波方程（S0≈Sf）来模拟。水流运动方程中的重力项占了主导地位，因而忽略局部加速度引起的能量损失的渐变流方程便可以用来描述水流的运动。
(3)泥石流由流域出口输送至公路的距离是比较短的。若以事先给定的流域出口处的泥石流作为模拟河段的来流条件，那么在非常短的时间内，来流的非恒定性便可传播到河段的下游端，结果，洪峰流量的变形是很小的，非恒定的泥石流便可以用一系列概化的恒定流过程来描述。
(4)根据当地的地形条件，设计过路涵洞拟由底坡不同的若干涵洞段组成，表现为由陡坡段与缓坡段相接的变化。水流在由陡坡段过渡到缓坡段时，在一定的流量下将可能以水跃的形式相衔接。因而，模型应具有在不间断的条件下模拟过渡水跃的功能。
考虑到以上所述研究流域与泥石流的特点及泥石流过路涵洞的设计要求，本研究拟将非恒定的泥石流简化为一维恒定渐变流进行模拟计算。为此选择了美国高速公路局研究开发的BRI-STARS模型（Molinas，1989；Molinas and Wu，2000）[5，6]作为基础，经过修改后增加了泥石流模拟功能，包括了适合于模拟粘性泥石流的粘塑技术和通用的综合阻力系数技术。BRI-STARS模型是一个准二维的渐变流通用流管模型，在河流模拟和桥渡设计中广泛采用，一个较为突出的特点是该模型既能够模拟缓流的回水水面线，又能够模拟急流水面线，当由急流过渡到缓流产生水跃时，模型自动以水跃的形式将急流和缓流衔接起来。模型之所以能够同时模拟急流和缓流以及相应的水跃，是因为模型的水流计算同时使用了能量和动量两个控制方程。
4 模型的验证
1997年9月在所研究区域内的编号为M的流域发生了一场泥石流，本次泥石流是由暴雨触发形成的。由于设计和研究的需要, H-P土木技术公司于1999年就该次泥石流形成的河槽进行了调查(Mock, 1999)[2]，并给出了若干横断面的最高水面线。根据所测量的冲积扇面积和厚度，估计本次泥石流堆积的总体积约为2 523m3，相应的径流过程约为15～20年一遇。图2是泥石流所塑造河槽的一个典型断面，相应的河床比降平均约为12%。以下我们将用该泥石流来对所建立的模型进行验证。
| 为了使用上述有限的观测资料来验证修改后的模型，按照图2所示断面形态构造了一个底坡为12%的概化河段。在验证计算中，采用了洪峰流量为1.67m3/s和2.46m3/s的两组洪水过程, 相应的径流过程分别为15年和25年一遇。取泥沙的体积比含量为50%（容积增大系数为2），相应的泥石流最大流量分别为3.34m3/s和4.92m3/s。泥石流的颗粒体积比含量一般在40%和50%之间，这里取容积增大系数取为2是比较保守的上限，对于工程设计是比较合适的。同时由径流洪水过程放大得到相应的泥石流过程，所计算出的泥石流体积与估算量2523m3(Mock, 1999)[2]大体吻合，由此说明采用容积增大系数取为2是可行的。对于综合糙率系数法，曼宁糙率系数选为0.16。对于粘塑技术法，τy和μ均取为622。4N/m2。图2给出了模拟计算的水面线与实测水面线的比较。相应于径流过程约为15～20年一遇的实测泥石流水面线，落在了模拟计算所给出的15年和25年一遇水面线之间。这一结果说明，采用综合糙率系数法和粘塑技术法所计算的泥石流水面线均与实测结果是一致的。 |
图2 横断面计算与实测水面线比较 |
5 涵洞设计
本文泥石流模拟的主要任务，是确定有关设计涵洞能否将给定的泥石流安全地通过涵洞送往公路下侧的谷地。图3给出了SH-82号高速公路扩建工程中所设计的泥石流过路涵洞的概化简图。其中图3(a)是高速公路穿过泥石流堆积扇时的典型平面布置情况，包括了泥石流自流域沟口以下的导流渠道和过路涵洞的衔接。图3(b)则是过路涵洞的纵剖面图，从中可以看到过路涵洞与公路路面的相对位置。一般来讲，过路涵洞由坡度明显不同的上段和下段所组成。由于受到下游河谷底部高程和公路建筑物的限制，下段涵洞的底坡通常要小于上段，但为了保持泥石流的有效流动，该段的底坡必须大于泥石流运动所需要的最小临界值。Bagnold（1956）[7]曾给出了临界坡降的如下计算公式
式中γs和γw分别为泥沙颗粒和水的容重。对于流动的沙粒，Chen(1988)[8]曾经使用Cv=0.53、tanφbed=0.32和γw=1793kg/m3,得出θ大约为3°到4°(tanθ≈0.05)。与此相似，Dietrich(1997)[9]认为当河道的坡降小于5% 时就不在具有输送泥石流的能力。因而，我们将过路涵洞底坡的最小临界值取为5%。 |
图3 SH 82号高速公路的典型泥石流过路涵洞设计布置简图 |
图4给出了流域M和I相应25年一遇径流的设计来流过程线，图中的虚线为洪水过程线，而实线则为泥石流过程线。需要指出的是，直接由产流产沙模型给出泥石流过程是十分困难的，因而设计的泥石流是由相应的设计洪水放大2倍而得到，即认为泥石流中泥沙颗粒的体积比含量为50%（容积增大系数为2）。
| 图5和6给出了使用不同模拟方法给出的水面线计算结果。所使用的方法包括：（a）一般洪水模拟；（b）使用放大洪水模拟；（c）使用粘塑技术的泥石流模拟；和（d）使用综合糙率技术的泥石流模拟。对于采用一般水流计算方法的（a）和（b）来说，曼宁糙率系数取为0.04,这一糙率系数与泥石流堆积扇的颗粒组成是相适应的。根据以上用野外观察资料对模型的验证计算结果，对于采用综合糙率技术方法的泥石流模拟来说，曼宁糙率系数取为0.16。对于使用粘塑技术的泥石流模拟，τy和μ均取为957.6N/m2。该值大于验证计算中所采用的相应值，这是为了设计的安全起见而考虑了τy和μ取值上的不确定性。并且，这里τy和μ的取值与相似情况下Jin 和Fread (1999)[4]所采用的数值大小是一致的。 图5所示流域M设计涵洞横断面的宽和高均为2.44m |
图4 洪水过程线与泥石流过程线：（a）流域M；（b）流域I |
| （8ft），图6所示流域I设计涵洞横断面的宽和高分别为2.44m（8ft）和1.83m（6ft）。由图5和6的模拟结果可以看出，如果按照一般水流的模拟方法来考虑，结果在发生泥石流的情况下，将会导致设计尺寸的严重不足。图中所示模拟结果还进一步表明，设计所采用的最终过路涵洞的底坡和断面尺寸，完全可以满足输送设计泥石流的要求。需要指出的是，上述涵洞设计是根据泥石流模拟结果经过多次修改后选定的；在没有考虑泥石流模拟之前所给出的最初概化设计的涵洞底坡和断面尺寸，与上述最终设计涵洞底坡和断面尺寸是很不一致的。另外，最终选用涵洞断面的宽度与H-P Geotechnical, Inc. (Mock, 1999)[2]所测量的自然条件下的泥石流河槽尺度是一致的。现有河网是以往泥石流活动的结果，一定程度上反映了河道的几何尺寸与泥石流之间存在的河相关系。 |
图5 流域M相应25年一遇泥石流水面线计算结果 |
图6 流域I相应25年一遇泥石流水面线计算结果 |
本文对过路涵洞设计中的泥石流进行了模拟计算。考虑到所研究的流域与泥石流的特点及泥石流过路涵洞的设计要求，本研究在美国联邦高速公路局现有的准二维河流模拟与桥渡冲刷流管模型（BRI-STARS）的基础上，增加了泥石流的模拟计算功能，并根据实测资料对修改后的模型进行了验证。然后采用修改后的模型，对位于美国科罗拉多州斯奴马西峡谷段的第SH-82号高速公路扩建工程中的泥石流涵洞设计，进行了不同涵洞设计方案的泥石流模拟计算，给出了可靠的涵洞与渠道水力设计参数。通过本文的研究可以得到以下结论：
(1) 传统的过路涵洞设计，是根据流域产流计算所得洪水过程，进行一般清水水流的有关水力学计算。对于输送泥石流的过路涵洞，若按照一般清水水流的模拟方法来考虑，将会导致设计过路涵洞尺寸的严重不足。
(2) 设计泥石流过程线必须在清水径流的基础上使用适当的容积增大系数，泥石流的峰值流量一般是相应清水洪峰流量的2倍，即容积增大系数为2。传统过路涵洞设计中常用的混凝土糙率系数（n=0.016～0.018），对于泥石流是不适用的。由于泥石流附加阻力的存在，应采用较大的曼宁糙率系数（n=0.14～0.16），或者采用考虑泥石流运动机理的有关模拟方法。
(3) 在泥石流过路涵洞的水力设计中,相应于泥石流设计峰值流量的河道坡降和最大水深是决定因素,而水深随时间的变化和洪峰流量的传播时间是相对次要的。修改后的BRI-STARS模型是过路涵洞设计中泥石流模拟的有效工具，可用于同类条件下过路涵洞设计的水力计算。
(4) 修改后的BRI-STARS模型目前包括了适合于模拟粘性泥石流的粘塑技术和通用的综合阻力系数技术。原则上不仅较适合于模拟粘性泥石流，也适合于模拟其它不同类型的泥石流。但对于颗粒泥石流，应增加适合模拟颗粒泥石流的颗粒滑动技术，以便于和综合阻力系数技术的结果比较。
参 考 文 献[1] Cotton, G. (1999).“Preliminary Drainage Report, Snowmass Canyon Segment, SH-82 Reconstruction." Chang Engineering Company, Inc., Denver, CO.
[2] Mock, R. G. (1999).“Phase I - Debris Flow Hazard and Mitigation Strategies Study, Colorado Sate Highway 82-Snowmass Canyon, Pitkin County, Colorado." Hepworth-Pawlak Geotechnical, Inc., 5020 County Road 154, Glenwood Springs, CO 81601.
[3] National Research Council (1982). “Selecting a Methodology for Delineating Mudslide Hazard Areas for the National Flood Insurance Program." National Academy of Science Report by the Advisory Board on the Built Environment, Washington, D。C.
[4] Jin, M. and Fread, D. L. (1999). “1D Modeling of Mud/Debris Unsteady Flows." Journal of Hydraulic Engineering, ASCE, Vol. 125, No. 8, pp. 827-834.
[5] Molinas, A. (1989). “BRI-STARS Model and Its Applications." U.S. Department of Transportation, Federal Highway Administration Bridge Scour Symposium, McLean, VA, Oct. 16-19.
[6] Molinas, A.and Wu, B.(2000).“User's Primer for BRI-STARS (BRIdge Stream Tube model for Alluvial River Simulation)." U.S. Department of Transportation, Federal Highway Administration, Pub. No. FHWA-RD-99-191, McLean, VA.
[7] Bagnold, R. A. (1956).“Flow of Cohesionless Grains in Fluids." Philos. Trans., Royal Society of London, Ser. B, 249, 235-297.
[8] Chen, C.-L. (1988).“Generalized Viscoplastic Modeling of Debris Flow." Journal of Hydraulic Engineering, ASCE, Vol. 114, No. 3, pp. 237-258.
[9] Dietrich, W. and Sitar, N. (1997). “Geoscience and Geotechnical Engineering Aspect of Debris-Flow Hazard Assessment." Debris-Flow Hazards Mitigation: Mechanics, Prediction, and Assessment, C.-L. Chen ed., ASCE, San Francisco, CA, pp. 656-676.
