跳到主要内容

HYDRUS-1D

HYDRUS模型和STANMOD模型发展概述

地下环境中流动、运移与生物地球化学过程建模

在地下环境中解决与水流、溶质运移及生物地球化学过程相关的问题,需要与具体应用相适应的建模工具。 对于某些问题,相对简单的解析或半解析模型即可求解;而其他更复杂的问题,则可能需要一维或多维的数值模型,以模拟水流、溶质运移及多种生物地球化学反应

为了在解决一般性和场地特定的环境问题时具备足够的灵活性,研究人员通常需要一个包含不同复杂程度的计算机模型工具箱。 在过去约 30 年间,美国农业部**盐渍化实验室(USSL)加州大学河滨分校(UCR)**合作开发了大量此类计算机程序,并向公众发布。本文旨在介绍其中最重要的一些程序,并讨论其若干应用。


HYDRUS 与相关模型

本文描述了 HYDRUS 模型及其相关软件包(见表 1)的发展历程、所涉及的主要过程以及典型应用。 我们首先重点介绍一系列数值型 HYDRUS 模型、其前身及后续改进或扩展版本(如 SWMS_2D、HYDRUS-1D、HYDRUS-2D、HYDRUS (2D/3D) 以及 HP1),这些模型由美国、捷克、以色列、荷兰和比利时等地的研究团队共同开发完成。

模型维度操作系统主要功能与过程简介参考文献
HYDRUS (DOS 3.0–6.0)1DDOS可变饱和多孔介质中的水流与溶质运移;采用 Richards 方程求解;后期版本加入根系吸水、线性吸附与降解项。Kool & van Genuchten (1991); Vogel et al. (1996); Šimůnek et al. (1998b)
SOILCO₂1DDOS模拟 CO₂ 在水气两相中的传输;考虑热传导、CO₂ 生成与根系吸收。Šimůnek & Suarez (1993c)
UNSATCHEM1DWindows模拟水流、热量、主要离子与 CO₂ 传输;含化学平衡与动力学反应;16 bit GUI。Šimůnek et al. (1996b)
HYDRUS-1D 2.0 - 4.01DWindows水、热、多溶质运移;根系吸水(含盐胁迫补偿);VG/MVG/BC 模型;可逆吸附、两区非平衡、温度依赖、逆问题;3.0 版增 Durner & Kosugi 水力函数、病毒/胶体/细菌运移、UNSATCHEM 模块。Šimůnek et al. (1998b, 2005)
SWMS_2D 1.0–2.02DDOS基于有限元的二维非饱和水流与溶质运移;线性溶质传输、根系吸水;2.0 版加入迭代求解器。Šimůnek et al. (1992, 1994)
CHAIN_2D2DDOS扩展 SWMS_2D;加入非线性溶质传输、串级一阶衰变、气体扩散、温度依赖、双位点吸附。Šimůnek & van Genuchten (1994)
UNSATCHEM-2D2DDOS多离子化学与 CO₂ 运移;含平衡/动力学反应、Pitzer 活动度模型。Šimůnek & Suarez (1993b)
SWMS_3D3DDOS三维变饱和水流与溶质运移;线性吸附;为 HYDRUS (2D/3D) 前身。Šimůnek et al. (1995)
HYDRUS-2D 1.0–2.02DWindows基于 SWMS_2D 与 CHAIN_2D;32 bit GUI;非平衡双孔隙流、滞后性、逆问题。Šimůnek et al. (1996a, 1999a)
HYDRUS (2D/3D)2D / 3DWindows统一 HYDRUS-2D 与 SWMS_3D;VG/MVG/Durner/Kosugi 函数;气体扩散、温度依赖、双孔隙流、根系补偿吸水、病毒/胶体/细菌运移、人工湿地模块;高端 3D GUI。Šimůnek et al. (2006b); Šejna & Šimůnek (2007)
DISK2D (轴对称)Windows反演张力盘入渗实验数据以估算土壤水力参数;32 bit GUI。Šimůnek & van Genuchten (2000)
STANMOD 2.01–3D (解析)Windows含 CFITM/CFITIM/CXTFIT/3DADE/N3DADE/CHAIN 等解析模型;顺序衰变、多维运移;32 bit GUI。Šimůnek et al. (1999b)
RETCNAWindows土壤水分特征与水力导度参数反演;VG、BC、Durner、Kosugi 模型;32 bit GUI。van Genuchten et al. (1991) / Released online
ROSETTA 1.0NAWindows分层神经网络 PTF 预测 van Genuchten–Mualem 参数(提供 5 级输入)。Schaap et al. (2001)
UNSODANAWindows不饱和土壤水力性质数据库(保留 ≈ 800 样品);用于 PTF 验证与模型输入。Leij et al. (1996)
HP1 1.01DWindows耦合 HYDRUS-1D 与 PHREEQC;多组分传输、平衡–动力学反应、热传导;用于重金属、放射性核素、爆炸物迁移等复杂反应。Jacques & Šimůnek (2005)

模型发展脉络

  1. UNSAT → SWMII → SWMS_2D → CHAIN_2D → HYDRUS-2D → HYDRUS (2D/3D):二维–三维演进线;逐步引入溶质链、热传导、化学反应与 GUI。
  2. HYDRUS (DOS) → HYDRUS-1D (Win) → HP1:一维方向;扩展至热–化学–生物–根吸多过程耦合。
  3. UNSATCHEM 系列:与 CO₂ 与主离子化学耦合。
  4. STANMOD、RETC、ROSETTA、UNSODA:为 HYDRUS 系列提供解析解、参数反演与数据库支撑。

模型简要介绍

UNSAT 与 SWMII

UNSAT 模型最早由 Neuman 等人 (1972–1975) 在以色列理工学院开发,是一个二维有限元模型,用于模拟变饱和介质中的水流(基于 Richards 方程),并考虑根系吸水和多种边界条件。 后来 Davis 和 Neuman (1983) 将其移植至个人电脑。此版本成为 Vogel (1987) 在瓦赫宁根大学开发的 SWMII 模型的基础。 SWMII 扩展了功能,引入 van Genuchten (1980a) 及其改进的土壤水力学模型(Vogel & Císlerová, 1988),整合根吸模型 (SWATRE),并引入尺度因子以处理非均质土壤和局部各向异性。

SWMS_2D

SWMS_2D (Šimůnek et al., 1992) 扩展了 SWMII 的能力,引入溶质运移模块。为提升计算效率,模型只在压力头变化区域进行迭代求解。溶质运移使用对流–弥散方程,允许线性吸附与一阶衰变过程Simunek et al., VZJ - 2008, Hyd…。


CHAIN_2D

CHAIN_2D 在 SWMS_2D 基础上进一步加入非线性溶质运移、串级一阶衰变、气体扩散、双位点吸附模型、温度依赖性等过程。可模拟挥发性污染物(如溴甲烷、1,3-二氯丙烯)的运移Simunek et al., VZJ - 2008, Hyd…。


UNSATCHEM-2D

Šimůnek 和 Suarez (1993b, 1994) 在 SWMS_2D 上扩展开发 UNSATCHEM-2D,用于模拟主要离子(Ca, Mg, Na, K, SO₄, Cl, NO₃ 等)迁移以及平衡与动力学化学反应,包括络合、阳离子交换、溶解–沉淀。 模型还考虑 CO₂ 的迁移与生成、温度与水分变化的影响,并采用 Debye–HückelPitzer 活度模型 计算离子活度Simunek et al., VZJ - 2008, Hyd…。


SWMS_3D

SWMS_3D (Šimůnek et al., 1995) 是 SWMS_2D 的三维扩展,使用线性四面体有限元求解 Richards 方程和对流–弥散方程。 由于三维问题计算量巨大,Hardelauf et al. (2007) 开发了并行版本 PARSWMS,利用 MPI、PETSc 和 PARMETIS 在多核系统上并行计算。

HYDRUS-1D

HYDRUS-1D 是一个基于 Windows 的一维有限元建模系统,用于模拟非饱和带水流、热传导和多溶质运移。 其图形界面(GUI)用于输入数据、剖分网格、分配参数、执行计算并可视化结果。用户可在界面中定义土层分布、根系吸水分布及初始条件,输出包含压力头、含水量、水/溶质通量、温度与浓度等随时间变化的剖面。

HYDRUS-2D

HYDRUS-2D 在 HYDRUS-1D 的基础上扩展到二维空间,包含大部分相同过程:根系吸水(考虑水分与盐分胁迫)多种水力学函数滞后性非线性及非平衡反应双孔隙运移参数反演。 V1.0 基于 SWMS_2D,V2.0 基于 CHAIN_2D,引入了 MeshGen2D 网格生成器,可定义任意不规则几何域或具有轴对称特征的三维区域

HYDRUS (2D/3D)

HYDRUS (2D/3D) 是 HYDRUS-2D v2.0 与 SWMS_3D完全重写与整合,支持二维与三维几何。 除继承前版功能外,还包括: (i) 双孔隙系统下的水流与溶质运移(优先生流)、 (ii) 补偿型根吸水与空间根分布、 (iii) Durner (1994)、Kosugi (1996) 的水力模型、 (iv) 病毒/胶体/细菌运移的吸附–脱附与过滤理论、 (v) 人工湿地模块(2D)、 (vi) 滞后模型与动态边界、 (vii) 粒子流可视化及累积通量计算等Simunek et al., VZJ - 2008, Hyd…。

界面特点: 新 GUI 采用高端三维图形库、支持多项目视图、多窗口结构、对象浏览、拖放操作、自定义显示样式与导出并行计算输入(PARSWMS)

HP1

HP1(HYDRUS + PHREEQC)将 HYDRUS-1D 的流动与运移模块与 PHREEQC 地球化学求解器耦合,模拟变饱和带中的多组分反应运移与热传递。 适用于重金属、放射性核素、有机污染物及爆炸物迁移等含复杂化学反应的场景,是 HYDRUS 家族中生物地球化学最完备的版本

STANMOD

STANMOD 是一组溶质运移解析模型工具包,包含 CFITM、CFITIM、CHAIN、CXTFIT、3DADE、N3DADE、SCREEN 等程序,用于实验数据拟合与解析解验证。 适用于一维至三维的溶质迁移与顺序衰变分析,具有独立 GUI 与可视化功能;在文献中常用于实验数据反演与模型对比。

土壤参数模型

用来估算土壤参数。

Rosetta

Rosetta 是基于分层神经网络的土壤转移函数 (PTFs) 模型。 提供五个层次的输入复杂度:

  • 1️⃣ 土壤分类平均参数;
  • 2️⃣ 砂/粉/黏比例;
  • 3️⃣ 加入体积密度;
  • 4️⃣ 增加 33 kPa 吸力的含水量;
  • 5️⃣ 增加 1500 kPa 吸力的含水量(永久萎蔫点)。

输出参数均附带不确定度估计,基于自助法 (bootstrap) 和神经网络联合训练。

UNSODA

UNSODA(UNsaturated SOil Hydraulic DAtabase)是一个包含约 800 个土样的数据库,存储实测的不饱和土壤水力性质(含水量曲线、水力导度、扩散率等)及相关属性(粒径组成、有机质含量等)。 数据库可用于查询、筛选、输出与分析数据,并支持 Brooks–Corey 与 van Genuchten 模型拟合;已广泛用于 PTF 推导与水力性质研究。

RETC

RETC 是专门用于分析土壤水力性质的独立程序,可依据实验数据或假设模型(van Genuchten、Brooks–Corey、Durner、Kosugi 等)拟合土壤水分特征曲线与水力导度函数。 在 HYDRUS 系列中被广泛用于土壤参数前处理及 PTF 验证。


一维有限元模型HYDRUS-1D介绍

HYDRUS 以有限元法求解 Richards 方程与对流-弥散方程,可在 1D 与 2D/3D 下模拟变饱和介质的水-热-溶质耦合过程,并通过反演同化观测数据。其模块体系(如 HPx、UNSATCHEM、Wetland、PFAS 等)支持从盐渍化、反应运移到湿地处理与边坡入渗的多场景应用,是当前变饱和带研究与工程设计的主流工具之一。

HYDRUS-1D(软件包)专指“一维”模型的完整包:计算内核 + 图形界面 + 示例/文档。用于 1D 的水流(Richards)、热传输、溶质运移及根吸等,由USDA-ARS U.S. 盐度实验室提供,长期作为免费/公共领域软件发布,目前版本为4.17仍免费开放,合并到HYDRUS5之后为收费版本。

HYDRUS-1D GUI(界面)Windows 图形化环境,负责建模前处理(几何/网格/参数/边界)→ 调用 1D 内核计算 → 结果可视化(剖面/时序/通量);官方文档常用“Windows-based modeling environment / interactive graphics-based interface”等表述。

Hydrus-2D 是一个基于 Microsoft Windows 的建模环境,用于分析可变饱和多孔介质中的水流和溶质传输。该软件包包括二维有限元模型HYDRUS2,用于模拟水、热和多种溶质在可变饱和介质中的运动。Hydrus-2D 已不再发行,并于 2007 年被 HYDRUS 2D/3D 完全取代**。

HYDRUS(2D/3D)是面向 变饱和多孔介质 的二维/三维数值模拟平台,采用有限元法求解 Richards 水流方程对流–弥散方程(溶质/热),带 Windows 图形界面用于建模、计算与可视化。与 1D 版本相比,它能处理平面、轴对称与完整 3D 的任意几何与非均质介质,并支持根系吸水、化学反应与多组分运移等过程。

HYDRUS 版本 5,合并了两个以前独立的软件包 HYDRUS-1D(版本 4,用于一维应用)和 **HYDRUS(2D/3D)(**版本 3,用于二维和三维应用),全部为收费。

image-20251025082728988

主要能力与物理过程

  • 方程体系:水流采用 Richards 方程;热传输采用基于 Fick 的对流—弥散(含传导与随水对流);溶质采用液相对流—弥散与气相扩散耦合方程,支持根系吸水(作为汇项)。ars.usda.gov
  • 反应/转化:含非线性与/或非平衡固—液反应液—气线性平衡零级生成,以及两类一级降解(含级联反应耦合)。可用于非饱和、部分饱和或饱和条件下的水—溶质联运分析。ars.usda.gov
  • 边界与介质:介质可非均质;流动方向可为垂直、水平或倾斜;支持规定水头/通量、大气控制、自由排水等边界,并允许土表在通量—水头边界间自动切换。ars.usda.gov
  • 土壤水力学:提供 van Genuchten (1980)Brooks & Corey (1964) 及改进型函数;滞后按 Scott (1983) 与 Kool & Parker (1987) 的经验模型实现;并可通过线性尺度变换表示剖面内的水力变异。ars.usda.gov
  • 作物与胁迫:根系生长按逻辑斯谛函数;水分/盐分胁迫函数可选 Feddes (1978)van Genuchten (1987)ars.usda.gov

数值与反演

  • 数值实现:采用 Galerkin 线性有限元隐式(反向)时间离散,提供自适应时间步长并控制 Courant/Péclet 数范围;水含量项按 Celia et al. (1990)质量守恒方法计算;并提供减小数值振荡的上风、人工弥散等选项。ars.usda.gov
  • 参数反演:内置 Marquardt–Levenberg 反演,可利用含水量、压力头、浓度、瞬时/累积通量等数据估计土壤水力与溶质反应参数,并可加入先验约束/惩罚函数(贝叶斯)ars.usda.gov

HYDRUS-1D GUI图形界面

Windows 图形界面(GUI)负责剖分与编辑、参数分配、任务执行、结果可视化。支持在任意深度设置观测点输出水分/压力头/温度/浓度时序;可绘制剖面分布图边界通量/累积量曲线;内置小型土壤水力性质目录一组教学算例(含正算与反演)ars.usda.gov

HYDRUS-1D 获取

页面提供软件下载与手册入口及HYDRUS-1D 教程;该页列示的外壳版本为 3.01(2005 年 4 月),模型内核为 HYDRUS v7.0;作者包括 J. Šimůnek, K. Huang, M. Šejna, M.Th. van Genuchten;页面“Last Modified”日期为 2020-09-09ars.usda.gov

参考文献:

HYDRUS 模型下载: USDA ARS

HYDRUS 模型:USDA ARS

HYDRUS-1D 教程:USDA ARS

土壤环境影响评价中的HYDRUS-1D

土壤环境影响评价中的HYDRUS-1D采用的版本为4.08版本,仅采用计算内核,可以满足土壤环境影响评价一维模拟要求,但不包含其所有功能。

DOS程序下载

说明文档下载


HYDRUS—1D界面版使用说明

1 HYDRUS—1D

1.1 介绍

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

θt=z[K(θ)(hz1)]S\frac{\partial\theta}{\partial t}=\frac{\partial}{\partial z}\left[K(\theta)\left(\frac{\partial h}{\partial z}-1\right)\right]-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版相比包括):

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

Hydrus-1D 3.0版的新功能(2005年4月发布;与2.0版相比,包括):

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

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

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

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

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

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

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

1.3.3 后期处理

后处理也在GUI中进行。

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

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

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

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

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

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

1.3.4 模型附带的示例

直接问题:

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

反问题:

  1. 一步式流出实验-来自Kool等的数据。(1987)
  2. 多步骤流出实验-来自Hopmans的数据
  3. 蒸发实验-来自Wendroth的数据
  4. 向上渗透
  5. 涉及磁滞的瞬态流
  6. 非线性阳离子吸附的溶质运移-赖和尤里纳克的数据
  7. 非线性阳离子吸附的溶质运移-Selim的数据
  8. 硝化链溶质运输
  9. 水平渗透-来自George Vachaud的数据
  10. 水平渗透和再分布-来自Vachaud的数据
  11. 沙柱中的排水-来自Vachaud的数据
  12. 草下田间土壤剖面中的水流-季节模拟
1.3.5 源代码

HYDRUS-1D直接计算模块的(开放)源代码已包含在GNU通用公共许可证中,其详细条件在此处列出。仅当您同意本文档中列出的条件时,才能下载HYDRUS-1D计算模块的源代码。我们欢迎您对本源代码提出任何意见或对其进行任何改进。每当您使用此源代码时,请参考HYDRUS-1D技术手册或以下手稿:[Simunek等,VZJ-2008,Hydrus和Stanmod models.pdf](https://www.pc-progress.com/Documents/Simunek et al., VZJ - 2008, Hydrus and 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中的模块相似(几乎相同)。可能的应用包括研究农业土壤的盐碱化/复垦,各种灌溉系统的可持续性以及采矿作业中盐水的处理。

image-20251024103529649

2 HYDRUS—1D GUI版使用说明

2.1 主过程

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

image-20251024103719941

  • 水流(water flow)

    • 蒸汽运输

    • 土壤表面积雪

  • 溶质运移(solute transport)

    • 标准溶质运移

    • 主要化学过程

    • PHREEQC生物地球化学模型

  • 热传导(heat transport)

  • 根系吸水(root water uptake)

  • 根系生长(root growth)

  • 反解(Inverse Solution)

2.2 基本信息

2.2.1 模型几何结构

(1)土壤结构输入参数(右侧)

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

(2)长度单位(Length Units)(左侧)

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

image-20251024104132660

2.2.2 时间信息

(1)时间单位(Time Units)

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

(2)时间离散(Time Discretization)

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

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

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

(3)时间变化边界状态(Time-Variable Boundary Conditions)

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

image-20251024104339684

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)。

image-20251024104524762

2.3 土壤水流模型

2.3.1 迭代规则(iteration criteria)

一般可以采用默认值。

1、迭代规则

  • 最大迭代次数(Maximum number of iteration),任意时间步长的最大迭代次数,同时用改进的Picard方法求解非线性Richards方程。建议的默认值是10。
  • 含水量迭代精度(water content tolerance),流区非饱和区节点的绝对含水量容差-。此参数表示在特定时间步长内连续两次迭代之间的含水量值的最大期望绝对值变化。建议的默认值是0.001。
  • 压力水头迭代精度(pressure head tolerance),流动区域饱和部分节点的绝对压头公差[L](推荐值为1 cm)。该参数表示在特定时间步长内连续两次迭代之间的压力头值的最大期望绝对值变化。建议的默认值是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/(即空气进口价值)。

image-20251024104928508

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)

image-20251024105313357

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、QsKs、l 六个参数。其中 l 一般取 0.5。其他参数与土质、土壤结构等有关系。参数取得是否合理和定解条件(初始条件和边界条件)是模拟是否能反映实际的关键。在土壤水力参数很难直接得知的情况下,需要通过其他方法推求参数,或是通过模拟结果与实际比较,来调节参数。

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

英文中文
QrResidual soil water content, θr\theta_{r} 土壤残余含水量
QsSaturated soil water content, θs\theta_{s} 饱和土壤含水量
AlphaParameter α\alpha in the soil water retention function [L1]\left[ \mathrm{L}^{-1} \right] 参数α\alpha 在土壤水分保持函数中
NParameter n in the soil water retention function参数n在土壤水分保持函数中
KsSaturated hydraulic conductivity, k_{s} $$\left[ \mathrm{L}^{-1} \right] 饱和导水率
lTortuosity parameter in the conductivity function [-]电导率函数中的弯曲度参数[-],一般取值0.5

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

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

英文中文
QaParameter θa\theta_{a} in the soil water retention function土壤水分保持函数中的参数**θa\theta_{a} **
QmParameter θm\theta_{m} in the soil water retention function土壤水分保持函数中的参数 θm\theta_{m}
KkMeasured value of the hydraulic conductivity, KkK_{k} , corresponding to θk\theta_{k} [L1]\left[ \mathrm{L}^{-1} \right] 对应**θk\theta_{k} ** [L1]\left[ \mathrm{L}^{-1} \right] 的水力传导率实测值KkK_{k}
QkSoil water content, θk\theta_{k} , corresponding to KkK_{k} 土壤含水量,θk\theta_{k} ,对应KkK_{k}

当Dumer模型使用时,

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

英文中文
wParameter w for material M [-]. Relative weighting factor for the subcurve of the second overlapping subregion.第二重叠子区域子曲线的相对权重因子,参数w
Alpha2Parameter a for material M [L1][ {L}^{-1} ] , for the second overlapping subregion.对于第二重叠子区域,分层的参数*\alpha,*M/L
n2Parameter n for material M [-], for the second overlapping subregion.对于第二重叠子区域,分层的参数n

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

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

英文中文
QrImParameter θr\theta_{r} for the immobile region of material M.
QsImParameter θs\theta_{s} for the immobile region of material M.
OmegaParameter ω\omega (mass transfer coefficient) for material M [T1][{T}^{-1} ] .

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

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

英文中文
QrImParameter θr\theta_{r} for the immobile region of material M.
QsImParameter θs\theta_{s} for the immobile region of material M.
AlphaIParameter α\alpha for the immobile region of material M [L1][ {L}^{-1} ]
nIParameter nn for the immobile region of material M
OmegaParameter ω\omega (mass transfer coefficient) for material M [T1][{T}^{-1} ] .

双重渗透率模型,

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

英文中文
QrFtParameter θr\theta_{r} for the fracture region of material M.
QsFrParameter θs\theta_{s} for the fracture region of material M.
AlphaFrParameter α\alpha for the fracture region of material M [L-1]
nFrParameter nn for the fracture region of material M
KsFrParameter ksk_{s} for the fracture region of material M [LT-1]
lFrTortuosity parameter in the conductivity function for the fracture region of material M [-]
wParameter ω\omega for material M [-]. w is the ratio of the volumes of the macropore or fracture domain and the total soil system.
BetaParameter β\beta (a shape factor that depends on the geometry) for material M [-].
GammaParameter γ\gamma (a scaling factor) for material M [-].
aParameter dd (an effective diffusion pathlength) for material M [L].
KsaParameter kak_{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)

image-20251024105710693

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 一般信息

image-20251024111023672

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

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

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

2、空间权重方案

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

在HYDRUS-2D中,上游加权是一个选项,它可以在模拟相对陡峭的浓度锋时最小化数值振荡的一些问题。为此,平流-色散方程的第二项(通量)并不是由常规线性基函数Φn\Phi_n ,而是使用非线性函数Φnu\Phi_{n}^{u} [Yeh和Tripathi, 1990]。赋权函数Φ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要打印到输出文件或显示在各种图形中的单元。采用的质量单位对计算没有影响。浓度单位一般应在[ML3][ML^{-3}]中给出,其中M为质量单位,L为几何信息对话框窗口中指定的长度单位。然而,由于浓度变量出现在每一项执政的溶质运移方程(Eq。技术手册的3.1和3.2),可以使用不同长度单位比用于定义几何和通量(例如,几何可能被指定在米浓度给出毫克/立方厘米)。在这种情况下,溶质通量(cq)将在单位MLc3LgT1ML_c^{-3}L_gT^{-1} ,LcL_c是长度单位(如厘米)用于定义浓度和LgL_g长度单元定义几何和通量(例如,米)。同样的溶质质量(cθVc{\theta}V)通过集成溶质在交通领域将MLc3Lg2ML_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 溶质运移参数

image-20251024112409208

土壤参数

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

溶质的具体参数

Diffus.W分子在自由水中的扩散系数,DwD_w [L2T1][L^2T^{-1}]
Diffus.G分子在土壤中的扩散系数 DaD_a* [L2T1][L^2T^{-1}]

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

Bulk.d.体积密度Bulk density, ρ[ML3]\rho [ML^{-3}]
Disp.M.矩阵中的纵向弥散性Longitudinal dispersivity, DLD_L [L]
Frac.M.在矩阵中分类为类型-1的吸附位点的无量纲分数,即,即考虑化学非平衡选择时具有瞬时吸附的位点[-]。
ThImob.M.基质中的固定含水量。当不考虑物理非平衡期权时,设为0。
Diffus.W在自由水中的分子扩散系数,DwD_w
Mass Tr.有效扩散系数[L2T-1]在基质与裂缝间溶质传递的传质系数中,Γs\Gamma_s [T1][T^{-1}]

地球化学模块[UNSATCHEM](mk:@MSITStore:C:\Program Files (x86)\PC-Progress\Hydrus-1D 4.xx\Hydrus1D.chm::/HYDRUS1D/UNSATCHEM_Module.htm)

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

image-20251024112500077

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

2.4.3 溶质反应参数

image-20251024112548459

类型描述
Kd吸附等温式系数,ks[M1 L3]k_{s} \cdot \left[\mathrm{M}^{-1} \mathrm{~L}^{3}\right]
Nu吸附等温线系数 v[M1 L3]v \cdot \left[\mathrm{M}^{-1} \mathrm{~L}^{3}\right]。 (手册中称h,朗缪尔系数)
Beta吸附等温线系数,β\beta (Freundlich指数)
Henry液气平衡分配常数,kgkg []
SinkWater1溶解相一阶速率常数,uwu_w [T1][T^{-1}]
SinkSolid1固相一阶速率常数,u_s$$[T^{-1}]
SinkGas1气相一阶速率常数,u_g$$[T^{-1}]
SinkWater1'溶解相的一阶速率常数,μw[T1]\mu_{w}^{\prime} \cdot \left[\mathrm{T}^{-1}\right]。,表示链式反应
SinkSolid1'固相一阶速率常数, μs[T1]\mu_{s}^{\prime} \cdot\left[\mathrm{T}^{-1}\right]。表示链式反应
SinkGas1'气相一阶速率常数,μg[T1]\mu_{g^{\prime}} \cdot \left[\mathrm{T}^{-1}\right],表示链式反应
SinkWater0溶解相零级速率常数, γw[ML3 T1]\gamma_{w} \cdot \left[\mathrm{ML}^{-3}\mathrm{~T}^{-1}\right]
SinkSolid0固相零级速率常数, γs[T1]\gamma_{s} \cdot \left[\mathrm{T}^{-1}\right]
SinkGas0固相零级速率常数, γg[ML3 T1]\gamma_{g} \cdot \left[\mathrm{ML}^{-3}\mathrm{~T}^{-1}\right]
Alpha单点或双点非平衡吸附的一阶速率系数,流动和不流动液相区间溶质交换的传质系数,为 ω[T1]\omega \cdot \left[\mathrm{T}^{-1}\right]

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

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

当\beta=1时,上述吸附方程成为 Langmuir equation, 当\eta=0, 方程成为 Freundlich equation, 当 \beta=1 和*\eta*=0, 它会导致一个线性吸等温线。无吸附的溶质运移用ksk_s =0来描述。

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

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

附件/分离模型:

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

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

未完待续…

2.4.4 温度依赖

image-20251024140841107

2.4.5 溶质运移边界条件

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

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

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

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

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

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

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

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

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

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

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

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

S=θc+ρs+avg=θc+ρksc+avkgc=c(θ+ρks+avkg)c=Sθ+ρks+avkg\begin{aligned} S &= \theta c + \rho s + a_{v} g = \theta c + \rho k_{s} c + a_{v} k_{g} c = c (\theta + \rho k_{s} + a_{v} k_{g}) \\ c &= \frac{S}{\theta + \rho k_{s} + a_{v} k_{g}} \end{aligned}

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

S=θc+ρs+avg=θc+ρkscβ1+ηcβ+avkgcS = \theta c + \rho s + a_{v} g = \theta c + \frac{\rho k_{s} c^{\beta}}{1 + \eta c^{\beta}} + a_{v} k_{g} c

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

ks=θkaρkdk_{s} = \frac{\theta k_{a}}{\rho k_{d}}

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

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

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

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

对于线性吸附:

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

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

s=θkaρkdcs = \frac{\theta k_{a}}{\rho k_{d}} c

式中,kak_akdk_d 分别为附着系数和脱离系数。

2.5 根系吸水模型(不做)

2.5.1 模型设定

image-20251024141134974

2.5.2 参数

image-20251024141146917