HYDRUS公开课-第一部分

主讲人：王珍Hydrus软件基础及应用培训第1部分：HYDRUS软件理论基础及技术手册12个人简介王珍，博士，中国水利水电科学研究院工程师教育经历：2004.9-2010.6西北农林科技大学农业水利工程本硕连读2010.9-2014.6中国水利水电科学研究院水利水电工程博士学术成果：已发表学术论文17篇，其中第一作者发表论文14篇，SCI论文3篇，EI论文4篇，其中与HYDRUS有关论文3篇，分别发表于AgriculturalWaterManagement（SCI），水利学报（EI），农业工程学报（EI）。汇报提纲2饱和/非饱和土壤水分运动1软件介绍3非稳定流中的溶质运移5物理问题描述6程序输入与输出文件的解读与说明34热量传导过程模拟41.软件介绍软件开发者软件开发者UniversityofCaliforniaRiverside教授CzechTechnicalUniversityandCzcchAcademyofSciences（布拉格，捷克）研究方向：土壤水溶质运移过程，平衡/非平衡物质转化迁移，田间尺度空间变异，土壤参数估算方法等已发表200余篇论文VadoseZoneJournal，JournalofHydrologicalSciences，JournalofHydrology等期刊副主编JiríŠimůnek51.软件介绍-应用方向农业领域降雨灌溉蒸发蒸腾根系吸水毛管上升深层排水施肥灌溉农药有机污染物病原体61.软件介绍-应用方向工业领域工业污染市政污染污水处理池污染物修复。。。71.软件介绍-应用方向生态领域碳的储存和释放热量交换与通量养分传输土壤呼吸微生物过程。。。81.软件介绍-结构91.软件介绍-结构适应于不同土壤饱和程度的Richards方程描述水分运动多种形式的土壤水力特性表达参数滞后效应模块根系吸水项优先流随水传递的热传导方程水分运动热量传导1.软件介绍-结构10液态形势下的对流-弥散方程气态形式下的扩散线性/非线性固态和液态转化液态和气态线性平衡放映0阶/一阶动力学反应过程非平衡溶质运移根系对肥料的吸收溶质运移1.软件介绍-历史发展111.软件介绍-历史发展(HYDRUS-1D)122.饱和/非饱和土壤水分运动132.1土壤水分运动基本方程2.2初始条件2.3边界条件2.3根系吸水模块2.1土壤水分运动基本方程14土壤中水流运动-Darcy’sLaw非饱和土壤中水流运动2.1土壤水分运动基本方程15土壤中水流运动控制方程-Richards方程θ—土壤体积含水率[L3L-3]h—土壤压力水头[L]S—根系吸水项或其它源汇项[T-1]K—非饱和土壤水力传导度[LT-1]KijA—各向异性张量的分量[-]xi—空间坐标[L]z—垂向坐标（向上为正）[L]t—时间[T]2.1土壤水分运动基本方程16土壤水分特征曲线2.1土壤水分运动基本方程17土壤水力传导度2.1土壤水分运动基本方程18土壤水分特征曲线方程BrooksandCorey（1964）：vanGenuchten（1980）：θs—土壤饱和含水率θr—土壤残余含水率α,n,h0,σ—经验参数Se—有效饱和度2.1土壤水分运动基本方程19土壤水力传导度BrooksandCorey（1964）：vanGenuchten（1980）：θs—土壤饱和含水率θr—土壤残余含水率α,n,h0,σ,l—经验参数Se—有效饱和度Ks—饱和水力传导度2.1土壤水分运动基本方程20水分特征曲线导水率变化曲线2.1土壤水分运动基本方程21RETC软件界面作用：由已测定的含水率与土壤水吸力的对应散点拟合土壤水分特征曲线参数拟合方程：BrooksandCorey或vanGenuchten模型2.1土壤水分运动基本方程22作用：通过神经网络对水分特征曲线和饱和导水率进行预测输入参数：砂粒、粉粒、粘粒、容重、33和1500kPa时的土壤含水率ROSETTA软件界面2.2初始条件230(,,,)(,,)0hxyzthxyzfort控制方程：2.3边界条件24压力水头边界（Dirichlet边界）：通量边界（Neumann边界）：水力梯度边界：2.3边界条件25大气边界（降雨和入渗过程）：2.3边界条件26大气边界（蒸发过程）：2.3边界条件27渗出面边界：排水管边界：2.3边界条件28随时间变化边界：降雪边界地下滴灌滴头出水边界地面滴灌湿润区变化边界2.4根系吸水模块29控制方程：b(z)—标准化根系分布系数Tp—潜在蒸腾量[LT-1]Sp—潜在蒸腾量[T-1]Lt—根区在土壤表面的宽度[L]2.4根系吸水模块30控制方程：b(z)—标准化根系分布系数Tp—潜在蒸腾量[LT-1]Sp—潜在蒸腾量[T-1]S—实际蒸腾量[T-1]α—根系吸水胁迫系数[-]4443413211121≤≤≤≤≤≤≤≥0hhh-hhhhh-h3α(x,z,h)=hhhh-hhhh2h-h0hh2.4根系吸水模块31标准化根系分布系数计算：Xm—x方向根系区分布最大长度x*,y*,z*—x,y,z方向根系密度最大处对应的横向坐标px,py,pz—表征x,y,z方向根系对称性的数值参数3.非稳定流中的溶质运移323.1溶质运移基本方程3.2初始条件3.3边界条件3.1溶质运移基本方程33一维溶质运移方程：θ—土壤体积含水率[L3L-3]c—溶液中溶质浓度[ML-3]s—吸附态溶质浓度[MM-1]ρ—土壤容重[ML-3]D—弥散系数[L2T-1]q—溶液通量[L]Φ—反应速率[ML-3T-1]t—时间[T]3.1溶质运移基本方程34物质反应转化链示意图w,s,g—下标分别表示液态、固态和气态c,s,g—分别表示液态、固态和气态溶质浓度3.1溶质运移基本方程35常见模拟溶质放射性元素农田氮素转化农药等有机污染物3.1溶质运移基本方程36w,s,g—下标分别表示液态、固态和气态c,s,g—分别表示液态、固态和气态溶质浓度3.1溶质运移基本方程37溶质扩散系数Dd—溶质在自由水中的分子扩散系数τ—曲率系数DL，DT—径向和横向弥散系数()AAarAAurETTRTTTra=aexp随温度变化的溶质参数aT—受温度影响后的溶质参数ar—未受温度影响时的参数参考值3.2初始条件380(,,,)(,,)0hxyzthxyzfort控制方程：3.3边界条件39第一类边界条件（Dirichlet边界，边界溶质浓度已知）：第三类边界条件（Cauchy边界，边界上溶质通量已知）：第二类边界条件（Neumann边界）：4.热量传导过程模拟404.1控制方程4.2初始条件和边界条件4.1热传导控制方程41控制方程（Sophocleous，1979）λij(θ)—土壤表观导热率C(θ)—多孔介质比热容Cw—液体比热容4.1热传导控制方程42导热率计算λ0(θ)—无水流交换条换条件下的土壤-水混合多孔介质导热率λL,λT—径向和横向热量传导系数4.2初始条件和边界条件43初始条件第一类边界条件（Dirichlet边界，边界温度已知）：(),,,()xyzt0DTx,y,z,t=Tforx,y,z第三类边界条件（Cauchy边界，边界上温度通量已知）：()CjTxijiwii0wii-n+Tcqn=TCqnforx,y,z第二类边界条件（Neumann边界）：()NjTxijin=0forx,y,z5.问题描述445.1时间控制5.2有限单元网格的划分5.3土壤分层和土壤分区的划分5.4边界条件的编码5.1时间控制45HYDRUS包括了三种不同的时间离散：（1）数值方法上的时间离散，（2）边界条件处理时的时间离散，（3）提供模拟成果输出时的时间离散。第二种和第三种离散是相互独立的，一般包括了在输入文件中指定的可变时间步长，而第一种时间离散则以给定的初始时间步长开始，在每一个时间层中这一时间步长将根据以下规则自动调整，并有如下要求：①离散1的时间必须能对应上离散2、3的时间点②时间步长必须介于预先给定的一个范围之内5.1时间控制46③在某一时间层，若达到收敛时的迭代次数少于3次，则下一时间层的时间步长σt必须乘以一个预先指定的大于1的常数（通常在1.1和1.5之间），相应地，如果迭代次数大于7，则下一时间层的时间步长σt必须乘以一个小于1的常数（通常在0.3和0.9之间）④在任何一个时间层如果迭代数大于一个给定的最大值（通常在10和50之间），该时间层的迭代将终止，并将步长设置为σt/3再重新开始进行该时间层的迭代过程5.2有限单元网格划分47在程序将计算区域剖分成有限个四边形或三角形（三维条件下为四面体、五面体或六面体）单元组合的网格系统，将每个单元顶点作为计算节点按逆时针进行编号，并在运行中自动将四边形处理成两个三角形单元单元边所形成的横向线应该沿着最短尺度的方向切分网格系统，这些横向线必须连续无交叉，但不必是直向。节点沿着每一条横向线向同一个方向连续地编号，从1到NumNp（节点总数），单元也按照类似的方法进行编号。在横向线上节点数的最大值IJ通常决定了有限单元矩阵的有效大小（例如矩阵的带宽），为了减小存储和运行时间，IJ应该尽可能地小5.2有限单元网格划分48有限元网格单元的大小应该根据实际问题来进行相应调整。在水力梯度比较大的区域单元的尺寸应划分的相对小一些。水力梯度较大的区域通常位于内部源汇项的附近、土壤表面。因此我们建议在土壤表面或附近划分较小的网格单元。一般随着深度的增加压力水头变化减小，单元大小可划分的相对大一些。此外，单元大小也决定于土壤水力参数，例如，粗糙的土壤比均匀质地的土壤有较高的值和较小的值，因此要求更细致的单元划分。HYDRUS软件建议使用近似大小的单元划分以减少数值误差。对于轴对称形式的三维流，对称轴必须是网格的左边界，或者位于左边界的左边。对于作物根系区土壤的单元划分，没有特别的要求5.3土壤分层和分区划分49软件用一套从1到NMat（土壤类型的总数）的整数编码来区分每一种土层，水流区域内每一个节点都对应了一个土层编码区域内土层之间的分界面不能与单元的边界重合。如果某个单元上各节点的土层编码不同，该单元的土壤参数将按照有限元运算法则在该单元上自动取平均，这样可以在某种程度上平滑土壤界面必须为每一种土壤指定一套水力特性参数和溶质运移参数，也必须为每一个单元指定一个各向异性张量的主方向和局部坐标与全局坐标轴的夹角5.3土壤分层和分区划分50采用了标定系数来简单地描述渗流区域内非饱和土壤水分运动参数的空间变异性，用户必须为水流区域内的每一个节点指定三个标定系数HYDRUS计算了每一个指定的分区内水流和溶质的质量平衡，指定的分区不一定就是土层区域HYDRUS用一套从1到NLay（土壤分区的总数）的整数编码来区分每一个分区，水流区域内每一个单元都对应了一个分区编码。5.4边界条件编码51程序代码中，每一个节点n都被赋予一个Kode(n)编码。总的编码原则是：如果该节点在一个时间步长内水头值是已知的，则在该时间步长内该节点的Kode(n)为正值，反之如果其流量是已知的话，则其Kode(n)为负值或为零归纳了水流边界条件有以下几类：稳定的边界条件；可变的边界条件；渗出面边界；排水管边界；根系吸水；深部排水；自由排水5.4边界条件编码52稳定的边界条件——这种边界的值由该节点的初值决定，包括两种情况：如果该节点属于第一类边界点，则在以后的计算过程中其水头值一直维持其初值不变；如果该节点为第二类边界节点，则其流量值在以后的计算过程中维持其初值不变。节点的编码Kode(n)、节点的流量Q(n)、水头值h(n)初始设定详情见下表：节点类型Kode(n)Q(n)h(n)内部节点（非内部源汇点）00.0初值内部源汇点（内部水头边界点