曹程鵬, 張 飛, 段建明
(1.綏化學院 農業與水利工程學院, 黑龍江 綏化 152061; 2.東北農業大學 水利與土木工程學院, 哈爾濱 150030)
我國東北黑土區犁底層普遍存在,導致降雨過程中土壤內部壤中流的產生,壤中流作為徑流不可或缺的重要組成部分,對農業生產及生態環境治理有著重要意義。基于不同的假設,國內外學者提出了許多壤中流產生的機理模型,但每種模型都有其應用局限性,不能將已有模型直接應用于實際田間生產活動。Richard模型作為經典非飽和土壤水分運動模擬模型在壤中流的模擬中有所應用[1-3],但其不能直接求出滲流區的側向出流量。貯水泄流模型從宏觀的水量平衡角度出發,雖然可以直接獲得坡面出口斷面處側向出流量[4-5],但其應用局限性較大。動力波模型假設不透水層或準不透水層邊界上飽和區域內流線平行于底板,且水力梯度等于基巖坡度[6-7],該模型早期應用較多。近年來,“雙超”模型及其改進模型逐漸應用于土壤壤中流過程的模擬中,由于其模擬參數變動較大,實際模擬效果并不十分理想[8-9]。
黑土區坡耕地犁底層對土壤水分入滲雖然存在阻滯作用,但犁底層中的垂向入滲仍然是土壤內部水分遷移的主要形式,這就造成了前人提出的貯水泄流模型及動力波等模型應用受到限制,能夠同時反映犁底層中土壤水分垂直遷移與犁底層上部壤中流水平遷移相互作用的模型更具合理性。為此,本文基于能夠連續模擬犁底層上部土壤及犁底層內部土壤水分垂直遷移的Richard模型,借鑒動力波模型及“雙超”模型原理,提出適合黑土區坡耕地壤中流水平遷移的機理模擬模型,并利用人工模擬降雨試驗對模型進行驗證,以期為黑土區坡耕地農田土壤水肥管理提供科學依據。
土壤水分入滲采用一維垂直土壤水分運動的Richard方程計算[10]:
(1)
式中:t表示時間(min);θ表示含水量(%);D(θ)表示土壤水分擴散率(cm2/min);k(θ)表示非飽和導水率(cm/min);θ(z,t)表示在深度為z處的土壤含水量(%)。方程(1)的初始和邊界條件為:

(2)
式中:R為入滲(cm/min);θi為土壤初始含水量(%);tp表示產流時間(min);i表示入滲率(cm/min)。
根據試驗條件,得到方程(1)的上邊界條件Parlange模型[11-12]:
(3)
采用Brooks-Corey[13]的水分特征曲線模型:
k(Θ)=ksΘM
(4)

(5)
(6)
式中:K(Θ)為非飽和導水率;D(Θ)為土壤水分擴散率;Θ為有效飽和度;N,M為形狀系數;M=3N+2;hd為土壤進氣吸力;l為彎曲度;θ0(t)為土壤表層含水量;Ks為土壤飽和導水率(cm/min);θs為土壤飽和含水量;θr為土壤滯留含水量。通過數值方法求解方程(3)得到土壤飽和前的上邊界條件,通過垂直和水平土柱(擾動土)入滲試驗得到土壤參數見表1。
由于土壤容重在犁底層處發生顯著變化,犁底層以上土壤與犁底層內部土壤進行土壤水分垂直遷移模擬時需要采用與土壤特性相對應的參數。2組參數需要分別測量獲得(表1),結合公式(1)將犁底層及其以上土壤作為一個整體,不同位置采用不同參數,利用數值法連續模擬濕潤峰下移過程。計算犁底層上邊界入滲通量時,采用公式(2)中入滲率i作為犁底層上邊界入滲通量計算公式,得到犁底層上邊界入滲能力。假設犁底層以上入滲水未觸及犁底層時段,土壤飽和層與集水槽擋板相接觸的飽和土壤中土壤水分無水平運移,結合動力波模型中假設犁底層區域內流線平行于底板(犁底層),水力梯度等于犁底層坡度,利用雙超模型原理,當犁底層上表面處土壤含水量超過土壤田間含水量,來水量大于犁底層下滲能力,在水力梯度作用下產生橫向流動,得到側向出流模型:
q=iAtanα
(7)
式中:q為側向出流量(ml);i為犁底層上表面入滲率;A為出流斷面面積(cm2);α為坡面坡度。公式(7)反映出壤中流側向出流流量變化過程與降雨強度及土壤初始含水量無關,與犁底層上表面入滲率有關。以上側向出流模型雖然忽略了側向出流土壤導水率與犁底層土壤導水率之間的差異,但在重力作用下的垂向遷移能力與相對疏松的犁底層以上土壤水平遷移能力相近,由此認為公式(7)是成立的。
利用RMSE[14-15]評價模型模擬值和試驗實測值之間的差異。公式RMSE表示為:
(8)
式中:n為總數據點數;pi為模擬值;oi為實測值。

表1 土壤基本參數
注:表中含水量為重量含水量。
利用自制坡度可調節的集雨槽(長100 cm,寬40 cm,高50 cm)裝土,坡度調節范圍0°~30°。土槽設計上下2個徑流收集口,上側收集口下邊緣低于土面0.5 cm,收集土表徑流,下側收集口位于上側收集口以下11.5 cm處,收集位于12 cm處犁底層以上壤中流。采用側噴式人工模擬降雨器(NLJY-09-2型,降雨均勻度系數>0.86,降雨高度4 m)模擬降雨。試驗裝置見圖1。設計坡面坡度7°和10°兩個水平,降雨強度0.12 cm/min,土壤初始含水量18.8%,根據取土樣地實測犁底層深度及土壤容重,設計犁底層位于土表以下12 cm深處,犁底層以上控制土壤容重1.23 g/cm3,犁底層土壤容重1.38 g/cm3,降雨歷時40 min。試驗2次重復。供試土壤理化性質見表2。將供試土壤過4 mm孔篩網除去碎石塊、植物根茬等雜質,經風干、均勻混合處理備用。利用烘干法測量土壤含水量,計算出達到設計含水量需加水量。將試驗用土、水按照重量(利用電子秤稱重)平均分成數份,將第一份土壤均勻攤鋪在塑料上,土層厚約2 cm,將水按每份重量均勻噴灑在土層上,按重量覆蓋上第二層土壤,再噴灑上設計水量的水,如此往復,最后用塑料將土包好,防止水分蒸發。24 h以后打開塑料,將土混合均勻。通過上述方法,可以獲得水分與土壤混合均勻的供試土壤,土壤含水量誤差可控制在士1%(重量含水量)范圍內。將混合好的土壤按設計容重分層裝入土槽中,每5 cm為一層,用秤稱量土重,裝土時用土錘輕輕地將土面砸平,為了得到平整的土面,用鋒利的刀將多余土量刮到土少的地方。填裝下一層土前抓毛下層土壤表面,以防土層之間出現分層現象。土壤填裝完成后用塑料將土槽蓋好。裝土深40 cm,預留10 cm土槽側壁防止雨水擊濺造成水量減少出現試驗誤差。約24 h以后開始降雨試驗。用容量為2 000 ml的塑料桶在土槽接水口接取地表徑流及壤中流,產流后每1 min收集1個水樣,停止降雨后在土槽中部沿坡面垂直方向取土壤剖面,收集土樣,每10 mm取1個土樣。飽和導水率的測定采用環刀法獲取土槽試驗后原狀土,環刀高70 mm,直徑50 mm,設3個重復,采用馬氏瓶定水頭供水測定飽和導水率,取平均值。
圖1 試驗裝置

表2 土壤基本農化狀況
東北黑土區犁底層普遍存在,犁底層界面以上當土壤含水量超過土壤田間持水量后,在水力梯度影響下產生橫向流動,本研究進行的2種坡面坡度試驗是模擬野外大田實際犁底層分布情況進行的,壤中流的產生與野外大田實際產流情況較為接近。
通過試驗獲取壤中流隨時間變化過程見圖2。隨著降雨時間的延續,前期(0~23 min左右)沒有壤中流出現,在這一階段,入滲水僅在重力及毛管吸力等作用下做垂向運動,壤中流產生以后流量隨時間推移首先呈增加趨勢,然后趨于平緩(圖2),與地面徑流變化規律相似。對比坡面坡度7°與10°條件下的壤中流流量變化圖,坡面坡度為7°時,壤中流穩定在1.5 ml/min附近,而坡面坡度為10°時,壤中流流量穩定在了2.4 ml/min,坡面坡度是壤中流的重要影響因素。2種處理條件下壤中流產流時間實測值分別為23.3,22.9 min,可見坡面坡度對壤中流產流時刻影響較小,對穩定流量影響顯著,2種坡度下壤中流流量均在37 min左右趨于穩定。從壤中流流量過程線可知,關于犁底層以上入滲水未觸及犁底層時段,土壤飽和層與集水槽擋板相接觸的飽和土壤中土壤水分無水平運移的假設接近實際情況。
對比壤中流模型模擬值與實測值,模型模擬值與實測值基本一致,2種坡度條件下RMSE分別為0.136 ml/min和0.138 ml/min,R2分別為0.89,0.92,模型能夠反映出存在犁底層條件下的壤中流產生情況,說明了壤中流水力坡度可以近似認為是地面坡度,同時將犁底層入滲能力作為壤中流計算依據的假設具有一定的合理性。2種處理條件下模擬值均小于實測值,產生這一現象的原因可能是試驗過程中出現了邊界效應,入滲水沿土槽壁面下滲快于在土壤中垂直下滲速度,使得犁底層以上靠近土槽側壁土壤含水量高于正常沿垂向土壤遷移條件下土壤含水量,壤中流流量實測值高于模擬值。

圖2 2種坡度處理條件下壤中流流量實測值與模擬值
圖3為降雨結束后土壤沿垂向含水量實測值與模擬值及犁底層分布位置。2種坡度處理條件下濕潤峰下移深度基本相同,可見坡度對土壤垂向入滲影響不顯著,2種處理濕潤峰均觸及犁底層并下移。壤中流的產生是土壤含水量超過田間持水量且入滲受到犁底層明顯阻礙后在犁底層表面產生的側向流動,其對土壤水分的垂直入滲影響很小,土壤水分的遷移仍然以垂直入滲為主,本研究試驗采用的土壤容重犁底層以上與農田耕層相同,飽和導水率為0.025 4 cm/min,而犁底層飽和導水率為0.012 0 cm/min,入滲能力減小顯著,壤中流流量與垂直入滲量相比較取決于犁底層坡面坡度,我國東北黑土區坡耕地大部分較為平坦,取土樣地平均坡度7°,壤中流相對于垂直入滲量較小。雖然圖3顯示出邊界效應的存在,但是模型模擬值與實測值仍然較為接近,2種坡度條件下RMSE分別為1.31%,1.24%(重量含水量),R2分別為0.94,0.97,說明Richard模型可以反映出東北黑土區坡耕地存在犁底層條件下的土壤水分入滲過程。

圖3 降雨結束后土壤剖面含水量實測值與模型模擬值
圖4為2種處理條件下降雨過程中坡面徑流過程,兩種處理條件下產流時刻較為接近,實測值分別為18.21,17.43 min,可見緩坡對于徑流產生時刻影響較小,這一研究成果與Dong等[10]研究成果較為接近。坡度對于隨著降雨時間的延續,流量顯著增加,2種處理條件下均接近250 ml/min趨于穩定,壤中流流量約為地表徑流量1/100,這一比例與降雨強度關系較大,降雨強度越大,比值越小。降雨強度也決定著土壤剖面含水量的變化,降雨強度越大,濕潤峰向下遷移越早,壤中流產生時刻越早,使得壤中流流量過程整體前移,形狀并未改變。雖然壤中流占徑流總量比例較小,但壤中流養分濃度可近似認為是相應深度的土壤水養分濃度,據相關研究顯示坡面徑流養分濃度較土壤水養分濃度小幾個數量級[10],所以土壤養分的徑流流失主要是通過壤中流這一途徑。同時,由于壤中流攜帶養分的水平遷移,也影響了植物根系對土壤養分的利用。
對模型的敏感性分析結果見圖5。利用模型分別模擬了不同坡面坡度、土壤初始含水量及降雨強度條件下的壤中流徑流過程。各條件下從降雨開始至產流時刻歷時見表3。坡面坡度對壤中流產流時刻影響較小(p>0.05),土壤初始含水量和降雨強度的變化引起壤中流出現時刻的變化,且變化規律與地表徑流相似。不同條件下壤中流產流時刻模擬結果較為合理,說明該模型對壤中流產生機理反映接近實際情況。圖5A為改變坡面坡度后壤中流過程,坡面坡度對產流時刻影響較小,但對穩定出流流量影響較大,隨著坡度增加顯著增加(p<0.05)。圖5B—5C顯示不同土壤初始含水量和降雨強度對壤中流過程穩定出流之前影響較大(p<0.05),土壤初始含水量越高,降雨強度越大壤中流出流越早趨于穩定,這與不同條件下壤中流出流過程理論分析基本一致,說明模型對不同影響因素較為敏感,同時說明模型能夠較為合理的反映壤中流出流過程。

圖4 2種坡度處理條件下降雨過程中坡面徑流過程


圖5 模型模擬不同坡度、初始含水量及降雨強度條件下的壤中流過程表3 不同條件下壤中流模型模擬產流時間
犁底層的存在導致壤中流的產生,土壤質地一定條件下,坡面坡度是決定壤中流強度的關鍵影響因素。隨著坡度的增加,壤中流強度增加顯著,坡度對于壤中流及坡面地表徑流產生時刻影響較小。初始含水量及降雨強度僅對產流時間影響較大,對壤中流流量大小影響較小。降雨強度增加、土壤初始含水量增高均能使壤中流產流時間縮短。壤中流流量顯著小于土壤表面徑流,但其對土壤養分水平遷移的貢獻可能高于土表徑流。基于連續模擬犁底層上部土壤及犁底層內部土壤水分垂直遷移的Richard模型,結合動力波模型及“雙超”模型原理構建的壤中流水平遷移機理模型能夠準確反映黑土區坡耕地壤中流徑流過程。這一模型的建立是基于黑土區坡耕地這一特定條件,其他下墊面條件下的壤中流水平遷移模擬適用性有待進一步驗證,同時,該模型適用于壤中流流量上升及穩定階段的模擬計算,對于下降段不適用。