基于HYDRUS-1D的秸稈覆蓋條件下黑土區(qū)玉米田土壤水鹽運移規(guī)律

劉繼龍,呂 航,曹曉強,李濟楨,張玲玲,奧列格

(1.東北農(nóng)業(yè)大學(xué)水利與土木工程學(xué)院,黑龍江 哈爾濱 150030;2. 農(nóng)業(yè)農(nóng)村部農(nóng)業(yè)水資源高效利用重點實驗室,黑龍江 哈爾濱 150030;3.東北農(nóng)業(yè)大學(xué)文理學(xué)院,黑龍江 哈爾濱 150030)

東北黑土區(qū)每年產(chǎn)出約占全國1/3的商品糧食,是中國重要的商品糧生產(chǎn)基地[1],然而由于不合理的人為和自然因素導(dǎo)致黑土區(qū)農(nóng)田土壤蓄水等能力降低,嚴重威脅了黑土區(qū)土壤資源安全[2]。秸稈覆蓋作為一種保護性耕作措施,能顯著抑制土壤水分無效蒸發(fā),提高土壤蓄水能力,抑制土壤表層鹽分積聚,是保護黑土區(qū)土壤資源的有效途徑[3-4]。探究秸稈覆蓋下黑土區(qū)土壤水鹽運移過程,揭示秸稈覆蓋對土壤水鹽運移機制,對黑土區(qū)土壤資源高效可持續(xù)利用具有重要意義。

目前國內(nèi)外學(xué)者已開展了許多關(guān)于秸稈覆蓋對土壤水鹽運移影響的研究,如孫博等[5]探究了不同秸稈覆蓋量下鹽漬土壤水鹽變化規(guī)律;紀文寧等[6]研究了秸稈帶狀單雙行覆蓋和全覆蓋下半干旱地區(qū)土壤水分分布情況;Peng等[7]分析了不同秸稈覆蓋時期對土壤水分含量及水分利用效率的影響。秸稈覆蓋條件下,秸稈與土壤具有相互作用,而且秸稈覆蓋量、覆蓋時間和覆蓋方式對農(nóng)田土壤水鹽運移過程的影響有所差異,同時作物生長會使土壤水鹽運移過程更加復(fù)雜。此外研究發(fā)現(xiàn)由于農(nóng)田中土壤水鹽動態(tài)信息不易連續(xù)監(jiān)測,系統(tǒng)表征土壤水鹽隨時間變化的連續(xù)動態(tài)難度較大。采用數(shù)值模擬方法是分析土壤水鹽連續(xù)動態(tài)變化的一種有效途徑,HYDRUS模型由于能夠綜合模擬大氣、作物共同作用下的田間土壤水鹽運移過程而受到廣泛關(guān)注[8-9]。目前國內(nèi)外已有許多學(xué)者利用HYDRUS模型對不同條件下土壤水鹽運移規(guī)律進行了研究,如Chen等[10]利用HYDRUS模型研究了不同地膜覆蓋下土壤水分運移過程;余根堅等[11]使用HYDRUS模型分析了不同灌水方式對農(nóng)田土壤水鹽運移的影響;Kanzari等[12]模擬出不同水質(zhì)灌溉下土壤鹽分累積情況;李亮等[13]利用HYDRUS模型探究了半干旱地區(qū)荒地水鹽運移機理。盡管國內(nèi)外學(xué)者在這一領(lǐng)域已經(jīng)取得了許多成果,關(guān)于秸稈覆蓋措施下作物冠層覆蓋與生理耗水作用的農(nóng)田土壤水鹽運移機理還有待深入研究。

本文以東北黑土區(qū)玉米地為例,利用HYDRUS-1D模型模擬不同秸稈覆蓋量與生理耗水復(fù)合機制驅(qū)動下的農(nóng)田土壤水鹽運移過程,并對模型進行適用性評價,揭示不同秸稈覆蓋量對黑土區(qū)農(nóng)田土壤水鹽運移的影響機制,以期為黑土區(qū)農(nóng)田土壤水鹽狀態(tài)監(jiān)測及土壤資源高效可持續(xù)利用提供理論參考。

1 材料與方法

1.1 試驗區(qū)概況

試驗于2022年5—9月在東北農(nóng)業(yè)大學(xué)阿城試驗基地(127°03′E, 45°31′N)進行,該試驗基地位于黑龍江省哈爾濱市阿城區(qū),海拔高程162 m,平均氣溫3.3 ℃,無霜期135～140 d,平均日照時數(shù)2 550 h,年平均降雨量580～600 mm。土壤類型為黑土,有機質(zhì)含量為18.42 g·kg-1,堿解氮量為141.6 mg·kg-1,有效磷量為22.3 mg·kg-1,速效鉀量為146.21 mg·kg-1。作物生育期降雨量與土壤蒸發(fā)量由試驗區(qū)內(nèi)自動氣象站采集,如圖1所示。

圖1 2022年玉米生育期降雨量與蒸發(fā)量Fig.1 Rainfall and evaporation during maize growth stages in 2022

1.2 試驗設(shè)計

試驗地種植玉米品種為‘翔禾88號’,設(shè)置對照(CK,無秸稈覆蓋)、秸稈半量覆蓋(CM1,0.33 kg·m-2)、秸稈全量覆蓋(CM2,0.65 kg·m-2)3個處理。所覆蓋秸稈的長度約為20 cm,秸稈全量覆蓋量為試驗前一年各小區(qū)平均秸稈產(chǎn)量。每個處理3次重復(fù),共9個小區(qū)。每個小區(qū)面積為36 m2(6 m×6 m),保護行、隔離帶均為6 m×2 m,種植密度為株距25 cm、壟距65 cm,其他種植管理方式與當(dāng)?shù)爻Ｒ?guī)種植方式一致。

在玉米苗期(6月4日—6月20日)、拔節(jié)期(6月21日—7月9日)、抽雄期(7月10日—7月31日)、灌漿期(8月1日—8月19日)、成熟期(8月20日—9月11日)內(nèi)每隔7～15 d分別采集試驗小區(qū)0～10、10～20、20～30、30～40 cm土樣,使用環(huán)刀法測定土壤容重,利用激光粒度分析儀(S3500 &Bluewave,大昌華嘉商業(yè)(中國)有限公司)測定土壤粒徑分布(采用國際制土壤分類標(biāo)準),如表1所示。土壤含水率用烘干法測定。土壤含鹽量采用土壤電導(dǎo)率表征,土壤電導(dǎo)率測定時,將風(fēng)干土樣過1 mm篩,以土水比1∶5配置待測液,提取上清液后采用電導(dǎo)率儀(DDS-307A,上海儀電科學(xué)儀器股份有限公司)測定。

表1 土壤容重與粒徑分布Table 1 Soil bulk density and particle size distribution

2 HYDRUS-1D模型建立

2.1 基本方程

2.1.1 水分運動基本方程 采用增加匯源項的一維Richards方程來描述土壤水分運移過程[14]。

(1)

式中,θ為土壤體積含水率(cm3·cm-3);t為時間(d);z為垂向空間坐標(biāo),向上為正;h為基質(zhì)勢(cm);K(θ)為非飽和導(dǎo)水率(cm·d-1);S(z)為匯源項,即單位時間內(nèi)的根系吸水量(cm3·m-3·d-1)。

其中,公式(1)中參數(shù)采用HYDRUS-1D中的Van Genuchten模型擬合[15]。

θ(h)=θr+(θs-θr)·[1+|ah|n]-m,h<0

(2)

(3)

(4)

式中,θ(h)為基質(zhì)勢h處的土壤體積含水率(cm3·cm-3);Se為有效飽和度(cm3·cm-3);θs為飽和含水率(cm3·cm-3);θr為殘余含水率(cm3·cm-3);Ks為飽和導(dǎo)水率(cm·d-1);a、m、n、l為經(jīng)驗系數(shù)。

2.1.2 根系吸水方程 根系吸水方程即匯源項S(z),可用Fedds模型進行計算[16]。

S(z)=α(h)·R(c)·β(z)·Tp

(5)

式中,α(h)為水分脅迫反應(yīng)函數(shù),無量綱;R(c)為鹽分脅迫反應(yīng)函數(shù),無量綱;β(z)為標(biāo)準化的根系吸水分布函數(shù),無量綱;Tp為作物潛在蒸騰速率(cm·d-1)。本文玉米根系吸水參數(shù)選用HYDRUS-1D軟件提供的玉米根系吸水參數(shù)如表2所示。

表2 玉米根系吸水參數(shù)Table 2 Water absorption parameters of maize root system

2.1.3 溶質(zhì)運移基本方程 采用對流-彌散方程來描述土壤溶質(zhì)的運移過程[17]。

(6)

式中,C為溶質(zhì)濃度(g·cm-3);D為水動力彌散系數(shù)(cm2·d-1);qz為水流通量(cm·d-1)。

2.2 初始條件與邊界條件

2.2.1 土壤水分運動的初始條件與邊界條件 土壤水分運動的初始條件、上邊界和下邊界條件分別如式(7)～(9)所示[18]:

h(z,t)=h0(z,0)

(7)

(8)

(9)

式中,ho為初始土壤剖面壓力水頭(cm);qo為初始土壤剖面水流通量(cm·d-1)。土壤水分運動初始條件基于土壤實測含水率設(shè)定。上邊界設(shè)為可隨時間變化的大氣邊界條件,需添加模擬時期內(nèi)逐日降雨量及蒸散量;由于試驗區(qū)平均地下水埋深大于8 m[19],不考慮下邊界水流通量,故設(shè)為自由排水邊界。

2.2.2 土壤鹽分運動的初始條件與邊界條件 土壤鹽分運動的初始條件、上邊界和下邊界條件分別如式(10)～(12)所示[18]:

C(z,t)=C0(z,0)

(10)

(11)

(12)

式中,Co為初始土壤溶質(zhì)濃度(g·cm-3)。土壤鹽分運動采用土壤溶液電導(dǎo)率表征,其初始條件基于土壤實測電導(dǎo)率設(shè)定,上邊界為濃度通量邊界,下邊界為零濃度梯度邊界。

2.3 單元劃分

模擬區(qū)域設(shè)為農(nóng)田0～40 cm土層,劃分為4層,分別于10、20、30、40 cm深度處設(shè)置觀測點。時間離散單位為 d,模擬時間從6月4日—9月11日,共100 d。設(shè)定初始時間步長、最小時間步長、最大時間步長分別為0.001、0.001、1 d。

2.4 模型參數(shù)

土壤水力參數(shù)與溶質(zhì)運移參數(shù)是土壤水鹽模擬的重要參數(shù),土壤水力參數(shù)包括飽和導(dǎo)水率(KS)、殘余含水率(θr)、飽和含水率(θs)、形狀參數(shù)(α)、經(jīng)驗參數(shù)(n)、孔隙連續(xù)度(l),根據(jù)實測土壤容重及粒徑分布使用Rosetta模型初步確定,土壤中溶質(zhì)運移主要與水動力彌散系數(shù)(D)有關(guān),其值基于文獻初步確定[20-21],并根據(jù)試驗實測數(shù)據(jù)對模型參數(shù)初值進行校準優(yōu)化,結(jié)果如表3所示。

表3 土壤水力參數(shù)與溶質(zhì)運移參數(shù)Table 3 Soil hydraulic parameters and solute transport parameters

3 結(jié)果與分析

3.1 模型檢驗

根據(jù)農(nóng)田實測土壤含水率與電導(dǎo)率對模擬結(jié)果進行驗證,采用決定系數(shù)(R2)、標(biāo)準誤差(RMSE)和平均相對誤差(MRE)評價模型模擬精度[21]。R2、RMSE和MRE閾值在0～1之間,其中R2值越大、RMSE和MRE值越小,表明模擬精度越高。模擬精度如表4所示,土壤含水率實測值與模擬值的R2介于0.903～0.940,RMSE介于0.010 ～0.014 cm3·cm-3,MRE介于3.46%～4.64%;土壤電導(dǎo)率實測值與模擬值的R2介于0.817～0.853,RMSE介于0.091～0.111 mS·cm-1,MRE介于4.06%～5.46%。模擬的模擬值與實測值吻合程度較高,模型模擬結(jié)果能較好地反映各處理土壤水鹽動態(tài)變化特征。其中,各處理土壤水分模擬精度均高于土壤鹽分模擬精度,究其原因是鹽分對流-彌散過程不易刻畫,造成其模擬誤差較大。而隨秸稈覆蓋量的增加模擬精度隨之降低,這是由于秸稈覆蓋減少了表層土壤與大氣的接觸面積,使其邊界條件更為復(fù)雜,致使模擬精度降低。

3.2 不同秸稈覆蓋量下土壤水鹽動態(tài)

3.2.1 不同秸稈覆蓋量下土壤水分動態(tài) 不同秸稈覆蓋量下玉米生育期各土層土壤水分動態(tài)變化如圖2所示,同一土層各處理土壤含水率隨時間的波動趨勢基本一致,波動幅度隨土層深度增加而降低,這是因為大氣降雨與土壤蒸發(fā)優(yōu)先作用于上層土壤。玉米各生育時期不同處理0～40 cm土層平均含水率如表5所示,各處理土壤平均含水率在成熟期達到最大值,原因是在此期間降雨較為豐富,玉米生長、葉面積增大較快,降低了土壤蒸發(fā)[22]。在玉米不同生育期內(nèi) CM1和CM2處理土壤平均含水率均顯著高于CK處理。與CK處理相比,苗期CM1和CM2處理土壤平均含水率分別增加7.11%和12.26%;拔節(jié)期CM1和CM2處理土壤平均含水率分別增加5.34%和8.68%;抽雄期CM1和CM2處理土壤平均含水率分別增加4.49%和7.16%;灌漿期CM1和CM2處理土壤平均含水率分別增加4.10%和6.36%;成熟期CM1和CM2處理土壤平均含水率分別增加3.51%和5.48%。

表5 各處理不同生育期0～40 cm土層平均含水率 /(cm3·cm-3)Table 5 Average soil water content of 0～40 cm soil layer in different growth stages of each treatment

圖2 玉米生育期間土壤含水率變化過程及模擬情況Fig.2 Change process and simulation of soil water content during maize growth

相同土層深度下,土壤水分含量隨秸稈覆蓋量增加而增加。在玉米全生育期內(nèi),CM1和CM2處理0～10 cm土層平均含水率較CK處理分別增加5.78%和9.62%,10～20 cm土層平均含水率較CK處理分別增加5.19%和8.87%,20～30 cm土層平均含水率較CK處理分別增加4.83%和7.43%,30～40 cm土層平均含水率較CK處理分別增加3.38%和5.26%。這是由于秸稈覆蓋阻斷了農(nóng)田土壤與大氣的直接接觸,有效降低大氣直接輻射到土壤表層,抑制了土壤水分蒸發(fā),提高了土壤含水率[23]。隨土層深度增大,CM1和CM2處理較CK處理土壤含水率增幅降低,這是因為淺層土壤受地面環(huán)境因素影響較大,使得秸稈覆蓋對淺土壤的影響效應(yīng)更加明顯。

3.2.2 不同秸稈覆蓋量下土壤鹽分動態(tài) 圖3反映了玉米生育期內(nèi)不同土層各處理土壤鹽分動態(tài)變化情況,由圖3可知,同一土層各處理土壤電導(dǎo)率隨時間的波動趨勢相似,波動幅度隨土層深度增加而降低。不同秸稈覆蓋量下玉米各生育時期0～40 cm土層平均電導(dǎo)率如表6所示,在玉米不同生育期內(nèi) CM1和CM2處理土壤平均電導(dǎo)率均顯著低于CK處理。與CK處理土壤平均電導(dǎo)率相比,苗期CM1和CM2處理分別降低7.76%和12.94%,拔節(jié)期CM1和CM2處理分別降低7.21%和11.77%,抽雄期CM1和CM2處理分別降低6.64%和10.51%,灌漿期CM1和CM2處理分別降低6.19%和9.53%,成熟期CM1和CM2處理分別降低5.42%和7.82%。

表6 各處理不同生育期0～40 cm土層平均電導(dǎo)率/(mS·m-1)Table 6 Average soil electrical conductivity of 0～40 cm soil layer in different growth stages of each treatment

圖3 玉米生育期間土壤電導(dǎo)率變化過程及模擬情況Fig.3 Change process and simulation of soil electrical conductivity during maize growth

相同土層深度下,土壤電導(dǎo)率隨秸稈覆蓋量增加而降低。在玉米全生育期內(nèi),CM1和CM2處理0～10 cm土層平均電導(dǎo)率較CK處理分別降低8.93%和12.41%,10～20 cm土層平均電導(dǎo)率較CK處理分別降低7.08%和10.50%,20～30 cm土層平均電導(dǎo)率較CK處理分別降低5.53%和9.29%,30～40 cm土層平均電導(dǎo)率較CK處理分別降低4.61%和8.87%。分析認為是秸稈覆蓋通過降低光照對土壤的輻射作用,減少土壤水分蒸發(fā),延長水分在土層停滯時間,使土壤中可溶性鹽加快溶解速率,從而提高土壤的抑鹽效果[24]。

4 結(jié) 論

HYDRUS-1D模型能夠模擬研究區(qū)不同秸稈覆蓋量條件下玉米生育期間農(nóng)田生態(tài)系統(tǒng)的土壤水鹽運移過程,模擬結(jié)果較好地反映出了土壤水鹽動態(tài)分布和隨生育階段的變化規(guī)律。土壤水分模擬精度高于土壤鹽分模擬精度,且隨秸稈覆蓋量的增加土壤水鹽模擬精度稍有降低。秸稈覆蓋可提高土壤水分含量,抑制土壤鹽分積累。隨秸稈覆蓋量增加,秸稈覆蓋對土壤水鹽運移的作用效果更加明顯。作物生育期內(nèi)土壤含水率與電導(dǎo)率的波動幅度隨土層深度增加而降低。0～40 cm土層中,秸稈半量覆蓋較裸地作物生育期土壤平均含水率提高3.38%～5.78%,土壤平均電導(dǎo)率降低4.61%～8.93%;秸稈全量覆蓋較裸地作物生育期土壤平均含水率提高5.26%～9.62%,平均電導(dǎo)率降低8.87%～12.41%。