Hydrus1D简明使用手册.docx

1、Hydrus1D简明使用手册用HYDRUS-1D模拟剖面变饱和度地下水流（简明手册）王旭升中国地质大学(北京)！HYDRUS-1D是一个共享专业软件，用于模拟一维变饱和度地下水流、根系吸水、溶质运移和热运移。本手册只介绍应用HYDRUS1D模拟垂向剖面水流和根系吸水的操作方法。1. 如何获取HYDRUS-1DHYDRUS-1D由位于欧盟捷克的PC-Progress工程软件开发公司发行，用户可以登录该公司首页: 。为了下载HYDRUS-1D，应先注册成为用户，然后下载Hydrus-1D的安装文件：。这个文件对应目前HYDRUS-1D的最高版本。2. 版权声明 HYDRUS-1D的作者为:(1)

2、J. Simunek, Department of Environmental Sciences, University of California Riverside, Riverside, California, USA.(2) M. Sejna, PC Progress, Prague, Czech Republic.(3) . van Genuchten, Department of Mechanical Engineering, Federal University of Rio de Janeiro, Rio de Janeiro, Brazil.感谢他们提供了一个如此精美而又免费

3、使用的专业软件，帮助我们从事有关的科学和教育工作。当你运行解压文件后，会在您的电脑中产生一个安装目录，其中包含可执行文件。运行这个文件即可安装HYDRUS-1D软件。当您安装HYDRUS-1D时，象安装其它软件一样，会出现一个许可协议，从中可知本共享软件也受到美国法规的保护。3. 参考资料 HYDRUS-1D安装之后，在软件运行目录下有HYDRS-1D 文件。从这个文件您可以了解到HYDRUS-1D的一些技术细节，如水流、溶质运移、热流的方程、一些处理专门问题的模型、输入输出文件等等。有一个Examples目录，包含大量的模拟算例可供参考。用户还可以参考以下文献： imnek, J., M.

4、Th. van Genuchten, and M. ejna, Development and applications of the HYDRUS and STANMOD software packages, and related codes, Vadose Zone Journal, doi:, Special Issue ”Vadose Zone Modeling”, 7(2), 587-600, 2008. Jacques, D., J. imnek, D. Mallants, and M. Th. van Genuchten, Modeling coupled hydrologic

5、al and chemical processes: Long-term uranium transport following mineral phosphorus fertilization, Vadose Zone Journal, doi:, Special Issue ”Vadose Zone Modeling”, 7(2), 698-711, 2008. | imnek, J. and M. Th. van Genuchten, Modeling nonequilibrium flow and transport with HYDRUS, Vadose Zone Journal,

6、doi:, Special Issue ”Vadose Zone Modeling”, 7(2), 782-797, 2008. 这些文献都可以从 下载。4. HYDRUS-1D的WINDOWS界面运行HYDRUS-1D，可以看到一个Windows的界面如下：图 1 所有的前后处理在界面中一目了然，左边是前处理工具，右边是后处理工具。其中前处理的各项功能如下图所示。图 25. 设计模型 在使用HYDRUS-1D之前，您需要对饱和-非饱和水流模拟的基本原理有所了解，并设计出自己想做的模型，准备好数据。一个剖面水流模型通常包含以下几个要素：(1)土壤剖面从地面算起的深度，准备模拟那个时间段的水分变

7、化。(2)土壤分几层，每层土壤的渗透性参数和水分特征曲线是怎样的。(3)根系是怎么分布的。(4)是否已经确定地面降雨入渗、蒸发蒸腾的信息，特别是它们随时间的变化。(5)是否已经确定剖面底部的状态属于哪种类型的边界条件。下面是一个参考模型的设计图:图 36. 使用HYDRUS-1D创建模型 打开HYDRUS-1D软件，选择”File/new”菜单，新建一个模型。在name一栏中输入本模型的名称”test”，更改模型存放的目录。￥图 4需要注意的是，HYDRUS-1D模型本身在计算机中就表现为一系列的输入输出文件，它们存放在与模型名称一致的目录中。本例中，软件会自动创建一个名称为”test”的目录

8、，而”CATOOLSHYDR1DProjects”中除了test目录之外，还有一个文件。这是一个模型项目(project)文件，告诉软件下次到哪里去寻找模型。 模型创建之后，会显示前处理和后处理窗口(图5)。由于是新模型，还没有任何模拟结果，所以后处理窗口是空白的。图 57. 输入模型控制信息 首先，在前处理窗口双击Main Processes，在弹出的对话框中输入模型的描述: a test model. 然后在Simulate 一栏中选中Root Water Uptake，表示想处理根系吸水问题。电击OK之后，前处理窗口将增加处理根系吸水的工具条。图 6 下一步，是输入模型的几何信息和土层划

9、分信息。在前处理窗口双击Geometry Information，在弹出的对话框中输入如图7所示的数据。图 7 接下来输入时间信息，在前处理窗口双击Time Information，会弹出一个对话框（图8）。图 8这个对话框中提供了一些灵活的选项来处理上边界条件的变化，下面简要加以说明： (1) 蒸腾量的每日周期变化 HYDRUS-1D 可以使用一个经验公式来处理每天24小时潜在蒸腾量的变化，设某天的潜在蒸腾量为（例如用Pemman公式获取的, cm/d），则其中Tp(t)是瞬时潜在蒸腾量，t为时间。模型假设早上6点之前以及晚上18点-24点的蒸腾量总和只占全天蒸腾量的1%。注意本例中蒸腾量的

10、单位是 cm/d。) (2) 降水量的周期变化如果在你的模型中降水量是周期性变化的，HYDRUS-1D也可以用一个公式来处理其中是周期 t内的平均降雨量。(3) 使用气象数据也可以在HYDRUS-1D中输入气象数据，它将自动利用这些数据计算潜在蒸散量ETp。可以选择FAO组织推荐的Penman-Monteith公式，也可以选择Hargreaves公式。这些公式需要辐射、气温、湿度之类的气象数据。模型的另一个控制信息是对模拟结果的输出如何进行设置。在前处理窗口双击Print Information工具条，弹出一个对话框。本例中确定输出30组模拟结果，每天输出1组。图 98. 水流模型迭代计算参数

11、HYDRUS-1D是采用迭代法来处理非线性Richards方程的。在前处理窗口双击Water Flow- Iteration Criteria工具条，弹出一个设置迭代参数的对话框（图10）。迭代控制参数的设置具有高度的专业技术性，除非特别了解，一般可以使用默认值。如果模拟结果出现不收敛的情况，需要对最大迭代次数、迭代精度等参数进行调整，但是在缺乏经验的情况下很难操作。图 10HYDRUS-1D 采用自动控制时间步长的方法来处理迭代的收敛性。对于每个时步，如果迭代次数太多，就缩小时间步长；如果没经过几次迭代就达到收敛精度，则适当增大时间步长。9. 水流模型土壤水力特性模型水分特征曲线是非饱和土壤

12、的重要物理性质，HYDRUS-1D提供了几种方法来处理与之有关的参数。在前处理窗口双击Water Flow- Soil Hydraulic Properties工具条，弹出一个设置水力特性模型的对话框（图11）。图 11在一般情况下，选择单孔介质模型，并选择用van Genuchten-Mualem公式处理土壤的水力特性就可以了。如果还要模拟溶质运移，可能需要考虑双重介质模型。双重介质在同一个点有两个孔隙度或两个渗透率，相当于两种介质的混杂。双重介质模型能够模拟这两种“介质”之间的水分和盐分交换。10. 水流模型土壤水分特征曲线在前处理窗口双击Water Flow- Soil Hydrauli

13、c Parameters工具条，弹出一个设置水分特征曲线参数的对话框（图12）。本例中选择van Genuchten-Mualem公式处理水分特征曲线，其中 , n, l均为控制因子。HYDRUS-1D软件中提供了一组土壤经验参数库，可供用户参考。本例中两层土壤的参数直接从数据库中调出：第1层对应Sandy loam，第二层对应sand。图 12&在输入参数时，请注意参数的单位。11. 水流模型边界条件在前处理窗口双击Water Flow- Boundary Conditions工具条，弹出一个设置边界条件的对话框（图13）。图 13上边界条件有6种类型，下边界条件有8种类型。边界类型的确定需

14、要考虑实际条件，在本算例中，上边界选择大气边界条件，在降雨量很大时地表可以产生积水。植被蒸腾量和土壤蒸发量分开处理，HYDRUS-1D推荐使用一个经验公式来把潜在蒸散量分割为蒸腾潜力和土壤蒸发潜力：其中ETp为潜在蒸散量（可以使用Penman-Monteith公式处理气象数据得到, cm/d），Tp为潜在蒸腾量(cm/d)，Ep为土壤潜在蒸发量(cm/d)，LAI是叶面积指数，k为消光系数，取决于太阳角度、植被类型及叶片空间分布特征。SCF是一个中间参数，即土壤覆盖度(Soil cover fraction)。在阔叶植被发育的情况下，消光系数的经验值为k=。12. 水流模型定水头或通量边界设置

15、如果边界条件中包含定水头或定通量的边界，则在前处理窗口双击Water Flow- Constant BC工具条，弹出一个设置边界数据的对话框。 本算例模型中，下边界为定流量边界，实际上就是隔水边界，因此直接输入0即可。13. 根系吸水吸水模型在前处理窗口双击Root Water Uptake- Models工具条，弹出一个处理根系吸水模型的对话框（图14）。图 14HYDRUS-1D使用水分胁迫和盐分胁迫模型处理根系的吸水。对于水分胁迫模型，计算公式为其中Tp是潜在蒸腾量(cm/d)，Ta是实际蒸腾量(cm/d)，S(x)是吸水强度函数(cm/，注意x坐标实际表示深度)， (h)是水分胁迫函数

16、，h为土壤压力水头(cm)，b(x)是根系吸水分配（密度）函数，LR为根系层的深度。水分胁迫函数有2种经验表示方法，即Feddes模型和S-Shape模型。Feddes模型是一个梯形函数，只需要知道h值。而S-Shape模型把水分胁迫和叶片气孔的压力水头联系起来，需要知道气孔压力水头的数值h 。土壤的湿润度可以表示为 但是如果直接用这种方法来计算实际蒸腾量有一定的问题。植被其实可以调节不同深度的水分胁迫响应特征；某个深度土壤干燥吸不上水，植被可以加大在比较湿润的土层的吸水量，以补偿不足。这种现象称为补偿吸水。为了模拟根系补偿吸水，HYDRUS-1D提供了一种简化的模型，即如果湿润度高于某个临界

17、值( c)，植被根系可以通过补偿机制充分吸水达到潜在蒸腾量。如果湿润度低于这个临界值，补偿机制受到抑制，发生整体的水分胁迫，根系吸水总量将低于潜在蒸腾量，并正比于湿润度。如果不考虑这种补偿吸水机制，可以令 c =1。14. 根系吸水水分胁迫参数在前处理窗口双击Root Water Uptake- Water Stress Reduction工具条，弹出一个处理水分胁迫参数的对话框（图15）。本算例中直接从数据库中调入Wheat的经验值。￥图 15 根系水分胁迫Feddes模型参数15. 输入可变边界条件的信息在前处理窗口双击Variable Boundary Conditions工具条，弹出一

18、个处理时间序列数据的对话框（图16）。本算例中，在步骤(11)中已经把地面处理大气边界，同时又选择使用消光系数法划分植被蒸腾和土面蒸发，因此需要输入每天的降水、潜在蒸散量、叶面积指数等数据。还有一个需要输入的数据是最小压力水头值，即地面土壤达到最干燥状态时的压力水头。从理论上讲，当土壤十分干燥时，吸力很大，而液态孔隙水的压强很小，与空气湿度保持平衡关系，因此有其中hA为最小压力水头，Hr为空气绝对湿度，RT/Mg为空气的摩尔气体常数。空气湿度虽然可以通过气象数据得到，但这里公式需要的是近地面的空气湿度。一般情况下，取饱和水汽湿度是可取的，因为2 cm深度以下土壤空气的湿度往往都是饱和的，只不过

19、随温度发生变化。因此，可以根据近地面气温的变化来推算地表土壤的空气湿度（饱和水汽湿度），再换算成压力水头。HYDRUS-1D中需要输入的是最小压力水头的绝对值，缺省值为hCritA=|hA|=106 cm=104 m这个数值只会对土壤蒸发起作用。HYDRUS-1D建议：hCritA 所对应的土壤含水量应该至少比残余含水量大，在模拟根系吸水的情况下，hA 还应该低于图15中的P3。否则（hAP3），当根系吸水的临界值压力水头(P3)和地面蒸发的最小压力水头 (hA) 满足时,会导致回流(inflow)现象，这是不合理的。除非存在特别干燥的情况，模型一般不需要仔细处理这些问题。 图 16表1时间(

20、d)降水量(cm/d)ETp (cm/d)hCritA (cm)LAI时间(d)降水量(cm/d)ETp (cm/d)hCritA (cm)LAI10100000160100000201000001701000003010000018010000040100000190100000【50100000200%10000060100000210！100000701000002201000008100000!230100000901000002410000010010000025100000110*10000026100000120$100000271000001331000002801000001

21、40100000290100000！15100000300100000数据可以先在Excel中准备好，如表1。这些数据可以拷贝到图16的电子表格中。这些数据显示的降水量、蒸散潜力和叶面积指数变化特征如图16右图所示。在第17日由于庄稼收割，叶面及指数大幅度下降。 16. 编辑土壤剖面使用图形界面在前处理窗口双击Soil Profile- Graphical Editor工具条，程序将弹出一个处理土壤剖面的软件（图17）。这个图形软件的使用比较简单，我们需要注意的是在Conditions菜单下面有很多子菜单，包含处理各种问题的功能菜单。 首先要做的事情，是确定把土壤剖面离散化为多少个节点。本算例

22、土壤模型深度为3 m，我们希望节点间距达到 1 cm，因此需要301个节点。选择菜单Conditions/ Profile Discretization，在下拉工具条中把Number修改为301。图 17 接下来，确定土层的分布，本模型有2个土层。缺省的土层编号为index=1，就是图12中的细砂壤土层。选择菜单Conditions/ Material Distribution，在下拉工具条中使用Edit condition，把下部土壤层设置index=2，这个土层编号为2，实际上就是图12中的中砂层。、图 18再就是确定根系随深度的分布，根系分布函数b(x)是一个很特殊的函数，它满足以下条件

23、其中LR为根系层厚度, x是深度； z是节点间距，bn是每个间距中的根系分布函数值，M是根系层占节点数。本算例中根系层的厚度为1 m，假设根系分布函数为线性，并有, x 1容易证明上式满足积分为1的条件。设置方法为，选择菜单Conditions/ Root Distribution，在在下拉工具条中使用Edit condition，划定根系层范围(100 cm)，把顶部数值设置为0，底部数值设置为2，让程序自动进行线性差值。图 19下一步处理初始条件，选择菜单Conditions/ Initial Condition。本算例中初始地下水位高于底板100 cm，假设土壤剖面初始状态是静力平衡态，

24、则模型底部的压力水头为+100 cm，地面的压力水头为 200 cm。使用Edit Condition工具条，分别设置顶部和底部的压力水头，并让程序自动插值形成初始条件。最后，选择菜单Conditions/ Observation Points，使用Insert工具条添加若干观察点，有必要在靠近模型底部的位置加一个观察点以判断地下水位的变化。关闭图形程序退回主程序。17. 编辑土壤剖面使用表格在前处理窗口双击Soil Profile- Summary工具条，程序将弹出一个表格（图20）。在这个表格中可以进一步修改土壤剖面数据。图 2018. 运行模型我们已经把所有需要的数据都输入了模型中，下面

25、就可以运行模拟计算程序了。选择菜单Calculation/Execute Hydrus，在弹出的对话框中选择OK，则会进入模拟计算模块。软件将调用HYDRUS-1D的核心程序进行计算，并输出有关的信息，见图21。图 21一般情况下，如果迭代计算的收敛性较好，模拟时期又不是很长，将很快计算完毕。但是，一旦出现难以收敛的情况，则可能等待很长时间才有结果，或者非正常中断。本算例模型在很短的时间内就能够完成计算。19. 察看结果模拟运算完毕之后，模拟结果可以通过后处理窗口中的各个工具条来察看，见图22。图 22 模型模拟的结果包括土壤含水量、压力水头、实际根系吸水（蒸腾）量、实际地表通量、土壤水分存贮

26、总量等的变化。20. 输出结果如果我们想把模拟结果保存为其它文件格式，用其它软件进行后处理，则可以做进一步的转换。HYDRUS-1D的结果图包含曲线信息，点击右键，会看到有一个Edit Chart Data菜单，用这个菜单可以通过电子表格输出模拟结果，见图23。这个电子表格的个别数据可以修改后重新显示图像。图像可以打印。图 23 然而，图23中的电子表格不具备输出数据的功能。为此，你可以打开模型目录（文件夹），找到以OUT为后缀的文件，图24。这些文件都可以用记事本打开察看数据。图 24有大量的计算结果保存在文件中，查阅HYDRUS-1D的用户手册，就知道文件保存了那些结果变量，例如根系吸水强度vROOT等等。附表 给出了变量缩写名称的含义。