基于不同灌水模式下HYDRUS模型土壤水盐运移模拟研究

王志坤

(辽宁省水利水电科学研究院有限责任公司，辽宁 沈阳 110003)

在我国东北部地区以辽宁省辽河流域为例，灌溉对于改善该区域生态环境和促进农业生产发展具有重要作用[1]。但是，灌溉可造成一系列的负面作用，如土壤出现次生盐渍化、地下水位显著升高等危害。辽河流域为解决以上问题分别采取了加大灌溉定额、建设排水系统以及利用秋浇淋洗盐分等措施并取得了良好的效果[2]。近年来，随着社会经济的发展和人口的急剧增加水资源短缺日益加剧，为此，辽河流域开始大规模地改造和实施节水灌溉工程。在干旱半干旱条件下辽河流域原有的地表水、地下水与土壤水之间的变化规律随着节水灌溉的推广而发生显著改变，并使生态环境各要素之间发生连锁反应直至达到新的平衡，而且新的平衡状态下也可能引发自然生态一系列问题的出现[3]。据此，研究节水灌溉可能造成的生态环境影响以及节水灌溉条件下区域水盐运动变化规律对于保障区域节水灌溉的持续健康发展具有重要意义，本文利用所构建的HYDRUS模型分别对不同灌水定额下沟灌和畦灌2种模式的土壤水盐变化规律进行模拟研究[4]。

1 区域概况及试验方法

1.1 自然概况

选取辽宁省辽河流域某试验站为试验地进行研究，试验区属于半干旱半湿润大陆性季风季候区，区域内雨热同季、冬季寒冷干燥、夏季高温多雨，年际降水量较少而蒸发量较大，全年无霜期较短，土壤处于较长的封冻期属于典型的季节性冻土区。多年平均降水量和蒸发量分别为350～450mm、700～900mm，其土壤主要以沙土壤和粉黏土为主，其中盐渍化土和辽河淤灌土为该区域主要土质类型，沙壤土主要位于浅水位以上属于非饱和带的土壤质地，而粉黏土主要位于不同厚度的夹层中[5]。试验区不同夹层土壤的主要理化性质指标见表1。

1.2 试验区布置

试验区域农作物主要以玉米、水稻和葵花为主，传统的地面灌溉模式是该区域田间灌水主要方式并以大规格的畦田灌溉为常用的地面灌溉方式[6]。在当前经济条件和作物种植下沟灌和小畦灌为辽河流域试验区主要田间节水灌溉方式。本文选择畦灌和沟灌2种灌水模式并以玉米作物为例分别采用3种灌水定额进行试验研究，其中畦灌和沟灌灌水定额分别为3.6、4.3、5.2m3/hm2和2.4、3.2、3.8m3/hm2，试验编号及处理结果见表2。灌水沟长约35m，玉米株距为20cm，大行和小行的行距分别为60、40cm，玉米作物的行列排序为大行和小行交叉排列，对大行开沟而不对小行开沟，两行玉米同时由一条沟控制，沟的尺寸如开口、深、收底分别50、30、20cm。依据辽河来水时间和实际状况对玉米试验田进行灌溉。分别对灌溉前、后进行土样信息采集，将沟尾30m、沟中15m和沟首5m处3个断面作为沟灌的取样断面，沟顶和沟底作为土壤取样具体位置。以40m×1.5m的尺寸设定畦灌的畦田规格，畦灌的取样断面位置分别与沟灌断面相对应并以此提高试验结果的可对比性[7]。为了对畦灌MB- 2位置的首部、中部和尾部进行土壤含水测定，分别在此3个部位设立TDR观测管进行测定，其他各取样点的深度分别与表1中土层深度一一对应，并进行土壤盐分和含水量的测定。

表1 辽河流域某试验站土壤理化性质主要指标

表2 玉米试验田灌溉试验处理

1.3 观测指标

8大离子、土壤pH值、土壤EC值以及土壤含水率等参数为土壤主要的观测指标，其中畦灌模式中各条沟的面积为40m×1.5m=60m2。

2 构建数学模型

2.1 基本方程

土壤水分运动在采用畦灌模式时其有水层的水分入渗符合一维垂直运动；而土壤水分在采用沟灌模式时沿垂直和水平2个方向的入渗属于二维土壤水分运动[8]。对沟灌和畦灌的土壤水分运动模型方程可利用修正的Richards方程进行描述，其中畦灌土壤水分运动模型可采用下式表征：

(1)

可采用下述公式代表沟灌模式下的土壤水分运动模型：

(2)

式中，θ—土壤体积含水量，%；t—时间，d；x、z—空间坐标，其中z以向下为正，cm；K(θ)—非饱和土壤导水率，cm/d。

对模型的初始条件、边界条件以及土壤水力函数等基本问题可通过对上述方程的求解确定。选择V G-Mualem公式作为土壤水力函数方程，函数方程表达式如下：

(3)

(4)

式中，Ks—土壤饱和导水率，cm/d；θe、θr、θs—土壤相对饱和度、剩余体积含水率和饱和体积含水率；α、n—经验拟合参数，并满足m=1-1/n，且α值与土壤物理性质相关；l—经验拟合参数，一般条件l值为0.5。

本研究选取可随水流移动的土壤可溶盐为研究对象，其主要参数指标为土壤水矿化度，然后利用多孔介质溶质运移理论分别建立畦灌和沟灌模式下的非饱和-饱和土壤溶质运移数学模型，其中畦灌溶质模型可采用下式进行表征：

(5)

沟灌溶质模型的表达式如下所示：

(6)

式中，Dij—弥散系数，cm2/d；qi—水流通量，cm/d；c—溶液质量浓度，g/cm3；xi—空间坐标,i=1,2,其中x1=x;x2=z;D11=Dxx,D12=Dxz。

2.2 初始与边界条件

畦灌模式的土壤水分运动初始条件可采用θ(z,0)=θ0(z);Z≤z≤0进行表征，其上、下边界条件分别为θ(0,t)=θs;z=0和θ(Z,t)=θ0(t)。

上述初始条件和边界条件表达式中，θ0、θs—土壤初始含水率和饱和含水率，%；qs—地表水分通量，cm/d，降水与灌溉入渗值为负而蒸散发为正；c0—初始剖面土壤水矿化度，g/cm3；cs、cb—上边界流量的矿化度和下边界浅水矿化度，g/cm3。若降水量或土壤水蒸散量作为边界流量则cs值为0；若灌溉水量作为边界流量则cs为灌水矿化度，g/cm3。

3 数值模拟与参数识别

3.1 HYDRUS模型简介

HYDRUS软件可对非饱和多孔隙媒介中水流和溶质的运移过程进行数值模拟，它是一种对非饱和土壤中水、热及溶质的一维或二维运动利用土壤物理参数进行模拟的有限元计算机模型。模型选择饱和-非饱和达西水流为其水流状态，并且不考虑土壤水流运动受空气的作用影响，采用修改后的Richards方程作为水流控制方恒。软件程序可对各类水流边界进行灵活处理，主要包括排水沟、自由排水边界、大气边界、变水头和定水头边界等。不规则的水流边界可作为水流区域本身的边界条件，特殊情况下还可以由非均质各向异性的土壤组成其边界[9]。

3.2 模型模拟

将地下0～100cm范围的土壤按照一定标准划分为5个层次并利用模型进行模拟分析，其模拟时长共98h并按照剖分方式划分变时间步长，对时间步长依据收敛迭代次数进行调整。设定的最小、最大步长以及初始时间步长分别为5、0.01和0.1d，其中压力水头和土壤含水量容许偏差分别为1cm和0.0005。

采用单孔隙模型中的V G-Mualen模型作为土壤水流模型，利用逆向求解法且不考虑水分滞后效应对水盐运动参数进行确定。选取开发大气边界作为盐分和水流模拟的熵边界条件，对土壤蒸发、降水、灌溉补给以及HYDRUS水流模拟分别赋值实测的蒸发量、灌溉量和降水量；对实测灌溉水矿化度赋值盐分模拟参数。选择已知变水头边界作为水流模拟的下边界条件，其压力水头赋于HYDRUS内并依据实测地下水埋深对其确定；选择已知浓度边界作为盐分模拟的下边界条件并赋值实测潜水矿化度[10]。

3.3 模型参数

依据实测土壤粒径构成土壤水力参数，参数初始值为Rosetta模型的初始值，然后利用2015年生育期试验区的实测数据利用模型进行拟合，对主要特征参数数值进行求解和确定，V G-Mualem公式各参数值经调整后的统计结果见表3。

表3 土壤水力特征各参数值计算结果

3.4 模型验证与效果评价

选取2015年试验区生育期土壤导电率(EC)和含水率实测值对模型的模拟结果进行验证分析，对所构建数值模型的合理性利用实测数据模拟数据进行对比分析[11]。为验证模型的精度分别选取实测值、电导率模拟值和土壤含水率进行两配对样本T检验见表4，其中畦灌和沟灌模式下参与检验的样本数分别为57和43。研究表明，在置信区间为0.05时显著性水平P值范围内未出现导电率和土壤含水率配对T检验结果，由此表明上述2参数实测值与模拟值无显著性差异，模型参数较为可靠，其模拟精度满足有关要求，相关参数的设定可应用于实际模型。

表4 模型模拟效果评价

4 模型应用与分析

分别对畦灌和沟灌模式下的3种灌水定额下的土壤水盐运移规律利用识别的田间水盐运移模型进行研究分析，对田间灌溉最佳模式进行确定，其中灌水等额分别为3.6、4.3、5.2m3/hm2和2.4、3.2、3.8m3/hm2。

4.1 土壤含水率影响分析

在沟灌和畦灌模式下的土壤含水率变化状况如图1所示，由图可知不同灌溉处理条件下的土壤含水率在玉米的整个生育期内均保持大致相同的变化趋势。根据不同土层深度的含水率变化可知，土壤含水率随土壤深度的增加而逐渐表现出增大的趋势，土壤含水率在0～20cm和60cm以下土层深度时约为30%和35%，土壤含水量在80～100cm处相对于0～20cm处有明显的提高约为30%，在80～100cm深度的土壤含水量相对于其上层土壤基本处于稳定状态，深处土壤含水量受灌溉作用的影响相对较低。

4.2 土壤盐分影响分析

不同时期不同灌溉模式下的各层土壤的盐分EC值分布状况如图2所示，由图可知土壤平均含盐量与灌水定额和灌水模式影响显著，而且土壤盐分含量降低的趋势可随灌水定额的增大而逐渐增加[12]。在相同灌水定额条件下沟灌模式相对于畦灌模式土壤平均含盐有所降低，在玉米作物生育期1m土体内的土壤盐分均值在沟灌模式下低于畦灌模式下，其原因主要为沟灌模式的入渗作用优于畦灌模式的入渗作用，特别是表层土壤水分入渗作用的增大可使得土壤盐分随水分运动的效果明显降低，并且盐分的减少幅度在0～40cm的表层土壤处更为明显，由此表明沟灌模式可有效降低土壤盐分[13]。

5 结论

(1)对不同灌水模式下土壤水盐运移时空变化规律以及土壤剖面水盐分布特征利用HYDRUS模型进行模拟分析，并将其模拟结果与田间试验实测数据进行对比验证分析。研究表明，土壤中水盐的分布随时间的变化规律可利用模型进行很好的模拟分析，所构建的模型表现出较好的适用性与可靠性，能较好地反映试验区水盐运移变化的实际状况。

(2)土壤含水率随土壤深度的增加而逐渐表现出增大的趋势，土壤含水率在0～20cm和60cm以下土层深度时约为30%和35%，土壤含水量在80～100cm处相对于0～20cm处有明显的提高约为30%，在80～100cm深度的土壤含水量相对于其上层土壤基本处于稳定状态，深处土壤含水量受灌溉作用的影响相对较低。

图1 不同灌溉模式下各层土壤含水量变化

图2 不同灌溉模式下各层土壤的EC值变化

(3)土壤平均含盐量与灌水定额和灌水模式影响显著，而且土壤盐分含量降低的趋势可随灌水定额的增大而逐渐增加，在相同灌水定额条件下沟灌模式相对于畦灌模式其土壤平均含盐有所降低。1m土体内畦灌和沟灌模式下的盐分变化量分别为70.66%、80.27%、58.42%、42.18%、45.45%、62.47%和52.88%、57.16%、47.62%、6.48%、21.65%、38.42%。

(4)对于表层土壤其水分入渗作用的增大可使得土壤盐分随水分运动的效果明显降低，并且盐分的减少幅度在0～40cm的表层土壤处更为明显，研究表明沟灌模式可有效降低土壤盐分；土壤盐分可通过采取合理的灌水定额以及灌溉模式进行控制，依据节水控盐综合标准可将其对应的灌溉制度和沟灌模式进行推广和应用。