曹程鵬, 張 飛, 段建明
(1.綏化學(xué)院 農(nóng)業(yè)與水利工程學(xué)院, 黑龍江 綏化 152061; 2.東北農(nóng)業(yè)大學(xué) 水利與土木工程學(xué)院, 哈爾濱 150030)
我國(guó)東北黑土區(qū)犁底層普遍存在,導(dǎo)致降雨過(guò)程中土壤內(nèi)部壤中流的產(chǎn)生,壤中流作為徑流不可或缺的重要組成部分,對(duì)農(nóng)業(yè)生產(chǎn)及生態(tài)環(huán)境治理有著重要意義?;诓煌募僭O(shè),國(guó)內(nèi)外學(xué)者提出了許多壤中流產(chǎn)生的機(jī)理模型,但每種模型都有其應(yīng)用局限性,不能將已有模型直接應(yīng)用于實(shí)際田間生產(chǎn)活動(dòng)。Richard模型作為經(jīng)典非飽和土壤水分運(yùn)動(dòng)模擬模型在壤中流的模擬中有所應(yīng)用[1-3],但其不能直接求出滲流區(qū)的側(cè)向出流量。貯水泄流模型從宏觀的水量平衡角度出發(fā),雖然可以直接獲得坡面出口斷面處側(cè)向出流量[4-5],但其應(yīng)用局限性較大。動(dòng)力波模型假設(shè)不透水層或準(zhǔn)不透水層邊界上飽和區(qū)域內(nèi)流線平行于底板,且水力梯度等于基巖坡度[6-7],該模型早期應(yīng)用較多。近年來(lái),“雙超”模型及其改進(jìn)模型逐漸應(yīng)用于土壤壤中流過(guò)程的模擬中,由于其模擬參數(shù)變動(dòng)較大,實(shí)際模擬效果并不十分理想[8-9]。
黑土區(qū)坡耕地犁底層對(duì)土壤水分入滲雖然存在阻滯作用,但犁底層中的垂向入滲仍然是土壤內(nèi)部水分遷移的主要形式,這就造成了前人提出的貯水泄流模型及動(dòng)力波等模型應(yīng)用受到限制,能夠同時(shí)反映犁底層中土壤水分垂直遷移與犁底層上部壤中流水平遷移相互作用的模型更具合理性。為此,本文基于能夠連續(xù)模擬犁底層上部土壤及犁底層內(nèi)部土壤水分垂直遷移的Richard模型,借鑒動(dòng)力波模型及“雙超”模型原理,提出適合黑土區(qū)坡耕地壤中流水平遷移的機(jī)理模擬模型,并利用人工模擬降雨試驗(yàn)對(duì)模型進(jìn)行驗(yàn)證,以期為黑土區(qū)坡耕地農(nóng)田土壤水肥管理提供科學(xué)依據(jù)。
土壤水分入滲采用一維垂直土壤水分運(yùn)動(dòng)的Richard方程計(jì)算[10]:
(1)
式中:t表示時(shí)間(min);θ表示含水量(%);D(θ)表示土壤水分?jǐn)U散率(cm2/min);k(θ)表示非飽和導(dǎo)水率(cm/min);θ(z,t)表示在深度為z處的土壤含水量(%)。方程(1)的初始和邊界條件為:
(2)
式中:R為入滲(cm/min);θi為土壤初始含水量(%);tp表示產(chǎn)流時(shí)間(min);i表示入滲率(cm/min)。
根據(jù)試驗(yàn)條件,得到方程(1)的上邊界條件Parlange模型[11-12]:
(3)
采用Brooks-Corey[13]的水分特征曲線模型:
k(Θ)=ksΘM
(4)
(5)
(6)
式中:K(Θ)為非飽和導(dǎo)水率;D(Θ)為土壤水分?jǐn)U散率;Θ為有效飽和度;N,M為形狀系數(shù);M=3N+2;hd為土壤進(jìn)氣吸力;l為彎曲度;θ0(t)為土壤表層含水量;Ks為土壤飽和導(dǎo)水率(cm/min);θs為土壤飽和含水量;θr為土壤滯留含水量。通過(guò)數(shù)值方法求解方程(3)得到土壤飽和前的上邊界條件,通過(guò)垂直和水平土柱(擾動(dòng)土)入滲試驗(yàn)得到土壤參數(shù)見表1。
由于土壤容重在犁底層處發(fā)生顯著變化,犁底層以上土壤與犁底層內(nèi)部土壤進(jìn)行土壤水分垂直遷移模擬時(shí)需要采用與土壤特性相對(duì)應(yīng)的參數(shù)。2組參數(shù)需要分別測(cè)量獲得(表1),結(jié)合公式(1)將犁底層及其以上土壤作為一個(gè)整體,不同位置采用不同參數(shù),利用數(shù)值法連續(xù)模擬濕潤(rùn)峰下移過(guò)程。計(jì)算犁底層上邊界入滲通量時(shí),采用公式(2)中入滲率i作為犁底層上邊界入滲通量計(jì)算公式,得到犁底層上邊界入滲能力。假設(shè)犁底層以上入滲水未觸及犁底層時(shí)段,土壤飽和層與集水槽擋板相接觸的飽和土壤中土壤水分無(wú)水平運(yùn)移,結(jié)合動(dòng)力波模型中假設(shè)犁底層區(qū)域內(nèi)流線平行于底板(犁底層),水力梯度等于犁底層坡度,利用雙超模型原理,當(dāng)犁底層上表面處土壤含水量超過(guò)土壤田間含水量,來(lái)水量大于犁底層下滲能力,在水力梯度作用下產(chǎn)生橫向流動(dòng),得到側(cè)向出流模型:
q=iAtanα
(7)
式中:q為側(cè)向出流量(ml);i為犁底層上表面入滲率;A為出流斷面面積(cm2);α為坡面坡度。公式(7)反映出壤中流側(cè)向出流流量變化過(guò)程與降雨強(qiáng)度及土壤初始含水量無(wú)關(guān),與犁底層上表面入滲率有關(guān)。以上側(cè)向出流模型雖然忽略了側(cè)向出流土壤導(dǎo)水率與犁底層土壤導(dǎo)水率之間的差異,但在重力作用下的垂向遷移能力與相對(duì)疏松的犁底層以上土壤水平遷移能力相近,由此認(rèn)為公式(7)是成立的。
利用RMSE[14-15]評(píng)價(jià)模型模擬值和試驗(yàn)實(shí)測(cè)值之間的差異。公式RMSE表示為:
(8)
式中:n為總數(shù)據(jù)點(diǎn)數(shù);pi為模擬值;oi為實(shí)測(cè)值。
表1 土壤基本參數(shù)
注:表中含水量為重量含水量。
利用自制坡度可調(diào)節(jié)的集雨槽(長(zhǎng)100 cm,寬40 cm,高50 cm)裝土,坡度調(diào)節(jié)范圍0°~30°。土槽設(shè)計(jì)上下2個(gè)徑流收集口,上側(cè)收集口下邊緣低于土面0.5 cm,收集土表徑流,下側(cè)收集口位于上側(cè)收集口以下11.5 cm處,收集位于12 cm處犁底層以上壤中流。采用側(cè)噴式人工模擬降雨器(NLJY-09-2型,降雨均勻度系數(shù)>0.86,降雨高度4 m)模擬降雨。試驗(yàn)裝置見圖1。設(shè)計(jì)坡面坡度7°和10°兩個(gè)水平,降雨強(qiáng)度0.12 cm/min,土壤初始含水量18.8%,根據(jù)取土樣地實(shí)測(cè)犁底層深度及土壤容重,設(shè)計(jì)犁底層位于土表以下12 cm深處,犁底層以上控制土壤容重1.23 g/cm3,犁底層土壤容重1.38 g/cm3,降雨歷時(shí)40 min。試驗(yàn)2次重復(fù)。供試土壤理化性質(zhì)見表2。將供試土壤過(guò)4 mm孔篩網(wǎng)除去碎石塊、植物根茬等雜質(zhì),經(jīng)風(fēng)干、均勻混合處理備用。利用烘干法測(cè)量土壤含水量,計(jì)算出達(dá)到設(shè)計(jì)含水量需加水量。將試驗(yàn)用土、水按照重量(利用電子秤稱重)平均分成數(shù)份,將第一份土壤均勻攤鋪在塑料上,土層厚約2 cm,將水按每份重量均勻噴灑在土層上,按重量覆蓋上第二層土壤,再噴灑上設(shè)計(jì)水量的水,如此往復(fù),最后用塑料將土包好,防止水分蒸發(fā)。24 h以后打開塑料,將土混合均勻。通過(guò)上述方法,可以獲得水分與土壤混合均勻的供試土壤,土壤含水量誤差可控制在士1%(重量含水量)范圍內(nèi)。將混合好的土壤按設(shè)計(jì)容重分層裝入土槽中,每5 cm為一層,用秤稱量土重,裝土?xí)r用土錘輕輕地將土面砸平,為了得到平整的土面,用鋒利的刀將多余土量刮到土少的地方。填裝下一層土前抓毛下層土壤表面,以防土層之間出現(xiàn)分層現(xiàn)象。土壤填裝完成后用塑料將土槽蓋好。裝土深40 cm,預(yù)留10 cm土槽側(cè)壁防止雨水擊濺造成水量減少出現(xiàn)試驗(yàn)誤差。約24 h以后開始降雨試驗(yàn)。用容量為2 000 ml的塑料桶在土槽接水口接取地表徑流及壤中流,產(chǎn)流后每1 min收集1個(gè)水樣,停止降雨后在土槽中部沿坡面垂直方向取土壤剖面,收集土樣,每10 mm取1個(gè)土樣。飽和導(dǎo)水率的測(cè)定采用環(huán)刀法獲取土槽試驗(yàn)后原狀土,環(huán)刀高70 mm,直徑50 mm,設(shè)3個(gè)重復(fù),采用馬氏瓶定水頭供水測(cè)定飽和導(dǎo)水率,取平均值。
圖1 試驗(yàn)裝置
表2 土壤基本農(nóng)化狀況
東北黑土區(qū)犁底層普遍存在,犁底層界面以上當(dāng)土壤含水量超過(guò)土壤田間持水量后,在水力梯度影響下產(chǎn)生橫向流動(dòng),本研究進(jìn)行的2種坡面坡度試驗(yàn)是模擬野外大田實(shí)際犁底層分布情況進(jìn)行的,壤中流的產(chǎn)生與野外大田實(shí)際產(chǎn)流情況較為接近。
通過(guò)試驗(yàn)獲取壤中流隨時(shí)間變化過(guò)程見圖2。隨著降雨時(shí)間的延續(xù),前期(0~23 min左右)沒(méi)有壤中流出現(xiàn),在這一階段,入滲水僅在重力及毛管吸力等作用下做垂向運(yùn)動(dòng),壤中流產(chǎn)生以后流量隨時(shí)間推移首先呈增加趨勢(shì),然后趨于平緩(圖2),與地面徑流變化規(guī)律相似。對(duì)比坡面坡度7°與10°條件下的壤中流流量變化圖,坡面坡度為7°時(shí),壤中流穩(wěn)定在1.5 ml/min附近,而坡面坡度為10°時(shí),壤中流流量穩(wěn)定在了2.4 ml/min,坡面坡度是壤中流的重要影響因素。2種處理?xiàng)l件下壤中流產(chǎn)流時(shí)間實(shí)測(cè)值分別為23.3,22.9 min,可見坡面坡度對(duì)壤中流產(chǎn)流時(shí)刻影響較小,對(duì)穩(wěn)定流量影響顯著,2種坡度下壤中流流量均在37 min左右趨于穩(wěn)定。從壤中流流量過(guò)程線可知,關(guān)于犁底層以上入滲水未觸及犁底層時(shí)段,土壤飽和層與集水槽擋板相接觸的飽和土壤中土壤水分無(wú)水平運(yùn)移的假設(shè)接近實(shí)際情況。
對(duì)比壤中流模型模擬值與實(shí)測(cè)值,模型模擬值與實(shí)測(cè)值基本一致,2種坡度條件下RMSE分別為0.136 ml/min和0.138 ml/min,R2分別為0.89,0.92,模型能夠反映出存在犁底層條件下的壤中流產(chǎn)生情況,說(shuō)明了壤中流水力坡度可以近似認(rèn)為是地面坡度,同時(shí)將犁底層入滲能力作為壤中流計(jì)算依據(jù)的假設(shè)具有一定的合理性。2種處理?xiàng)l件下模擬值均小于實(shí)測(cè)值,產(chǎn)生這一現(xiàn)象的原因可能是試驗(yàn)過(guò)程中出現(xiàn)了邊界效應(yīng),入滲水沿土槽壁面下滲快于在土壤中垂直下滲速度,使得犁底層以上靠近土槽側(cè)壁土壤含水量高于正常沿垂向土壤遷移條件下土壤含水量,壤中流流量實(shí)測(cè)值高于模擬值。
圖2 2種坡度處理?xiàng)l件下壤中流流量實(shí)測(cè)值與模擬值
圖3為降雨結(jié)束后土壤沿垂向含水量實(shí)測(cè)值與模擬值及犁底層分布位置。2種坡度處理?xiàng)l件下濕潤(rùn)峰下移深度基本相同,可見坡度對(duì)土壤垂向入滲影響不顯著,2種處理濕潤(rùn)峰均觸及犁底層并下移。壤中流的產(chǎn)生是土壤含水量超過(guò)田間持水量且入滲受到犁底層明顯阻礙后在犁底層表面產(chǎn)生的側(cè)向流動(dòng),其對(duì)土壤水分的垂直入滲影響很小,土壤水分的遷移仍然以垂直入滲為主,本研究試驗(yàn)采用的土壤容重犁底層以上與農(nóng)田耕層相同,飽和導(dǎo)水率為0.025 4 cm/min,而犁底層飽和導(dǎo)水率為0.012 0 cm/min,入滲能力減小顯著,壤中流流量與垂直入滲量相比較取決于犁底層坡面坡度,我國(guó)東北黑土區(qū)坡耕地大部分較為平坦,取土樣地平均坡度7°,壤中流相對(duì)于垂直入滲量較小。雖然圖3顯示出邊界效應(yīng)的存在,但是模型模擬值與實(shí)測(cè)值仍然較為接近,2種坡度條件下RMSE分別為1.31%,1.24%(重量含水量),R2分別為0.94,0.97,說(shuō)明Richard模型可以反映出東北黑土區(qū)坡耕地存在犁底層條件下的土壤水分入滲過(guò)程。
圖3 降雨結(jié)束后土壤剖面含水量實(shí)測(cè)值與模型模擬值
圖4為2種處理?xiàng)l件下降雨過(guò)程中坡面徑流過(guò)程,兩種處理?xiàng)l件下產(chǎn)流時(shí)刻較為接近,實(shí)測(cè)值分別為18.21,17.43 min,可見緩坡對(duì)于徑流產(chǎn)生時(shí)刻影響較小,這一研究成果與Dong等[10]研究成果較為接近。坡度對(duì)于隨著降雨時(shí)間的延續(xù),流量顯著增加,2種處理?xiàng)l件下均接近250 ml/min趨于穩(wěn)定,壤中流流量約為地表徑流量1/100,這一比例與降雨強(qiáng)度關(guān)系較大,降雨強(qiáng)度越大,比值越小。降雨強(qiáng)度也決定著土壤剖面含水量的變化,降雨強(qiáng)度越大,濕潤(rùn)峰向下遷移越早,壤中流產(chǎn)生時(shí)刻越早,使得壤中流流量過(guò)程整體前移,形狀并未改變。雖然壤中流占徑流總量比例較小,但壤中流養(yǎng)分濃度可近似認(rèn)為是相應(yīng)深度的土壤水養(yǎng)分濃度,據(jù)相關(guān)研究顯示坡面徑流養(yǎng)分濃度較土壤水養(yǎng)分濃度小幾個(gè)數(shù)量級(jí)[10],所以土壤養(yǎng)分的徑流流失主要是通過(guò)壤中流這一途徑。同時(shí),由于壤中流攜帶養(yǎng)分的水平遷移,也影響了植物根系對(duì)土壤養(yǎng)分的利用。
對(duì)模型的敏感性分析結(jié)果見圖5。利用模型分別模擬了不同坡面坡度、土壤初始含水量及降雨強(qiáng)度條件下的壤中流徑流過(guò)程。各條件下從降雨開始至產(chǎn)流時(shí)刻歷時(shí)見表3。坡面坡度對(duì)壤中流產(chǎn)流時(shí)刻影響較小(p>0.05),土壤初始含水量和降雨強(qiáng)度的變化引起壤中流出現(xiàn)時(shí)刻的變化,且變化規(guī)律與地表徑流相似。不同條件下壤中流產(chǎn)流時(shí)刻模擬結(jié)果較為合理,說(shuō)明該模型對(duì)壤中流產(chǎn)生機(jī)理反映接近實(shí)際情況。圖5A為改變坡面坡度后壤中流過(guò)程,坡面坡度對(duì)產(chǎn)流時(shí)刻影響較小,但對(duì)穩(wěn)定出流流量影響較大,隨著坡度增加顯著增加(p<0.05)。圖5B—5C顯示不同土壤初始含水量和降雨強(qiáng)度對(duì)壤中流過(guò)程穩(wěn)定出流之前影響較大(p<0.05),土壤初始含水量越高,降雨強(qiáng)度越大壤中流出流越早趨于穩(wěn)定,這與不同條件下壤中流出流過(guò)程理論分析基本一致,說(shuō)明模型對(duì)不同影響因素較為敏感,同時(shí)說(shuō)明模型能夠較為合理的反映壤中流出流過(guò)程。
圖4 2種坡度處理?xiàng)l件下降雨過(guò)程中坡面徑流過(guò)程
圖5 模型模擬不同坡度、初始含水量及降雨強(qiáng)度條件下的壤中流過(guò)程表3 不同條件下壤中流模型模擬產(chǎn)流時(shí)間
犁底層的存在導(dǎo)致壤中流的產(chǎn)生,土壤質(zhì)地一定條件下,坡面坡度是決定壤中流強(qiáng)度的關(guān)鍵影響因素。隨著坡度的增加,壤中流強(qiáng)度增加顯著,坡度對(duì)于壤中流及坡面地表徑流產(chǎn)生時(shí)刻影響較小。初始含水量及降雨強(qiáng)度僅對(duì)產(chǎn)流時(shí)間影響較大,對(duì)壤中流流量大小影響較小。降雨強(qiáng)度增加、土壤初始含水量增高均能使壤中流產(chǎn)流時(shí)間縮短。壤中流流量顯著小于土壤表面徑流,但其對(duì)土壤養(yǎng)分水平遷移的貢獻(xiàn)可能高于土表徑流?;谶B續(xù)模擬犁底層上部土壤及犁底層內(nèi)部土壤水分垂直遷移的Richard模型,結(jié)合動(dòng)力波模型及“雙超”模型原理構(gòu)建的壤中流水平遷移機(jī)理模型能夠準(zhǔn)確反映黑土區(qū)坡耕地壤中流徑流過(guò)程。這一模型的建立是基于黑土區(qū)坡耕地這一特定條件,其他下墊面條件下的壤中流水平遷移模擬適用性有待進(jìn)一步驗(yàn)證,同時(shí),該模型適用于壤中流流量上升及穩(wěn)定階段的模擬計(jì)算,對(duì)于下降段不適用。