黃土區(qū)裸露坡地徑流養(yǎng)分流失模型的建立與驗(yàn)證*

邵凡凡，吳軍虎，李玉晨

（西安理工大學(xué)省部共建西北旱區(qū)生態(tài)水利國家重點(diǎn)實(shí)驗(yàn)室，西安 710048）

溶質(zhì)從表層土壤進(jìn)入地表徑流的過程非常復(fù)雜。雨滴擊濺、土壤中的養(yǎng)分濃度以及對流、擴(kuò)散作用和土壤顆粒吸附均會對這一過程產(chǎn)生影響[1-3]。Ahuja等[4]1981年通過在飽和土層不同深度位置放置一定量的32P，發(fā)現(xiàn)表層土壤的溶質(zhì)進(jìn)入徑流的可能性最高，隨著土層深度的增加，進(jìn)入徑流的概率呈指數(shù)遞減。然后，Ahuja和Lehman[5]在1983年發(fā)現(xiàn)，觀察到的交換層深度較通過擬合模型和實(shí)測數(shù)據(jù)獲得的深度小得多。結(jié)合交換層理論和Rose土壤侵蝕模型，Gao等[6]在2003年建立了基于物理過程并考慮雨滴擊濺和擴(kuò)散作用的溶質(zhì)運(yùn)移模型；該模型的所有參數(shù)均可通過實(shí)驗(yàn)測量獲得。鑒于黃土高原地區(qū)養(yǎng)分流失的特點(diǎn)，王全九等[7]提出了一種新的方法來改進(jìn)等效對流傳質(zhì)模型。在此基礎(chǔ)上，Dong等[8]假設(shè)交換層被混合層代替，且交換率被雨滴誘導(dǎo)水分轉(zhuǎn)移率所代替。Yang等[9]結(jié)合質(zhì)量守恒方程和降雨對表層土壤的剝蝕過程，建立了可預(yù)測黃土高原地表養(yǎng)分隨徑流流失的數(shù)學(xué)模型。但是，這些模型無法描述徑流發(fā)生之前土壤溶質(zhì)濃度的變化，并且僅能通過擬合曲線來獲得初始土壤溶質(zhì)濃度。Tong等[10]基于質(zhì)量守恒方程和水平衡方程建立了一維兩層的溶質(zhì)運(yùn)移模型。該模型結(jié)合了入滲和擴(kuò)散作用，并由不完全混合參數(shù)來描述。基于該模型，Tong等[11]使用集合卡爾曼濾波數(shù)據(jù)同化方法（EnKF）來校準(zhǔn)參數(shù)并更新可溶性化學(xué)物質(zhì)從土壤至地表徑流的轉(zhuǎn)移過程，并消除了實(shí)驗(yàn)觀測數(shù)據(jù)的誤差。因此，基于混合層理論的模型由于具有明晰的物理意義而被廣泛用于預(yù)測斜坡上的溶質(zhì)運(yùn)移。

養(yǎng)分流失的模擬是在徑流模擬的基礎(chǔ)上進(jìn)行的。通常使用 Saint-Venant方程（即連續(xù)性方程和動量方程）來描述地表徑流過程[12]。然而，由于Saint-Venant方程是高度非線性的，很難獲得解析解，這意味著僅能使用數(shù)值方法對其進(jìn)行求解。但當(dāng)忽略 Saint-Venant方程的加速度項(xiàng)時，可使用擴(kuò)散波方程來對其進(jìn)行簡化[13]。若同時忽略Saint-Venant方程的加速度和壓力項(xiàng)時，Saint-Venant方程可表示為運(yùn)動波方程。Luce和Cundy[14]通過使用菲利普（Philip）入滲方程修改了運(yùn)動波方程來預(yù)測超滲降雨條件下的產(chǎn)流過程。Yang等[15]通過假設(shè)水深與入滲率之間的線性關(guān)系簡化運(yùn)動波模型中的水深項(xiàng)，并結(jié)合Philip入滲方程得到了運(yùn)動波模型的近似解析解。該模型因其參數(shù)簡單易獲取而被廣泛使用于坡面徑流的模擬。在養(yǎng)分隨徑流流失的研究中，Gao等[16]在2004年提出了一個基于溶質(zhì)守恒方程的模型，該模型考慮了雨滴飛濺和徑流沖刷作用。但是，該模型僅用于模擬積水條件下飽和土壤的養(yǎng)分流失過程，這與黃土區(qū)初始非飽和土壤條件下的流失過程存在較大差異。因此，本文以 Yang等[15]建立的坡面徑流近似解析解為基礎(chǔ)，進(jìn)一步延伸于養(yǎng)分隨地表徑流流失過程的模擬中，并修改了Gao等[16]的模型以適應(yīng)本文的測試條件。通過模型參數(shù)分析揭示了不同因素對養(yǎng)分流失的貢獻(xiàn)作用，提出了防止養(yǎng)分流失的有效措施。該研究可為防治農(nóng)田退化和農(nóng)業(yè)面源污染提供有力基礎(chǔ)。

1 材料與方法

1.1 理論與模型

1.1.1 坡面徑流運(yùn)動過程 采用運(yùn)動波模型來描述次降雨條件下的坡面水流流動過程[15]，其中超滲凈雨可用降雨強(qiáng)度與入滲率的差值來表示，見式（1）：

式中，h為徑流深，cm；t為徑流時間，min；q為單寬流量，cm2·min–1；x為坡面任一位置距離入流口的長度，cm；p為降雨強(qiáng)度，mm·h–1；i為土壤入滲率，cm·min–1。

由于坡面水深與入滲率之間存在關(guān)聯(lián)關(guān)系[9]，徑流水量為超滲凈雨所產(chǎn)生的，Yang等[15]用Philip公式表示降雨條件下的入滲過程，進(jìn)一步求解了單寬流量和坡面水深，見式（2）和式（3）：

式中，c為入滲率參數(shù)[15]；S為土壤吸滲率，cm·min–0.5；Δt=3S2/（16p2）。

式中，n為曼寧糙率系數(shù)，；S0為水力梯度，本文中坡度為15°，故S0為Sin 15°。

1.1.2 徑流養(yǎng)分流失過程 降雨條件下土壤表層養(yǎng)分在雨滴擊濺和水分入滲的作用下隨徑流遷移并在土壤中重新分配，因此土壤剖面的水和養(yǎng)分運(yùn)移系統(tǒng)自上而下可分為3層：徑流積水層、養(yǎng)分交換層和交換層以下土壤，如圖1所示。

養(yǎng)分交換層是徑流積水層與土壤剖面交界面以下厚度較薄的土層。交換層中化學(xué)物質(zhì)的傳輸主要受入滲、水動力彌散和雨滴飛濺侵蝕控制[16-17]。

式中，de為交換層深度，cm；Ce為交換層中溶質(zhì)濃度，mg·L–1；Cw為徑流中溶質(zhì)濃度，mg·L–1；is為徑流水進(jìn)入交換層的入滲速率，cm·min–1；ix為交換層水分進(jìn)入更深土層的入滲速率，cm·min–1；er為雨滴誘導(dǎo)水分轉(zhuǎn)移速率，cm·min–1；λCw為徑流層進(jìn)入交換層的溶質(zhì)濃度（0 ≤λ≤ 1，Gao等[6]研究表明計(jì)算模型對該參數(shù)不敏感，取λ=0），mg·L–1；J為較深土壤層與交換層內(nèi)部的溶質(zhì)擴(kuò)散通量，mg·cm–2·min–1。

為了簡化計(jì)算過程，Gao等[16]對J進(jìn)行近似求解如下：

式中，Ds為養(yǎng)分在土壤中的擴(kuò)散性，cm2·min–1；Cs為更深土層的溶質(zhì)濃度，mg·g–1；γ為土壤容重，g·cm–3；K為土壤吸附系數(shù)，mL·g–1；β=er/（αde）。

從降雨開始，將整個降雨過程劃分為3個階段。

第一階段：從降雨開始t0至交換層完全飽和tsa。在這一階段，土壤入滲率為降雨強(qiáng)度，土壤表層未產(chǎn)生徑流，故i=p，q=0，交換層完全飽和的時間tsa可以表示為：

式中，tsa為交換層完全飽和所需時間，min；θs為飽和含水率，cm3·cm–3；θ0為初始含水率，cm3·cm–3。

第二階段：從交換層完全飽和tsa至土壤表層出現(xiàn)徑流tp。在這個階段，徑流層溶質(zhì)濃度Cw和雨滴誘導(dǎo)水分轉(zhuǎn)移速率er的取值為0，ix = p。交換層中溶質(zhì)濃度見式（9）：

式中，C0為初始溶質(zhì)濃度，mg·L–1。

當(dāng)式（9）中t=tp時，開始產(chǎn)流時刻交換層中的溶質(zhì)濃度見式（10）：

式中，A=ix/（αde）。

第三階段：從開始產(chǎn)流至降雨結(jié)束。這一過程中，徑流中的養(yǎng)分濃度遠(yuǎn)低于交換層中的養(yǎng)分濃度，因此為了簡化計(jì)算過程，忽略了徑流養(yǎng)分在入滲作用下對交換層養(yǎng)分的微小補(bǔ)給作用，ix=0.01 cm·min–1。結(jié)合這一階段的起始產(chǎn)流時間，即：t=tp，可求解得到交換層中溶質(zhì)濃度表示如下：

式中，B=（ix+er）/（αde）。

產(chǎn)流過程中，徑流中化學(xué)溶質(zhì)的質(zhì)量守恒關(guān)系可表示為：

結(jié)合式（1）和式（12），可得到：

為了簡化計(jì)算，忽略了入滲和擴(kuò)散作用，式（13）可表示為：

對式（14）進(jìn)行積分得到徑流中化學(xué)溶質(zhì)的濃度，見式（15）：

徑流中化學(xué)溶質(zhì)的濃度與徑流量的乘積便為徑流中化學(xué)溶質(zhì)的流失速率，結(jié)合式（2）、式（3）和式（15），徑流中溶質(zhì)流失速率可表示為：

式中，Mw為養(yǎng)分隨徑流流失的速率，mg·min–1。

1.2 試驗(yàn)區(qū)概況

試驗(yàn)于2019年5月在中國科學(xué)院水利部水土保持研究所長武黃土高原農(nóng)業(yè)生態(tài)試驗(yàn)站擔(dān)水溝流域野外模擬降雨小區(qū)（35°12′N，107°10′E）進(jìn)行，試驗(yàn)區(qū)平均海拔為1 200 m，氣候?qū)倥瘻貛О霛駶櫞箨懶约撅L(fēng)氣候，年平均氣溫為 9.1℃，年平均降水量580 mm，地下水位50～80 m，無灌溉條件，屬典型的旱作雨養(yǎng)農(nóng)業(yè)區(qū)。該流域內(nèi)典型土壤為粉砂質(zhì)壤土，母質(zhì)為深厚的中壤質(zhì)馬蘭黃土，具體土壤物理特性見表1，流域塬面面積占35%，梁坡占35.6%，溝谷占29.4%，各約占1/3；流域地貌屬典型的黃土高原溝壑區(qū)[18-20]。

表1 試驗(yàn)區(qū)土壤物理化學(xué)特性Table 1 Physical and chemical properties of the soil tested

1.3 試驗(yàn)方法

試驗(yàn)用地為3年閑置坡耕地，模擬降雨試驗(yàn)的小區(qū)設(shè)置尺寸為1.0 m×1.0 m，并根據(jù)當(dāng)?shù)氐湫推赂仄露群推旅媲治g的臨界坡度，設(shè)置小區(qū)坡度為15°；使用針孔式人工模擬降雨裝置進(jìn)行降雨試驗(yàn)（圖2），其主要由主體支架（可調(diào)節(jié)高度）、底板布設(shè)有針孔的水槽（可根據(jù)降雨強(qiáng)度更換不同孔徑的針孔）和供水裝置 3部分組成，有效降雨面積為1.0 m2。經(jīng)過測試：該套人工降雨器的平均雨滴直徑為2 mm，降雨均勻度在80%以上，雨滴終速符合天然降雨特征[3]。

為了消除土壤前期含水率對試驗(yàn)結(jié)果的影響，每次開始降雨24 h前以25 mm·h–1的降雨強(qiáng)度在試驗(yàn)小區(qū)進(jìn)行預(yù)降雨，直至開始產(chǎn)流時停止降雨。開始正式降雨試驗(yàn)前測定小區(qū)內(nèi)土壤表層 0～20 cm剖面的初始含水率，采用烘干法測定質(zhì)量含水率為0.11± 0.003 7 g·g–1（ 即 體 積 含 水 率 為 0.15±0.005 cm3·cm–3）時開始試驗(yàn)。為了提高土壤初始養(yǎng)分濃度值，使養(yǎng)分在土壤中均勻分布，預(yù)降雨結(jié)束后，在小區(qū)土壤表面均勻噴灑氯化銨和硝酸鉀混合溶液，其具體操作方式為：將預(yù)先配置好的5.0 g·L–1的氯化銨溶液和10.0 g·L–1的硝酸鉀溶液各取1 L混合均勻（為消除溶液噴施次序?qū)ζ浞植嫉挠绊懀⒂脡毫妷卦诿總€小區(qū)分別定量噴灑2 L氯化銨和硝酸鉀的混合溶液，為盡可能減小噴壺壓力對表層土壤造成的壓實(shí)作用，將噴嘴調(diào)節(jié)至霧化度最強(qiáng)位置處，噴嘴霧化半徑為5 cm，采用左右往復(fù)的方式將混合溶液均勻噴灑在小區(qū)表層，使其在2 cm范圍內(nèi)均勻分布。并在小區(qū)坡面上、中、下三個部位分別取土樣來測定表層 2 cm土壤中的溶質(zhì)濃度作為模型計(jì)算所需的初始土壤養(yǎng)分濃度（表2）。根據(jù)研究區(qū)暴雨實(shí)測資料及降雨分級標(biāo)準(zhǔn)[20]，設(shè)計(jì)30、45、60、75、90 mm·h–1的 5種降雨強(qiáng)度，設(shè)計(jì)總降雨歷時為120 min，按照0～10 min之間，每隔2 min承接1次徑流，10～120 min之間，每隔5 min承接1次徑流的頻率用量杯承接出口處徑流，并用量筒進(jìn)行精確測量，通過沉淀過濾除去徑流中的泥沙，用50 mL的塑料瓶收集徑流水樣并存放于實(shí)驗(yàn)室冰箱中，用全自動高通量間斷分析儀（SmartChem450，AMS Allinace公司，意大利）測定徑流和土壤中的養(yǎng)分濃度。

1.4 模型基本參數(shù)

通過現(xiàn)場測定和文獻(xiàn)查閱等方式獲取了模型計(jì)算所需參數(shù)[20-21]（表2）：

表2 模型基本參數(shù)Table 2 Basic parameters of the calculation model

1.5 數(shù)據(jù)處理

所有試驗(yàn)實(shí)測數(shù)據(jù)均為 3次重復(fù)試驗(yàn)的平均值，使用Matlab 2015b進(jìn)行參數(shù)求解和模型模擬；使用SPSS 24.0進(jìn)行數(shù)據(jù)分析，使用Origin 2018進(jìn)行圖表繪制和函數(shù)擬合。

2 結(jié)果與討論

2.1 產(chǎn)流過程分析及模擬

2.1.1 產(chǎn)流過程分析 土壤表層的養(yǎng)分通常會隨著地表徑流而流失。因此，探究產(chǎn)流規(guī)律是模型準(zhǔn)確預(yù)測養(yǎng)分隨徑流流失過程的基礎(chǔ)。在5種降雨強(qiáng)度下，開始產(chǎn)流的時間點(diǎn)分別為20.5、8.5、4.8、3.0和 1.8 min（圖3）；90 mm·h?1較 30 mm·h?1提前19 min產(chǎn)流；說明隨降雨強(qiáng)度的增大，開始產(chǎn)流時間開始顯著縮短。降雨強(qiáng)度與產(chǎn)流時間的關(guān)系可用冪函數(shù)來描述，R2=0.997 6。單寬流量開始產(chǎn)流后快速增大，而后進(jìn)入穩(wěn)定產(chǎn)流階段（圖4）。這是由于表層土體中的黏粒分散堵塞了土壤的孔隙，并伴隨著雨滴的飛濺使表層土壤變得密實(shí)，降低了土壤的入滲能力[13]。5種降雨強(qiáng)度下實(shí)測單寬流量的標(biāo)準(zhǔn)差分別為 0.01～0.16、0.01～0.16、0.10～0.41、0.13～0.71和0.23～0.69。這可能是由于土壤的非均質(zhì)性、蟻穴和植物根系對小區(qū)土壤入滲過程的影響，以及雨滴飛濺和徑流發(fā)育過程中微地形的形成，可能導(dǎo)致徑流滯后。進(jìn)一步分析發(fā)現(xiàn)，90 mm·h?1在穩(wěn)定產(chǎn)流階段的單寬流量分別較其他降雨強(qiáng)度依次增加6.3倍、2.7倍、1.6倍和1.2倍，這表明降雨強(qiáng)度的增加顯著增大了坡面徑流率。

2.1.2 產(chǎn)流過程模擬 通過將已知參數(shù)S、p（表2）和實(shí)測單寬流量代入式（2）中推求入滲率參數(shù)c，采用R2、均方根誤差（RMSE）和納什效率系數(shù)（NSE）對模擬結(jié)果進(jìn)行評價。從表3可看出，入滲率參數(shù)c隨著降雨強(qiáng)度的增大呈減小趨勢，且分布在0.003 1～0.006 0之間；R2均在0.89以上，隨著降雨強(qiáng)度的增大，RMSE也隨之增大，取值分布在0.406～1.052之間，NSE均大于0.397，而當(dāng)降雨強(qiáng)度大于等于 60 mm·h–1時，NSE 則進(jìn)一步增大至0.783以上；說明降雨強(qiáng)度越大，模型計(jì)算值和實(shí)測值的匹配度也隨之提高。指數(shù)函數(shù)可很好地?cái)M合參數(shù)c與降雨強(qiáng)度之間的關(guān)系（圖5），決定系數(shù)R2為0.976 5，表達(dá)式為：

表3 入滲率參數(shù)c的最佳擬合值Table 3 Optimal fitting values of c，R2，RMSE and NSE relative to rainfall intensity

圖6分別顯示了5種降雨強(qiáng)度下的單寬流量模擬過程。可以看出，產(chǎn)流模型能夠較好地模擬地表徑流過程，且隨著降雨強(qiáng)度的增大，模擬趨勢變得更加準(zhǔn)確。在產(chǎn)流初期，計(jì)算值的上升趨勢均慢于實(shí)測值；在穩(wěn)定產(chǎn)流階段，30 和45 mm·h?1下的模擬值均大于實(shí)測值，而在 60、75、90 mm·h?1下，實(shí)測值與計(jì)算值的匹配程度較好。這可能是由于降雨強(qiáng)度較小時，雨滴動能的濺蝕增加了前期土壤表面粗糙度；雨滴擊濺形成的微地形和洼地?cái)r截部分地面徑流，從而削減了連續(xù)徑流的沖刷作用。這表明產(chǎn)流模型可準(zhǔn)確模擬大于等于 60 mm·h?1的產(chǎn)流過程。此外，雨滴的擊濺使表層土壤被壓實(shí)，容重增大，土壤表層形成密封層，降低了土壤的入滲能力。但該模型未考慮地表土壤容重和孔隙率的變化，導(dǎo)致模型計(jì)算的土壤入滲能力明顯大于實(shí)測值。

2.2 養(yǎng)分流失過程分析及模擬

2.2.1 養(yǎng)分隨徑流流失過程分析 養(yǎng)分從土壤至徑流的傳輸是通過雨滴擊濺作用和徑流溶解作用來完成的[16]。不同降雨強(qiáng)度下硝態(tài)氮和銨態(tài)氮隨時間的流失過程可用單峰形式來描述；即：徑流初期養(yǎng)分流失速率迅速增大，到達(dá)峰值后開始減少，最后進(jìn)入穩(wěn)定流失階段的趨勢（圖7）。這是由于土壤水和土壤顆粒表面吸附的養(yǎng)分在雨滴擊濺作用下進(jìn)入徑流所引起的；隨著土壤結(jié)皮厚度和徑流深度的增大，土壤表層形成“堅(jiān)實(shí)的保護(hù)殼”削弱了雨滴動能，延緩?fù)寥浪蛷搅鞯慕粨Q作用，使得進(jìn)入徑流的土壤水和溶解態(tài)氮顯著減少。此外，徑流中氮素濃度的降低也是由于隨著降雨過程的推移，表層土壤中氮含量的逐漸減少所引起的。同時，降雨強(qiáng)度、入滲能力、養(yǎng)分濃度和徑流率對養(yǎng)分的峰值流失速率及其發(fā)生時間均有一定的影響。簡言之，硝態(tài)氮和銨態(tài)氮的峰值損失率隨降雨強(qiáng)度的增大而增大。以硝態(tài)氮流失過程為例，當(dāng)降雨強(qiáng)度為 30 mm·h–1時，硝態(tài)氮流失速率在 25 min左右達(dá)到峰值5.74 mg·min–1，而當(dāng)降雨強(qiáng)度為 45、60、75 和90 mm·h–1時，分別在 13、8、6和 5 min達(dá)到硝態(tài)氮流失速率的峰值：35.21、121.3、280.4和468.4 mg·min–1。因此，降雨強(qiáng)度對硝態(tài)氮的峰值流失速率具有較大貢獻(xiàn)。通過對比養(yǎng)分的峰值流失速率出現(xiàn)時間和穩(wěn)定產(chǎn)流時間，各降雨強(qiáng)度下養(yǎng)分流失速率峰現(xiàn)時間分別為25、13、8、6、5 min，而穩(wěn)定產(chǎn)流時間分別為40、28、18、16、14 min；由此看出，硝態(tài)氮流失速率的峰值出現(xiàn)時間要早于徑流速率達(dá)到穩(wěn)定階段所需的時間，這可能是隨著產(chǎn)流時間的推移，交換層土壤中硝態(tài)氮濃度的不斷減小和徑流量的增大共同作用所造成的。5種降雨強(qiáng)度下硝態(tài)氮損失率的標(biāo)準(zhǔn)誤差分別分布在0.03～0.53、0.05～5.13、0.03～6.9、0.1～21.0 和 1.1～31.7。

2.2.2 養(yǎng)分流失過程模擬 交換層深度de和雨滴誘導(dǎo)水分轉(zhuǎn)移速率er是養(yǎng)分流失模型中的兩個重要參數(shù)。由于受室外實(shí)驗(yàn)條件的限制，交換層深度de很難通過實(shí)地測量得到，因此需借助模型擬合實(shí)測的養(yǎng)分流失速率來反推交換層深度de。研究[4-5]發(fā)現(xiàn)交換層的深度de在2～3 mm的范圍內(nèi)。Tong等[11]指出，交換層深度隨入滲率的增加而減小。有研究[21-23]指出，交換層深度隨著初始含水量的增大而增大。關(guān)于雨滴誘導(dǎo)水分轉(zhuǎn)移速率er，Gao等[16]在2004年提出了適用于初始飽和土壤的雨滴誘導(dǎo)水分轉(zhuǎn)移率er的計(jì)算方法，由于黃土區(qū)坡耕地在降雨前為非飽和土壤，因此，借助 Matlab非線性擬合的方法，將式（3）所計(jì)算出的坡面徑流深度代入式（16）來計(jì)算養(yǎng)分流失速率，并結(jié)合硝態(tài)氮流失速率的實(shí)測值推求出了式（16）中的交換層深度de和雨滴誘導(dǎo)水分轉(zhuǎn)移速率er（表4），并進(jìn)一步模擬了銨態(tài)氮流失過程。可以看出：de和er均隨著降雨強(qiáng)度的增大而增大，其分別從0.68增至1.32、從0.006增至0.023。這與前述單寬流量隨著降雨強(qiáng)度的增大而增加是一致的[24]。硝態(tài)氮和銨態(tài)氮流失速率的R2值分別分布于 0.834～0.922和 0.800～0.921之間，RMSE值分別分布在1.188～58.50和0.974～58.37之間，NSE值分別分布在0.653～0.881和0.546～0.775之間。從圖8中可以看出，計(jì)算出的曲線可很好地模擬養(yǎng)分流失過程。當(dāng)降雨強(qiáng)度為30 mm·h–1時，初始增大階段的測量值與計(jì)算值之間的差異較大，但隨著降雨強(qiáng)度的增大，差異逐漸減小。而在養(yǎng)分流失的穩(wěn)定階段，實(shí)測數(shù)據(jù)大于計(jì)算值。這可能是由于本文建立的養(yǎng)分流失模型近似求解了對流彌散項(xiàng)所造成的，這使得穩(wěn)定減小階段養(yǎng)分流失速率衰減得過快[9]。同時模型忽略了徑流層養(yǎng)分對交換層的微弱補(bǔ)給作用，然而在降雨開始時可能存在從交換層至徑流層的擴(kuò)散過程[16]。從交換層完全飽和（tsa）到坡面開始產(chǎn)流（tp）的時間段內(nèi)，模型假設(shè)交換層中的養(yǎng)分隨入滲水向土壤深處遷移的速率大小即為對應(yīng)的降雨強(qiáng)度，從而使得計(jì)算出的開始產(chǎn)流時土壤表層的濃度Ce（tp）小于理論值。

表4 不同降雨強(qiáng)度下參數(shù)de和er的最佳擬合值Table 4 Optimal fitting values of de，er，R2，RMSE and NSE relative to rainfall intensity

以上分析表明本文建立的模型可以很好地模擬裸露坡面徑流和養(yǎng)分隨徑流遷移過程。但該模型未用于模擬不同坡度、坡長和土壤初始含水量條件下的養(yǎng)分流失過程，本研究所獲得的參數(shù)是否具有普遍適用性，需要在以后的研究中加以驗(yàn)證?？梢灶A(yù)見的是，坡長的增大將顯著增加徑流量和泥沙量，初始含水量的增大將提前產(chǎn)流時間并增加養(yǎng)分的峰值流失速率[25，26]。同時，土壤中植物根系的生長和土壤生物活性可能形成連通的土壤孔隙結(jié)構(gòu)，導(dǎo)致優(yōu)先流的出現(xiàn)，這將對模型的模擬精度產(chǎn)生較大影響。泥沙顆粒中通常吸附有大量養(yǎng)分，但本研究的模型并未將泥沙考慮在內(nèi)，使得模型并不能完整模擬徑流過程所帶走的養(yǎng)分總量。簡言之，該近似解析模型充分考慮了非飽和土壤水分入滲對交換層中養(yǎng)分運(yùn)移過程的影響，因此，該模型可用于預(yù)測干旱和半干旱氣候條件下裸露坡耕地的養(yǎng)分流失過程。但是，徑流過程的精確計(jì)算是進(jìn)行養(yǎng)分流失模擬的基礎(chǔ)，應(yīng)根據(jù)土壤質(zhì)地、養(yǎng)分類型和雨水中養(yǎng)分濃度選擇合適的入滲公式和溶質(zhì)吸附系數(shù)。

3 結(jié) 論

本研究以交換層理論為基礎(chǔ)，根據(jù)黃土區(qū)降雨量少，降雨前土壤通常為非飽和狀態(tài)，其產(chǎn)流需要較長時間的實(shí)際情況對降雨過程進(jìn)行劃分，建立了基于坡面徑流養(yǎng)分遷移理論的機(jī)理模型，并通過 5個降雨強(qiáng)度的模擬降雨試驗(yàn)對模型進(jìn)行了驗(yàn)證。3組實(shí)測重復(fù)試驗(yàn)間的標(biāo)準(zhǔn)誤差均較小，試驗(yàn)結(jié)果具有可靠性。模型驗(yàn)證的結(jié)果表明，本文建立的養(yǎng)分流失近似解析模型能夠準(zhǔn)確描述不同降雨強(qiáng)度下的坡面流和養(yǎng)分流失特征（R2> 0.8，NSE > 0.347）。參數(shù)c（入滲率參數(shù)）、de（交換層深度）、er（雨滴誘導(dǎo)水分轉(zhuǎn)移速率）均隨降雨強(qiáng)度的增大而增大。養(yǎng)分流失模型對交換層的深度de較雨滴誘導(dǎo)水分轉(zhuǎn)移速率er更敏感，de可顯著影響可交換溶質(zhì)的量。因此，在施肥過程中應(yīng)采取一些措施，如施肥后覆蓋坡面土壤、暴雨前避免施肥等，以達(dá)到防止土壤貧瘠化并控制農(nóng)業(yè)面源污染的目的。