HYDRUS—1D界面版使用说明

2025-11-01-Hydrus-1d界面版使用说明

HYDRUS—1D界面版使用说明

1 HYDRUS—1D

1.1 介绍

HYDRUS—1D是国际地下水模型中心公布的，计算包气带水分、盐分运移规律的软件，用它可以解算在不同边界条件制约下的数学模型。若将坐标原点选在地面，取z轴向下为正，则一维饱和—非饱和带水分运移基本方程为：

\frac{\partial θ}{\partial t} = \frac{\partial}{\partial z} [K (θ) (\frac{\partial h}{\partial z} - 1)] - S

HYDRUS-1D是美国农业部盐渍土实验室开发的模拟非饱和土壤中的水、热、溶质运移的软件，它在模拟土壤中水分运动、盐分、污染物(如农药)和养分(如土壤氮素)运移方面得到广泛应用。软件分为1维、2维、3维三种，分别命名维Hydrus-1D、Hydrus-2D、Hydrus-3D。由于非饱和土壤水主要是1维垂向运动的形式，因此Hydrus-1D的应用非常广泛。

1.2 HYDRUS—1D模型说明

该程序用数值方法求解了饱和-不饱和水流的Richards方程和热和溶质传输的基于Fickian的对流扩散方程。

流动方程采用了汇项，以考虑水的吸收植物根系。

热传输方程认为传导以及与流动的水的对流。

所述溶质输运方程考虑了溶质在液相中的对流和分散和在气相中的扩散。

迁移方程还包括以下方面的规定：

固相和液相之间的非线性和/或非平衡反应，
液相和气相之间的线性平衡反应，
零阶生产
两个一级降解反应：
- 独立降解：一种独立于其他溶质的物质
- 耦合降解：它提供了涉及顺序一阶衰变反应的溶质之间的偶联

该程序可用于分析不饱和，部分饱和或完全饱和的多孔介质中的水和溶质运动。

流动区域本身可能由不均匀的土壤组成。流动和运输可以在垂直，水平或大致倾斜的方向上进行。模型的水流部分可以处理（恒定或随时间变化）规定的水头和通量边界，受大气条件控制的边界以及自由排水的边界条件。在模拟过程中，土壤表面边界条件可能会从规定的通量变为规定的水头类型条件（反之亦然）。

对于溶质传输，该代码支持（恒定和变化）规定的浓度（狄利克雷或第一类Dirichlet or first-type）和浓度通量（柯基或第三类Cauchy or third-type)）边界条件。分散系数包括反映分子扩散和曲折度影响的项。

使用van Genuchten [1980]，Brooks and Corey [1964]和改进的van Genuchten类型的分析函数描述了非饱和土壤的水力特性。进行了修改以改善对接近饱和的水力特性的描述。HYDRUS代码通过使用Scott等人引入的经验模型来合并滞后现象。[1983]和Kool and Parker [1987]。该模型假设干燥扫描曲线是从主干燥曲线定标的，而润湿扫描曲线是从主润湿曲线定标的。

HYDRUS还通过一组线性比例转换将相应的土壤水力特性与参考土壤的水力特性联系起来，从而执行定标程序来近似给定土壤剖面中的水力变化。

根的生长通过逻辑生长函数进行模拟。水分和盐分应力响应函数可以根据Feddes等人提出的函数定义。[1978]或van Genuchten [1987]。

使用Galerkin型线性有限元方案对控制流和输运方程进行数值求解。对于饱和和非饱和条件，都使用隐式（向后）有限差分方案实现了时间积分。采取了其他措施来提高解决瞬态问题的效率，包括自动调整时间步长以及遵守Courant和Peclet编号的预设范围。使用Celia等人提出的质量保守方法评估含水量。[1990]。最小化运输解决方案中的数字振荡的可能选项包括上游称重，人工分散和/或性能指标。

HYDRUS实施了Marquardt-Levenberg类型的参数估计技术，用于根据测得的瞬态或稳态流量和/或传输数据对选定的土壤水力和/或溶质运移和反应参数进行反推。该程序允许根据观察到的水含量，压头，浓度和/或瞬时或累积边界通量（例如，渗透或流出数据）来估算几个未知参数。可以在参数估计过程中选择性地包含其他保留或水力传导率数据，以及用于约束优化参数以保留在某个可行区域中的罚函数（贝叶斯估计）。

1.3 Hydrus-1D软件版本

1.3.1 版本更新

Hydrus-1D 4.0版中的新功能（2007年发布；与3.0版相比包括）：

蒸气流量
水，蒸汽和能量传输耦合
双渗透式水流与溶质运移
在流动区域中具有溶质迁移和两点吸附的双孔隙水流，
Penman-Monteith组合和Hargreaves方程计算潜在的蒸散量
蒸发，蒸腾和沉淀的每日变化
通过将HYDRUS与PHREEQC生物地球化学代码耦合获得的HP1代码支持。
一种考虑根源溶质吸收的选择，包括被动吸收和主动吸收[Šimunekand Hopmans，2009]。另请参阅FAQ18。

Hydrus-1D 3.0版的新功能（2005年4月发布；与2.0版相比，包括）：

补偿摄取根水
Kosugi [1996]（对数正态模型）和Durner [1994]（双孔隙率模型）建议的土壤水力特性的附加分析模型
双孔系统中的水流
具有附着/分离系数的溶质迁移，可以模拟胶体，病毒和细菌的迁移
两个动力学吸附位点（例如一个可以用于空气-水界面）
基于过滤理论的附着系数评估
二氧化碳生产和运输模块（Unsatchem模块）
考虑主要离子的迁移，沉淀/溶解，阳离子交换和络合反应的地球化学碳酸盐化学模块（Unsatchem模块）
新模型比旧模型快大约3倍。

1.3.2 用户界面

基于Microsoft Windows的图形用户界面（GUI）管理运行HYDRUS以及节点离散化和编辑，参数分配，问题执行以及结果可视化所需的输入数据。

在图形环境中指定了所有空间分布的参数，例如土壤类型/层，根系水分吸收分布以及水，热和溶质运动的初始条件。

离散化节点的位置可以由用户以图形方式进行编辑，以优化不同元素的厚度。

该程序包括允许用户构建特定于应用程序的流和传输模型以及动态执行图形分析的控件。

输入和输出都可以使用图形工具进行检查。

Hydrus-1D Shell程序将所有几何和参数数据转换为HYDRUS输入格式。

1.3.3 后期处理

后处理也在GUI中进行。

Hydrus-1D提供了有关预选时间的压头分布，水含量，水和溶质通量，根系吸水量，温度以及土壤剖面中浓度的图表。

输出还包括随时间变化的图，例如跨边界或离开根区域的实际，潜在和累积通量。

可以在配置文件中的任何位置添加观察点，以获得水含量，压头，温度和/或浓度的图形输出。

支持的外围设备包括最流行的打印机和绘图仪类型。

该程序包括一小部分土壤水力学特性。

广泛的上下文相关的联机帮助是该界面的一部分。

1.3.4 模型附带的示例

直接问题：

草场下土壤剖面中的水流和溶质运移-季节性模拟
大型沉箱中的渗透和排水
涉及磁滞的瞬态流
色谱柱渗透测试-来自Skaggs的数据
非线性阳离子吸附的溶质运移-赖和尤里纳克的数据
带有非线性阳离子吸附的溶质运移-来自Selim等的数据。
硝化链溶质运输
非平衡阳离子吸附的溶质运移
大气条件波动下的热传递

反问题：

一步式流出实验-来自Kool等的数据。（1987）
多步骤流出实验-来自Hopmans的数据
蒸发实验-来自Wendroth的数据
向上渗透
涉及磁滞的瞬态流
非线性阳离子吸附的溶质运移-赖和尤里纳克的数据
非线性阳离子吸附的溶质运移-Selim的数据
硝化链溶质运输
水平渗透-来自George Vachaud的数据
水平渗透和再分布-来自Vachaud的数据
沙柱中的排水-来自Vachaud的数据
草下田间土壤剖面中的水流-季节模拟

1.3.5 源代码

HYDRUS-1D直接计算模块的（开放）源代码已包含在GNU通用公共许可证中，其详细条件在此处列出。仅当您同意本文档中列出的条件时，才能下载HYDRUS-1D计算模块的源代码。我们欢迎您对本源代码提出任何意见或对其进行任何改进。每当您使用此源代码时，请参考HYDRUS-1D技术手册或以下手稿：Simunek等，VZJ-2008，Hydrus和Stanmod models.pdf。

HYDRUS-1D计算直接模块的源代码可以在这里下载。可以使用FORTRAN编译器编译和链接此源代码。生成的可执行代码应命名为h1d_calc.exe。然后可以将其直接放置在HYDRUS-1D安装文件夹中，并直接从HYDRUS-1D GUI执行。

请注意，当使用不来自Microsoft PowerStation的编译器（例如，在工作站上）编译此源代码时，您需要在对函数GetDat和GeTim（获取日期和获取时间）的主文件调用中删除（或注释掉），并可能替换它们具有相应的功能。此外，您需要删除对MSFLIB库的调用（使用MSFLIB，“ interface”和“ end interface”（包括）之间的代码以及Signal handler例程h_sig和hand_fpe）。这部分代码对程序的其余部分没有影响。它仅捕获异常，例如被零除，并在这种情况下关闭程序和输出文件。

1.4 地球化学UNSATCHEM模块

地球化学UNSATCHEM模块（Šimůnek和Suarez，1994年；Šimůnek等人（1996）已分别在HYDRUS-1D和HYDRUS（2D / 3D）软件包的一维和二维计算模块中实现。地球化学UNSATCHEM模块模拟主要离子在可变饱和多孔介质中的传输，包括主要离子平衡和动力学非平衡化学。生成的代码旨在预测主要离子化学以及瞬态流动过程中土壤中的水和溶质通量。UNSATCHEM-2D中化学系统的主要变量是Ca，Mg，Na，K，SO4，Cl，NO3，H4SiO4，碱度和CO2。该模型考虑了这些组分之间的各种平衡化学反应，例如络合，阳离子交换和沉淀-溶解。对于方解石的沉淀溶解和白云石的溶解，可以使用平衡或多组分动力学表达式，包括正向和反向反应。考虑的其他溶解-沉淀反应包括石膏，菱镁矿，菱镁矿和海泡石。由于土壤溶液的离子强度会在时间和空间上发生很大变化，并且经常达到很高的值，因此修改后的Debye-Hückel表达式和Pitzer表达式都被纳入模型中以计算单个离子的活性。还考虑了溶液化学性质对水力传导率的影响。水流和热传输模块与常规HYDRUS中的模块相似（几乎相同）。可能的应用包括研究农业土壤的盐碱化/复垦，各种灌溉系统的可持续性以及采矿作业中盐水的处理。

2 HYDRUS—1D GUI版使用说明

2.1 主过程

用来选择计算模块，水流计算模块时其他计算模块的基础，为必选模块。为默认勾选。

水流（water flow）
- 蒸汽运输
- 土壤表面积雪
溶质运移（solute transport）
- 标准溶质运移
- 主要化学过程
- PHREEQC生物地球化学模型
热传导（heat transport）
根系吸水（root water uptake）
根系生长（root growth）
反解（Inverse Solution）

2.2 基本信息

2.2.1 模型几何结构

（１）土壤结构输入参数（右侧）

土壤质地的种类；（Number of Soil Materials）
土壤剖面分区的数目；(Numbe of Layers for Mass)
土壤表层的倾斜程度，取值在 0 到 1 之间。（土壤表层与水平面夹角的余弦函数值。例如表层是水平时，取值为 1。）；(Decline from Vertical Axes)
所模拟土壤的厚度。(depth of the soil)

（２）长度单位（Length Units）（左侧）

长度单位可以选择 mm， cm 或 m。注意:此处选择的长度单位将是以后所有涉及长度量纲的量的默认长度单位。

2.2.2 时间信息

（1）时间单位（Time Units）

时间单位可选择秒、分、小时、天、年。注意:此处选择的时间单位将是以后所有涉及时间量纲的量的默认长度单位。

（2）时间离散（Time Discretization）

模拟起始时间
模拟结束时间
初始时间步长
最小时间步长
最大时间步长

前两项由模拟的需要而定。后三项是控制迭代计算的时间步长参数，若填入的时间步长过大，很有可能导致迭代结果不收敛；而填入的时间步长过小，会使得模拟耗费过多的时间。根据土壤质地的复杂程度填入时间步长，若计算时出现不收敛或运算时间过长等情况，再进行调整。

初始时间步长可设置0.0001天，最小时间步长=初始时间步长/10，最大时间步长=初始时间步长*1000

（3）时间变化边界状态（Time-Variable Boundary Conditions）

有一类边界条件是随时间变化的。例如要模拟田间蒸发蒸腾条件下的土壤含水量，就必须输入模拟时间段每天的蒸发量和蒸腾量。若是时变边界条件，需要将 Time-Variable Boundary Condition 选中，并填入一共有多少时变的边界条件记录。
增加了利用气象资料计算蒸腾量的模块。可以选择蒸腾量的计算公式（彭曼-蒙特斯公式和哈格里夫斯公式及表层能量平衡）；气象数据可以输入，也可以由 Hydrus-1D 生成。

2.2.3 输出信息

T-level Information 一般是默认选中，在代码计算的结果中(后缀名为.out)，有专门存储如平均水头等信息的文件。
- Every n time steps: 每隔几个时间步长输出一组这样的信息存在文件中。
print at Regular Time Interval 是专门指定按固定时间间隔输出信息的。
Screen Out 如果选中，则在代码执行时， DOS 程序框中会快速闪现计算的过程。
Print Fluxes for Observation Nodes：对设定的观察点输出流量，代替输出温度。运行过程中打印这些信息
Hit Enter at End：选中后，代码执行完毕，按 Enter 键才会关闭 DOS 计算框。
Number of Print Times：在指定时间专门输出水头、含水量等详细信息的个数。具体的时间点可以选择。默认的有均匀分布时间点（ default）和对数分布时间点（ default log）。

2.3 土壤水流模型

2.3.1 迭代规则（iteration criteria）

一般可以采用默认值。

1、迭代规则

最大迭代次数（Maximum number of iteration），任意时间步长的最大迭代次数，同时用改进的Picard方法求解非线性Richards方程。建议的默认值是10。
含水量迭代精度（water content tolerance），流区非饱和区节点的绝对含水量容差–。此参数表示在特定时间步长内连续两次迭代之间的含水量值的最大期望绝对值变化。建议的默认值是0.001。
压力水头迭代精度（pressure head tolerance），流动区域饱和部分节点的绝对压头公差L。该参数表示在特定时间步长内连续两次迭代之间的压力头值的最大期望绝对值变化。建议的默认值是0.1

2、时间步长控制

HYDRUS采用了三种不同的时间离散控制:(1)时间离散的数值解,(2)时间相关的离散边界条件的实现,和(3)时间离散的提供打印输出仿真结果(如节点因变量的值,水和溶质质量平衡组件,和其他流态的信息)。

离散化2和离散化3是相互独立的;它们通常涉及输入数据文件中描述的可变时间步长。离散化1从规定的初始时间增量Δt开始。此时间增量按以下规则在每个时间级别自动调整:

（1）离散化1必须与离散化2和3得到的时间值一致。

（2）时间增量不能小于预先选定的最小时间步长Δtmin，也不能超过最大时间步长Δtmax（即Δtmin≤Δt≤Δtmin）。

（3）如果在特定的时间步长中，达到收敛所需的迭代次数为≤3，则通过将Δt乘以预定的常数>1(通常在1.1到1.5之间)来增加下一个时间步长的时间增量。如果迭代次数≥7, 则下一个时间步长Δt下乘以一个常数< 1(通常在0.3和0.9之间)。

（4）如果在特定的时间步长中，任意时间级的迭代次数大于规定的最大值(通常在10到50之间)，则该时间级的迭代过程终止。时间步长随后重置为Δt≤/3，迭代过程重新开始。

下最优迭代范围（lower optimal iteration range），当水流达到收敛所需的迭代次数小于该次数时，将时间步长乘以较高的时间步长乘数因子(时间步长增加)。建议的默认值是3
上最优迭代范围（upper optimal iteration range），当水流达到收敛所需的迭代次数大于该次数时，将时间步长乘以较低的时间步长倍增因子(时间步长减小)。建议的默认值是7。
较低的时间步长倍增因子（lower time step multiplication factor），如果水流收敛所需的迭代次数小于较低的最优迭代范围，则将时间步长乘以该次数(时间步长增加)。推荐的默认值是1.3。
较高的时间步长乘数因子（upper time step multiplication factor），如果水流达到收敛所需的迭代次数大于上最优迭代范围，则乘以时间步长(时间步长减小)。建议默认值为0.7。

3、内部插值表（internal interpolation tables）

在数值模拟开始时，HYDRUS根据给定的一组水力参数，为流域中的每种土壤类型生成一个表，表中包括含水量、水力传导性和特定的水容量。在迭代求解过程中，利用表中各条目之间的线性插值计算水力特性值。如果参数h超出了规定的区间(ha, hb)，则直接从水力函数(即,没有插值。在计算上，上述插值技术比直接计算整个压力头范围内的水力函数要快得多。可以通过将两个值都指定为零来避免表中的插值表达式。因此，通常直接从水力函数来评价土壤的水力性质。,没有插值。

张力区间下限（low limit of the tension interval），压力头区间的下限值[L]的绝对值，每一土壤层的内部都将生成一个水力属性表。默认值1e-006。
张力区间上限（low limit of the tension interval），压力头区间的上限[L]的绝对值，每一土壤层的内部都将生成一个水力属性表。默认值10000。

将给出区间(ha, hb)土壤水力特性的输出图。

建议用户使用张力区间的上下限的默认值。张力区间的下限通常被选择为一个非常小的数(例如，1.e – 3比1.e-6厘米)。上限可以修改，使它包括在模拟过程中遇到的大多数压力头(例如，hCritA或萎蔫点)。为了避免在使用Brooks和Corey土壤水动力特性模型时数值求解的问题，其下限(绝对值)应大于1/(即空气进口价值)。

2.3.2 土壤非饱和水力模型（Soil Hydraulic model）

在土壤水力模型对话框中，用户选择用于土壤水力特性的水力模型，并指定在计算过程中是否考虑Hysteresis滞后。

水力模型:

程序允许用户选择六种模型对土壤水力性质:

the van Genuchten-Mualem model [van Genuchten, 1980]
the van Genuchten-Mualem model with an air-entry value of -2 cm
modified van Genuchten type equations [Vogel and Cislerova, 1988]
Brooks and Corey 方程[1964]
Kosugi对数正态分布模型[1996]
双孔隙度模型 [Durner, 1994].

这些模型概述如下。此外，用户还可以选择两个双重孔隙度非平衡流动模型，假定流动区和静止区之间的传质与f)含水率或g)压头和双重渗透率模型成正比。有关这些型号的详细说明，请参阅HYDRUS技术手册。另外两种方法(双渗透模型和查找表)在当前版本的HYDRUS中不可用。

由于 Richards 方程中涉及含水量q 、负压 h ，非饱和水力传导度 K 这三个未知量，所以必须补充两个关系式。由于负压 h 和饱和水力传导度 K 都是含水量q的函数，有很多的模型来描述这两个关系。模型分为两大类，单孔模型和多孔模型。最常用的是单孔模型中的 vanGenuchten-Mualem 模型。

Hysteresis 表示滞后现象，即吸湿和疏干并不是的负压—含水量曲线不一致。一般选择不考虑滞后。

一、水力模型

单孔隙度模型
- van Genuchten-Mualem模型(van Genuchten, 1980)
- van Genuchten-Mualem模型的考虑2cm的空气进入
- 修正van Genuchten类型方程(Vogel和Cislerova, 1988),
- Brooks and Corey 方程[1964]
- Kosugi对数正态分布模型[1996]
双孔隙度模型(Durner,1994)
- 双孔隙度模型（Durnner，dual van Genuchten ）
- 双孔隙度模型（mobile-immobile，water c mass transfer ）
- 双孔隙度模型（mobile-immobile，head mass transfer ）
- 双渗透率模型（kinematic wave equation ）
- 双渗透率模型（Gerk and Genuchten，1993）
查找表

二、滞后现象（Hysteresis）

没有滞后现象（No Hysteresis）
保留曲线滞后（Hysteresis in retention curve）
保留曲线电导率滞后（Hysteresis in retention curve conductivity）
保留曲线滞后（无泵送，Bob lenhard） Hysteresis in retention curve (no pumping ,Bob lenhard)
- 初始干燥曲线（initially drying curve）
- 初始湿润曲线（initially wetting curve）

2.3.3 土壤水力参数（water flow parameters）

HYDRUS规范中的非饱和土水力性质是用一组与van Genuchten[1980]相似的封闭方程来描述的，van Genuchten利用Mualem[1976]的统计孔隙大小分布模型来获得非饱和水电导率函数的预测方程。原始的van Genuchten方程可以进行修改，以便在描述接近饱和的水力特性时增加额外的灵活性[Vogel和Cislerova, 1988]。用户还可以选择Brooks和Corey[1964]、Kosugi[1996]和Durner[1994]提出的模型。

不同的水力模型就有不同的水力参数。

对于 van Genuchten-Mualem 模型，有α 、n、Qr、Qs、Ks、l 六个参数。其中 l 一般取 0.5。其他参数与土质、土壤结构等有关系。参数取得是否合理和定解条件（初始条件和边界条件）是模拟是否能反映实际的关键。在土壤水力参数很难直接得知的情况下，需要通过其他方法推求参数，或是通过模拟结果与实际比较，来调节参数。

表 2‑1 van Genuchten-Mualem 模型参数表

	英文	中文
Qr	$\theta_{r}$	土壤残余含水量
Qs	$\theta_{s}$	饱和土壤含水量
Alpha	$\alpha$ $\left[ \mathrm{L}^{-1} \right]$	$\alpha$ 在土壤水分保持函数中
N	Parameter n in the soil water retention function	参数n在土壤水分保持函数中
Ks	$k_{s}$ $\left[ \mathrm{L}^{-1} \right]$	饱和导水率
l	Tortuosity parameter in the conductivity function [-]	电导率函数中的弯曲度参数[-]，一般取值0.5

当修正van Genuchten-Mualem模型（[Vogel and Cislerova, 1988]）使用时，需要以下参数

表 2‑2 van Genuchten-Mualem 模型参数表

	英文	中文
Qa	Parameter $\theta_{a}$ in the soil water retention function	土壤水分保持函数中的参数 $\theta_{a}$
Qm	$\theta_{m}$ in the soil water retention function	$\theta_{m}$
Kk	$K_{k}$ $\theta_{k}$ $\left[ \mathrm{L}^{-1} \right]$	$\theta_{k}$ $\left[ \mathrm{L}^{-1} \right]$ $K_{k}$
Qk	$\theta_{k}$ $K_{k}$	$\theta_{k}$ $K_{k}$

当Dumer模型使用时，

表 2‑3 van Genuchten-Mualem 模型参数表

	英文	中文
w	Parameter w for material M [-]. Relative weighting factor for the subcurve of the second overlapping subregion.	第二重叠子区域子曲线的相对权重因子，参数w
Alpha2	Parameter aM $[ {L}^{-1} ]$ , for the second overlapping subregion.	对于第二重叠子区域，分层的参数α，M/L
n2	Parameter n for material M [-], for the second overlapping subregion.	对于第二重叠子区域，分层的参数n

双孔隙度模型，质量转换为含水量水力梯度驱动

表 2‑4 van Genuchten-Mualem 模型参数表

	英文	中文
QrIm	Parameter $\theta_{r}$ for the immobile region of material M.
QsIm	Parameter $\theta_{s}$ for the immobile region of material M.
Omega	Parameter $\omega$ M $[{T}^{-1} ]$ .

双孔隙度模型，质量转换为含压力水头驱动

表 2‑5 van Genuchten-Mualem 模型参数表

	英文	中文
QrIm	Parameter $\theta_{r}$ for the immobile region of material M.
QsIm	Parameter $\theta_{s}$ for the immobile region of material M.
AlphaI	Parameter $\alpha$ M $[ {L}^{-1} ]$
nI	Parameter $n$ for the immobile region of material M
Omega	Parameter $\omega$ M $[{T}^{-1} ]$ .

双重渗透率模型，

表 2‑6 van Genuchten-Mualem 模型参数表

	英文	中文
QrFt	Parameter $\theta_{r}$ for the fracture region of material M.
QsFr	Parameter $\theta_{s}$ for the fracture region of material M.
AlphaFr	$\alpha$ for the fracture region of material M [L-1]
nFr	Parameter $n$ for the fracture region of material M
KsFr	Parameter $k_{s}$ for the fracture region of material M [LT-1]
lFr	Tortuosity parameter in the conductivity function for the fracture region of material M [-]
w	Parameter $\omega$ for material M [-]. w is the ratio of the volumes of the macropore or fracture domain and the total soil system.
Beta	Parameter $\beta$ (a shape factor that depends on the geometry) for material M [-].
Gamma	$\gamma$ (a scaling factor) for material M [-].
a	$d$ (an effective diffusion pathlength) for material M [L].
Ksa	Parameter $k_{a}$ (the effective hydraulic conductivity Ka [LT-1] of the fracture-matrix interface) for material M [LT-1].

2.3.4 土壤水力边界状况（water flow boundary conditions）

用户需要定义上边界和下边界条件

1、以下上边界条件是可用的：

上边界
- 定水头（Constant pressure head）
- 定流量（Constant flux）
- 大气边界（Atmospheric boundary condition with surface layer）
- 带径流大气边界（Atmospheric boundary condition with surface run off）
- 变水头（Variable Pressure Head）
- 变水头/流量（Variable Pressure Head/Flux）

2、以下下边界条件是可用的：

下边界
- 定水头（Constant pressure head）
- 定水分通量（Constant flux），需要输入水分通量。如果是底部隔板边界，则水分通量为零。
- 变压力水头（Variable pressure head）
- 变水分流量（Variable flux）
- 自由下泄/排水边界（Free drainage）适用于地下水埋深的情况，包气带为自有排水的情况。
- 深层排水（Deep drainage）边界排泄与地下水位相关，这个需要给定地下水位，经验公式在hydrus里面。
- 渗漏面（Seepage face）适用于实验室土柱
- 水平排水沟（Horizontal drains）；均匀土壤剖面

3、初始条件

初始条件可以根据压力水头或含水量来指定。

需要注意的是，当初始条件以含水率的形式给出时，对于双重孔隙度模型，规定的初始含水率将在可动域和不可动域之间重新分配，并与各域的饱和含水率成比例。对于双渗透模型，规定的初始含水率将被解释为矩阵域的初始含水率，并根据该值计算矩阵域的初始压头。然后将裂缝域分配给与矩阵域相同的压头。

压力水头（in pressuer heads）
含水量（in water contents）

有些边界条件是时变边界，需要输入时变的记录。

大气边界条件Atmospheric BC

输入PET和LAI:代替输入大气通量(即，也可以根据LAI和消光系数分别输入潜在蒸散发和独立潜在蒸散发通量的组合值。也可以规定地表径流开始前地表水层的最大厚度(土壤表面的最大h)。

2.4 溶质运移模型

2.4.1 一般信息

1、时间权重方案，选择时间权重系数

= 0.0 显式方案
= 0.5 Crank-Nicolson隐式方案
= 1.0 全隐式方案

从求解精度的角度出发，建议采用Crank-Nicholson克兰克-尼克尔森隐式格式。完全隐式格式可能导致数值色散过大，但避免了数值不稳定性。显式格式最容易出现非理想振荡的数值不稳定性。

2、空间权重方案

Galerkin方案
上游加权方案
带人工分散的Galerkin方案

$\Phi_n$ $\Phi_{n}^{u}$ $\Phi^u$ 确保相对较多的权重放在位于单元上游一侧的节点的流速上。

在Galerkin有限元结果中，可以加入额外的人工扩散来稳定数值解，避免不必要的振荡。在此基础上加入人工扩散，从而实现了一个涉及Pe.Cr (Peclet数* Curant数)的稳定性准则[Perrochet and Berod, 1993]。Pe.Cr的推荐值是2.0。

Number of Solutes	同时模拟或参与衰变链反应的溶质数量。
Pulse Duration	浓度脉冲持续时间[T]。所有浓度边界条件(没有指定时变边界条件)在大于“脉冲持续时间”几倍时等于零。当使用时变边界条件(如大气BCs或时变磁头或通量)时，忽略脉冲持续时间值，必须使用时变边界条件窗口指定边界浓度。
Mass Units	$[ML^{-3}]$ $ML_c^{-3}L_gT^{-1}$ $L_c$ $L_g$ $c{\theta}V$ $ML_c^{-3}L_g^2$ 单位。需要对其他涉及浓度和长度单位的变量进行类似的单位调整。
Stability Criterion	无量纲Peclet和Curant数的乘积(Pe.Cr)。该准则既可用于在具有人工离散格式的Galerkin有限元中加入人工离散，也可用于限定Galerkin有限元格式的时间步长(给定的Peclet数导致更低的Courant数)。
Use Tortuosity Factor	当分子在水和气相中的扩散系数要乘以一个弯曲系数时，请选中此框。
Temperature Dependence	如果假定溶质输送和反应参数与温度有关，则勾选此框。
Water Content Dependence	如果假定溶质输送和反应参数与含水量有关，则勾选此框。

3、溶质输运模型。这些模型包括:

均衡模型
单点吸附模型(化学非平衡)
双点吸附模型(化学非平衡)
两个动力学位点模型(使用附着/分离的粒子运输，化学非平衡)。该模型常用于病毒、胶体或细菌的运输。
两个动力学位点模型(基于过滤理论，化学非平衡)
双重孔隙度(移动-不移动水)模型(物理非平衡)
流动带双点吸附的双重孔隙模型(物理和化学非平衡)
双渗透模型(物理非平衡)
在基质中加入不动水或动态吸附的双渗透模型（物理和化学非平衡）

4、收敛标准

当考虑非线性吸附时，平流-弥散溶质输运方程变为非线性。与Richards方程相似，在每一个新的时间步上，必须使用迭代过程来获得全局矩阵方程的解。在每次迭代过程中，用高斯消去法或共轭梯度法导出并求解线性化的代数方程组。反演后，利用初始解对系数重新求值，重新求解新的方程。这个迭代过程一直进行，直到获得一个令人满意的收敛程度，即，直到在所有节点上，连续两次迭代之间浓度的绝对变化小于某个浓度容差(在HYDRUS中定义为绝对浓度容差与浓度与相对浓度容差的乘积之和(推荐的默认值为0.001))。需要指定某个时间步长的最大迭代次数(建议值为10)。当达到最大迭代次数时，数值解可以(a)终止涉及瞬态水流的问题，或(b)以减少稳态水流问题的时间步长重新启动。

绝对收敛误差|
相对收敛误差
最大迭代次数
弯曲度
- 使用弯曲度因子
- Millington & Quirk
- Moldrup
溶质数量
浓度脉冲持续时间[T]

2.4.2 溶质运移参数

土壤参数


Bulk.d.	$\rho [ML^{-3}]$
Disp	$D_L$ [L]
Frac	无量纲分数的吸附位点分类为-1型，即，即考虑化学非平衡选择时具有瞬时吸附的位点[-]。如果考虑平衡迁移，则设为1。考虑物理非平衡时，与流动水接触的吸附点的无量纲分数[-]。如果所有吸附点都与流动水接触，则设置为1。
ThImob	束缚水含量。当不考虑物理非平衡时设为0。 ThImob(当>0)被解释为一种含水量，当使用两个动力学位点模型时，其中的颗粒(如病毒、胶体、细菌)被排除在外。

溶质的具体参数

Diffus.W	$D_w$ $[L^2T^{-1}]$
Diffus.G	$D_a$ $[L^2T^{-1}]$

基质中存在不动水或动力学吸附的双渗透模型(物理和化学非平衡)


Bulk.d.	$\rho [ML^{-3}]$
Disp.M.	$D_L$ [L]
Frac.M.	在矩阵中分类为类型-1的吸附位点的无量纲分数，即，即考虑化学非平衡选择时具有瞬时吸附的位点[-]。
ThImob.M.	基质中的固定含水量。当不考虑物理非平衡期权时，设为0。
Diffus.W	$D_w$
Mass Tr.	$\Gamma_s$ $[T^{-1}]$

地球化学模块UNSATCHEM

指定了一组溶液、吸附物和沉淀浓度组合。溶液、吸附物和沉淀浓度组合的数量在普通溶质传输信息窗口中指定。这种(溶液、吸附和沉淀)浓度组合可用于指定初始和边界条件。

溶液浓度必须以meq/L为单位(即， mmol_c / L土壤溶液)。吸附和沉淀(固相)浓度必须以meq/kg为单位(即， mmol_c / kg土壤基质;例如mmol_c中的方解石(每公斤土壤基质中的钙)。详见技术手册。注意，在Unsatchem模块中，“溶质传输-一般信息”窗口中指定的“质量单位”被忽略。

2.4.3 溶质反应参数

类型	描述
Kd	$k_{s} \cdot \left[\mathrm{M}^{-1} \mathrm{~L}^{3}\right]$ 。
Nu	$v \cdot \left[\mathrm{M}^{-1} \mathrm{~L}^{3}\right]$ 。 (手册中称h，朗缪尔系数)
Beta	$\beta$ (Freundlich指数)
Henry	$kg$ []
SinkWater1	$u_w$ $[T^{-1}]$
SinkSolid1	$u_s$ $[T^{-1}]$
SinkGas1	$u_g$ $[T^{-1}]$ 。
SinkWater1'	$\mu_{w}^{\prime} \cdot \left[\mathrm{T}^{-1}\right]$ 。,表示链式反应
SinkSolid1'	$\mu_{s}^{\prime} \cdot\left[\mathrm{T}^{-1}\right]$ 。表示链式反应
SinkGas1'	$\mu_{g^{\prime}} \cdot \left[\mathrm{T}^{-1}\right]$ ,表示链式反应
SinkWater0	$\gamma_{w} \cdot \left[\mathrm{ML}^{-3}\mathrm{~T}^{-1}\right]$ 。
SinkSolid0	$\gamma_{s} \cdot \left[\mathrm{T}^{-1}\right]$ 。
SinkGas0	$\gamma_{g} \cdot \left[\mathrm{ML}^{-3}\mathrm{~T}^{-1}\right]$ 。
Alpha	$\omega \cdot \left[\mathrm{T}^{-1}\right]$ 。

吸附等温线(上表中前三个参数)如下:

s = \frac{k_{s} c^{β}}{1 + η c^{β}}

当β=1时，上述吸附方程成为 Langmuir equation, 当η=0, 方程成为 Freundlich equation, 当 βη $k_s$ =0来描述。

流动带双点吸附的双重孔隙模型(物理和化学非平衡)

Henry	=0
Frac_M	1型吸附位点在移动区内的比例，考虑化学非平衡时的瞬时吸附位点[-]。
Omega	非平衡吸附一阶速率系数w [T-1]

附件/分离模型:

D_Soil	砂粒直径，dc [L]。
iPsi2	阻塞类型，iPsi，在第二个吸收地点。 =0:无阻塞。 = 1: Langmuirian动力学。 = 2: ripening.。 = 3:随机顺序吸附模型。 = 4:深度相关的阻塞系数。
iPsi1	第一个吸附点也是如此
SMax2	(1)、(2)、(3)、(4)中β的smax模型。
AttachSolid2	$k_{a} \cdot [\mathrm{T}^{-1}]$ 。
DetachSolid2	$k_{d} \cdot \left[\mathrm{T}^{-1}\right]$
SMax1	参数为第一个吸附点的阻塞函数。
AttachSolid1	$k_{a} \cdot \left[\mathrm{T}^{-1}\right]$ ,
DetachSolid1	$k_{d} \cdot \left[\mathrm{T}^{-1}\right]$ 。

过滤理论用于计算附着系数

未完待续…

2.4.4 温度依赖

2.4.5 溶质运移边界条件

使用此命令指定上下边界条件的类型。

1、下列类型的上边界条件是可用的:

浓度边界条件(用户指定边界处液相浓度)
浓度通量边界条件(用户指定渗入水的液相浓度)
易挥发溶质的停滞边界层

当用户为易挥发的溶质选择一个不流动的边界层时，他/她必须指定该边界层的厚度和上面的浓度。

2、下边界条件的种类如下:

浓度边界条件(用户指定边界处液相浓度)
浓度通量边界条件(用户指定渗入00水的液相浓度)
自由排水(零坡度)

3、溶质通量与溶质浓度的边界条件

第三类(Cauchy)柯西)边界方程规定了在边界处的溶质通量(不是浓度)。然后，由于与水(和溶质)混合最初出现在剖面，一个人不会立即获得边界上的c0浓度。由于规定了溶质流量，就可以完全控制质量平衡和溶质进入输运域的量。

如果使用第一类型(狄利克雷)BC，则规定边界处的浓度(而不是流入域的通量)。由于溶质通量由对流通量和弥散通量组成，所以这两种通量在边界上都是活动的。初始浓度梯度大，色散通量大。在这种情况下，进入定义域的通量要比使用第三类溶质通量边界条件时大得多。

我们总是建议使用第三种类型的BC，因为这是一种更实际的BC。狄利克雷BC不是物理的，由于上述原因，它不能保持质量(例如，van Genuchten和Parker, 1984;Leij等，1991)。只有在，例如，有大量的污染物与输运区域接触的情况下，人们才能使用它，因此，人们可以假定边界浓度是固定的。

4、初始条件是总浓度(而不是液体)

溶质运移的初始条件是根据总浓度[溶质质量/土壤体积]而不是液体浓度[溶质质量/水体积]来规定的。线性吸附时液相浓度计算如下:

\begin{aligned} S & = θ c + ρ s + a_{v} g = θ c + ρ k_{s} c + a_{v} k_{g} c = c (θ + ρ k_{s} + a_{v} k_{g}) \\ c & = \frac{S}{θ + ρ k_{s} + a_{v} k_{g}} \end{aligned}

对于非线性吸附，则求以下非线性方程的根:

S = θ c + ρ s + a_{v} g = θ c + \frac{ρ k_{s} c^{β}}{1 + η c^{β}} + a_{v} k_{g} c

对于两个动力学吸附点模型(仅实现线性情况，即假设分布系数ks定义为:

k_{s} = \frac{θ k_{a}}{ρ k_{d}}

5、非平衡相浓度最初与平衡相浓度处于平衡状态

对于双重孔隙度模型，不动含水率浓度将等于动含水率浓度。

动力学吸附点处的吸附浓度设为:

s = (1 - f) k_{s} c

对于线性吸附:

s = (1 - f) \frac{ρ k_{s} c^{β}}{1 + η c^{β}}

对于非线性吸附。对于两个动力学吸附点模型，吸附浓度设为:

s = \frac{θ k_{a}}{ρ k_{d}} c

$k_a$ $k_d$ 分别为附着系数和脱离系数。

目录