袁成福,冯绍元,季泉毅,霍再林
(1.扬州大学水利与能源动力工程学院,江苏 扬州 225009;2.江西水利职业学院,江西 南昌 330013;3.中国农业大学中国农业水问题研究中心,北京 100083)
甘肃省石羊河流域地处西北干旱内陆区,该地区降雨稀少、气候干燥、蒸发强烈,水资源短缺成为制约该地区农业生产的主要因素之一[1]。石羊河流域水资源总量为18.17亿m3,地表水资源量为15.52亿m3,地下水资源量为2.65亿m3[2]。石羊河流域地下水资源存在着不同矿化度的咸水,其中上游地区地下水矿化度为0.5~1.0 g·L-1,中游地区为1.0~3.0 g·L-1,下游地区为3.0~9.0 g·L-1[3-4]。在有限的水资源情况下,为了满足农业生产用水需要和充分利用地下水资源,可以将地表淡水资源和地下咸水资源相联合。国内外有关学者对咸水灌溉技术及利用进行了大量研究,实践证明咸水或微咸水灌溉可使部分农作物产量接近或达到淡水灌溉的产量[5-8]。咸水灌溉的原则是控制土壤中盐分不能超过作物的耐盐度。目前,国内外利用咸水进行农田灌溉,主要有直接利用咸水灌溉,咸水与淡水混合灌溉,咸水与淡水轮灌。不同的咸水灌溉方式,其效果不同。然而,由于野外田间试验受各种因素的制约,很难全面开展各种咸水灌溉的试验,也不易得到不同作物的咸水灌溉利用模式。在田间试验的基础上,应用数学模型来模拟和预测不同咸水灌溉方式下的土壤水盐运移规律及对作物产量的影响,是一种可行的科学研究方法[9-11]。其中,国内外广泛利用SWAP模型来模拟不同咸水灌溉方式下的土壤水盐分布及作物产量,并预测长时期采用咸水灌溉对农田土壤环境的影响。杨树青等[12]在内蒙古河套灌区利用SWAP模型与Visual ModFlow相结合,模拟了不同灌溉定额(淋洗和正常)的微咸水灌溉对土壤盐分累积效应以及作物产量的影响,并预测了长时期微咸水灌溉后土壤根系层盐分分布与平衡。王相平等[13]在苏北地区利用SWAP模型模拟分析了水稻生育期土壤水盐运移规律和水稻水分利用效率,并预测了长时期采用1.5 g·L-1的微咸水灌溉对土壤盐分分布的影响。唐秀楠等[14]在内蒙古河套灌区利用SWAP模型模拟了枸杞不同咸淡水轮灌模式下的土壤盐分运移规律及盐分平衡,并预测了长时期咸淡水轮灌方式对土壤环境的影响。Kumar等[15]在印度新德里利用SWAP模型模拟了小麦不同咸水灌溉条件下根区土壤盐分动态及小麦的相对产量,并预测了小麦长时期咸水灌溉下的相对产量。Jiang等[16]在石羊河流域利用SWAP模型对春玉米咸水非充分灌溉下土壤水盐分布及春玉米产量进行了模拟,并预测了10年间不同土层的土壤水盐分布和春玉米相对产量。石羊河流域由于具有独特的地理气候条件及丰富的光热资源,是我国重要的玉米制种生产基地[17]。在石羊河流域有关咸水灌溉的利用方式主要是直接利用咸水灌溉,而关于咸水与淡水联合灌溉的研究报道较少。因此,本研究在制种玉米咸水与淡水灌溉田间试验的基础上,利用SWAP模型模拟不同咸淡水轮灌下的土壤水盐平衡,并预测较长时期采用咸淡水轮灌模式下的土壤含盐量及制种玉米相对产量,为研究区合理利用地下咸水资源及农业生产实践提供理论依据。
田间试验于2014年4—9月在甘肃省武威市的中国农业大学石羊河试验站进行,所在经纬度为E102°52′, N37°52′,海拔为1 581 m。该试验站位于甘肃省武威市的石羊河流域中上游区域。该研究区地处西北干旱内陆区,降雨稀少,蒸发强度大,多年平均降雨量为164.4 mm,多年平均蒸发量为2 000 mm,研究区地下水位埋深为48 m。咸水和淡水灌溉试验在试验站内测坑中进行,试验站内共有9个测坑,每个测坑的面积为6.66 m2(3.33 m×2 m),深度为3 m,每个测坑之间用混凝土分隔,可防止侧渗。测坑内0~100 cm土层的土壤理化性质见表1。
根据石羊河流域上、中、下游典型区域的地下水矿化度情况试验设置3个处理,分别为s0(灌溉水矿化度为0.71 g·L-1,淡水)、s3(灌溉水矿化度为3.0 g·L-1)和s6(灌溉水矿化度为6.0 g·L-1),每个处理设3次重复,共9个试验小区,采用随机排列方式布置。试验所用淡水为当地井水,咸水根据当地地下水化学组成,采用质量比为2∶2∶1的NaCl、MgSO4和CaSO4混合地下水配制而成。利用管道对试验地进行灌溉,由供水干管引水到田间,再由支管分配到每个测坑,每个支管安装水表,用来控制每次的灌溉水量。
试验作物为当地制种玉米富农963号,于2014年4月19日播种,9月19日收获,全生育期153 d。制种玉米按父本和母本1∶7的比例方式进行种植,制种密度为每个小区56株,制种玉米株距为25 cm,行距为35 cm。各处理制种玉米生育期内共灌溉5次,灌溉水量与当地实际情况保持一致,其中第5次灌溉是为了使制种玉米在成熟期籽粒饱满而进行的灌溉,灌溉制度见表2。其他各种农艺措施均与当地实际情况保持一致。
试验期间在制种玉米播种前、每次灌溉前后和收获后通过土钻田间分层取土的方法获取土样,每个小区每次取1个取样点,每个取样点均分为5层,取土深度分别为0~20、20~40、40~60、60~80 cm和80~100 cm。采用烘干法测定土壤含水量;田间取完土样后预留部分土样带回实验室,将土样风干,进行研磨和过1 mm筛后,采用SG-3型电导率仪(SG3-ELK742,Mettler-toledo international Inc.,Switzerland)测定土壤饱和浸提液的电导率EC1∶5(单位为mS·m-1),并根据已有换算公式(S=0.0275EC1∶5+0.1366)将EC1∶5转化为土壤含盐量,其中土壤含盐量S的单位为mg·cm-3[18]。
在制种玉米出苗后每隔7~10 d获取制种玉米不同生育期的株高、叶面积指数、根长分布等资料。制种玉米收获时玉米脱粒和晒干后称重得到每个处理的产量,再折算单位为kg·hm-2的产量。土壤水分特征曲线参数采用高速离心机测定,利用van Genuchten-Mulaem 模型对土壤水分特征曲线的参数进行拟合。饱和导水率采用渗透仪(TST-55,China),按常水头法测定,为消除不同温度对饱和导水率大小的影响,换算成10℃时的饱和导水率值。气象数据从试验站内的自动气象站采集获取,2014年制种玉米生育期内有效降雨量为141.2 mm。
表1 土壤理化性质
表2 各处理灌溉制度
SWAP(Soil-water-atmosphere-plant)模型是由荷兰Wageningen大学开发的一种用于模拟农田尺度下土壤水分、溶质和热量在土壤-植物-大气-作物系统(SPAC系统)中运移及作物生长过程的综合模型。该模型的上边界位于植物冠层上方,下边界位于地下水系统(饱和层)的上部,分别考虑大气环境因素和区域地下水动态变化的影响。在非饱和带中,SWAP模型假定水流运动的主方向是垂直的,水流运动主要按垂直一维运动考虑,在垂直方向上,SWAP模型将土层分为不同的单元,在每个单元上,耦合求解水分及溶质运动方程和热量传输方程。该模型在国内外干旱地区或半干旱地区模拟土壤水盐运移及作物生长方面得到了较广泛的应用。
土壤水流采用Richards方程:
(1)
式中,θ为体积含水率(cm3·cm-3);K为土壤饱和导水率(cm·d-1);h为土壤水头(cm);Z为垂向坐标(cm),向上为正;t为时间(d);C为容水度(cm-1);S为作物根系吸水项(cm3·cm-3·d-1)。
溶质运移采用对流弥散方程:
(2)
式中,J为总溶质通量浓度(g·cm-2·d-1);q为在边界处的垂向水流通量(cm·d-1);Ddif为溶质扩散系数(cm2·d-1);Ddis为溶质弥散系数(cm2·d-1);∂c/∂z为溶质浓度梯度;c为溶质浓度(g·cm-3)。
SWAP模型模拟的作物生长过程是采用WOFOST作物生长模型,其中包括详细作物模型和简单作物模型,本研究采用简单作物模型。简单作物模型是静态模型,只描述作物最终产量与水分的关系。简单作物模型计算作物的实际产量与潜在产量的比值为相对产量,运用各生育阶段相对产量连乘的数学模型表示整个生育阶段的相对产量。其计算公式如下:
(3)
式中,Ya,k为各生育阶段作物实际产量,Yp,k为各生育阶段作物最大产量,Ta,k、Tp,k为各生育阶段实际蒸腾量和最大蒸腾量,Ky,k为各生育阶段产量反应系数,k为作物不同生育阶段。
(4)
式中,Ya为整个生育阶段累积作物实际产量,Yp为整个生育期作物累积最大产量,n为作物不同生育期阶段的数量。
水盐联合胁迫作用下,SWAP模型描述的作物根系吸水过程是根据作物潜在根系吸水速率,引入水分胁迫修正系数和盐分胁迫修正系数相乘进行计算。计算公式如下:
(5)
Sa(z)=arwarsSp(z)
(6)
式中,Sp(z)是指作物潜在根系吸水速率(d-1);Tp是指潜在腾发速率(cm·d-1);Droot是指作物根系深度(cm);Sa(z)是指作物的实际根系吸水速率(d-1);arw是指水分胁迫修正系数;ars是指盐分胁迫修正系数。
SWAP模型需要输入气象数据、灌溉资料、作物生长资料、土壤理化参数和水力特性参数、初始和边界条件、初始压力水头和溶质浓度等资料。输入气象数据包括每天的太阳辐射量、最高温度、最低温度、平均风速、平均湿度和降雨量,气象数据由自动气象站获取;灌溉资料、作物生长资料、土壤理化参数和水力特性参数均采用田间试验实测数据;土壤剖面的初始上边界为气象因素决定的降雨、蒸发、植物蒸腾和灌溉,由于地下水位埋深较大,土壤剖面下边界条件为自由排水边界;初始压力水头由初始土壤含水量通过水分特征曲线换算得到;初始溶质浓度由土壤初始含盐量换算得到。具体有关SWAP模型的详细介绍参见SWAP模型手册[19]。
模型模拟值与实测值吻合度采用均方误差(RMSE)、平均相对误差(MRE)2个指标进行评价。
(7)
(8)
式中,N为观测值的个数,Pi表示第i个模拟值,Oi表示第i个观测值,其中RMSE和MRE值越小,模型模拟效果越好。
利用研究区实测土壤水盐数据、土壤水力特性参数、制种玉米生长资料、灌溉资料以及气象数据等对SWAP模型进行率定和验证。其中以s3处理为模型的率定处理,s6处理为模型的验证处理。图1为模型率定和验证时不同土层土壤含水量的模拟值与实测值的比较。从图1可以看出,土壤含水量的模拟值与实测值吻合较好,模拟值较好地反映了实测值的变化趋势。表3为土壤含水量率定与验证时模拟值与实测值的判别指标。由表3可知土壤含水量率定与验证过程中土壤含水量RMSE值在0.05 cm3·cm-3以下,MRE值在15%以下,在允许的误差精度范围25%之内。率定后得到的土壤水力特性参数见表4。
图1 土壤含水量模拟值与实测值的比较Fig.1 Comparison of the simulated and measured soil water content in calibration and validation
表3 土壤含水量模拟值与实测值的RMSE和MRE
图2为模型率定和验证时不同时期土壤含盐量的模拟值与实测值的比较。从图2可以看出,不同时期的土壤含盐量的模拟值与实测值吻合较好,模拟值基本上反映了实测值的变化趋势。表5为土壤含盐量率定与验证时模拟值与实测值的判别指标。由表5可知土壤含盐量率定与验证过程中土壤含盐量RMSE值均在4.2 mg·cm-3以下,MRE值均在25%以下,在允许的误差精度范围25%之内。率定后得到的分子扩散系数和弥散度见表4。
表4 率定与验证后的土壤水力特性与溶质运移参数
SWAP模型模拟得出的产量结果为相对产量,本研究假定2014年试验s0处理得到的制种玉米产量(6 303.36 kg·hm-2)为最大实际产量,根据模拟的相对产量与最大实际产量之间的换算可得到各处理的模拟产量。以s0处理为模型的率定处理,s3、s6处理为模型的验证处理,制种玉米产量的率定与检验结果如图3所示,制种玉米的产量模拟值与产量实测值基本一致。制种玉米产量的率定与检验过程中,RMSE值均在380 kg·hm-2以内,MRE值均在10%以下,符合误差精度要求。制种玉米产量率定后得到的最小冠层阻力为60 s·m-1,产生盐分胁迫时土壤含盐量的临界值为1.7 dS·m-1,盐分胁迫引起的根系吸水系数的变化比率为4.0%。
上述对SWAP模型的率定与验证结果表明,率定参数后的SWAP模型能够较好地模拟土壤水盐运动规律,可以用于研究区咸水与淡水灌溉的模拟与预测。
图2 土壤含盐量模拟值与实测值的比较Fig.2 Comparison of the simulated and measured soil salinity in calibration and validation
表5 土壤含盐量模拟值与实测值的RMSE和MRE
图3 制种玉米产量实测值与模拟值比较Fig.3 Comparison of the simulated and measured values for seed maize yield
模拟研究采用咸淡水轮灌方式,设置不同的咸淡水轮灌灌溉方案,其中咸淡水轮灌模式是以淡水灌溉制度而拟定。根据当地灌溉经验和前人的研究结果,制种玉米生育期内一般灌溉5次或者4次,灌溉5次主要是为了使制种玉米收获时籽粒饱满,获得更高的产量和经济收入,在制种玉米成熟期增加了1次灌溉[20]。本次模拟不考虑第5次灌溉,按照制种玉米生育期内需水量进行灌溉,即生育期内共灌溉4次,灌溉定额为465 mm,其中苗期(6月10日)灌溉120 mm、拔节期(7月1日)灌溉120 mm、抽穗期(7月25日)灌溉120 mm、灌浆期(8月15日)灌溉105 mm;制种玉米第1次灌溉时期为制种玉米苗期~拔节期阶段,该阶段制种玉米抗盐能力较差,为了使制种玉米在生育前期顺利生长,第1次灌溉各咸淡水轮灌溉方案均拟定为淡水灌溉,咸淡水轮灌溉只是在第2、第3和第4次这3次灌水期间进行。因此,根据这3次灌水的咸淡水轮灌方式以及考虑3.0、6.0 g·L-1两种灌溉水矿化度,设置的咸淡水轮灌方案见表6,共有12种。由于石羊河流域地处干旱荒漠地带,降雨稀少,枯水年份、平水年份和丰水年份降雨量年际之间差别不大,在模拟过程中,气象数据是采用平均年法(P=50%),根据1960-2014年的降雨资料进行降雨频率分析,得出P=50%的典型年份,对应的年份是2011年,制种玉米生育期降雨量为118.8 mm;初始土壤含水量与土壤含盐量、制种玉米生长资料与2014年田间试验实测数据一致;模拟的土层深度为0~100 cm。根据上述条件,利用率定好的SWAP模型分别对12种咸淡水轮灌灌溉方案进行模拟,表6为模拟不同咸淡水轮灌方案的土壤水盐平衡计算结果。以土体盐分增加量最少和制种玉米产量最高为原则,进行较优咸淡水轮灌模式的筛选。从表6可知,3.0 g·L-1的微咸水条件下,2淡1咸轮灌方案中,淡淡咸模式下的土体盐分增加量最少,产量最高;2咸1淡轮灌方案中,淡咸咸模式下的土体盐分增加量较少,产量最高。6.0 g·L-1的咸水条件下,2淡1咸轮灌方案中,淡淡咸模式下的土体盐分增加量最少,产量最高;2咸1淡方案下各轮灌模式的土体盐分增加量均较大,不考虑。因此,3.0 g·L-1微咸水条件下,采用淡淡咸和淡咸咸两种咸淡水轮灌模式为较优轮灌模式,6.0 g·L-1咸水条件下,采用淡淡咸的咸淡水轮灌模式为较优轮灌模式。
表6 模拟不同咸淡水轮灌方案的土壤水盐平衡计算结果
注:土壤水分变化量“-”表示土壤水分被消耗,水分底部通量为负表示土壤水分向下渗漏,盐分底部通量为负表示土壤盐分向下运动,土体盐分增加量为负表示土壤盐分被淋洗。
Notes: The “-” of soil water change indicates that soil water content was consumed. The negative sign of bottom soil water flux indicates that soil water was leaking downward. The bottom soil salt flux was negative, which means that soil salinity moved downward. The negative increase of soil salinity means that soil salinity was leached.
通过SWAP模型对不同咸淡水轮灌方案进行模拟,已经筛选出了研究区制种玉米较优咸淡水轮灌模式,即3.0 g·L-1微咸水灌溉条件下,采用淡淡咸和淡咸咸两种咸淡水轮灌模式,6.0 g·L-1咸水灌溉条件下,采用淡淡咸的咸淡水轮灌模式。但这3种较优咸淡水轮灌模式是否可以在研究区进行长时期利用,需要利用SWAP模型对长时期利用咸淡水轮灌模式下的土壤盐分及制种玉米产量进行预测。模拟预测时灌溉制度和气象资料不变,以每一年末的土壤水分和盐分模拟结果作为下一年度的初始条件,将上述3种咸淡水轮灌模式连续计算10年。图4为不同咸淡水轮灌模式模拟10年间制种玉米生育期结束后的土壤含盐量变化规律。从图4可以看出,随着运行年份的增加,3种轮灌模式下的土壤含盐量均有不同程度的增加,模拟早期土壤含盐量增幅较大,模拟后期土壤含盐量增幅较小。以3.0 g·L-1微咸水条件下的淡咸咸轮灌模式为例,第2年较第1年增加了3.535 mg·cm-3,第3年比第2年增加了1.739 mg·cm-3,第4年比第3年增加了0.921 mg·cm-3,第5年比第4年增加了0.436 mg·cm-3,第6年比第5年增加了0.226 mg·cm-3,从第7年后土壤含盐量比前一年增加量在0.1 mg·cm-3以内,说明随着咸淡水轮灌的长时期进行,土壤含盐量增幅呈现逐渐减小的趋势,土壤含盐量能够达到平稳,不会造成土壤盐渍化。图5为不同咸淡水轮灌模式模拟10年间制种玉米的产量。由图5可以看出,随着运行年份的增加,制种玉米的产量有一定程度的下降,但减产程度较小,以3.0 g·L-1微咸水条件下的淡咸咸轮灌模式为例,第1年制种玉米的产量为4 664.5 kg·hm-2,第2年制种玉米的产量为4 412.4 kg·hm-2,第3年、第4年制种玉米的产量均为4 349.3 kg·hm-2,第5年后制种玉米的产量均为4 286.3 kg·hm-2,制种玉米产量在模拟期第5年达到平稳,说明随着咸淡水轮灌的长时期进行,制种玉米的相对产量能够保持平稳,不会造成制种玉米大量减产。
图4 模拟10年内土壤盐分动态变化Fig.4 Simulation of soil salinity dynamic for 10 years
图5 模拟10年内制种玉米产量Fig.5 Simulation of seed maize yield for 10 years
王卫光等[21]提出了适合内蒙古河套灌区条件的咸淡水轮灌方式作为指导灌溉方案,灌溉次数为3次,即分别在春小麦分蘖期、拔节期和灌浆期进行灌溉。米迎宾等[22]研究表明在河南东北部地区利用3.0 g·L-1的微咸水作为小麦和玉米的灌溉用水,最好采用咸淡水交替的方式,综合土壤积盐状况和对作物产量的影响,淡淡咸组合灌溉顺序为最优轮灌方案。杨树青等[12]研究表明在内蒙古河套地区利用3.84 g·L-1的微咸水灌溉,春小麦生育期采用淡咸咸的咸淡水轮灌方案,作物根层土壤积盐量较其他咸淡水轮灌方案更少,对研究区的土壤环境和春小麦产量的影响较小。唐秀楠等[14]在内蒙古河套地区利用3.84 g·L-1的微咸水作为枸杞的灌溉用水,采用2淡1咸轮灌方案中,淡淡咸模式在枸杞盛果期盐分较低,对枸杞的产量影响较小;2咸1淡轮灌方案中,淡咸咸模式对枸杞的产量影响较小,并且利用SWAP模型预测了淡咸咸轮灌方式5 年后,土壤盐分增幅不大,能够达到盐分的进出平衡,且5年后盐分依然适宜枸杞的生长,不会对枸杞产量造成影响。本研究得到的3种较优的咸淡水轮灌模式,在预测10 年后,土壤含盐量增加较少,能够达到平稳,不会造成土壤盐渍化,制种玉米的减产幅度较小,这与前人研究结果基本一致。
(1)SWAP模型参数率定与验证结果表明:土壤含水量、土壤含盐量和制种玉米产量的模拟值与实测值吻合度较好,率定参数后的SWAP模型能够较好的模拟土壤水盐运动规律和制种玉米产量,可以用于研究区咸水与淡水灌溉的模拟与预测。
(2)不同咸淡水轮灌模式下的土壤水盐平衡分析表明:3.0 g·L-1微咸水条件下,采用淡淡咸和淡咸咸两种咸淡水轮灌方案为较优轮灌模式,即1水120 mm(淡水)、2水120 mm(淡水)、3水120 mm(淡水)、4水105 mm(3.0 g·L-1微咸水)和1水120 mm(淡水)、2水120 mm(淡水)、3水120 mm(3.0 g·L-1微咸水)、4水105 mm(3.0 g·L-1微咸水);6.0 g·L-1咸水条件下,采用淡淡咸的咸淡水轮灌方案为较优轮灌模式,即1水120 mm(淡水)、2水120 mm(淡水)、3水120 mm(淡水)、4水105 mm(6.0 g·L-1咸水)。
(3)较长时期土壤盐分动态变化及制种玉米产量预测结果表明:3种较优的咸淡水轮灌模式在模拟期内随着时间的增加土壤含盐量增加量较少,能够达到平稳,不会造成土壤盐渍化;制种玉米的减产幅度较小,产量能够保持平稳,不会造成制种玉米大量减产。本研究所得到的3种较优的咸淡水轮灌模式有待在生产实践中进行进一步验证和补充完善。