基于HYDRUS對稻田不同階段土壤水分的模擬與分析

馬黎華，胡笑濤，蔣先軍，王文娥*

(1. 西南大學(xué)資源環(huán)境學(xué)院，重慶 400715； 2. 西北農(nóng)林科技大學(xué)旱區(qū)農(nóng)業(yè)水土工程教育部重點實驗室，陜西 楊凌 712100)

在水稻生育后期，一般要進行排水干田或自然落干的過程，稻田土壤水分狀況實現(xiàn)由飽和態(tài)到非飽和態(tài)的轉(zhuǎn)變.在這一轉(zhuǎn)變過程中，土壤水分運移方式以及涉及土壤水分各平衡要素的轉(zhuǎn)化均存在差異性.加深對稻田不同水分狀況的認識和理解，對于農(nóng)業(yè)水管理、水分高效利用以及減緩農(nóng)業(yè)面源污染等方面具有實際意義.在水稻生產(chǎn)過程中的水量平衡要素如地表蒸發(fā)量、作物蒸騰量以及底層滲漏量等難以實現(xiàn)連續(xù)和有效地監(jiān)測.土壤水分的數(shù)值模擬可以利用模擬結(jié)果為一些較難實現(xiàn)連續(xù)監(jiān)測的指標(biāo)提供更多的數(shù)據(jù)信息，加強對土壤水分運移過程以及環(huán)境因素影響的認識和理解.稻田的水分運移與消耗包括地表蒸發(fā)、植物蒸騰、深層滲漏和其他邊界水分損失等.蒸發(fā)量與蒸騰量發(fā)生在稻田的土壤上邊界，受降雨與日照氣象因素的影響較為直接[1].底層滲漏量取決于土壤的水力特性、壓力水頭以及田塊的狀態(tài)[2]等多種因素.底層滲漏量占稻田水分損失總量的50%～80%[3].稻田采用干濕交替和節(jié)水控灌等不同水分管理措施會對稻田水分運移規(guī)律和稻田水分利用效率產(chǎn)生影響[1,4].

HYDRUS模型基于Richards方程在飽和-非飽和不同狀態(tài)下的對流-彌散方程開展數(shù)值運算分析[5-6]，可以用于稻田土壤水分狀況、水稻的水分利用效率和稻田氮淋失等方面的研究[1,4,7].HYDRUS在水分運移模擬過程中對土壤水力特性參數(shù)的選取至關(guān)重要，考慮到田間實際水分變化的復(fù)雜性，往往進行求逆解(inverse solution)模塊對模型輸入的參數(shù)進行校驗，校驗過程也是對土壤水分運移規(guī)律的再認識過程[1].

文中充分考慮稻田土壤水分的狀態(tài)因素和作物生長情況，利用HYDRUS-1D模型，對常規(guī)水旱輪作(conventional paddy-upland rotation tillage, CT)水稻生育后期和收獲后的田間水分狀況進行定量分析.① 選擇水稻種植后期淹水、排水干田及收獲后3種不同土壤水分狀況的典型階段(階段Ⅰ——淹水、階段Ⅱ——排水干田、階段Ⅲ——收獲后)，使用HYDRUS-1D模型對不同土壤深度的土壤含水量進行模擬，通過試驗數(shù)據(jù)率定土壤水力參數(shù)，確定模擬的有效性；② 根據(jù)水量平衡分析，估算3個不同階段底層滲漏量的差異性；③ 根據(jù)3個不同階段的模擬結(jié)果，分析降雨和日照時數(shù)對水稻土壤上下邊界流的影響，闡明氣象因素影響下土壤水分的變化規(guī)律.

1 材料和方法

1.1 試驗設(shè)計

試驗于2017年在重慶市北碚區(qū)西南大學(xué)稻田長期定位試驗站(106°26′E，30°26′N)開展.試驗站海拔230 m，年均溫度為18.3 ℃，年日照時數(shù)為1 277 h，年均降水量為1 105 mm，無霜期約為334 d；土壤類型為灰棕紫色沙泥巖母質(zhì)上發(fā)育的中性紫色水稻土，屬于重慶地區(qū)典型丘陵紫色土稻田區(qū).水稻品種為炳優(yōu)6號(炳優(yōu)6028)，為早秈三系雜交稻品種.處理CT在夏季種植水稻，收獲后粉碎秸稈翻耕20 cm深；冬春季節(jié)種植油菜.

2017年試驗水稻全生育期101天(5月5日移栽，8月14日收獲)，第101天進行水稻收割.在水稻生長季，處理CT使表層土壤的水深保持在5 cm左右，從成熟期后期進行田間排水直至收獲.降水量P、田間水分管理以及所選擇的3個階段如圖1所示，圖中Hw為水層深度，ds為播種后天數(shù).選取的3個階段分別為階段Ⅰ(7月14—28日)，持續(xù)14 d(336 h)，稻田淹水，土壤水分為飽和狀態(tài)；階段Ⅱ(7月30日—8月13日)，持續(xù)14 d(336 h)，稻田排水干田，土壤水分從飽和狀態(tài)到非飽和狀態(tài)；階段Ⅲ(8月17—24日)，持續(xù)7 d(168 h)，稻田收獲后，該階段為持續(xù)無降雨階段，稻田土壤為非飽和狀態(tài).

在水稻生長季節(jié)和收獲后收集氣象和土壤含水量等相關(guān)數(shù)據(jù).氣象數(shù)據(jù)采用HOBO小型氣象站(H21-USB, ONSET, America)和HOBO U30 USB氣象站數(shù)據(jù)記錄器(U30-NRC, ONSET, America)采集，包括降水量(mm)、日最高氣溫(℃)，日最低氣溫(℃)，日照時數(shù)(h)，地面2 m處風(fēng)速(m/s).土壤含水量(體積，下同)(cm3/cm3)采用EC5土壤水分智能傳感器(S-SMC-M005, ONSET, America)監(jiān)測，安裝在10，20，30，40和50 cm土層深度處.

圖1 常規(guī)水旱輪作(CT)降雨、稻田水分管理以及三階段選取示意圖

1.2 HYDRUS模型描述

HYDRUS用于描述飽和-非飽和態(tài)土壤水分運動的方程為Richards方程，即

(1)

式中：θ為土壤體積含水量，cm3/cm3；t為時間，h；z為空間坐標(biāo)，cm；K為土壤導(dǎo)水率，cm/h；h為壓力水頭，cm；S為根系吸水量，cm3/(cm3·h).

土壤體積含水量θ(h)和土壤導(dǎo)水率K(h)由Van Genuchten(1980)公式計算，即

(2)

(3)

式中：θs為飽和含水量，cm3/cm3；θr為滯留含水量，cm3/cm3；m，α和n為決定水分特征曲線形狀的經(jīng)驗參數(shù)，其中m=1-1/n；Ks為飽和導(dǎo)水率，cm/h；l為土壤孔隙連通性參數(shù)，一般在HYDRUS模擬中選擇l=0.5；Se為相對飽和度，定義為

(4)

1.2.1 土壤機械組成與水力參數(shù)確定

對于水力參數(shù)的估計采用了Neural Network Prediction模塊的Rosetta Lite v.1.1[8]轉(zhuǎn)換實現(xiàn)，主要由土壤質(zhì)地即土壤的機械組成決定[6].土壤的機械組成和容重轉(zhuǎn)化為6個土壤水力參數(shù)，分別為飽和含水量(θs)，滯留含水量(θr)，飽和導(dǎo)水率(Ks)，進氣吸力倒數(shù)(α)，孔隙尺寸分布指數(shù)(n)(n大于1)，土壤孔隙連通性參數(shù)(l)取常數(shù)0.5.

分層由種植前測得的土壤分層與土壤組成，見表1，表中Hs，τ，ρ分別為土層深度、顆粒體積分數(shù)、土壤容重.表1表明，根據(jù)Rosetta Lite v.1.1得到的第2層的土壤水力參數(shù)Ks比第1層減小了27.3%.

表1 水稻土的土壤機械組成

考慮到實際情況的復(fù)雜性，采用田間土壤含水量的觀測數(shù)據(jù)和求逆解(inverse solution)模塊，對土壤的水力參數(shù)進行參數(shù)擬合，結(jié)果見表2.

表2 土壤水分特征參數(shù)初始值以及Hydrus的優(yōu)化值

在土壤水處于飽和狀態(tài)的階段Ⅰ，擬合值表明θs的擬合結(jié)果為測量數(shù)據(jù)的最大值，其他值保持初始設(shè)定不變.在非飽和階段Ⅲ，HYDRUS的擬合值主要是關(guān)于α，n，Ks的擬合值變化，其中17～33 cm的Ks與初始值的差異最大.具體有(0，17] cm，Ks擬合值的變化范圍為1.46～1.48 cm/h，在3個階段的對比中，階段Ⅰ略小于階段Ⅱ和階段Ⅲ；(17，33] cm，Ks擬合值的變化范圍為1.04～1.71 cm/h，在3個階段的比較中，階段Ⅰ(1.04 cm/h)小于階段Ⅱ(1.11 cm/h)，小于階段Ⅲ(1.71 cm/h)；(33，63] cm，Ks擬合值的變化范圍為1.23～1.32 cm/h，階段Ⅱ大于階段Ⅰ和Ⅲ.

1.2.2 HYDRUS邊界條件設(shè)定

HYDRUS-1D的模擬設(shè)定為1 m深度的土壤剖面上，邊界條件不考慮土壤水分的側(cè)向流動.考慮到該區(qū)域的地下水位深度大約為3 m，底部的邊界條件設(shè)定為自由排水，土壤剖面的上邊界條件設(shè)定為大氣邊界.初始條件為每個階段初始時刻的測量土壤含水量(cm3/cm3).

1.2.3 模型誤差分析

通過計算決定系數(shù)R2和均方根誤差RMSE將不同土壤深度土壤含水量的HYDRUS模型模擬值與觀測值進行比較.

R2的值在0到1之間，越接近于1表明模型的模擬效果越好；RMSE與土壤含水量觀測值具有相同的量綱，便于對不同階段之間的模擬誤差進行比較.

1.3 潛在地表蒸發(fā)量Ep和潛在植物蒸騰量Tp的估算方法

HYDRUS模型中的潛在蒸發(fā)量(Ep，cm/h)與潛在蒸騰量(Tp，cm/h)之和為作物的蒸發(fā)蒸騰量，以作物系數(shù)Kc[9]結(jié)合Penman-Monteith公式[10]組合計算得到[11].Ep和Tp因作物在土壤表面覆蓋度差異而不同，與太陽輻射值Kgr[12]以及葉面積指數(shù)LAI密切相關(guān)，計算公式為

ETC=ET0·KC，

(5)

EP+TP=ETC，

(6)

Ep=ETC·eKgr·LAI.

(7)

水稻的葉面積指數(shù)LAI根據(jù)生育期的劃分[6]實測得到：階段Ⅰ取6.2，階段Ⅱ取4.7，階段Ⅲ為無植株，故LAI為0.由此得到的地表蒸發(fā)量(Ep)和植物蒸騰量(Tp)如圖2所示.

圖2 常規(guī)水旱輪作(CT)稻田3個階段的潛在地表蒸發(fā)量(Ep)與植物蒸騰量 (Tp)

1.4 根系吸水模型

在HYDRUS-1D根系吸水模塊的計算中，基于利用FEDDES等[13]提出的水分脅迫響應(yīng)函數(shù).

采用SINGH等[14]優(yōu)化的水稻作物[h1=100 cm，h2=55 cm，h3(高)=-160 cm，h3(低)=-250 cm，h4=-15 000 cm]的參數(shù)值作為水稻根系吸水模型的設(shè)定[6].以不同生育期實際測得的實際根系參數(shù)[15]為根系分布的輸入項.

1.5 水量平衡分析

水量平衡分析計算式為

P+I=E+T+D+ΔV，

(8)

式中：P為累積降雨量，cm；I為累積灌溉量，cm(選取的階段Ⅰ和Ⅱ無灌溉)；E為實際和累積蒸發(fā)量，cm；T為實際和累積根系吸水量，cm；D為累積底層滲漏量；ΔV為土壤蓄水量、土壤耗水或補充水量的變化，用于反映從初始時間含水量與結(jié)束時間含水量的差異.

2 結(jié)果與分析

2.1 實測土壤含水量

稻田土壤不同深度Hs的實測土壤含水量在3個選取時段的最大值、最小值和均值見表3.

表3 不同階段不同深度的土壤含水量

在選取的階段Ⅰ，不同深度的土壤含水量均處于飽和含水量的狀態(tài)，各個土層之間的最大值、最小值和平均值的差異也很小.在選取的階段Ⅱ，除了10 cm深度的最大值和最小值的差異不大，深度20，30，40和50 cm土層含水量的最大值接近飽和含水量0.5 cm3/cm3，最小值為0.1～0.2 cm3/cm3.在水稻收獲后階段Ⅲ中除50 cm，其他土層的含水量明顯處于非飽和狀態(tài)，最大值和最小值均為0.1～0.2 cm3/cm3.階段Ⅰ，Ⅱ，Ⅲ的0～50 cm平均土壤含水量分別為0.521，0.316，0.169 cm3/cm3.選取的3個階段基本代表了稻田土壤處于淹水飽和狀態(tài)(階段Ⅰ)、排水干田由飽和狀態(tài)到非飽和態(tài)的轉(zhuǎn)換過程(階段Ⅱ)和收獲后完全非飽和狀態(tài)(階段Ⅲ).

2.2 模型模擬與誤差分析

稻田土壤10，30和50 cm深度土壤含水量的實測值θt和HYDRUS模擬值θm如圖3所示.

在階段Ⅰ，不同深度稻田土壤含水量的實測值均接近其飽和含水量，即0.5 cm3/cm3左右，其中10 cm深度的實測值存在微小的波動變化，各土層的HYDRUS的模擬值為一條直線，基本等于該土層的飽和含水量值；30 cm深度土壤含水量的模擬值略低于實測值，模擬值和實測值的平均值分別為0.52和0.53 cm3/cm3；50 cm深度土壤含水量的模擬值略高于實測值，模擬值和實測值的平均值分別為0.51和0.49 cm3/cm3.

在階段Ⅱ，土壤含水量的變化范圍為0.20～0.50 cm3/cm3.在變化節(jié)點附近，實際的土壤含水量的變化比模擬值減小得更為迅速.10 cm深度的土壤含水量為0.20 cm3/cm3左右，由于降水呈現(xiàn)2個峰值的變化，模擬值的變化比實測值的變化更加明顯.30 cm在這一階段部分仍處于飽和狀態(tài)，存在飽和態(tài)與非飽和狀態(tài)的轉(zhuǎn)換過程，這可能與犁底層的存在或者水稻土特殊的土層結(jié)構(gòu)有關(guān).在模擬中調(diào)整壓力水頭可以取得較好的模擬效果，由圖3b部分下降趨勢中，實測值比模擬值的變化更陡.

在階段Ⅲ，HYDRUS能夠較好地模擬非飽和狀態(tài)下無降雨土壤水分持續(xù)降低的過程，30 cm深度處的模擬值略高于實測值，模擬值和實測值的平均值分別為0.17和0.15 cm3/cm3.

圖3 稻田10，30和50 cm深度階段Ⅰ，Ⅱ和Ⅲ土壤含水量

表4為CT在不同土壤深度、不同階段的實測和模擬含水量的誤差分析結(jié)果.除了階段Ⅰ，R2均符合模擬值與觀測值存在顯著(P<0.05)線性關(guān)系.RMSE在3個階段均為0.00～0.10 cm3/cm3.階段Ⅰ中，除了在50 cm處，其他土壤深度處于飽和土壤含水量(0.51～0.53 cm3/cm3)中，RMSE最小(0.02).階段Ⅱ中的R2最高(0.72)，RMSE較低(0.06)；階段Ⅲ中的R2較高(0.65)，RMSE最低(0.05).可以認為HYDRUS可以對稻田的水分進行有效模擬，對階段Ⅲ土壤水分非飽和狀態(tài)下的模擬效果好于階段Ⅱ和Ⅰ，即土壤水分飽和或部分飽和狀態(tài)下的模擬效果.

表4 土壤含水量HYDRUS模擬誤差分析

2.3 水量平衡分析與底層滲漏量估算

由式(12)，考慮不同階段的降雨P(guān)，E，T和土壤水分變化的水量平衡結(jié)果見表5，表中d1，θa，Tp，Ep，Dc，Ds分別為天數(shù)、0～50 cm含水量均值(cm)、實際蒸騰量、實際蒸發(fā)量、底層滲漏量計算值、底層滲漏量HYDRUS模擬值.在階段Ⅰ，水分輸入總量來自降水量5.84 cm，其他補充量為3.56 cm，土壤水分損失提供為1 cm；水分損耗中，植株蒸騰量占57%，地表蒸發(fā)量占43%.在階段Ⅱ，水分輸入總量來自降水量5.26 cm，土壤水分損失提供的水分補給為5 cm；水分損耗中，植株蒸騰量占36%，地表蒸發(fā)量占12%，底層滲漏量占52%.在階段Ⅲ，水分輸入總量來自降水量0 cm，土壤水分損失提供為12 cm；水分損耗中，植株蒸騰量占0%，地表蒸發(fā)量占5%，底層滲漏量占95%.

表5 常規(guī)水旱輪作(CT)稻田3個階段的水量平衡分析

在階段Ⅱ和階段Ⅲ，底層滲漏量是稻田水分損失的主要來源；階段Ⅰ以水分補充為主，階段Ⅱ以底層滲漏為主.在階段Ⅰ和階段Ⅱ，土壤含水量變化不大，土壤含水量基本維持和接近土壤飽和含水量.在階段Ⅲ，由于作物已經(jīng)被移除，導(dǎo)致土壤含水量的變化主要有地表蒸發(fā)和底層滲漏導(dǎo)致，底層滲漏量是水分減少的主要方式.

2.4 HYDRUS模擬降雨、日照對土壤水分的影響

2.4.1 降雨對HYDRUS模擬上、下邊界流的影響

根據(jù)HYDRUS模擬結(jié)果中對于上邊界流的計算值，降雨對于上、下邊界流q的影響如圖4所示.降雨量的大小影響了邊界流的大小.階段Ⅰ中，3次降雨的累積降雨量分別為1.98，3.14和0.06 cm；對應(yīng)邊界流的變化值為-0.010～0，-0.030～-0.007和-0.007～-0.003 cm/h.

圖4 降水量對上、下邊界流的影響

在階段Ⅱ，上邊界流的變化依然受到降水量的影響較大，下邊界流的變化受降雨量影響較小，下邊界流先迅速減少然后較為穩(wěn)定.

在持續(xù)無降水的階段Ⅲ，由于缺少降水量的補給，上邊界流為0，下邊界流先迅速下降隨后趨于穩(wěn)定.

2.4.2 日照時數(shù)對HYDRUS模擬實際蒸散發(fā)的影響

根據(jù)HYDRUS的模擬結(jié)果，日照時數(shù)tsun對于實際作物蒸騰量和實際地表蒸發(fā)量產(chǎn)生直接的影響如圖5 所示.

圖5 日照時數(shù)對于實際作物蒸騰量和實際地表蒸發(fā)量的影響

由圖可以看出，日照時數(shù)的日變化規(guī)律導(dǎo)致了實際作物蒸騰量和實際地表蒸發(fā)量也呈現(xiàn)了明顯的日變化規(guī)律.總體上，階段Ⅱ的蒸發(fā)蒸騰量低于階段Ⅰ.在階段Ⅰ和Ⅱ，對應(yīng)日照時數(shù)最大的正午，實際作物蒸騰量和實際地表蒸發(fā)量最大值分別為0.006和0.003 cm/h；對應(yīng)日照時數(shù)為0的夜間，實際作物蒸騰量和實際地表蒸發(fā)量最小值分別為0.003和0.001 cm/h.在階段Ⅲ，實際作物蒸騰量均為0，對應(yīng)日照時數(shù)最大的正午，實際地表蒸發(fā)量最大值為0.005 cm/h；對應(yīng)日照時數(shù)為0的夜間，實際地表蒸發(fā)量最小值為0.003 cm/h.階段Ⅲ的地表實際蒸發(fā)介于階段Ⅰ和Ⅱ之間.

3 討 論

不同階段稻田土壤水分的差異來源之一體現(xiàn)在土壤水力特性參數(shù)的差異上.文獻[16-17]證實了HYDRUS反推模塊在實際使用中的有效性.文中(見表2)以Ks的差異為例，0～17和17～33 cm 2個土層的Ks均為階段Ⅰ的小于階段Ⅱ和Ⅲ的，即在飽和態(tài)到非飽和態(tài)的轉(zhuǎn)變過程中，Ks是逐漸變大的，這與土壤自身飽和-非飽和入滲的規(guī)律具有一致性[18].Ks的變化由階段Ⅰ的1.04 cm/h和階段Ⅱ的1.11 cm/h，變?yōu)殡A段Ⅲ的1.71 cm/h.在圖3的結(jié)果中也可以看到30 cm的土壤水分在階段Ⅱ中處于飽和與非飽和態(tài)交替變化中.在分層交界面，當(dāng)分層土壤水力特性差別較大時，土壤入滲規(guī)律存在變化上的復(fù)雜性[17].

稻田水量平衡要素的差異很大程度上取決于降水量和灌溉量的差異[6].文中研究中選取水稻生育后期，不考慮灌溉因素，水量平衡要素的影響主要來自于降水量與土壤水分變化，以及底層滲漏量的差異上.由圖4可知，在假設(shè)不產(chǎn)生地表徑流的前提下，降雨過程與上邊界流和下邊界流的變化密切相關(guān)；而階段Ⅱ和Ⅲ中的變化更接近非飽和土壤的變化，即先迅速降低然后趨于平衡的變化[17].日照對于3個階段的影響也存在一定的差異性.由圖5所示，對地表實際蒸發(fā)量而言，階段Ⅲ中，日照時數(shù)的變化與地表實際蒸發(fā)量關(guān)系的決定系數(shù)為0.83，階段Ⅰ和Ⅱ的該值分別為0.64和0.61.根據(jù)HYDRUS這一模擬結(jié)果，可以認為在水稻收獲后的階段Ⅲ，稻田土壤的地表實際蒸發(fā)量與日照時數(shù)的關(guān)系更密切.

4 結(jié) 論

1) HYDRUS可以實現(xiàn)對稻田不同水分狀態(tài)3個階段土壤含水量的有效模擬.

2) 根據(jù)水量平衡結(jié)果，階段Ⅰ稻田土壤水分以水分補充為主；階段Ⅱ和Ⅲ的土壤水分以水分損失為主，底層滲漏量是水分損失的主要途徑，分別占總水分損失的52%和95%.

3) 降雨是土壤上邊界流發(fā)生變化的直接原因，日照時數(shù)為蒸發(fā)蒸騰量帶來了日變化規(guī)律.