热湿循环作用下节理砂岩重复剪切力学特性及本构模型

张景昱; 詹润禾; 邓华锋; 李建林; 王文东; 宛良朋

doi:10.11779/CJGE20230652

热湿循环作用下节理砂岩重复剪切力学特性及本构模型 English Version

1.
三峡大学三峡库区地质灾害教育部重点实验室, 湖北宜昌 443002
2.
中国长江电力股份有限公司, 湖北武汉 430014
3.
中国三峡建设有限管理公司, 四川成都 610041

基金项目:

国家自然科学基金项目 52009067

国家自然科学基金项目 U22A20600

三峡库区地质灾害教育部重点实验室开放基金项目 2022KDZ23

详细信息

作者简介:
张景昱(1987—)，男，副教授，主要从事水工岩土工程方面的研究工作。E-mail：zhangjingyu@ctgu.edu.cn

通讯作者:
邓华锋, E-mail: dhf8010@ctgu.edu.cn

中图分类号: TU470
计量
- 文章访问数: 373
- HTML全文浏览量: 51
- PDF下载量: 73
出版历程
- 收稿日期: 2023-07-10
- 网络出版日期: 2024-03-24
- 刊出日期: 2024-09-30

Repeated shear mechanical properties and constitutive model of jointed sandstone under heat-wet cycles

1.
Key Laboratory of Geological Hazards on Three Gorges Reservoir Area, China Three Gorges University, Ministry of Education, Yichang 443002, China
2.
China Yangtze Power Co., Ltd., Beijing 100033, China
3.
China Three Gorges Construction Co., Ltd., Chengdu 610041, China

摘要

摘要: 乌东德、白鹤滩等大型水库处于金沙江典型干热河谷内，其岸坡消落带岩土体除经受库水位每年大幅升降以及汛期小幅频繁变化的影响外，还包括干热河谷高温低湿气候条件以及区域高频地震和水库地震的影响。为研究干热河谷地区消落带节理岩体的损伤劣化特性，在前期水-岩作用的研究基础上，以乌东德库区典型岸坡节理砂岩为研究对象，开展了热湿循环作用下节理砂岩重复剪切试验，分析了节理岩体的剪切力学特性劣化规律和劣化机理，并分段建立了考虑热湿循环作用损伤的节理砂岩剪切本构模型。研究表明：①在热湿循环以及重复剪切作用下，节理砂岩的抗剪强度、剪切刚度、剪切重点剪胀值和屈服点剪切位移均呈现“快-缓-稳”的劣化趋势，其中在前3个循环周期内劣化幅度显著，劣化占比近70%，3~8个循环周期内劣化速率放缓，8个循环周期后劣化速率趋稳；②相比于干湿循环，热湿循环对岩石孔隙度和纵波波速的影响更大，且岩石颗粒间的连接作用在较高温度的影响下更易遭受破坏；③根据剪切应力-剪切位移曲线特征，建立了以Goodman模型和Clough-Duncan双曲线模型分段修正描述的本构模型，模型曲线与试验结果能够较好吻合。研究成果及方法可为干热河谷水库岩质岸坡消落带稳定分析提供参考。
- 干热河谷 /
- 水库消落带 /
- 热湿循环 /
- 节理岩体 /
- 劣化效应
Abstract: The large reservoirs of Wudongde and Baihetan are located in the typical dry and hot valley of the Jinsha River. The rock and soil masses in the fluctuating zone of the bank slope are affected by the annual large rise and fall of the reservoir water level and the small frequent changes during the flood season, as well as the high temperature and low humidity climate conditions in the dry and hot valley, regional high-frequency earthquakes and reservoir earthquakes. In order to study the damage and degradation characteristics of jointed rock masses in the dry hot valley region, based on the previous researches on water-rock interaction, a repeated shear test on the jointed sandstone under heat- wet cycles is conducted in a typical bank slope in Wudongde reservoir area. The degradation law and mechanism of shear mechanical properties of jointed rock masses are analyzed, and a shear constitutive model for the jointed sandstone considering heat-wet cycling damage in sections is established. The results show that: (1) Under the action of heat-wet cycles and repeated shear, the shear strength, shear stiffness, shear expansion value and yield-point shear displacement of the jointed sandstone exhibit a "fast-slow-stable" deterioration trend. Among them, the degradation amplitude is significant in the first three cycles, accounting for nearly 70% of the degradation rate. The degradation rate slows down within 3 to 8 cycles, and stabilizes after 8 cycles. (2) Compared to the dry-wet cycle, the heat-wet cycle has a greater impact on rock porosity and longitudinal wave velocity, and the connection between rock particles is more susceptible to damage at higher temperatures. (3) According to the characteristics of shear stress-shear displacement curve, the modified constitutive model described by Goodman model and Clough-Duncan hyperbola model is established, and the model curve is in good agreement with the test results. The research results and methods can provide reference for the stability analysis of the rocky bank of reservoirs in dry hot valley.
- dry hot valley /
- reservoir drawdown zone /
- heat-wet cycle /
- jointed rock mass /
- deterioration effect

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. Typical dry and hot valley area of Jinsha River

下载: 全尺寸图片幻灯片

图 2 试验岩样

Figure 2. Test rock samples

下载: 全尺寸图片幻灯片

图 3 YRK-2岩石浸泡-风干试验仪

Figure 3. YRK-2 rock soaking air drying tester

下载: 全尺寸图片幻灯片

图 4 节理面剪切应力-剪切位移曲线

Figure 4. Shear stress-shear displacement curves of joint plane

下载: 全尺寸图片幻灯片

图 5 节理面法向位移-剪切位移曲线

Figure 5. Normal displacement-shear displacement curves of joint plane

下载: 全尺寸图片幻灯片

图 6 节理面抗剪强度劣化规律曲线

Figure 6. Curves of shear strength degradation of joint plane

下载: 全尺寸图片幻灯片

图 7 黏聚力、内摩擦角劣化规律曲线

Figure 7. Curves of degradation law of cohesion and internal friction angle

下载: 全尺寸图片幻灯片

图 8 砂岩纵波波速和孔隙度变化规律图

Figure 8. Variation patterns of longitudinal wave velocity and porosity in sandstone

下载: 全尺寸图片幻灯片

图 9 不同循环周期典型试样SEM照片

Figure 9. SEM photos of typical samples with different cycles

下载: 全尺寸图片幻灯片

图 10 参数r对曲线弹性阶段形态的影响

Figure 10. Effects of parameter on shape of elastic stage of curves

下载: 全尺寸图片幻灯片

图 11 参数r随试验周期N的变化规律曲线

Figure 11. Curves of variation of parameter r with test period N

下载: 全尺寸图片幻灯片

图 12 屈服点处剪切刚度K_sρ随试验周期的劣化规律曲线

Figure 12. Deterioration law curves of shear stiffness K_sρ at yield point ρ with test period N

下载: 全尺寸图片幻灯片

图 13 参数n对曲线弹性阶段形态的影响

Figure 13. Effects of parameter n on shape of elastic stage of curves

下载: 全尺寸图片幻灯片

图 14 剪切本构模型理论曲线与试验值对比图

Figure 14. Comparison between theoretical curve and experimental value of shear constitutive model

下载: 全尺寸图片幻灯片

表 1 节理面变形参数统计表

Table 1 Parametric statistics of deformation of joint plane

法向应力/MPa	周期	剪切刚度/ (GPa·m^-1)	剪切终点剪胀值/mm	屈服点剪切位移/mm
0.5	0	1.10	0.45	1.25
	1	0.90	0.26	1.40
	3	0.76	0.19	1.50
	5	0.67	0.18	1.59
	8	0.58	0.17	1.64
	12	0.56	0.16	1.68
1.0	0	1.55	0.40	1.31
	1	1.24	0.23	1.47
	3	1.03	0.17	1.59
	5	0.88	0.15	1.71
	8	0.75	0.14	1.77
	12	0.71	0.14	1.82
1.5	0	1.90	0.35	1.33
	1	1.50	0.19	1.52
	3	1.20	0.14	1.65
	5	0.98	0.13	1.77
	8	0.82	0.12	1.85
	12	0.77	0.12	1.88
2.0	0	2.34	0.47	1.36
	1	1.82	0.26	1.55
	3	1.44	0.19	1.69
	5	1.16	0.17	1.83
	8	0.95	0.16	1.92
	12	0.88	0.15	1.96

下载: 导出CSV

表 2 本构模型相关参数

Table 2 Relevant parameters of constitutive model

周期	0	1	3	5	8	12
u_ρ/mm	1.25	1.4	1.5	1.59	1.64	1.68
u_m/mm	14.98	14.98	14.98	14.99	14.96	14.97
τ_ρ/MPa	0.59	0.51	0.46	0.41	0.38	0.37
τ_m/MPa	0.78	0.7	0.6	0.57	0.55	0.54
K_sρ/(GPa·m^-1)	0.34	0.3	0.26	0.24	0.23	0.23
r	0.73	0.8	0.86	0.91	0.99	1.05
n	0.65	0.57	0.65	0.72	0.74	0.64

下载: 导出CSV

参考文献(20)

[1]	刘新荣, 景瑞, 缪露莉, 等. 巫山段消落带岸坡库岸再造模式及典型案例分析[J]. 岩石力学与工程学报, 2020, 39(7): 1321-1332. LIU Xinrong, JING Rui, MIAO Luli, et al. Reconstruction models and typical case analysis of the fluctuation belt of reservoir bank slopes in Wushan[J]. Chinese Journal of Rock Mechanics and Engineering, 2020, 39(7): 1321-1332. (in Chinese)
[2]	黄波林, 殷跃平, 李滨, 等. 三峡工程库区岩溶岸坡岩体劣化及其灾变效应[J]. 水文地质工程地质, 2020, 47(4): 51-61. HUANG Bolin, YIN Yueping, LI Bin, et al. Rock mass deterioration and its catastrophic effect of karst bank slope in the Three Gorges Project Reservoir area[J]. Hydrogeology and Engineering Geology, 2020, 47(4): 51-61. (in Chinese)
[3]	刘新荣, 傅晏, 王永新, 等. (库)水-岩相互作用下砂岩抗剪强度劣化规律的试验研究[J]. 岩土工程学报, 2008, 30(9): 1298-1302. doi: 10.3321/j.issn:1000-4548.2008.09.006 LIU Xinrong, FU Yan, WANG Yongxin, et al. Deterioration rules of shear strength of sand rock under water-rock interaction of reservoir[J]. Chinese Journal of Geotechnical Engineering, 2008, 30(9): 1298-1302. (in Chinese) doi: 10.3321/j.issn:1000-4548.2008.09.006
[4]	傅晏, 王子娟, 刘新荣, 等. 干湿循环作用下砂岩细观损伤演化及宏观劣化研究[J]. 岩土工程学报, 2017, 39(9): 1653-1661. FU Yan, WANG Zijuan, LIU Xinrong, et al. Meso damage evolution characteristics and macro degradation of sandstone under wetting-drying cycles[J]. Chinese Journal of Geotechnical Engineering, 2017, 39(9): 1653-1661. (in Chinese)
[5]	姚华彦, 张振华, 朱朝辉, 等. 干湿交替对砂岩力学特性影响的试验研究[J]. 岩土力学, 2010, 31(12): 3704-3708, 3714. doi: 10.3969/j.issn.1000-7598.2010.12.002 YAO Huayan, ZHANG Zhenhua, ZHU Chaohui. Experimental study of mechanical properties of sandstone under cyclic drying and wetting[J]. Rock and Soil Mechanics, 2010, 31(12): 3704-3708, 3714. (in Chinese) doi: 10.3969/j.issn.1000-7598.2010.12.002
[6]	姜永东, 阎宗岭, 刘元雪, 等. 干湿循环作用下岩石力学性质的实验研究[J]. 中国矿业, 2011, 20(5): 104-106, 110. doi: 10.3969/j.issn.1004-4051.2011.05.030 JIANG Yongdong, YAN Zongling, LIU Yuanxue, et al. Experimental study on mechanical properties of rock under the conditions of wet and dry cycles[J]. China Mining Magazine, 2011, 20(5): 104-106, 110. (in Chinese) doi: 10.3969/j.issn.1004-4051.2011.05.030
[7]	崔凯, 吴国鹏, 王秀丽, 等. 不同水岩作用下板岩物理力学性质劣化实验研究[J]. 工程地质学报, 2015, 23(6): 1045-1052. CUI Kai, WU Guopeng, WANG Xiuli, et al. Drying-wetting-saturating experiments for deteriora-tion of physical and mechanical properties of slate[J]. Journal of Engineering Geology, 2015, 23(6): 1045-1052. (in Chinese)
[8]	杨和平, 刘艳强, 李晗峰. 干湿循环条件下碾压膨胀土的裂隙发展规律[J]. 交通科学与工程, 2012, 28(1): 1-5. YANG Heping, LIU Yanqiang, LI Hanfeng. The development of cracks of compacted expansive soil under dry-wet cycling[J]. Journal of Transport Science and Engineering, 2012, 28(1): 1-5. (in Chinese)
[9]	申林方, 冯夏庭, 潘鹏志, 等. 单裂隙花岗岩在应力-渗流-化学耦合作用下的试验研究[J]. 岩石力学与工程学报, 2010, 29(7): 1379-1388. SHEN Linfang, FENG Xiating, PAN Pengzhi, et al. Experimental research on mechano-hydro-chemical coupling of granite with single fracture[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(7): 1379-1388. (in Chinese)
[10]	邓华锋, 段玲玲, 支永艳, 等. 干湿循环作用下节理面剪切力学特性演化规律[J]. 岩石力学与工程学报, 2018, 37(增刊2): 3958-3967. DENG Huafeng, DUAN Lingling, ZHI Yongyan, et al. Evolution of shear mechanical properties of jointed surface under dry-wet cycle[J]. Chinese Journal of Rock Mechanics and Engineering, 2018, 37(S2): 3958-3967. (in Chinese)
[11]	邓华锋, 齐豫, 李建林, 等. 水-岩作用下断续节理砂岩力学特性劣化机理[J]. 岩土工程学报, 2021, 43(4): 634-643. DENG Huafeng, QI Yu, LI Jianlin, et al. Degradation mechanism of intermittent jointed sandstone under water-rock interaction[J]. Chinese Journal of Geotechnical Engineering, 2021, 43(4): 634-643. (in Chinese)
[12]	刘新荣, 缪露莉, 袁文, 等. 化学溶液和干湿循环作用下砂岩Ⅰ/Ⅱ复合型断裂韧性劣化特性研究[J/OL]. 岩土工程学报: 1-8[2023-07-07]. LIU Xinrong, MIAO Luli, YUAN Wen, et al. The degradation characteristics of Ⅰ/Ⅱ mixed-mode fracture toughness of sandstone under the action of chemical solution and drying-wetting cycles[J/OL]. Chinese Journal of Geotechnical Engineering: 1-8[2023-07-07]. (in Chinese)
[13]	王乐华, 金晶, 张冰祎, 等. 热湿循环作用下砂岩加卸荷力学特性研究[J]. 岩石力学与工程学报, 2018, 37(3): 699-708. WANG Lehua, JIN Jing, ZHANG Bingyi, et al. Experimental study on loading and unloading mechanical properties of sandstone under heat and wet cycles[J]. Chinese Journal of Rock Mechanics and Engineering, 2018, 37(3): 699-708. (in Chinese)
[14]	孙杰豪, 郭保华, 田世轩, 等. 峰前循环剪切作用下岩石节理剪切力学特性[J]. 岩土力学, 2022, 43(增刊2): 52-73. SUN Jiehao, GUO Baohua, TIAN Shixuan, et al. Shear mechanical properties of rock joints under pre-peak cyclic shearing condition[J]. Rock and Soil Mechanics, 2022, 43(S2): 52-73. (in Chinese)
[15]	刘新荣, 邓志云, 刘永权, 等. 峰前循环剪切作用下岩石节理损伤特征与剪切特性试验研究[J]. 岩石力学与工程学报, 2018, 37(12): 2664-2675. LIU Xinrong, DENG Zhiyun, LIU Yongquan, et al. An experimental study on the cumulative damage and shear properties of rock joints under pre-peak cyclic shear loading[J]. Chinese Journal of Rock Mechanics and Engineering, 2018, 37(12): 2664-2675. (in Chinese)
[16]	刘日成, 尹乾, 杨瀚清, 等. 恒定法向刚度边界条件下三维粗糙节理面循环剪切力学特性[J]. 岩石力学与工程学报, 2021, 40(6): 1092-1109. LIU Richeng, YIN Qian, YANG Hanqing, et al. Cyclic shear mechanical properties of 3D rough joint surfaces under constant normal stiffness(CNS)boundary conditions[J]. Chinese Journal of Rock Mechanics and Engineering, 2021, 40(6): 1092-1109. (in Chinese)
[17]	李超, 齐豫, 李涛, 等. 库岸边坡水-岩作用对节理砂岩力学特性影响研究[J]. 水利水电技术, 2019, 50(6): 163-168. LI Chao, QI Yu, LI Tao, et al. Study on influence of water-rock interaction of reservoir bank slope on mechanical characteristics of jointed sandstone[J]. Water Resources and Hydropower Engineering, 2019, 50(6): 163-168. (in Chinese)
[18]	夏才初, 宋英龙, 唐志成, 等. 反复直剪试验节理强度与粗糙度变化的研究[J]. 中南大学学报(自然科学版), 2012, 43(9): 3589-3594. XIA Caichu, SONG Yinglong, TANG Zhicheng, et al. Shear strength and morphology characteristic evolution of joint surface under cyclic loads[J]. Journal of Central South University (Science and Technology), 2012, 43(9): 3589-3594. (in Chinese)
[19]	邓华锋, 肖瑶, 李建林, 等. 重复剪切作用下节理强度和形貌特征劣化规律[J]. 岩土工程学报, 2018, 40(增刊2): 183-188. doi: 10.11779/CJGE2018S2037 DENG Huafeng, XIAO Yao, LI Jianlin, et al. Degradation laws of joint strength and micro- morphologyunder repeated shear tests[J]. Chinese Jounal ofGeotechnical Engineering, 2018, 40(S2): 183-188. (in Chinese) doi: 10.11779/CJGE2018S2037
[20]	沈明荣, 张清照. 规则齿型结构面剪切特性的模型试验研究[J]. 岩石力学与工程学报, 2010, 29(4): 713-719. SHEN Mingrong, ZHANG Qingzhao. Experimental study of shear deformation characteristics of rock mass discontinuities[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(4): 713-719. (in Chinese)