應(yīng)用CNMM-DNDC模擬小流域土壤水力侵蝕和顆粒態(tài)碳氮磷遷移*

李思琪, 李 勇**, 張 偉, 鄭循華, 胡鵬程, 范繼輝, 汪 濤, 朱 波**

(1.大氣邊界層物理與大氣化學(xué)國家重點(diǎn)實(shí)驗(yàn)室/中國科學(xué)院大氣物理研究所 北京 100029; 2.中國科學(xué)院水利部成都山地災(zāi)害與環(huán)境研究所 成都 610041)

土壤侵蝕及其導(dǎo)致的顆粒態(tài)碳氮磷流失是引起流域或區(qū)域土壤肥力下降、土地退化、河道淤積和河流湖泊富營養(yǎng)化等生態(tài)環(huán)境問題的主要因素之一。其中, 水力侵蝕(簡稱“水蝕”)總面積約占國土總面積的17%, 是我國分布最廣、危害最嚴(yán)重的水土流失類型。因此, 定量評估區(qū)域或流域土壤水力侵蝕強(qiáng)度并識別其空間分布特征, 對制定水土流失防治措施具有重要意義, 是國內(nèi)外學(xué)者一直以來的重要研究課題。某區(qū)域或流域單位時(shí)間內(nèi)侵蝕強(qiáng)度的大小一般用土壤侵蝕模數(shù)來表征。土壤侵蝕模數(shù)的大小與氣候條件、地形、土壤屬性、土地利用和管理措施等密切相關(guān)。我國不同地區(qū)的水土流失程度不同, 土壤侵蝕模數(shù)差異較大。趙明松等研究結(jié)果表明, 2010年安徽省土壤侵蝕模數(shù)平均值為257 t?km?a; 而高海東等指出2010年黃土高原土壤侵蝕模數(shù)平均值高達(dá)3355 t?km?a。

用徑流小區(qū)試驗(yàn)測定土壤水力侵蝕強(qiáng)度的現(xiàn)場測定方法人力和物力成本高, 期望通過這種直接測定方法量化各個(gè)地區(qū)不同條件下的水力侵蝕強(qiáng)度是不現(xiàn)實(shí)的, 更無法實(shí)現(xiàn)實(shí)時(shí)監(jiān)測和預(yù)測。利用土壤水力侵蝕模型和地理空間信息相結(jié)合的方式有助于彌補(bǔ)現(xiàn)場直接測量方法的不足。目前, 國內(nèi)外廣泛使用的通用土壤流失方程(universal soil loss equation,USLE)及其改進(jìn)版RUSLE (revised universal soil loss equation)是量化水蝕強(qiáng)度的一種經(jīng)驗(yàn)統(tǒng)計(jì)模型, 該模型利用降水、土壤性質(zhì)、地形、植被覆蓋和管理信息等因子直接估算水蝕強(qiáng)度。USLE僅對影響水蝕的各種因素進(jìn)行量化, 不涉及具體的地表徑流和泥沙輸移的過程和機(jī)理。近年來發(fā)展起來的水文物理模型—水文評價(jià)模型(soil and water assessment tool,SWAT), 耦合RUSLE來預(yù)測水文響應(yīng)單元的土壤水力侵蝕, 但SWAT是以水文響應(yīng)單元為基本模擬單元的半分布式水文模型, 且土壤水力侵蝕模塊仍是基于經(jīng)驗(yàn)統(tǒng)計(jì)方法計(jì)算, 且不考慮泥沙遷移路徑。美國農(nóng)業(yè)部推出的土壤侵蝕預(yù)報(bào)模型(water esrosion prediction project, WEPP)是基于物理過程模擬土壤水力侵蝕的機(jī)理模型, 但模型需要的細(xì)溝間可蝕性、細(xì)溝可蝕性和土壤臨界剪切力3個(gè)關(guān)鍵參數(shù)很難獲取, 同時(shí)WEPP模型的生物地球化學(xué)元素循環(huán)的建模理念相對簡單。ROSE模型(model of soil erosion and depostion processes)是一個(gè)基于物理過程的單純的土壤水力侵蝕模型, 通過計(jì)算降雨剝離、地表徑流夾帶和泥沙沉降等物理過程來模擬土壤水力侵蝕過程, 該模型在對美國亞利桑那州一個(gè)1.3 hm地塊的土壤水力侵蝕流失量的估算中表現(xiàn)出較好的模擬效果, 模擬值與觀測值之間的誤差僅為-4%～15%。但ROSE模型只考慮與水力侵蝕直接相關(guān)的物理過程, 不涉及碳氮循環(huán)過程的模擬。目前, 由地表徑流引起的顆粒態(tài)碳氮磷遷移的模擬和評估都是依托在土壤水力侵蝕模擬的基礎(chǔ)上, 應(yīng)用元素富集系數(shù)進(jìn)行估算。綜上所述, 目前用于估算土壤水力侵蝕以及隨之的顆粒態(tài)碳氮磷流失的物理模型,要么是半分布式, 要么其碳氮循環(huán)過程建模理念過于簡化。值得強(qiáng)調(diào)的是, 水蝕及其顆粒態(tài)碳氮磷營養(yǎng)物遷移在很大程度上受到流域生態(tài)系統(tǒng)的碳氮循環(huán)過程的影響。因此, 具備復(fù)雜碳氮循環(huán)過程的分布式水文生物地球化學(xué)模型與土壤水力侵蝕物理模型之間的耦合, 對水蝕及其導(dǎo)致的顆粒態(tài)碳氮磷遷移的準(zhǔn)確預(yù)測具有重要意義。

紫色土集中分布在我國的川中丘陵, 其養(yǎng)分儲(chǔ)量豐富, 多作為耕作土壤, 不合理的管理措施和頻繁耕作, 加之多以大到暴雨形式出現(xiàn)的夏季降水, 導(dǎo)致其水土流失問題不容忽視, 評估其土壤水力侵蝕強(qiáng)度具有重要意義。目前, 我國學(xué)者對紫色土區(qū)域的水力侵蝕模擬進(jìn)行了一定程度的研究, 但這些研究大多停留在對現(xiàn)有模型的評估和應(yīng)用上。例如, 繆馳遠(yuǎn)和嚴(yán)冬春等對紫色土水力侵蝕以及隨之的顆粒態(tài)碳氮磷流失的模擬研究結(jié)果表明, WEPP比USLE表現(xiàn)出更好的模擬效果。薛菲等應(yīng)用SWAT模型成功模擬了紫色土丘陵區(qū)典型農(nóng)業(yè)小流域的非點(diǎn)源氮磷月流失規(guī)律與負(fù)荷。Deng等嘗試把徑流曲線方程和RUSLE引入到DNDC (DeNitrification-DeComposition)模型, 改進(jìn)模型對川中紫色土典型農(nóng)林復(fù)合小流域的土壤水力侵蝕和顆粒態(tài)氮素?fù)p失均表現(xiàn)出良好模擬效果, 但該模型的地表徑流和土壤水力侵蝕過程仍是基于經(jīng)驗(yàn)統(tǒng)計(jì)方法建模, 應(yīng)用到其他流域仍需大量的評估驗(yàn)證工作且存在較大不確定性。

本研究使用的分布式流域水碳氮耦合模型(CNMM-DNDC)的構(gòu)建核心是將生物地球化學(xué)模型DNDC復(fù)雜碳氮循環(huán)等核心動(dòng)力學(xué)過程(如反硝化、硝化、分解和發(fā)酵過程)耦合到流域生源要素管理模型(catchment nutrients management model, CNMM)的分布式水文模塊DHSVM框架內(nèi), 實(shí)現(xiàn)基于水文過程和碳氮循環(huán)機(jī)理模擬流域水碳氮循環(huán)。該模型具有高時(shí)間分辨率、高空間分辨率和高過程分辨率, 能夠用于模擬小流域生態(tài)系統(tǒng)的能量平衡、水文過程、植物生長和碳氮遷移轉(zhuǎn)化過程,進(jìn)而獲取流域尺度生態(tài)系統(tǒng)生產(chǎn)力、溫室氣體凈排放、污染氣體排放和流域水氮流失等關(guān)鍵碳氮變量。目前, CNMM-DNDC模型已經(jīng)實(shí)現(xiàn)對亞熱帶紫色土農(nóng)林復(fù)合小流域(流域徑流、作物產(chǎn)量以及NH、NO、NO和CH排放)、亞熱帶茶園(NO和NO排放)、東北黑土(NO—流失)和全國農(nóng)田63個(gè)施肥事件氨揮發(fā)的模擬驗(yàn)證工作, 并在點(diǎn)位驗(yàn)證的基礎(chǔ)上實(shí)現(xiàn)了對亞熱帶紫色土農(nóng)林復(fù)合小流域生產(chǎn)力、NO—流失以及NH、NO、NO和CH排放和東北黑土NO—流失的區(qū)域模擬。然而, CNMM-DNDC模型目前還缺乏對地表徑流侵蝕引起的泥沙遷移的模擬功能, 進(jìn)而尚不具備模擬顆粒態(tài)碳氮磷隨泥沙遷移的能力。因此, 本研究為CNMM-DNDC開發(fā)了一個(gè)水土流失模塊并對其進(jìn)行模擬效果評估和應(yīng)用。從模型過程和代碼以及驗(yàn)證數(shù)據(jù)的可獲取性考慮,本研究以鹽亭截流紫色土農(nóng)林復(fù)合小流域(以下簡稱“截流小流域”)為例, 把單純的土壤水力侵蝕物理模型ROSE和顆粒態(tài)碳氮磷元素富集模塊耦合到CNMM-DNDC的水文過程中, 實(shí)現(xiàn)在地塊和流域尺度土壤水力侵蝕和碳氮磷隨侵蝕泥沙運(yùn)移的高時(shí)空分辨率模擬。

1 材料和方法

1.1 研究區(qū)域概況

截流小流域位于四川鹽亭農(nóng)田生態(tài)系統(tǒng)國家野外科學(xué)觀測研究站(以下簡稱“鹽亭站”, 31°16′N,105°28′E, 海拔400～600 m)所在的研究區(qū)(圖1), 地處四川盆地中北部涪江支流—彌江和湍江的分水嶺上。截流小流域總面積約為35 hm, 農(nóng)田面積約19 hm, 林地面積約11 hm, 其他用地約5 hm, 是紫色土丘陵區(qū)典型農(nóng)林復(fù)合小流域。該區(qū)屬于亞熱帶濕潤季風(fēng)氣候區(qū), 2005—2018年年平均氣溫為16.6 ℃,年平均降水量為920 mm, 約75%降水集中分布在6—9月。徑流小區(qū)的土壤類型為石灰性紫色土, 土壤容重為1.33 g?cm, pH為8.83, 有機(jī)質(zhì)含量為8.75 g?kg, 總氮含量為0.62 g?kg。該流域旱地采用冬小麥(Triticum aestivum)-夏玉米(Zea mays)輪作制度,兩季田采用冬小麥-水稻(Oryza sativa)輪作和油菜(Brassica napus)-水稻兩種輪作制度, 灌溉水田采用水稻-休耕制度。用于模型驗(yàn)證的徑流小區(qū)內(nèi)種植制度為冬小麥-夏玉米輪作, 施肥、耕作等田間管理措施均為當(dāng)?shù)爻Ｒ?guī)管理方式。小麥季氮肥(碳酸氫銨)和磷肥(過磷酸鈣)施用量分別為130 kg(N)?hm和90 kg(PO)?hm; 玉米季氮肥(碳酸氫銨)和磷肥(過磷酸鈣)施用量分別為150 kg(N)?hm和90 kg(PO)?hm。

圖1 四川鹽亭截流小流域位置示意圖Fig.1 Location diagram of the Jieliu catchment in the Yanting, Sichuan Province

1.2 CNMM-DNDC模型和土壤水力侵蝕模塊

本研究通過在CNMM-DNDC模型的水文框架中耦合土壤水力侵蝕過程模塊, 使CNMM-DNDC模型具備模擬由地表徑流引起的侵蝕產(chǎn)沙的能力。土壤水力侵蝕產(chǎn)沙量采用簡化的ROSE模型, 根據(jù)地表徑流量和坡度、植被覆蓋度等因子計(jì)算由單次降水事件引起的侵蝕量, 表達(dá)式如公式(1):

式中: S是產(chǎn)沙量, 即單位面積由地表徑流引起的泥沙量, 單位為kg?hm; Q是由降水引起的地表徑流量, 單位為m, 采用CNMM-DNDC模型水文模塊計(jì)算地表徑流; G表示坡度(弧度值)的正弦值的絕對值, 無量綱; P表示沒有植被覆蓋的地表面積百分比;η是地表徑流中的泥沙攜帶效率, 無量綱, 由公式(2)計(jì)算。模型中地表徑流在以下兩種情況下發(fā)生: 第一, 降水量、灌溉量等水分輸入大于土壤最大入滲量; 第二, 土壤表層水分含量超過其飽和含水量時(shí)仍出現(xiàn)降水或灌溉。地表徑流采用逐網(wǎng)格方式匯流,根據(jù)網(wǎng)格高程計(jì)算的其上、下、左、右4個(gè)方向坡度, 求算流向4個(gè)方向的地表徑流權(quán)重。

式中: c是裸地地表徑流的泥沙攜帶率, 該值一般需要用研究區(qū)的徑流小區(qū)試驗(yàn)觀測到的泥沙流失數(shù)據(jù)進(jìn)行率定, 例如, Loch和Donnollan報(bào)道的黏土、黏壤土和砂壤土的泥沙攜帶率分別為8.3%、4.8%和1.6%; c是模型系數(shù), 缺省值一般為0.15, 無量綱;e是自然常數(shù), 其值為2.718。本文用模型植物生長模塊計(jì)算的生育期指數(shù)(f, 無量綱)表征植被對地表的覆蓋面積百分比, 沒有植被覆蓋的地表面積百分比P由公式(3)計(jì)算:

研究表明, 被地表徑流侵蝕的泥沙中碳氮磷含量通常比地表土壤中相應(yīng)元素的含量更高。這種被泥沙侵蝕引起元素富集作用, 通常用富集比(R)表征, 即被侵蝕泥沙中元素濃度和原地表土壤中相應(yīng)元素濃度的比例。顆粒態(tài)碳氮磷隨地表徑流遷移的計(jì)算方法參照McElroy等和Williams和Hann的研究, 計(jì)算公式為:

式中: S是由地表徑流侵蝕引起的顆粒態(tài)碳氮磷的流失量, 單位為g?hm; C是表層土壤中顆粒態(tài)碳氮磷的含量, 單位為g?kg; R是由地表徑流引起的顆粒態(tài)碳氮的富集比, 無量綱。一般情況下, 土壤水力侵蝕量越大, 富集比越小。地表徑流引起的顆粒態(tài)碳氮的富集比(R)由公式(5)計(jì)算, 其中指前因子(k)需要利用研究區(qū)的徑流小區(qū)試驗(yàn)觀測到的地表徑流引起的顆粒態(tài)碳氮數(shù)據(jù)進(jìn)行率定。

地表徑流引起的顆粒態(tài)磷的富集比(R)參照Sharpley由公式(6)計(jì)算:

土壤水力侵蝕模塊計(jì)算的產(chǎn)沙量和元素富集模塊計(jì)算的顆粒態(tài)碳氮磷流失量隨著地表徑流遷移進(jìn)入河道, 在河道中運(yùn)移最終流出所研究的流域。同時(shí), 元素富集模塊計(jì)算的顆粒態(tài)碳氮磷流失量也會(huì)在表層土壤的相應(yīng)組分中扣除。

1.3 數(shù)據(jù)資料準(zhǔn)備

測試CNMM-DNDC模型對流域土壤水力侵蝕和顆粒態(tài)物質(zhì)流失的模擬能力需要準(zhǔn)備兩個(gè)數(shù)據(jù)庫:一個(gè)是用于驅(qū)動(dòng)模型運(yùn)轉(zhuǎn)的輸入數(shù)據(jù)庫, 另一個(gè)是用于模型校準(zhǔn)和驗(yàn)證的觀測數(shù)據(jù)庫。模型輸入庫包括地形、氣象、土壤和管理措施數(shù)據(jù)庫。地形數(shù)據(jù)庫是根據(jù)原始空間分辨率為5 m×5 m的數(shù)字高程圖,在兼顧水文模塊計(jì)算效率和過程描述準(zhǔn)確度的前提下, 將其分辨率重采樣為15 m×15 m, 并以此作為模型運(yùn)行的柵格空間分辨率。氣象數(shù)據(jù)庫采用鹽亭站自動(dòng)氣象站的3 h平均氣溫、風(fēng)速、空氣濕度和累積降水量、短波輻射和長波輻射的觀測值。土壤數(shù)據(jù)包括土壤容重、pH、黏粒含量、有機(jī)質(zhì)含量、總氮含量和水力學(xué)特征(如飽和導(dǎo)水率、飽和含水量、田間持水量和萎蔫點(diǎn)等)等, 是根據(jù)鹽亭站不同土壤類型分層觀測獲取的。管理措施數(shù)據(jù)庫包括從鹽亭站獲取的作物種類、種植制度和播種和收獲時(shí)間、耕作、施肥和灌溉等信息。土地利用數(shù)據(jù)來源于鹽亭站提供的斑塊化土地利用重新分類圖。土壤、管理措施和土地利用數(shù)據(jù)經(jīng)重采樣處理后與地形具有相同的空間分辨率。用于本研究新引進(jìn)的土壤水力侵蝕和碳氮磷元素富集模塊的模型率定和檢驗(yàn)的觀測資料包括2004—2006年觀測的徑流小區(qū)的地表徑流量、水力侵蝕產(chǎn)沙量、由地表徑流引起顆粒態(tài)氮和總氮流失量。其中, 徑流小區(qū)具備極端降水事件的2005—2006年觀測數(shù)據(jù)被用作模型率定, 2004年觀測數(shù)據(jù)被用作模型驗(yàn)證。由于缺乏地表徑流引起的顆粒態(tài)磷流失的觀測數(shù)據(jù), 本研究未涉及顆粒態(tài)磷流失的模型驗(yàn)證。此前, Zhang等已經(jīng)利用鹽亭站多個(gè)典型生態(tài)系統(tǒng)進(jìn)行多年度同步連續(xù)觀測獲取的多變量數(shù)據(jù), 對CNMM-DNDC模擬作物產(chǎn)量、出水口流量、NH、NO、NO和CH排放的性能進(jìn)行了系統(tǒng)全面的測試和檢驗(yàn), 本研究在此基礎(chǔ)上繼續(xù)開展模型優(yōu)化和驗(yàn)證。

1.4 模型評價(jià)指標(biāo)

本文采用相對均方根誤差(RMSE)和線性回歸擬合方程兩種方法, 定量評估改進(jìn)的CNMM-DNDC模型對土壤水力侵蝕以及隨之的顆粒態(tài)氮流失的模擬能力。其中, RMSE是用來評判模擬值和觀測值之間的一致性, RMSE越趨近于0表明模型模擬效果越好; 線性回歸擬合方程的斜率和決定系數(shù)(R)是評價(jià)模擬值與觀測值之間的系統(tǒng)偏差和動(dòng)態(tài)吻合程度的重要指標(biāo), 二者越接近1表明模型模擬效果越好。

2 結(jié)果與分析

2.1 地表徑流、土壤水力侵蝕和氮流失的模型驗(yàn)證評估

2004—2006年地表徑流、土壤水力侵蝕、隨地表徑流的顆粒態(tài)氮和總氮流失的觀測值和模擬值如圖2所示。經(jīng)2005—2006年觀測數(shù)據(jù)的率定, 公式(2)中的c和公式(5)中的k分別取值0.85%和1.2。改進(jìn)的CNMM-DNDC模型能夠較好地模擬截流小流域典型玉米-小麥輪作農(nóng)田徑流小區(qū)的地表徑流季節(jié)動(dòng)態(tài)變化趨勢(圖2a)。在模型驗(yàn)證階段(2004年),地表徑流觀測值范圍為1.3～12.0 mm (均值為6.2 mm),同期模型模擬的地表徑流范圍為1.8～9.6 mm(均值為6.0 mm)。統(tǒng)計(jì)評價(jià)指標(biāo)顯示, 模擬值和觀測值具有顯著線性擬合相關(guān)關(guān)系(R=0.96, P＜0.001,n=7), RMSE為15.2% (表1)。圖2b展示了由地表徑流引起的侵蝕產(chǎn)沙量的觀測值和模擬值的動(dòng)態(tài)變化。在模型驗(yàn)證階段(2004年), 由地表徑流引起的侵蝕泥沙產(chǎn)量觀測值范圍為0.03～0.6 t?hm(均值為0.2 t?hm), 相應(yīng)的模型模擬值范圍0.01～0.6 t?hm(均值為0.2 t?hm)。模型驗(yàn)證階段, 與觀測值相比, 模型能夠很好地捕獲徑流小區(qū)觀測的土壤水力侵蝕產(chǎn)沙季節(jié)變化特征(R=0.83, P＜0.01, n=7, RMSE=32.0%)。盡管如此, 模型對2006年9月4日由一個(gè)極端暴雨事件引起的侵蝕產(chǎn)沙量極端高值低估了約20%。由此可見, 模型在模擬由極端降水事件引起的侵蝕產(chǎn)沙時(shí)仍存在一定的不確定性, 需更多觀測數(shù)據(jù)支撐模型率定和驗(yàn)證工作。改進(jìn)的CNMM-DNDC模型能夠較好地模擬截流小流域典型農(nóng)田徑流小區(qū)的顆粒態(tài)氮和總氮流失的季節(jié)變化規(guī)律(圖2c-d)。在模型驗(yàn)證階段, 觀測到的顆粒態(tài)氮和總氮流失范圍分別為0.03～0.5 kg(N)?hm和0.1～0.6 kg(N)?hm, 均值分別為0.1 kg(N)?hm和0.3 kg(N)?hm; 對應(yīng)的模型模擬值范圍分別為0.02～0.6 kg(N)?hm和0.1～0.8 kg(N)?hm, 均值分別為0.3 kg(N)?hm和0.4 kg(N)?hm。模型驗(yàn)證階段, 顆粒態(tài)氮流失的觀測值與模擬值之間存在顯著線性回歸關(guān)系(R=0.85), RMSE為88.0%。而總氮流失的觀測值和模擬值之間的線性回歸關(guān)系不顯著。

表1 模型率定和驗(yàn)證階段觀測值和模擬值的統(tǒng)計(jì)評價(jià)指標(biāo)Table 1 Model evaluation indices for comparing observations and simulations during the calibration and verification periods

圖2 2004—2006年徑流小區(qū)地表徑流(a)、侵蝕產(chǎn)沙(b)、地表徑流引起的顆粒氮流失(c)和總氮流失(d)模擬結(jié)果Fig.2 Simulated runoff (a), sediment (b), particulate nitrogen loss (c) and total nitrogen loss (d) of the runoff plot from 2004 to 2006

2.2 流域土壤水力侵蝕空間分布模擬

截流小流域由地表徑流引起的侵蝕泥沙年總量空間差異較大(圖3), 2004—2006年不同柵格由地表徑流引起的土壤侵蝕模數(shù)范圍分別為0～11 045 t?km?a(平均值為400 t?km?a)、0～22 618 t?km?a(平均值為701 t?km?a)和0～31 581 t?km?a(平 均值為1550 t?km?a)。其中, 2004年和2005年水力侵蝕強(qiáng)度為微度和輕度的水文響應(yīng)單元總數(shù)目占比分別為99.6%和97.8%, 而2006年該數(shù)值為74.3%。綜合土壤侵蝕模數(shù)均值和侵蝕強(qiáng)度占比數(shù)據(jù), 根據(jù)中華人民共和國水利部發(fā)布的中華人民共和國水利行業(yè)標(biāo)準(zhǔn)《土壤侵蝕分類分級標(biāo)準(zhǔn)》(SL190—2007)規(guī)定的侵蝕強(qiáng)度分級, 截流小流域水力侵蝕強(qiáng)度分級為輕度。

圖3 截流小流域土地利用、坡度和2004—2006年地表徑流和侵蝕泥沙產(chǎn)量空間分布Fig.3 Spatial distributions of land use type, slope and the simulated runoff and sediment yield in the Jieliu catchment from 2004 to 2006

2.3 流域顆粒態(tài)碳氮磷流失空間分布模擬

圖4展示了截流小流域2004—2006年由地表徑流引起的碳氮磷流失量的空間分布。2004—2006年地表徑流引起的顆粒態(tài)碳流失年總量范圍分別為0～835.1 kg(C)?hm?a[平均值為63.9 kg(C)?hm?a]、0～1325.3 kg(C)?hm?a[ 平均值為107.2 kg(C)?hm?a] 和0～1800.4 kg(C)?hm?a[平均值為200.2 kg(C)?hm?a];同期顆粒態(tài)氮流失年總量范圍分別為0～80.8 kg(N)?hm?a[平均值為6.3 kg(N)?hm?a]、0～127.4 kg(N)?hm?a[平均值為10.5 kg(N)?hm?a]和0～172.8 kg(N)?hm?a[平均值為19.5 kg(N)?hm?a]; 相 應(yīng)的顆粒態(tài)磷流失年總量范圍分別為0～11.4 kg(P)?hm?a[平均值為0.9 kg(P)?hm?a]、0～18.0 kg(P)?hm?a[平均值為1.5 kg(P)?hm?a]和0～24.4 kg(P)?hm?a[平均值為2.8 kg(P)?hm?a]。

圖4 2004—2006年截流小流域由地表徑流引起的碳氮磷流失量的空間分布Fig.4 Spatial distributions of the simulated particulate carbon, nitrogen and phosphorus losses in the Jieliu catchment from 2004 to 2006

2.4 流域出水口泥沙、顆粒碳氮磷流出量模擬

模型模擬的2004—2006年截流小流域出水口泥沙年總流出量分別為138 t、97 t和110 t。2004—2006年截流小流域出水口顆粒態(tài)碳的年流出總量分別為1763 kg(C)、1101 kg(C)和1120 kg(C), 顆粒態(tài)氮的年流出總量分別為175 kg(N)、109 kg(N)和110 kg(N), 顆粒態(tài)磷的年流出總量分別為21 kg(P)、13 kg(P)和14 kg(P)。

3 討論

3.1 改進(jìn)CNMM-DNDC模型對土壤水力侵蝕和顆粒態(tài)氮流失模擬的適用性

本研究拓展了CNMM-DNDC模型對土壤水力侵蝕的模擬功能, 并在紫色土農(nóng)林復(fù)合截流小流域進(jìn)行了改進(jìn)模型的適用性驗(yàn)證評估。從驗(yàn)證結(jié)果可以看出, 增加了土壤水力侵蝕和元素富集模塊的CNMM-DNDC模型能夠成功模擬截流小流域的土壤水力侵蝕產(chǎn)沙和隨之的顆粒態(tài)碳氮磷流失, 尤其對小和中度降水事件引起的侵蝕產(chǎn)沙和顆粒態(tài)碳氮磷流失模擬精度較高。但遺憾的是, 盡管改進(jìn)模型成功捕捉2006年9月4日由一個(gè)極端暴雨事件(單日降水量高達(dá)107 mm)引起的地表徑流峰值, 但改進(jìn)模型對這一極端暴雨事件引起的侵蝕產(chǎn)沙極端高值低估了約20%。慶幸的是, 改進(jìn)模型對這一極端暴雨事件引起顆粒態(tài)氮流失僅低估7%, 且模擬值處在測量誤差范圍內(nèi)。由于可獲取觀測數(shù)據(jù)的匱乏, 目前本研究只對紫色土典型玉米-小麥輪作農(nóng)田徑流小區(qū)(地塊尺度)的水力侵蝕產(chǎn)沙量和顆粒態(tài)氮流失進(jìn)行了模擬驗(yàn)證, 未來亟需地塊尺度顆粒態(tài)碳和顆粒態(tài)磷流失的可靠觀測數(shù)據(jù)支撐改進(jìn)模塊的模擬驗(yàn)證工作的進(jìn)一步開展。另外, 亟需流域出水口流量、泥沙和顆粒態(tài)碳氮磷流出量等觀測資料對改進(jìn)模塊的流域尺度驗(yàn)證工作的支持。

3.2 土壤水力侵蝕和顆粒態(tài)氮流失的時(shí)空分布特征

表2以降水強(qiáng)度較大的2006年為例, 展示了截流小流域土壤水力侵蝕、顆粒態(tài)碳氮磷流失強(qiáng)度與土地利用、坡度的關(guān)系。截流小流域由地表徑流引起的侵蝕產(chǎn)沙、顆粒態(tài)碳氮磷流失強(qiáng)度與土地利用密切相關(guān), 與前人在該流域的研究結(jié)果一致。由于水田淹水期養(yǎng)分溢出引起的營養(yǎng)物流失不屬于土壤水力侵蝕范疇, 本研究暫不予討論, 但水田非淹水期仍存在地表徑流、水力侵蝕產(chǎn)沙和顆粒態(tài)碳氮磷損失。2006年截流小流域各種土地利用類型中, 坡耕旱地和居民區(qū)的土壤水力侵蝕強(qiáng)度最高, 分別為2550 t?km?a和 2178 t?km?a(表2)。坡耕旱地由于長期耕作導(dǎo)致表層土壤結(jié)構(gòu)被破壞, 易被雨水沖刷流失。而居民區(qū)主要是由于水泥和建筑物的覆蓋降低了雨水的下滲, 形成更多的地表徑流, 水力侵蝕引起的泥沙也就更多。因此, 坡耕旱地和居民區(qū)由地表徑流引起的顆粒態(tài)碳氮磷流失強(qiáng)度也最高, 成為截流小流域顆粒態(tài)碳氮磷流失的主要來源。而人為干擾較少, 植被覆蓋度較高的林地, 土壤水力侵蝕程度較弱, 2006年產(chǎn)沙量僅為498 t?km?a。這可能是由于林地地表植被的覆蓋對降水有截留作用, 有效降低了雨滴對土壤顆粒剪切破壞, 同時(shí)減小地表徑流流速, 大大降低了地表徑流對土壤顆粒的攜帶作用。另外, 林地植被根系發(fā)達(dá)能夠使土壤有機(jī)質(zhì)含量增加、組成結(jié)構(gòu)改變, 提高土壤孔隙度, 增加土壤水分的入滲量, 增加水穩(wěn)性團(tuán)聚體數(shù)量, 進(jìn)而提高土壤的抗侵蝕能力。而土壤水力侵蝕產(chǎn)沙量僅為坡耕旱地1/5的林地的顆粒態(tài)碳氮磷的流失量卻為坡耕旱地的60%～70%。這可能是由于樹木凋落物導(dǎo)致林地土壤有機(jī)質(zhì)含量高, 可供侵蝕的顆粒態(tài)碳氮磷含量較高。另外, 除土地利用外, 坡度也顯著影響水力侵蝕產(chǎn)沙和顆粒態(tài)碳氮磷流失強(qiáng)度。對2006年林地、坡耕地和居民區(qū)3種土地利用類型分別進(jìn)行的產(chǎn)沙量、顆粒態(tài)碳氮磷流失與坡度的相關(guān)性分析結(jié)果表明, 在同一種土地類型上, 由地表徑流引起的泥沙、顆粒態(tài)碳氮磷的流失強(qiáng)度與坡度正弦值的絕對值均呈顯著正相關(guān)關(guān)系(表2), 水力侵蝕產(chǎn)沙和顆粒態(tài)碳氮磷流失強(qiáng)度隨坡度的增加而增加, 皮爾森相關(guān)系數(shù)均不低于0.89。由于坡度的影響, 在灌溉水田和雙季稻田極陡地形的地塊柵格也出現(xiàn)了相當(dāng)大的水土流失和顆粒態(tài)碳氮磷損失(圖3和圖4)。

表2 2006年截流小流域土壤水力侵蝕產(chǎn)沙、顆粒態(tài)碳氮磷流失強(qiáng)度與土地利用、坡度的關(guān)系Table 2 Relationships between simulated sediment yield, particulate carbon/nitrogen/phosphorus loss and land use type, slope in the Jieliu catchment in 2006

截流小流域產(chǎn)沙量和顆粒態(tài)碳氮磷等養(yǎng)分流失量的觀測值和模型模擬值均呈現(xiàn)旱季(冬季和春季)低而雨季高的特征, 峰值均出現(xiàn)在極端暴雨事件之后。雖然截流小流域2006年年降水量低于2004年和2005年(2004年、2005年和2006年降水量分別為860 mm、835 mm和806 mm), 但2006年出水口泥沙以及碳氮磷流出量與2005年相當(dāng)。這可能與中雨和暴雨事件發(fā)生頻率有關(guān)。2004年、2005年和2006年中雨(24 h降水量大于10 mm)發(fā)生次數(shù)分別為33、21和13, 但中雨發(fā)生次數(shù)較少的2006年的出水口泥沙流出量甚至比2005年還高, 這可能是由于2006年7月3日和9月4日分別發(fā)生一次暴雨(24 h降水量大于50 mm)和特大暴雨(24 h降水量大于100 mm)導(dǎo)致的。大量研究也表明, 在雨強(qiáng)較小時(shí), 碳氮磷等養(yǎng)分主要以可溶態(tài)形式隨徑流遷移; 而在雨強(qiáng)較大時(shí), 碳氮磷等養(yǎng)分則主要以顆粒態(tài)形式隨泥沙遷移。土壤水力侵蝕以及隨之的碳氮磷流失隨雨強(qiáng)的增加而增大。這可能是由于大雨強(qiáng)雨滴在打擊地表時(shí)對土壤結(jié)構(gòu)的剪切破壞作用大, 使土粒更加分散, 從而使得地表徑流攜帶土壤團(tuán)聚體的能力增強(qiáng), 隨之的顆粒態(tài)養(yǎng)分流失增加。另外,大雨強(qiáng)雨滴的擊濺作用使吸附在土壤團(tuán)聚體上的碳氮磷等養(yǎng)分大量解吸, 使得地表徑流中的碳氮磷等養(yǎng)分濃度升高, 地表徑流中泥沙顆粒對碳氮磷養(yǎng)分的富集量增加, 導(dǎo)致顆粒態(tài)碳氮磷的流失量增加。

4 結(jié)論

本研究為分布式流域水碳氮模型CNMM-DNDC增加了土壤水力侵蝕和顆粒態(tài)碳氮磷元素富集模塊,改進(jìn)模型不僅能夠模擬地塊尺度土壤水力侵蝕強(qiáng)度和顆粒態(tài)氮磷流失動(dòng)態(tài), 而且能夠評估整個(gè)流域的水土流失強(qiáng)度。增加了土壤水力侵蝕和顆粒態(tài)碳氮磷元素富集模塊的CNMM-DNDC, 對截流小流域玉米-小麥輪作徑流小區(qū)的泥沙產(chǎn)量和顆粒態(tài)氮流失量的季節(jié)動(dòng)態(tài)表現(xiàn)出良好的模擬性能, 模型驗(yàn)證階段模擬值和觀測值具有顯著線性擬合關(guān)系(R分別為0.83和0.85)。改進(jìn)的CNMM-DNDC模型能夠成功模擬流域土壤水力侵蝕產(chǎn)流產(chǎn)沙強(qiáng)度和隨之的顆粒態(tài)碳氮磷流失強(qiáng)度的空間分布。模擬結(jié)果顯示,2004—2006年截流小流域柵格尺度土壤侵蝕模數(shù)的平均值分別為400 t?hm?a、701 t?hm?a和1550 t?hm?a, 根據(jù)《土壤侵蝕分類分級標(biāo)準(zhǔn)》(SL190—2007)規(guī)定的侵蝕強(qiáng)度分級, 屬于輕度侵蝕區(qū)。2004—2006年由水力侵蝕引起的顆粒態(tài)碳年流失量的平均值分別為63.9 kg(C)?hm?a、107.2 kg(C)?hm?a和200.2 kg(C)?hm?a, 顆粒態(tài)氮年流失量的平均值分別為6.3 kg(N)?hm?a、10.5 kg(N)?hm?a和19.5 kg(N)?hm?a, 顆粒態(tài)磷年流失量平均值分別為0.9 kg(P)?hm?a、1.5 kg(P)?hm?a和2.8 kg(P)?hm?a。截流小流域由地表徑流引起的侵蝕產(chǎn)沙、顆粒態(tài)碳氮磷流失強(qiáng)度與土地利用和地形密切相關(guān)。其中,坡耕旱地和居民區(qū)土壤水力侵蝕最為嚴(yán)重, 是截流小流域顆粒態(tài)碳氮磷流失的主要來源。在同一種土地利用類型上, 土壤水力侵蝕產(chǎn)沙和顆粒態(tài)碳氮磷流失強(qiáng)度隨坡度的增加而增加。本研究拓展了CNMM-DNDC模型對土壤水力侵蝕和面源污染的模擬功能, 在一定程度上為流域水土流失和碳氮磷面源污染的研究提供先進(jìn)的模型工具。