换填及排水改造对北疆输水渠道稳定性的影响

邓铭江; 蔡正银; 郭万里; 黄英豪; 张晨

doi:10.11779/CJGE202105001

换填及排水改造对北疆输水渠道稳定性的影响 English Version

邓铭江^1,,
蔡正银^{2, 3},
郭万里^3, ,,
黄英豪^{2, 3},
张晨³

1.
西北旱区生态水利国家重点实验室,陕西　西安　710048
2.
南京水利科学研究院土石坝破坏机理与防控技术水利部重点实验室,江苏　南京　210024
3.
南京水利科学研究院岩土工程研究所,江苏　南京　210024

基金项目:

国家重点研发计划重点专项项目 2017YFC0405100

国家自然科学基金项目 51909170

中央级公益性科研院所基本科研业务费专项资金项目 Y319010

中国博士后科学基金项目 2018M640500

详细信息

作者简介:
邓铭江(1960—),男,博士,湖南耒阳人,中国工程院院士,教授级高级工程师,主要从事干旱区水资源研究与水利工程建设管理工作。E-mail: xjdmj@163.com

通讯作者:
郭万里, E-mail: guowljs@163.com

中图分类号: TU43
计量
- 文章访问数: 349
- HTML全文浏览量: 19
- PDF下载量: 430
出版历程
- 收稿日期: 2020-06-04
- 网络出版日期: 2022-12-04
- 刊出日期: 2021-04-30

Influences of filling replacement and drainage modification on stability of water conveyance canals in North Xinjiang

1.
State Key Laboratory of Eco-Hydraulics in Northwest Arid Region of China, Xi'an 710048, China
2.
Key Laboratory of Failure Mechanism and Safety Control Techniques of Earth-rock Dam of the Ministry of Water Resources, Nanjing 210024, China
3.
Geotechnical Engineering Department, Nanjing Hydraulic Research Institute, Nanjing 210024, China

摘要

摘要: 北疆渠道渗漏严重,膨胀性渠基土性质易劣化,使得渠道滑坡频发。为解决这一难题,建管部门对渠道进行了渠坡换填改造,同时在渠底设置了渗水抽排系统。通过室内试验和数值计算,分析了改造后渠坡的稳定性,得到如下结论：①换填料白砂岩的性质较为稳定,抗剪强度指标c,φ值基本不受干湿循环、冻融循环的影响；②滑弧由渠底贯通到渠顶,主要穿过内部的膨胀土层,而膨胀土的性质会随着“湿干”循环而劣化,仅在渠坡表层设置0.5 m厚的碎石和1.5 m厚的白砂岩换填层起到的抗滑作用较低；③渠道改造后滑坡发生量减小、渠坡稳定性提高,主要是排水体系的作用,将渗水及时抽排,降低膜后水位,保障了渠道内部膨胀土性质不再劣化,显著提高了渠道的稳定性；④为了保障渠道的长期稳定性,建议对渠道渗水进行持续抽排作业,将膜后水位控制在2 m以内。
- 北疆 /
- 输水渠道 /
- 换填 /
- 排水改造 /
- 稳定性
Abstract: The leakage of canals in North Xinjiang is serious, and the behavior of expansive canal soil deteriorates rapidly, leading to frequent canal slides. In order to solve this problem, the canals are refilled and rebuilt, and a drainage system is set up at the bottom of the canals. Laboratory tests and numerical calculation are conducted to analyze the influences of reconstruction measures of canal slopes on their stability. The following conclusions are drawn: (1) The properties of white sandstone for filling replacement are relatively stable, and the shear strength indexes c and φ are basically not affected by the drying-wetting cycle and the freeze-thaw cycle. (2) The slide circle extends from the bottom to the top of the canals, mainly through the inner expansive soil layer. However, the properties of expansive soil will deteriorate with the drying-wetting cycle, and only 0.5 m- thick gravel and 1.5 m-thick white sandstone replacement layer on the surface of the canal slope shows a low resistance effect. (3) After canal reconstruction, the amount of landslide decreases, and the stability of canal slopes is improved, which is mainly due to the function of drainage system. The drainage system reduces the post-membrane water level, which ensures that the expansive soil inside the canals will not deteriorate, and significantly improves the stability of the canals. (4) In order to ensure the long-term stability of the canals, it is suggested to carry out continuous drainage operation for water seepage in the canals, and control the post-membrane water level within 2 m.
- North Xinjiang /
- water conveyance canal /
- filling replacement /
- drainage modification /
- stability

HTML全文

0. 引言

优势流会显著影响土体内水流与溶质迁移特征^[1]，究其原因在于土体的裂隙一般会为优势流提供必要通道。这些裂隙可能源于干旱气候作用下的土体干缩^[2]、土层植物根系扩展^[3]等大气-土体-植物相互作用过程。要分析雨水在裂隙土体内优势入渗特征，一般可将土体分成裂隙和团粒两个区域（图 1），从而建立双区域模型^[4]。现有模型要么仅考虑了雨水在裂隙内流动^[4]，要么认为雨水在裂隙土体内的优势入渗所占比例不随土体初始条件和时空变化^[5]。然而，对含裂隙土体内优势流所占比例的试验研究表明：优势流对雨水在含裂隙土体内入渗量的贡献主要受降雨强度、降雨持续时间^[6]、降雨前的土体初始湿度状态^[7]、裂隙与土体表面及相邻裂隙的连通程度^[8]等因素影响。

图 1 土体内团粒-裂隙双区域的特征示意

Figure 1. Schematic diagram of characteristics of aggregate and crack domains in soils

下载: 全尺寸图片幻灯片

由此可见，基于双区域模型推求优势流的解析时，如何考虑降雨前土体的初始条件、团粒性质、降雨强度及其持续时间的影响需要深入探讨。值得注意：在模拟雨水入渗土体过程的各类模型中，Green-Ampt入渗模型具有明确的物理意义，便于建立其特征参数与土物理性质之间的关系，虽已广泛用于模拟均质土体内水流的运动特征^[9]，却鲜见于模拟含裂隙土体内水流的运动特征^[10]。为此，将Green-Ampt模型用于具有团粒-裂隙双区域特征的土体（图 1）中。通过区分雨水在团粒与裂隙内的流动特征，使推得的解析表达式参数较少、求取简便，能描述裂隙内的优势入渗量所占比例及其相对于团粒的湿润锋深度。

1. 解析模型的推导

1.1 雨水入渗裂隙域与团粒域的时段分区

为模拟雨水在团粒-裂隙双区域土体内的流动特征，构建考虑开敞裂隙（裂隙直接贯通至土体表面，如图 1和图 2（a）所示）的解析模型（对封闭裂隙的解析模型构建将另文论述）。依据文献[11]，当在土体表面选取横、纵两个剖面对裂隙长度密度（即土体表面单位面积内的裂隙长度）的估算误差在-22%~10%变化时，则可将土体中多条裂隙按照等效面积和体积简化为单个形状规则的裂隙域（图 2（b））。关于土体裂隙的剖面形状，文献[11，12]中给出了长方形、三角形和抛物线形3种类型，而且当裂隙深度远大于其宽度时，这3种形状裂隙剖面的边界长度（即剖面内土-大气交界面的长度）趋于重合，可见，当裂隙较窄且较深时，裂隙剖面形状对雨水入渗计算参数的影响很小^[12]。为计算简便起见，本文选取长方形裂隙剖面为分析对象，如所示。在该解析模型中，降雨强度 $I$ 为常数，仅考虑水流在土体内的竖向入渗，不考虑水流在团粒域与裂隙域之间因水平向运动而产生的交换。假定裂隙域的饱和渗透系数（ ${k_{{\text{s(c)}}}}$ ）远高于团粒域的饱和渗透系数（ ${k_{{\text{s(a)}}}}$ ）（下标 $({\text{c}})$ 和 $({\text{a}})$ 分别表示裂隙域和团粒域，下文的下标定义相同），故团粒域先出现积水，而裂隙域后出现积水。这里将团粒表面出现积水的时间点定义为团粒域积水点（ ${t_{{\text{p(a)}}}}$ ），将裂隙被水充满的时间点定义为裂隙域积水点（ ${t_{{\text{p(c)}}}}$ ），显然 ${t_{{\text{p(a)}}}}$ ＜ ${t_{{\text{p(c)}}}}$ 。由此可将雨水入渗含裂隙土体的过程分为3个时段，如图 2（c）所示。

图 2 雨水在含裂隙土体中的入渗过程示意图

Figure 2. Schematic diagram of rainfall infiltration into soils with cracks

下载: 全尺寸图片幻灯片

（1）团粒表面出现积水前的时段（ $t$ ＜ ${t_{{\text{p(a)}}}}$ ）：团粒域和裂隙域的入渗能力高于降雨强度，因此这两个区域内无积水，相应的入渗率与按照其表面积分配的降雨强度相同。

（2）团粒表面出现积水后且裂隙域未被水充满的时段（ ${t_{{\text{p(a)}}}}$ ≤ $t$ ＜ ${t_{{\text{p(c)}}}}$ ）：当团粒表面出现积水后，该区域的入渗率小于该区域按其表面积所分配的降雨强度。如图 2（c）所示，团粒表面的积水会沿着裂隙边壁快速地流入裂隙区域的底部，从而使裂隙底面的积水深度增大。该假定与张家铭等^[13]、罗易等^[14]开展降雨条件下含裂隙土坡足尺模型试验所观测到的积水现象一致。团粒域在该时段的入渗率可采用Selker等^[15]提出的解析表达式计算，对其积分可计算雨水在团粒域内的累计入渗量。雨水在团粒域内的入渗率与降雨强度之差由入渗率更大的裂隙补充。裂隙底面在该时段的入渗率可采用Green-Ampt入渗模型计算。当裂隙底面的入渗率与流入裂隙的雨水流速相同时，裂隙域内不会出现积水，反之则会出现积水。

（3）裂隙域被水充满后的时段（ $t$ ≥ ${t_{{\text{p(c)}}}}$ ）：此时裂隙域处于饱和状态，土体的总入渗率小于降雨强度。裂隙域内的最大入渗率受限于其饱和渗透系数，团粒域的入渗率仍可采用文献[15]所述解析表达式计算。

依据本节提出的模拟雨水在含裂隙土体中入渗过程的时段分区方法，可确定雨水在团粒和裂隙两个区域内的入渗率，进而在不考虑两个区域间水流交换的假定下，根据这两个区域内的雨水入渗量来确定湿润锋的位置。

1.2 表征裂隙优势入渗过程的解析模型推导

依据文献[]，雨水在团粒和裂隙这两个区域内的总入渗率 $v$ 可表示为

$v = ({Q_{{\text{(a)}}}} + {Q_{{\text{(c)}}}})/({A_{{\text{(a)}}}} + {A_{{\text{(c)}}}}) = \xi {v_{{\text{(c)}}}} + (1 - \xi ){v_{{\text{(a)}}}} 。$

(1)

式中： $Q$ 为体积流速（cm³/min）， $A$ 为横截面积（cm²）， $v$ 为入渗率（cm/min）； $\xi$ 为裂隙体积分数， $\xi = {V_{{\text{(c)}}}}/V$ ，即裂隙体积（cm³）与土体总体积（cm³）之比。

由图 3可知，Green和Ampt开展了均质土层的薄层积水入渗试验，发现随着入渗时间的推移，整个土层的体积含水量分布呈阶梯状^[9]（中， ${t_1}$ ~ ${t_3}$ 为入渗过程的不同时刻； $z$ 为地表以下深度（cm）； ${z_{\text{f}}}$ 为某一时刻概化的湿润锋深度（cm）； ${h_0}$ 为土层表面的积水深度（cm）； ${\theta _{\text{s}}}$ 为土层湿润区的饱和体积含水量； ${\theta _{\text{i}}}$ 为土层初始含水率； ${\psi _{\text{f}}}$ 为湿润锋处的平均基质吸力（kPa））。关于雨水在裂隙内的运动状态，依据文献[]可知，裂隙底面的积水现象多出现于较窄裂隙内，雨水在这类裂隙内的雷诺数一般低于10，故可假定这类裂隙内的雨水呈层流运动，认为裂隙和团粒这两个区域内入渗水流运动均服从Darcy定律。因此，可将Darcy定律表达式（ $v = {k_{\text{s}}}({\text{d}}h/{\text{d}}z)$ ）代入式（1）得到

图 3 Green-Ampt入渗模型示意

Figure 3. Schematic diagram of Green-Ampt infiltration model

下载: 全尺寸图片幻灯片

${k_{\text{s}}}({\text{d}}h/{\text{d}}z) = \xi {k_{{\text{s(c)}}}}({\text{d}}{h_{{\text{(c)}}}}/{\text{d}}z) + (1 - \xi ){k_{{\text{s(a)}}}}({\text{d}}{h_{{\text{(a)}}}}/{\text{d}}z) 。$

(2)

式中： ${k_{\text{s}}}$ 为饱和渗透系数（cm/min）； ${\text{d}}h/{\text{d}}z$ 为土体的竖向平均水力坡降。此外，由文献[]所述Green- Ampt入渗模型可知，当土层表面无积水时（ ${h_0} = 0$ ），其竖向平均水力坡降可表示为

${\text{d}}h/{\text{d}}z = 1 + {\psi _{\text{f}}}/({\gamma _{\text{w}}}{z_{\text{f}}}) = 1 + {\psi _{\text{f}}}(\Delta \theta )/({\gamma _{\text{w}}}D) 。$

(3)

式中： $\Delta \theta$ 为土层的体积含水量增量， $\Delta \theta = {\theta _{\text{s}}} - {\theta _{\text{i}}}$ ； $D$ 为雨水在土层内的累计入渗量（cm）， $D =$ ${z_{\text{f}}}(\Delta \theta )$ ； ${\gamma _{\text{w}}}$ 为水的重度（kN/m³）。

为概化解析模型，依据文献[]所述方法可将时间 $t$ 进行无量纲化：

$T = {\gamma _{\text{w}}}{k_{{\text{s(a)}}}}t/[{\psi _{{\text{f(a)}}}}(\Delta {\theta _{{\text{(a)}}}})] 。$

(4)

式中： $T$ 为无量纲时间； ${\psi _{{\text{f(a)}}}}$ 为团粒域内湿润锋处的平均基质吸力（kPa）； $\Delta {\theta _{{\text{(a)}}}}$ 为团粒域的体积含水量增量，即团粒域的饱和体积含水量 ${\theta _{{\text{s(a)}}}}$ 与其初始体积含水量 ${\theta _{{\text{i(a)}}}}$ 之差。

由式（4）可知，团粒域的无量纲积水点为 ${T_{{\text{p(a)}}}} =$ ${\gamma _{\text{w}}}{k_{{\text{s(a)}}}}{t_{{\text{p(a)}}}}/[{\psi _{{\text{f(a)}}}}(\Delta {\theta _{{\text{(a)}}}})]$ ，而裂隙域的无量纲积水点为 ${T_{{\text{p(c)}}}} = {\gamma _{\text{w}}}{k_{{\text{s(a)}}}}{t_{{\text{p(c)}}}}/[{\psi _{{\text{f(a)}}}}(\Delta {\theta _{{\text{(a)}}}})]$ 。需注意：由于团粒域的水力特性参数较裂隙域变化范围更小，且本文重点关注雨水在裂隙域内的优势入渗占比，故采用团粒域水力特性参数对两个区域的积水点进行无量纲化。

由图 2（c）可知，在团粒域出现积水前，团粒域与裂隙域的入渗率之比与这两个区域的截面积或体积之比相同（假定降雨过程中不考虑土层的竖向变形）：

${v_{{\text{(a)}}}} = \frac{{{A_{{\text{(a)}}}}I}}{{{A_{{\text{(a)}}}} + {A_{{\text{(c)}}}}}} = \frac{{(V - {V_{{\text{(c)}}}})I}}{V} = (1 - \xi )I （ t ＜ {t_{{\text{p(a)}}}} ）\text{，}$

(5)

${v_{{\text{(c)}}}} = \frac{{{A_{{\text{(c)}}}}I}}{{{A_{{\text{(a)}}}} + {A_{{\text{(c)}}}}}} = \frac{{{V_{{\text{(c)}}}}I}}{V} = \xi I （ t ＜ {t_{{\text{p(a)}}}} ）。$

(6)

如1.1节所述，当团粒域出现积水时，可采用文献[]所述Green-Ampt入渗模型关于时间 $t$ 的显式近似解计算此时团粒域的入渗率：

$\begin{aligned} v_{\text {(a) }} & =(1-\xi)\left\{k_{\mathrm{s}(\mathrm{a})}+\frac{\beta k_{\mathrm{s}(\mathrm{a})}+\sqrt{\frac{\psi_{\mathrm{f}(\mathrm{a})} k_{\mathrm{s}(\mathrm{a})}\left(\Delta \theta_{(\mathrm{a})}\right)}{2 \gamma_{\mathrm{w}} t}}}{1+\frac{\beta \gamma_{\mathrm{w}} k_{\mathrm{s}(\mathrm{a})} t}{\psi_{\mathrm{f}(\mathrm{a})}\left(\Delta \theta_{\text {(a) }}\right)}+\sqrt{\frac{2 \beta \gamma_{\mathrm{w}} k_{\mathrm{s}(\mathrm{a})} t}{\psi_{\mathrm{f}(\mathrm{a})}\left(\Delta \theta_{\text {(a) }}\right)}}}\right\} \\ & =(1-\xi) k_{\mathrm{s}(\mathrm{a})}\left[1+\frac{\beta+\sqrt{1 /(2 T)}}{1+\beta T+\sqrt{2 T}}\right] 。 \end{aligned}$

(7)

式中： $\beta$ 为常数，依据文献[]取2/3。为简化公式表述，令 $m = 1 + \beta T + \sqrt {2T}$ ，则式（7）可简化为

${v_{{\text{(a)}}}} = (1 - \xi ){k_{{\text{s(a)}}}}\left\{ {1 + \left[ {\beta + \sqrt {1/(2T)} } \right]/m} \right\} 。$

(8)

再令式（5）与式（8）相等可得团粒域积水点 ${T_{{\text{p(a)}}}}$ 的隐式表达式：

$I/{k_{{\text{s(a)}}}} = 1 + \left[ {\beta + \sqrt {1/(2{T_{{\text{p(a)}}}})} } \right]/{m_{{p(a)}}} （ I ＞ {k_{{\text{s(a)}}}} ）。$

(9)

由式（2）可知，若假定裂隙域内的雨水在单位水力坡降下所得的最大入渗率为 ${v_{{\text{(c)}}}} = \xi {k_{{\text{s(c)}}}}$ ，则裂隙表面出现积水需满足的条件为

$\xi {k_{{\text{s(c)}}}} = I - {v_{{\text{(a)}}}} \quad \text { 。 }$

(10)

将式（8）代入式（10）可得

$\xi k_{\mathrm{s}(\mathrm{c})}=I-(1-\xi) k_{\mathrm{s}(\mathrm{a})}\left\{1+\left[\beta+\sqrt{1 /\left(2 T_{\mathrm{p}(\mathrm{c})}\right)}\right] / m_{\mathrm{p}(\mathrm{c})}\right\} \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(I>(1-\xi) k_{\mathrm{s}(\mathrm{a})}+\xi k_{\mathrm{s}(\mathrm{c})}\right) \quad \text { 。 }$

(11)

重新整理式（11）可得关于 ${T_{{\text{p(c)}}}}$ 的隐式表达式：

$\begin{gathered} \left(I-\xi k_{\mathrm{s}(\mathrm{c})}\right) /\left[(1-\xi) k_{\mathrm{s}(\mathrm{a})}\right]=1+\left[\beta+\sqrt{1 /\left(2 T_{\mathrm{p}(\mathrm{c})}\right)}\right] / m_{\mathrm{p}(\mathrm{c})} \\ \ \ \ \ \left(I>(1-\xi) k_{\mathrm{s}(\mathrm{a})}+\xi k_{\mathrm{s}(\mathrm{c})}\right) \quad \text { 。 } \end{gathered}$

(12)

为定量描述雨水优势入渗量的相对比例，引入雨水在裂隙域与团粒域内的累计入渗量之比 $f$ （下文简称累计入渗量相对占比）：

$f=D_{\text {(c) }} / D_{\text {(a) }} \quad \text { 。 }$

(13)

式中： ${D_{{\text{(c)}}}}$ 为雨水在裂隙域内的累计优势入渗量； ${D_{{\text{(a)}}}}$ 为雨水在团粒域内的累计入渗量。值得注意： $f$ 值通过雨水的累计入渗量反映了裂隙域与团粒域之间的几何尺度差异性。

当裂隙未被水充满时（即 $T$ ＜ ${T_{{\text{p(c)}}}}$ ），雨水在裂隙域内的累计优势入渗量可通过对累计降雨量与其在团粒域内的累计入渗量求差而得

${D_{{\text{(c)}}}} = It - {D_{{\text{(a)}}}} = \frac{{(\Delta {\theta _{{\text{(a)}}}}){\psi _{{\text{f(a)}}}}IT}}{{{\gamma _{\text{w}}}{k_{{\text{s(a)}}}}}} - {D_{{\text{(a)}}}} （ T ＜ {T_{{\text{p(c)}}}} ）。$

(14)

雨水在团粒域内的累计入渗量 ${D_{{\text{(a)}}}}$ 会随其表面出现积水前、后而发生变化。以下按照雨水入渗裂隙域与团粒域的3个时段分区进行计算：

（1）当团粒域表面出现积水前（ $T$ ＜ ${T_{{\text{p(a)}}}}$ ）：团粒域与裂隙域按各自体积所占比例吸收雨水：

${D_{{\text{(a)}}}} = (1 - \xi )It = \frac{{(1 - \xi )(\Delta {\theta _{{\text{(a)}}}}){\psi _{{\text{f(a)}}}}IT}}{{{k_{{\text{s(a)}}}}{\gamma _{\text{w}}}}} （ T ＜ {T_{{\text{p(a)}}}} ）\text{，}$

(15)

${D_{{\text{(c)}}}} = \xi {\psi _{{\text{f(a)}}}}(\Delta {\theta _{{\text{(a)}}}})IT/({k_{{\text{s(a)}}}}{\gamma _{\text{w}}}) \ \ \ \ \ （ T ＜ {T_{{\text{p(a)}}}} ）。$

(16)

因此，该时段内雨水在裂隙域与团粒域内的累计入渗量之比为

$f = {D_{{\text{(c)}}}}/{D_{{\text{(a)}}}} = \xi /(1 - \xi ) （ T ＜ {T_{{\text{p(a)}}}} ）。$

(17)

（2）当团粒表面出现积水后且裂隙域未被水充满时（ ${T_{{\text{p(a)}}}}$ ≤ $T$ ＜ ${T_{{\text{p(c)}}}}$ ）：当团粒表面出现积水后（ $T$ ≥ ${T_{{\text{p(a)}}}}$ ），其累计入渗量为

${D_{{\text{(a)}}}} = [(\Delta {\theta _{{\text{(a)}}}}){\psi _{{\text{f(a)}}}}/({\gamma _{\text{w}}}{k_{{\text{s(a)}}}})][(1 - \xi )IT + \int_{{T_{{\text{p(a)}}}}}^T {{v_{{\text{(a)}}}}({T^ * }){\text{d}}{T^ * }} ] \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ （ T ≥ {T_{{\text{p(a)}}}} ）。$

(18)

式中： ${T^ * }$ 为积分变量； ${v_{{\text{(a)}}}}({T^ * })$ 为关于 ${T^ * }$ 的被积函数。将式（8）代入式（18）并积分得

$D_{\text {(a) }}=\frac{(1-\xi)\left(\Delta \theta_{(\mathrm{a})}\right) \psi_{\mathrm{f}(\mathrm{a})}}{\gamma_{\mathrm{w}}}\left\{\left(\frac{I}{k_{\mathrm{s}(\mathrm{a})}}-1\right) T_{\mathrm{p}(\mathrm{a})}+T+\ln \left(\frac{m}{m_{\mathrm{p}(\mathrm{a})}}\right)\right\} \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(T \geqslant T_{\mathrm{p}(\mathrm{a})}\right) 。$

(19)

由此可算得 ${T_{{\text{p(a)}}}}$ ≤ $T$ ＜ ${T_{{\text{p(c)}}}}$ 时段内雨水在裂隙域与团粒域的累计入渗量之比 $f$ 可表示为

$\begin{aligned} f+1 & =D_{\text {(c) }} / D_{(\mathrm{a})}+1 \\ & =\frac{I / k_{\mathrm{s}(\mathrm{a})}}{(1-\xi)\left\{\left(\frac{I}{k_{\mathrm{s}(\mathrm{a})}}-1\right)\left(\frac{T_{\mathrm{p}(\mathrm{a})}}{T}\right)+1+\frac{1}{T} \ln \left(\frac{m}{m_{\mathrm{p}(\mathrm{a})}}\right)\right\}} 。 \end{aligned}$

(20)

（3）当裂隙域被水充满后（ $T$ ≥ ${T_{{\text{p(c)}}}}$ ）：假定裂隙域内的雨水在单位水力坡降作用下呈层流运动，则可算得该时段裂隙域内的累计入渗量：

${D_{{\text{(c)}}}} = [{\psi _{{\text{f(a)}}}}(\Delta {\theta _{{\text{(a)}}}})/{\gamma _{\text{w}}}]\left[ {I{T_{{\text{p(c)}}}}/{k_{{\text{s(a)}}}} + } \right.\xi {k_{{\text{s(c)}}}}(T - {T_{{\text{p(c)}}}})/{k_{{\text{s(a)}}}} - \left. \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ {T{D_{{\text{(a)}}}}(T = {T_{{\text{p(c)}}}})} \right] （ T ≥ {T_{{\text{p(c)}}}} ）。$

(21)

式中： ${D_{{\text{(a)}}}}(T = {T_{{\text{p(c)}}}})$ 为 $T = {T_{{\text{p(c)}}}}$ 时雨水在团粒域内的累计入渗量，可表示为

${D_{{\text{(a)}}}}(T = {T_{{\text{p(c)}}}}) = [(1 - \xi )(\Delta {\theta _{{\text{(a)}}}}){\psi _{{\text{f(a)}}}}/{\gamma _{\text{w}}}]\left[ {(I/{k_{{\text{s(a)}}}} - 1){T_{{\text{p(a)}}}}} \right. + \\ \ \ \ \ \ \ \ \ \ \ \left.T_{\mathrm{p}(\mathrm{c})}+\ln \left(m_{\mathrm{p}(\mathrm{c})} / m_{\mathrm{p}(\mathrm{a})}\right)\right] \text { 。 }$

(22)

联立式（19）和式（21），可推得 $f$ 为

$\begin{aligned} & f= D_{(\mathrm{c})} / D_{(\mathrm{a})} \\ & \ \ = \frac{{\frac{{I{T_{{\text{p(c)}}}}}}{{{k_{{\text{s(a)}}}}T}} + \frac{{\xi {k_{{\text{s(c)}}}}}}{{{k_{{\text{s(a)}}}}}}\left( {1 - \frac{{{T_{{\text{p(c)}}}}}}{T}} \right) - \frac{{{\gamma _{\text{w}}}{D_{{\text{(a)}}}}(T = {T_{{\text{p(c)}}}})}}{{{\psi _{{\text{f(a)}}}}(\Delta {\theta _{{\text{(a)}}}})T}}}}{{(1 - \xi )\left[ {\left( {\frac{I}{{{k_{{\text{s(a)}}}}}} - 1} \right)\frac{{{T_{{\text{p(a)}}}}}}{T} + 1 + \frac{1}{T}\ln \left( {\frac{m}{{{m_{{p(a)}}}}}} \right)} \right]}} 。 \end{aligned}$

(23)

综上所述，联立式（17），（20）和式（23）不仅可计算任意时刻雨水在裂隙域与团粒域之间入渗的相对分区，而且可定量描述裂隙域内累计入渗量在总入渗量中所占的比例（下文简称裂隙入渗量总占比）：

$F = {D_{{\text{(c)}}}}/D = {D_{{\text{(c)}}}}/({D_{{\text{(c)}}}} + {D_{{\text{(a)}}}}) = f/(f + 1) 。$

(24)

不仅如此，概化湿润锋深度在研究团粒吸附溶质诱发的污染物迁移规律时也具有重要意义，因此这里引入裂隙域与团粒域内的概化湿润锋深度比 $\lambda$ （下文简称湿润锋深度比）：

$\lambda = \frac{{{z_{\text{f}}}_{{\text{(c)}}}}}{{{z_{\text{f}}}_{{\text{(a)}}}}} = \frac{{{D_{{\text{(c)}}}}/[\xi (\Delta {\theta _{{\text{(c)}}}})]}}{{{D_{{\text{(a)}}}}/[(1 - \xi )(\Delta {\theta _{{\text{(a)}}}})]}} = \frac{{(1 - \xi )(\Delta {\theta _{{\text{(a)}}}})}}{{\xi (\Delta {\theta _{{\text{(c)}}}})}}f 。$

(25)

式中： ${z_{{\text{f(c)}}}}$ 为裂隙域内的概化湿润锋深度（cm）； ${z_{{\text{f(a)}}}}$ 为团粒域内的概化湿润锋深度（cm）； $\Delta {\theta _{{\text{(c)}}}}$ 为裂隙域的饱和体积含水量 ${\theta _{{\text{s(c)}}}}$ 与初始体积含水量 ${\theta _{{\text{i(c)}}}}$ 之差。值得注意： $\lambda$ 值通过湿润锋深度反映了裂隙域与团粒域之间的几何尺度差异性，与 $f$ 值的含义类似。

总体而言，本节提出的解析模型（式（17），（20），（23）和式（25））主要包括以下7个参数：团粒域和裂隙域的饱和渗透系数（ ${k_{{\text{s(a)}}}}$ 和 ${k_{{\text{s(c)}}}}$ ）及其相应的体积含水量增量（ $\Delta {\theta _{{\text{(a)}}}}$ 和 $\Delta {\theta _{{\text{(c)}}}}$ ）、裂隙体积分数 $\xi$ 、降雨强度 $I$ 、团粒域内湿润锋处的平均基质吸力 ${\psi _{{\text{f(a)}}}}$ 。其中，前5个参数可通过裂隙土的室内物理力学性质试验测定， $I$ 可通过国家气象中心发布的平均降雨量算得， ${\psi _{{\text{f(a)}}}}$ 可依据文献[18]所述近似计算公式算得

${\psi _{{\text{f(a)}}}} = \left( {\frac{1}{{{\alpha _{{\text{(a)}}}}}}} \right)\left( {\frac{{0.046{n_{{\text{(a)}}}} + 2.07n_{{\text{(a)}}}^2 + 19.5n_{{\text{(a)}}}^3}}{{1 + 4.7{n_{{\text{(a)}}}} + 16n_{{\text{(a)}}}^2}}} \right) 。$

(26)

式中： ${\alpha _{{\text{(a)}}}}$ 和 ${n_{{\text{(a)}}}}$ 为van Genuchten持水曲线模型参数。

2. 解析模型的计算结果分析与讨论

2.1 不同影响因素下解析模型的计算结果分析

采用式（17），（20）和式（23）可计算不同影响因素下累计入渗量相对占比的变化过程。在计算过程中，取裂隙体积分数 $\xi$ =0.05，0.25；对每种 $\xi$ 值采用3种降雨强度： $I/{k_{{\text{s(a)}}}}$ =1.5，3，12；令 ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ 取20（当 $\xi$ =0.05时）或60（当 $\xi$ =0.05，0.25时）；如1.2节所述，参数 $\beta$ 取2/3。由式（25）可知，湿润锋深度相对占比 $\lambda$ 是 $(\Delta {\theta _{{\text{(a)}}}})/(\Delta {\theta _{{\text{(c)}}}})$ 的函数，而联立不同的 $I/{k_{{\text{s(a)}}}}$ ， ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ ， $\xi$ 和 $\beta$ 参数组合可通过 $\lambda {\text{/}}f$ 算得 $\lambda$ 。

在不同的裂隙体积分数 $\xi$ 和裂隙域与团粒域的饱和渗透系数比（ ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ ）条件下，整理累计入渗量相对占比 $f$ 和裂隙入渗量总占比 $F$ 随无量纲时间 $T$ 的变化曲线，如图 4，所示。由可知， $f$ 起初为常数，直至团粒表面出现积水（ $T$ = ${T_{{\text{p(a)}}}}$ ）后， $f$ 随时间和相对降雨强度（ $I/{k_{{\text{s(a)}}}}$ ）均呈递增趋势，而且 $\xi$ 越大，则 $f$ 越大。当降雨强度较高时（ $I/{k_{{\text{s(a)}}}}$ =12）， $\xi$ =0.05且 ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ =20（）对应的 $f$ ≈1， $\xi$ =0.25且 ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ =60（）对应的 $f$ ＞10。当降雨强度较低时（ $I/{k_{{\text{s(a)}}}}$ =1.5）， $\xi$ 对 $f$ 的影响很小。当 $I/{k_{{\text{s(a)}}}}$ =12且 $\xi$ =0.05时，裂隙域积水点 ${T_{{\text{p(c)}}}}$ 远小于1，而该区域被水充满后， $f$ 继续递增。然而，饱和渗透系数比（ ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ ）控制着优势入渗量所占比例的上限：当 $I/{k_{{\text{s(a)}}}}$ =12时， ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ =60（）对应的 $f$ 值较 ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ =20时（）超过约3倍；当饱和渗透系数比较低时（ ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ =20），裂隙域在中等降雨强度下（ $I/{k_{{\text{s(a)}}}}$ =3）亦会被水充满，因此，当 $T$ ＞2时，中（ $I/{k_{{\text{s(a)}}}}$ =3）、高（ $I/{k_{{\text{s(a)}}}}$ =12）降雨强度的 $f$ 值会出现重合；反之，当 $\xi$ =0.25且 ${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ =60时，裂隙域在不同降雨强度下均未被水充满（图 4（c））。

图 4 不同的裂隙体积分数

$\xi$ 、饱和渗透系数比（

${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ ）与降雨强度（

$I/{k_{{\text{s(a)}}}}$ ）下累计入渗量相对占比与无量纲时间的关系

Figure 4. Relative accumulative infiltration versus dimensionless time under different volume fractions of crack

$\xi$ , ratios of saturated permeability

${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ and rainfall intensities

$I/{k_{{\text{s(a)}}}}$

下载: 全尺寸图片幻灯片

图 5 不同的裂隙体积分数

$\xi$ 、饱和渗透系数比（

${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ ）与降雨强度（

$I/{k_{{\text{s(a)}}}}$ ）下裂隙入渗量总占比

$F$ 与无量纲时间

$T$ 的关系

Figure 5. Proportion of infiltration from crack versus dimensionless time under different volume fractions of crack

$\xi$ , ratios of saturated permeability

${k_{{\text{s(c)}}}}/{k_{{\text{s(a)}}}}$ and rainfall intensities

$I/{k_{{\text{s(a)}}}}$

下载: 全尺寸图片幻灯片

裂隙入渗量总占比（ $F = {D_{{\text{(c)}}}}/D$ ）亦存在与 $f$ 类似的趋势（），且 $F$ 在0.05~0.92变化（当 $\xi$ =0.05且 $T$ ＜ ${T_{{\text{p(a)}}}}$ 时F=0.05，如，；当 $\xi$ =0.25且 $I/{k_{{\text{s(a)}}}}$ =12时F=0.92，如图 5（c））。

在不同的裂隙体积分数 $\xi$ 、降雨强度（ $I/{k_{{\text{s(a)}}}}$ ）和无量纲时间 $T$ 下，整理湿润锋深度比 $\lambda$ 随团粒域和裂隙域体积含水量增量比（ $(\Delta {\theta _{{\text{(a)}}}})/(\Delta {\theta _{{\text{(c)}}}})$ ）的变化曲线，如所示， $\lambda$ 随 $(\Delta {\theta _{{\text{(a)}}}})/(\Delta {\theta _{{\text{(c)}}}})$ 和 $I/{k_{{\text{s(a)}}}}$ 均呈单调递增趋势，随 $\xi$ 减小呈递增趋势。不仅如此，在持续降雨过程中，即 $T$ 从0.1（）增至1（）的过程中，由于雨水在团粒域内的入渗率逐渐减小，故 $\lambda$ 逐渐增大。此外，整理不同 $\xi$ 下 $\lambda /f$ 与 $(\Delta {\theta _{{\text{(a)}}}})/(\Delta {\theta _{{\text{(c)}}}})$ 的关系（）可知， $\lambda /f$ 随 $\xi$ 减小呈递增趋势。

图 6 不同的裂隙体积分数

$\xi$ 、降雨强度（

$I/{k_{{\text{s(a)}}}}$ ）与无量纲时间

$T$ 下湿润锋深度比与体积含水量增量比的关系

Figure 6. Relative wetting depth ratio versus increment of volumetric water content under different dimensionless times, volume fractions of crack

$\xi$ and rainfall intensities

$I/{k_{{\text{s(a)}}}}$

下载: 全尺寸图片幻灯片

图 7 不同

$\xi$ 下

$\lambda /f$ 与

$(\Delta {\theta _{{\text{(a)}}}})/(\Delta {\theta _{{\text{(c)}}}})$ 的关系

Figure 7.

$\lambda /f$ versus

$(\Delta {\theta _{{\text{(a)}}}})/(\Delta {\theta _{{\text{(c)}}}})$ under different values of

$\xi$

下载: 全尺寸图片幻灯片

2.2 解析模型与HYDRUS-1D数值模拟计算结果的对比分析

为进一步验证本文解析模型的有效性。选取HYDRUS-1D软件^[19]提供的壤土、壤质砂土和粉质黏壤土为研究对象，采用该软件对Richards方程进行数值求解，以对比分析该数值模拟方法与本文解析模型对 $f$ 和 $\lambda$ 的计算结果。在HYDRUS-1D数值模拟过程中，对团粒域和裂隙域均采用van Genuchten-Mualem渗透系数函数进行模拟，该软件所提供的水力特性参数如，所示。需注意：在数值模拟过程中亦不考虑这两个区域之间水流质量的转换，即令渗透系数（ ${k_{{\text{s(a - c)}}}}$ ）取0。此外，将参数 ${\alpha _{{\text{(a)}}}}$ 和 ${n_{{\text{(a)}}}}$ 代入式（26）即可估算 ${\psi _{{\text{f(a)}}}}$ 。

表 1 HYDRUS-1D提供的团粒域水力特性参数

Table 1. Hydraulic parameters in aggregate domain from HYDRUS-1D

参数	${\theta _{{\text{r(a)}}}}$	${\theta _{{\text{s(a)}}}}$	${\alpha _{{\text{(a)}}}}$ / kPa^-1	${n_{{\text{(a)}}}}$	${k_{{\text{s(a)}}}}$ / (cm·min^-1)	${L_{{\text{(a)}}}}$
壤土	0	0.40	0.50	1.49	0.010	0.5
壤质砂土	0	0.40	1.25	2.50	0.020	0.5
粉质黏壤土	0.09	0.43	0.10	1.23	0.0024	0.5

下载: 导出CSV

| 显示表格

表 2 HYDRUS-1D数值模拟所用的裂隙域水力特性参数

Table 2. Hydraulic parameters in crack domain from HYDRUS-1D

参数	${\theta _{{\text{r(c)}}}}$	${\theta _{{\text{s(c)}}}}$	${\alpha _{{\text{(c)}}}}$ / kPa^-1	${n_{{\text{(c)}}}}$	${k_{{\text{s(c)}}}}$ / (cm·min^-1)	$\xi$	${L_{{\text{(c)}}}}$
壤土	0	1	60	2.0	0.6	0.10	1
壤质砂土	0	1	40	2.5	1.0	0.075	1
粉质黏壤土	0	1	20	2.5	0.1	0.12	1
注：，中： ${\theta _{\text{r}}}$ 和 ${\theta _{\text{s}}}$ 分别为残余、饱和体积含水量； $L$ 为孔隙连通-迂曲参数；下标(a)和(c)同1.1节定义。

下载: 导出CSV

| 显示表格

在进行HYDRUS-1D数值模拟时，以均质土层为研究对象，其深度均为200 cm。对每种土采用两种工况模拟。在这两种工况的数值模型上端设置与大气连通的恒定降雨强度边界（其中第一种工况的降雨强度 $I$ = $3{k_{{\text{s(a)}}}}$ ，第二种工况的降雨强度 $I$ = ${7.5k_{{\text{s(a)}}}}$ ）。土层的初始基质吸力 $\psi_{\mathrm{i}}=1000 \mathrm{kPa}$ ，壤土、壤质砂土和粉质黏壤土的团粒域初始体积含水量 $\theta_{\mathrm{i}(\mathrm{a})}$ 分别为0.02，0和0.21，相应的裂隙初始体积含水量 $\theta_{\mathrm{i}(\mathrm{c})}$ 均为0。每种工况模拟的运行时间为120 min。采用本文的解析模型（对前两种工况采用式（17），（20），（23））每分钟计算一次 $f$ 值，而基于HYDRUS-1D中对团粒域与裂隙域的累计入渗量亦可算得 $f_{\text {HYDRUS }}$ 值。采用均方根误差RMSE对比分析每分钟计算的 $f$ 与 $f_{\text {HYDRUS }}$ 值的吻合程度：

$\mathrm{RMSE}=\sqrt{(1 / N) \sum\limits_{j=1}^N\left(f_{\text {HYDRUS }-j}-f_j\right)^2} \text { 。 }$

(27)

式中： $N$ 为数据点数（ $N$ =120）。采用式（25）计算这两个区域在最后一个时间步（ $t$ =120 min）的概化湿润锋深度即可算得 $\lambda$ 值。

采用解析模型与HYDRUS-1D数值模拟对累计入渗量相对占比 $f$ 计算结果进行对比，如所示，采用解析模型对 $f$ 的计算曲线与HYDRUS-1D模型的趋势基本吻合：对壤土的RMSE在0.081~0.15范围内变化，对壤质砂土的RMSE在0.047~0.12范围内变化，对粉质黏壤土的RMSE在0.0095~0.077范围内变化。在 $I / k_{\mathrm{s}(\mathrm{a})}=3$ （）的工况中，采用式（9）计算的无量纲团粒域积水点为 $T_{\mathrm{p}(\mathrm{a})}=0.093$ ，由此可算得壤土的团粒域积水点为15.8 min，而HYDRUS-1D数值模拟计算的团粒域积水点为23.6 min，故解析模型对 $Q$ 的计算值在t=16~49 min内高于HYDRUS-1D数值模拟值，而在t＞49 min以后低于HYDRUS-1D数值模拟值。类似的，采用式（9）可算得壤质砂土和粉质黏壤土的团粒域积水点分别为7.8，73.8 min，而HYDRUS-1D数值模拟计算的对应团粒域积水点分别为12.6，98 min。类似的，在降雨强度为 $I / k_{\mathrm{s}(\mathrm{a})}=7.5$ （图 8（b））的工况中，采用式（9）亦可算得壤土、壤质砂土和粉质黏壤土的积水点分别为1.5，0.8，8 min范围内，而HYDRUS-1D数值模拟计算的相应积水点分别为3.5，1.6，18 min。

图 8 本文的解析模型与HYDRUS-1D模型对雨水在不同种土体内累计入渗量相对占比随时间变化曲线的计算结果对比

Figure 8. Simulation of relative accumulative infiltration as a function of time for proposed analytical model versus HYDRUS-1D

下载: 全尺寸图片幻灯片

不仅如此，对比本文解析模型（式（25））与HYDRUS-1D数值模拟对不同种土的湿润锋深度比 $\lambda$ 计算值（）可知，本文解析模型在 $t$ =120 min时对湿润锋深度比的计算结果与HYDRUS-1D数值模拟较为接近。HYDRUS-1D数值模拟对壤土在 $I / k_{\mathrm{s}(\mathrm{a})}=3$ 条件下计算的裂隙域概化湿润锋深度约为17 cm，且团粒域概化湿润锋深度约为6.5 cm，从而可得 $\lambda$ =2.6，而式（25）算得 $\lambda$ =2.1，相对误差为19.2%；HYDRUS- 1D数值模拟对壤土在 $I / k_{\mathrm{s}(\mathrm{a})}=7.5$ 条件下算得 $\lambda$ =9.6，而式（25）算得 $\lambda$ =9.1，相对误差为5.2%。由亦可知，式（25）对壤质砂土与粉质黏壤土的 $\lambda$ 值计算结果与HYDRUS-1D数值模拟也较为吻合。两种计算方法对壤质砂土在 $I / k_{\mathrm{s}(\mathrm{a})}=7.5$ 条件下的计算结果相差最小（0.6%），对粉质黏壤土在 $I / k_{\mathrm{s}(\mathrm{a})}=3$ 条件下的计算结果相差最大（31.5%）。不仅如此，后续将开展含裂隙土柱的降雨入渗试验，采用土柱试验结果进一步验证、评价本文解析模型的适用性。

图 9 式（25）与HYDRUS-1D数值模拟对3种土

$\lambda$ 的计算对比

Figure 9. Relative wetting depths estimated by Eq. (25) versus HYDRUS-1D for three soils

下载: 全尺寸图片幻灯片

2.3 关于解析模型计算结果及其假定的讨论

为描述雨水在团粒域内与裂隙域内入渗的动态分区，本文将Green-Ampt入渗模型与团粒-裂隙双区域模型相结合，推得表征裂隙域相对于团粒域的优势入渗量及其湿润锋深度的解析模型。通过计算分析表明：裂隙体积分数 $\xi$ 、时间 $t$ 、降雨强度与团粒域饱和渗透系数比（ $I / k_{\mathrm{s}(\mathrm{a})}$ ）、裂隙域与团粒域饱和渗透系数比（ $k_{\mathrm{s}(\mathrm{c})} / k_{\mathrm{s}(\mathrm{a})}$ ）主要影响雨水优势入渗土体裂隙的过程。需注意：时间 $t$ 已无量纲化为 $T$ ，旨在考虑团粒域的水力特性（即 $k_{\mathrm{s}(\mathrm{a})}=7.5$ 与 $\psi_{\mathrm{f}(\mathrm{a})}$ ）与体积含水量增量（即 $\Delta \theta_{(\mathrm{a})}$ ）。

上述这些因素对雨水优势入渗的贡献比例受降雨强度及其持时影响。当降雨强度较低时，雨水的优势入渗在初始阶段仅取决于裂隙所占比例 $\xi$ 。一旦团粒表面出现积水时，降雨强度就开始显著影响优势入渗，尤其当降雨强度很高时，雨水入渗率越大，则越易出现优势入渗。不仅如此，通过计算分析表明：若不考虑裂隙域与团粒域之间的水流交换，则当降雨强度不变时，随着降雨持时增长，优势入渗比例也会随之增大。此外，通过对比图 4（a）和图 4（b）亦可知，当降雨强度较大或持时较长时，裂隙域与团粒域的饱和渗透系数比也会显著影响雨水的优势入渗过程。

不仅如此，裂隙体积对雨水优势入渗也有显著影响。当其他因素一定时， $\xi$ 值越小，则裂隙域内的雨水入渗量随之减小（即和中的 $f$ 和 $F$ 值越小），而裂隙域内的湿润锋深度较团粒域更大（即 $\lambda$ 值越大）。例如，当无量纲时间 $T$ =0.1时， $I / k_{\mathrm{s}(\mathrm{a})}$ =12，从而使 $\xi$ =0.05土体的 $f$ =0.953，而 $\xi$ =0.25土体的 $f$ =2.40，可见土体的裂隙越多，则雨水在裂隙域内的优势入渗量会随之增大两倍以上。此外， $\xi$ =0.05土体中优势入渗产生的湿润锋深度比 $\lambda$ =9.05（若 $\left(\Delta \theta_{(\mathrm{a})}\right) /\left(\Delta \theta_{(\mathrm{c})}\right)$ =0.5），而 $\xi$ =0.25土体的 $\lambda$ =3.60，表明雨水在裂隙域内优势入渗的湿润锋深度会随土体的裂隙增多而减小两倍以上。该计算结果可解释文献[7]所述试验中存在的普遍现象，即在裂隙或者大孔隙不易被观察到的土体中，优势流仍会使雨水快速入渗至土体深部。

本文的解析模型虽假定裂隙域内的雨水在单位水力坡降作用下呈层流运动，但依据Green-Ampt入渗模型描述了水在团粒域与裂隙域内的流动特征。只有裂隙域内的毛细作用可不考虑时，裂隙域内取单位水力坡降的假定才是合理的。由于湿润锋处的平均基质吸力产生的毛细作用与土的进气值及其孔径分布密切相关^[18]，而大部分裂隙的毛细作用较小，因此该假定是合理的。虽然采用Darcy定律和Richards方程描述裂隙内的水流运动特征所得解的合理性仍在讨论^[20]，但裂隙域内水流取单位水力坡降的假定表明本文采用的水流描述方法类似于文献[]所述饱和流动条件下的运动波近似方法。然而，假定单位水力坡降作用下的层流运动会使HYDRUS-1D软件与本文解析模型在模拟积水点时存在差异，主要原因在于HYDRUS-1D数值模拟过程中考虑裂隙内水流的水力坡降大于1（即 $\mathrm{d} h_{(\mathrm{c})} / \mathrm{dz}>1$ ）。

关于裂隙域内的饱和流动假定，文献[]所述研究认为裂隙域内的水流通常以薄膜形式存在。此外，在裂隙域被雨水充满前（ $t<t_{\mathrm{p}(\mathrm{c})}$ ），入渗分区则由团粒域水力特性（ $k_{\mathrm{s}(\mathrm{a})}, \Delta \theta_{(\mathrm{a})}$ 和 $\psi_{\mathrm{f}(\mathrm{a})}$ ）与裂隙体积分数 $\xi$ 控制，致使裂隙域的渗透系数与饱和度不协调。此外，随着雨水在裂隙内逐渐积聚，饱和流动假定也趋于合理，并与文献[]所述现场试验结果一致。尽管如此，本文假定裂隙域内水流服从层流也是HYDRUS-1D数值模拟与本文解析模型之间存在差异的另一原因。需注意：裂隙的湿润锋深度可能会受裂隙内薄膜形式的水流影响，这与本文解析模型推导过程中的层流假定不一致。该问题的解决方法有待进一步研究，既可通过修正裂隙域的体积含水量增量 $\Delta \theta_{(\mathrm{c})}$ 来描述非饱和薄膜流动对裂隙湿润锋深度的影响，也可在本文的解析模型推导中引入Li等^[24]提出的裂隙域非饱和渗透系数函数。

此外，本文解析模型在推导过程中还包括以下假定：①团粒和裂隙这两个区域的性质与初始条件是均匀的；②降雨强度可近似为常数；③团粒和裂隙这两个区域之间水流的转换可忽略，即不考虑雨水从裂隙域侧渗入团粒域的过程。假定①不适用于大空间尺度，而假定②~③不适用于长时间尺度。然而，在解析模型的推导过程中，引入文献[]所述方法，即可考虑降雨强度变化的条件。文献[]研究表明地下水的迁移对雨水入渗分区影响较小，因此对雨水优势入渗量的相对占比 $f$ 及其总占比 $F$ 的估算结果影响较小。然而，这两个区域之间的水流转换会影响概化湿润锋深度比 $\lambda$ 。由此可见，该解析模型仍需进一步改进，使之能够考虑这两个区域之间的水流转换，例如，可引入一阶水流转换函数^[27]；也可假定裂隙形如圆柱，采用水平向Green-Ampt入渗模型^[28]估算水流转换的速率与体积。不仅如此，本文在推导解析模型时未考虑土层收缩-膨胀产生的裂隙尺寸变化，若要考虑该因素，则在考虑团粒域与裂隙域的渗透系数动态变化时可采用文献[16]所述方法进行模拟。

3. 结论

为定量描述雨水在土体内的优势入渗过程，将Green-Ampt入渗模型与团粒-裂隙双区域模型相结合，以区分雨水在这两个区域中的运动特征。主要得到以下两点结论。

（1）采用无量纲参数对解析模型进行概化，使之可预测降雨强度不变时的优势入渗量与湿润锋深度。模型的计算结果表明：雨水的流动分区会因团粒与裂隙这两个区域的积水点而发生变化；优势入渗比例随降雨强度和裂隙所占相对体积均呈递增趋势；裂隙内的湿润锋深度随其所占相对体积的增大呈递减趋势。

（2）对比分析解析模型与基于Richards方程的HYDRUS-1D数值模拟对典型含裂隙土体的计算结果发现：解析模型与数值模拟的优势入渗量计算曲线趋势较为吻合，相应的均方根误差≤0.15，湿润锋深度的相对误差为0.6%~31.5%。

本文的解析模型适用于裂隙较窄且较深、土体表面横、纵两个剖面对裂隙长度密度的估算误差在-22%~10%的条件，不仅能定量证明优势流在团粒表面出现积水前、在裂隙域内未被水充满时、在裂隙域内被水充满后存在显著分区，而且能够显式描述雨水在裂隙域与团粒域内的入渗率及其累计入渗量，从而为定量描述雨水在土体裂隙内的优势入渗特征提供了一条新途径。

图 1 北疆输水渠道断面图

Figure 1. Sectional view of canal in North Xinjiang

下载: 全尺寸图片幻灯片

图 2 膨胀土强度参数的衰减规律

Figure 2. Attenuation law of strength parameters of expansive soil

下载: 全尺寸图片幻灯片

图 3 渠道断面有限元网格

Figure 3. Mesh of canal section

下载: 全尺寸图片幻灯片

图 4 典型膨胀土渠道滑坡

Figure 4. Typical landslide of expansive soil canal

下载: 全尺寸图片幻灯片

图 5 膜后水位对F_s的影响

Figure 5. Influences of post-membrane water level on F_s

下载: 全尺寸图片幻灯片

图 6 换填前后渠坡稳定性系数

Figure 6. Stability coefficient of canals before and after replacement

下载: 全尺寸图片幻灯片

图 7 2017年和2018年膜后水位监测结果

Figure 7. Monitoring results of post-membrane water level in 2017 and 2018

下载: 全尺寸图片幻灯片

图 8 换填及渗水抽排后渠坡稳定性系数

Figure 8. Coefficient of slope stability after filling replacement and drainage modification

下载: 全尺寸图片幻灯片

表 1 湿干冻融边界条件

Table 1 Conditions of wetting-drying and freeze-thaw cycle

状态	温度	时间
湿润	常温	抽气饱和
干燥	60℃	至含水率w_D为止
冻结	-20℃	24 h
融化	20℃	24 h

下载: 导出CSV

表 2 ^#1和^#2白砂岩的粒组含量

Table 2 Granular content of white sandstones No.1 and No. 2

土料	细砾含量(5~2 mm)	砂粒含量(2~0.075 mm)	细粒含量(<0.075 mm)
^#1白砂岩	4.7	65.2	30.1
^#2白砂岩	4.0	75.9	20.1

下载: 导出CSV

表 3 ^#1白砂岩的抗剪强度值

Table 3 Strength values of white sandstone No. 1

试样状态	制样干密度（压实度）/%	c/kPa	φ/(°)
普通状态	88.4	12	31.8
普通状态	90.9	18	32.0
普通状态	93.4	30	31.9
普通状态	96.0	41	32.0

下载: 导出CSV

表 4 干湿循环、冻融循环条件下^#2白砂岩的强度值

Table 4 Strength values under wetting-drying cycle and freeze-thaw cycle

试样状态	试样编号	c/ kPa	φ/(°)
普通状态	Y1-1	3.9	33.2
干湿循环1次	W1-1	4.7	33.6
干湿循环1次	W1-2	2.0	31.9
干湿循环3次	W3-1	3.9	31.8
干湿循环3次	W3-2	9.2	32.3
干湿循环7次	W7-1	5.1	33.2
干湿循环7次	W7-2	3.2	32.5
冻融循环1次	F1-1	3.7	33.1
冻融循环1次	F1-2	6.0	33.4
冻融循环3次	F3-1	0.5	32.8
冻融循环3次	F3-2	8.9	31.4
冻融循环7次	F7-1	1.0	33.0
冻融循环7次	F7-2	5.4	32.9
注：进行了平行试验,试样编号-1和-2表示的平行试样。

下载: 导出CSV

表 5 渠道换填前各工况下渠坡稳定性系数

Table 5 Stability of canals under various conditions before reconstruction

N/次	膜后水位
N/次	0 m	1 m	2 m	3 m	4 m	5 m
0	1.78	1.73	1.67	1.61	1.56	1.50
1	1.78	1.49	1.43	1.38	1.33	1.28
2	1.78	1.42	1.37	1.31	1.27	1.22
3	1.78	1.39	1.34	1.29	1.24	1.19
4	1.78	1.37	1.32	1.28	1.23	1.18
5	1.78	1.35	1.31	1.26	1.22	1.17
注：N表示循环次数,一年表示渠基土经历一次“湿干”、“冻融”或“湿干冻融”循环。

下载: 导出CSV

表 6 渠道换填后各工况下渠坡稳定性系数

Table 6 Stability of canals under various conditions after reconstruction

N/次	膜后水位
N/次	0 m	1 m	2 m	3 m	4 m	5 m
0	1.82	1.76	1.71	1.65	1.60	1.53
1	1.82	1.56	1.51	1.46	1.39	1.35
2	1.82	1.51	1.46	1.41	1.36	1.31
3	1.82	1.48	1.44	1.39	1.34	1.29
4	1.82	1.47	1.42	1.38	1.33	1.28
5	1.82	1.46	1.41	1.37	1.32	1.27

下载: 导出CSV

表 7 边坡稳定性状态划分

Table 7 Division of slope stability states

稳定性系数F_s	F_s<1.0	1≤F_s<1.05	1.05≤F_s<F_st	F_s>F_st
稳定性状态	不稳定	欠稳定	基本稳定	稳定
注：F_st为边坡稳定安全系数,一级、二级和三级安全等级的永久边坡,一般工况对应的F_st分别为1.35,1.30,1.25。本计算边坡取F_st=1.35。

下载: 导出CSV

参考文献(9)

[1]	王浩, 游进军. 中国水资源配置30年[J]. 水利学报, 2016, 47(3): 265-271. https://www.cnki.com.cn/Article/CJFDTOTAL-SLXB201603003.htm WANG Hao, YOU Jin-jun. Progress of water resources allocation during the past 30 years in China[J]. Journal of Hydraulic Engineering, 2016, 47(3): 265-271. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-SLXB201603003.htm
[2]	蔡正银, 黄英豪. 咸寒区渠道冻害评估与处治技术[M]. 北京: 科学出版社, 2015. CAI Zheng-yin, HUANG Ying-hao. Evaluation and Treatment Technology of Canal Frost Damage in Salt and Cold Region[M]. Beijing: Science Press, 2015. (in Chinese)
[3]	朱洵, 蔡正银, 黄英豪, 等. 湿干冻融耦合循环作用下膨胀土力学特性及损伤演化规律研究[J]. 岩石力学与工程学报, 2019, 38(6): 1233-1241. https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX201906014.htm ZHU Xun, CAI Zheng-yin, HUANG Ying-hao, et al. Research on mechanical properties and damage evolution law of expensive soils under the cyclic action of coupling wetting-drying and freeze-thaw[J]. Chinese Journal of Rock Mechanics and Engineering, 2019, 38(6): 1233-1241. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX201906014.htm
[4]	蔡正银, 朱洵, 黄英豪, 等. 冻融过程对膨胀土裂隙演化特征的影响[J]. 岩土力学, 2019, 40(12): 4555-4563. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201912001.htm CAI Zheng-yin, ZHU Xun, HUANG Ying-hao, et al. Influences of freeze-thaw process on evolution characteristics of fissures in expensive soil[J]. Rock and Soil Mechanics, 2019, 40(12): 4555-4563. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201912001.htm
[5]	唐朝生, 施斌. 干湿循环过程中膨胀土的胀缩变形特征[J]. 岩土工程学报, 2011, 33(9): 1376-1384. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201109013.htm TANG Chao-sheng, SHI Bin. Swelling and shrinkage behaviour of expansive soil during wetting-drying cycles[J]. Chinese Journal of Geotechnical Engineering, 2011, 33(9): 1376-1384. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201109013.htm
[6]	张晨, 蔡正银, 黄英豪, 等. 输水渠道冻胀离心模拟试验[J]. 岩土工程学报, 2016, 38(1): 109-117. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201601013.htm ZHANG Chen, CAI Zheng-yin, HUANG Ying-hao, et al. Centrifuge modelling of frost-heave of canals[J]. Chinese Journal of Geotechnical Engineering, 2016, 38(1): 109-117. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201601013.htm
[7]	王正中, 刘少军, 王羿, 等. 寒区弧底梯形衬砌渠道冻胀破坏的尺寸效应研究[J]. 水利学报, 2018, 49(7): 803-813. https://www.cnki.com.cn/Article/CJFDTOTAL-SLXB201807004.htm WANG Zheng-zhong, LIU Shao-jun, WANG Yi, et al. Size effect on frost heave damage for lining trapezoidal canal with arc-bottom in cold regions[J]. Journal of Hydraulic Engineering, 2018, 49(7): 803-813. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-SLXB201807004.htm
[8]	郑颖人, 赵尚毅. 有限元强度折减法在土坡与岩坡中的应用[J]. 岩石力学与工程学报, 2004, 23(19): 3381-3388. https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX200419037.htm ZHENG Ying-ren, ZHAO Shang-yi. Application of strength reduction FEM in soil and rock slope[J]. Chinese Journal of Rock Mechanics and Engineering, 2004, 23(19): 3381-3388. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX200419037.htm
[9]	建筑边坡工程技术规范：GB 50330[S]. 2013. Technical Code for Building Slope Engineering: GB 50330[S]. 2013. (in Chinese)