1. 引言

作为全球半干旱至干旱地区的代表，我国西北黄土地区长期面临水分短缺的挑战[1]，水分匮缺成为限制该区植被生长的关键因素。近年来，通过实施生态综合治理工程，该区植被覆盖和生态状况得到了显著改善，但整个生态系统的脆弱性依然显著[2]。甚至区域人工植被建设更加速了土壤水分消耗，引发了以土壤水分失衡和土壤干燥化为特点的新环境问题[3] [4]。特别是在黄土灌区，土壤水分时空变动更加剧烈，相关研究需求愈发迫切[5] [6]。而土壤水分作为农田管理的重要参数[7]，其动态变化与地下水之间存在密切联系[8]。地下水对土壤水的补充与涵养作用对维持土壤水量稳定至关重要，同时，地下水与土壤水之间的持续水分交换对表层土壤生态功能和作物正常生长具有决定性作用。

鉴于水循环过程的复杂性和介质空间的多样性，土壤水和地下水的运动过程与模拟研究在各自领域内独立发展，形成了一系列数值模拟软件模型[9] [10]。近年来，杨锋等[11]依据稻田土壤水转化的实测数据，构建了基于Hydrus-1D的不同灌水下限和地下水埋深情景下稻田地下水补给过程的动态模拟模型，探究控制灌排对稻田地下水补给过程的影响机制。张贺[12]以辽宁大洋河岫岩流域为例，构建了基于Hydrus-Modflow的地表水和地下水耦合模型，验证了模型模拟的高精度和广泛适应性；代锋刚等[13]则建立了Modflow-Hydrus耦合模型，分析了不同情境下水平井的排水效果。尽管学者们基于不同视角提出了多种土壤水和地下水耦合模型，但在实际应用中由于不同地区土壤水和地下水特性存在显著差异，选择适合特定研究区的模型至关重要。

本文提出了一种适用于典型黄土灌区、考虑不同降雨条件的Hydrus-Modflow模型，并将其应用于宝鸡峡灌区，旨在为该区的土壤水和地下水关系研究提供基础。

2. 研究区概况

宝鸡峡灌区地处陕西省关中地区中西部。西起渭河宝鸡峡口，东到泾河河畔，南临渭水，北抵渭北高原，与冯家山、羊毛湾灌区接壤，灌溉面积19.77万公顷(图1)。灌区属大陆性半干旱气候区，多年平均年降水量约600 mm，主要集中在7~9月，多年平均年蒸发量约1100 mm，集中在12~2月。同时受季风作用影响，降水年际变化大，且年内分配不均。灌区主要作物为小麦和玉米，冬小麦生育期11月~次年5月，夏玉米生育期为6~9月。

图1. 宝鸡峡灌区地层钻孔分布图

3. 模型构建

3.1. 模型简介

Hydrus模型由美国国家盐改中心研发，目前被广泛应用于分析水流和溶质在非饱和多孔隙介质中的运移过程，以及干井入渗地下水运动、生态屋顶土壤水分运动模拟等方面[14]，是模拟非饱和土壤水运动最有力工具之一。Hydrus模拟土壤水分运移的方程主要是基于Richards方程进行的，一维Richards方程如下[15] [16]：

$\frac{\partial \theta }{\partial t}=\frac{\partial }{\partial z}\left[K\left(\theta \right)\frac{\partial h\left(\theta \right)}{\partial z}\right]-\frac{\partial K\left(\theta \right)}{\partial z}$ (1)

式中：$K\left(\theta \right)$ 表示非饱和渗透性函数；$h\left(\theta \right)$ 表示压力水头，cm；$\theta$ 为土壤体积含水率，cm³/cm³；t为入渗时间，h。

Modflow模型[17]是由美国地质调查局开发的三维地下水模拟模型，属于分布式物理模型，能真实地反映地下水的空间特性及其运动[18]。该模型基于有限差分法，通过假设特定的边界条件、渗透系数以及其他含水层参数，建立地下水控制偏微分方程，通过迭代运算进行求解，基本方程[19]如下：

$S_x\frac{\partial h}{\partial t}=\frac{\partial }{\partial x}{K}_{xx}\frac{\partial h}{\partial x}+\frac{\partial }{\partial y}{K}_{yy}\frac{\partial h}{\partial y}+\frac{\partial }{\partial x}{K}_{zz}\frac{\partial h}{\partial z}\pm W$ (2)

式中，$K_{xx}$ ，$K_{yy}$ 和$K_{zz}$ 分别为在x轴，y轴和z轴方向上的渗透系数(LT⁻¹)；h为地下水位(L)；W为源汇项(T⁻¹)；S_x为孔隙介质的贮水率(L⁻¹)；t为时间(T)。

Hydrus模型在非饱和带土壤水模拟方面具有显著优势，而Modflow模拟软件在区域尺度地下水计算中表现出色。鉴于此，要探索一个区域的土壤水–地下水关联性的水分运动特征，应首先探索Hydrus土壤水模拟与Modflow地下水模型之间的有效结合方式，随后基于Modflow地下水模型，构建一个能够同时描述土壤水和地下水运动特征的集成模型，实现两类模型在时空尺度上的整合，从而提高模拟的精度和全面性。

3.2. Hydrus模型建立

1) 灌区地层结构

通过分析黄土台塬84号钻孔的资料，发现研究区域地层可分为4层，具体分层见表1。较为特殊的一点是，在土层深度205~300 m的范围内，主要土层质地为含砾砂土，其间仅夹一层较薄的砂质壤土。由于砂质壤土层厚度远小于含砾砂土层，所以在深度205~300 m范围内以砂土作为代表土层。

2) 土壤水分特征参数

本研究选取粉壤土(黄土)、壤砂土、砂土和粘壤土，各自的土壤水力参数见表2。

表1. 典型地层结构表

表2. 土壤水力特征参数

注：${\theta }_r$ 土壤剩余体积含水率：${\theta }_s$ 土壤饱和体积含水率；α和n是拟合参数；Ks饱和水力传导度。

3) 概念模型

宝鸡峡灌区土壤的物理性质在垂向上表现出显著分层现象，而水平方向上变化不明显。鉴于水分运动以垂向为主，本研究将宝鸡峡灌区的土壤水分运移模型概化为垂向一维非饱和水流模型。

4) 参数确定

a) 几何信息：设定土壤模拟深度为100 cm，采用0.5 cm的步长进行划分，共设置201个节点。根据钻孔资料，按等比例将土壤剖面细分为4层。

b) 时间信息：模拟时段总长为72 h，时间步长1 h。初始时间步长设置为0.001 h，最小时间步长保持为0.001 h，最大时间步长则设定为120 h。

c) 初始与边界条件：在模型设置中，上边界设置为大气边界，下边界为定水头。可变边界条件则根据黄土台塬84号钻孔附近礼泉县在典型平水年1992年的三日降水序列，选择频率为75%、50%和20%时的相应降水量作为输入。

3.3. Modflow模型建立

1) 概念模型

水文地质概念模型是在对地下水系统进行全面分析和对水文地质条件进行科学总结与处理的基础上构建的，它为数学模型的建立提供了基础框架[20]。通过GMS软件的Solids模块，利用钻井数据构建Solids 3D可视化地层模型，并据此确定地层的空间结构。首先，通过solids的数量确定最小层数，进行三维网格建模并初始化Modflow。其次，将Solids记录的地层空间信息转换为Modflow中的层(layer)顶板标高[21]。在此基础上设定了模型参数，并给出了各层的水文地质参数。地层类型主要分为四类：完整型、尖灭型、出露型和截断型。

2) 边界条件概化

研究区南部毗邻渭河，故将其视为定水头边界。鉴于河流水位未知，则根据两端点的高程设定了相应的定水头值。具体而言，西南角定水头设为480 m，东北角设为470 m。该定水头边界作用于模型1~3层。

3) 网格划分

研究范围总体为长方形，故选择水平方向上的网格划分为100 (横向) × 50 (竖向)个单元。在垂直方向上，地层被分为三层，因此三维网格也相应地分为三层，确保每层网格与地层一一对应。

4) 参数设定

模拟旨在研究日降雨量分别为10 mm (中雨)、30 mm (大雨)和60 mm (暴雨)时灌区地层中的水头变化。由于降雨为天然降水，可视为24小时内恒定水头。在最大降雨量60 mm的情况下，土表固定水头仅2.5 mm，因此垂向渗透系数的影响相对较小，可近似等于水平渗透系数。由于固定水头值较小，水分入渗主要发生在第一层，因此模拟参数均采用第一层的值。参考相关文献，本次模拟中三层黄土渗透系数均设定为0.015 m/24h。

3.4. 模型验证

本文选择以含水率实测值与模拟值的相关关系为主要验证指标，通过计算R²绝对系数来评估模拟效果。结果如表3所示，各土层土壤水分模拟精度较高，R²均在在0.8以上，表明模型可靠。

表3. 土壤含水率模拟结果评价

4. 研究结果分析

4.1. 降雨量对土壤水分运移的影响

为探究降雨量对土壤水分运移的影响，需要在多层土壤入渗模型中调节降雨量这一时间变量边界条件，其他参数、输入条件和边界条件保持不变。具体设置为：采用黄土台塬84号钻孔临近礼泉县典型平水年1992年三日降水序列频率为75% (1.53 mm)、50% (4.83 mm)、25% (16.8 mm)的降水量。

通过模拟，获得不同埋深点包气带剖面含水率随时间变化的曲线，并将数据导出至Excel表格，进而绘制成图表。图2展示了75%、50%和25%降雨量条件下的土壤含水量随时间变化的关系。

(a) 频率为75%时 (b) 频率为50%时 (c) 频率为25%时

图2. 不同降雨频率下的土壤含水量与时间变化图

对比分析不同降雨量条件下的土壤含水量与时间的变化关系图，发现：95 cm深度的土壤含水量几乎不随时间变化；67.5 cm深度的土壤含水量在一段时间内保持稳定后急剧下降了约0.275，并最终趋于稳定；32.5 cm深度的土壤含水量从0.41迅速减少至0.2，随后缓慢下降至稳定；10 cm深度的土壤含水量则缓慢减少至约0.325后保持稳定。稳定后，各层土壤含水量呈现为95 cm > 10 cm > 32.5 cm > 67.5 cm。这表明降雨入渗难以到达深层土壤，对深层土壤水分的影响微乎其微。

进一步分析模拟数据发现，入渗稳定后，25%降雨量条件下的各层土壤含水量均高于50%和75%降雨量条件下的土壤含水量。由此可知土壤水分变化量与降雨量的大小成正比，即随着降雨量的增大，土壤水分入渗量增大。

4.2. 降补给条件下的地下水运动特征

本研究设定了暴雨(60 mm/24h)、大雨(30 mm/24h)和中雨(10 mm/24h)三种降水补给条件，并将相应补给参数设置为0.0025、0.00125和0.00042，同时保持其他参数恒定。模拟后，得到不同降雨补给条件下的模型水头平面图和正视图(图3~5)。

由图3可以看出，在暴雨补给条件下，宝鸡峡灌区的地下水流向主要呈西北至东南趋势，与地表的起伏态势相吻合，地下三层中的水分运动状态也表现出高度的一致性。然而，当水流经过灌区中部时，由于下伏基底岩层的顶托作用，地下水面发生了显著的上翘和阻滞，这一现象可能阻碍地下水的自然流动，并可能加剧该区域的地下水侵蚀作用，严重时甚至导致地下空洞和地面塌陷，因此需引起高度重视并采取相应的预防措施。

图4显示，在大雨补给条件下，宝鸡峡灌区的地下水运动特征与暴雨时相似，即由西北向东南流动，并与地表起伏保持一致。同样地，地下三层中的水分运动状态也保持了一致性。流经灌区中部时，地下水面也出现了明显的上翘阻滞现象，对地下岩土结构的稳定性和地面安全构成威胁。但由于补给量较小，各层地下水的运动速度相对较慢，高水头区范围略小于暴雨补给条件。

图5表明中雨补给条件下的地下水运动特征与暴雨和大雨时基本一致，且高水头区范围与大雨补给条件相近。这一结果说明，在宝鸡峡灌区的特定水文地质条件下，暴雨以下的天然降水对地下水运动状态的影响并不显著，与Hydrus模拟结果吻合。

(a) 平视图 (b) 正视图

图3. 60 mm降水量时模型水头变化

(a) 平视图 (b) 正视图

图4. 30 mm降水量时模型水头变化

(a) 平视图 (b) 正视图

图5. 10 mm降水量时模型水头平面视图

4.3. 合理性分析

1) Hydrus模型选择以水分入渗点垂直方向上不同深度的含水率实测值与模拟值的相关关系对模型进行验证。运用计算R²绝对系数的方法，计算模拟值与其对应的实测值的相关系数。具体计算结果见表3，R²在0.8以上，说明该模型可以较好地模拟研究区的土壤水分运移特征。

2) Modflow模型所得到的地下水模拟结果，与Hydrus模型模拟的土壤水分运动特征相互佐证，模拟精度较高，同时，模拟所得的地下水流向也与实际情况相符。本文研究结果与多篇相关参考文献研究成果进行对比后也显示，该模型在黄土灌区地下水运动特征模拟上的精度可以满足生产应用要求。

5. 结论

本文采用Hydrus-Modflow模型，根据实际情况设置了不同降雨，模拟了不同条件下多层黄土土壤水分运移特征以及地下水横向运动特征。获得的主要结论有：

1) 各层土壤含水量稳定后，95 cm层土壤含水量最高，其次分别为10 cm、32.5 cm和67.5 cm。表明降雨入渗主要影响浅层土壤，对深层土壤水分的影响微乎其微。在频率为25%的降雨量条件下目标土壤含水量最高，其次是50%和75%，可知土壤水分变化量与降雨量的大小成正比。

2) Modflow模拟表明，在暴雨、大雨和中雨补给条件下，宝鸡峡灌区地下三层中的地下水流向均为西北至东南，与地面起伏趋势一致。但当水流经灌区中部时，由于下伏基底岩层的顶托作用，地下水面出现显著的上翘阻滞现象，不仅影响地下水的顺畅流动，还可能加剧地下水侵蚀作用，严重时甚至导致地下空洞和地面塌陷。因此，这一问题应引起高度重视，并采取相应的预防与应对措施。

3) Modflow模拟表明，针对宝鸡峡灌区的水文地质条件，暴雨以下的天然降水对地下水运动状态的影响几乎可以忽略不计，地下水在垂向上的运动相对缓慢。该结论与Hydrus软件模拟的土壤水分运动特征结果相吻合。

4) 基于Hydrus软件建立的典型黄土灌区多层土壤水分运移模型的决定参数R²大于0.8，验证了模型能够准确反映典型黄土灌区土壤水分运移的特征。同时，Modflow在模拟黄土灌区地下水运动特征方面也具有良好可行性。

6. 讨论

目前，关于宝鸡峡灌区地下水分的研究主要采用统计水文学和GIS相结合的方式，一般常用统计水文学方法探究地下水动态变化趋势及其驱动力因素，进而利用GIS技术研究地下水位动态的空间分布特征。本文采用的方法既可弥补统计水文学方法中实际观测的不足，又能结合实测值实现对地下水运动特征的科学表征。另外，本次研究的模拟值与实测值的R²均在0.8以上，证明了本文建立的Hydrus-Modflow模型对典型黄土灌区的土壤水分运移特征及其地下水运动特征模拟效果较好。但在实际情况中，由于缺少对灌区主要种植的植物作用的考虑，且地下水位的影响因素较多，尤其是灌溉水的影响，所以通过模型预测的地下水位的变化也存在一定的偏差。未来可综合考虑植物蒸腾、地表水灌溉等各项条件对地下水的影响，不断提高模型模拟的精度。

基金项目

国家自然科学基金项目“变化环境下水文非平稳强波动序列随机建模与预测方法研究”(52079110)。

NOTES

作者简介：李佳玉(2000.5-)：女，黑龙江宾县人，在读硕士研究生，研究方向：水资源利用与保护，Email: nancyjj9288@163.com

^*通讯作者Email: yanbaowen2000@163.com

参考文献