穆興民, 杜敏, 邵祎婷, 孫文義, 趙廣舉, 高鵬

(1.中國(guó)科學(xué)院 水土保持研究所,陜西 楊凌 712100; 2.西北農(nóng)林科技大學(xué) 黃土高原土壤侵蝕與旱地農(nóng)業(yè)國(guó)家重點(diǎn)實(shí)驗(yàn)室,陜西 楊凌 712100; 3.黃河水利委員會(huì) 黃河上中游管理局,陜西 西安 710021; 4.浙江水利水電學(xué)院水利與環(huán)境工程學(xué)院,浙江 杭州 310018)

黃土高塬溝壑區(qū)作為黃土高原的主要土壤侵蝕區(qū)和入黃泥沙的重要來(lái)源區(qū),是黃河流域生態(tài)保護(hù)和高質(zhì)量發(fā)展的重點(diǎn)關(guān)注區(qū)域。黃土塬作為黃土高塬溝壑區(qū)特有的地貌類(lèi)型,受極端暴雨和不合理的人為活動(dòng)影響較大,黃土高塬溝壑區(qū)劇烈的土壤侵蝕導(dǎo)致塬面破碎、溝壑密度增加、塬面面積銳減,嚴(yán)重制約了區(qū)域經(jīng)濟(jì)社會(huì)的發(fā)展。以“固溝保塬”為重點(diǎn)的水土流失治理成效顯著,但黃土高塬溝壑區(qū)依然侵蝕嚴(yán)重[1]。因此,認(rèn)識(shí)黃土高塬溝壑區(qū)的土壤侵蝕產(chǎn)沙格局演變特征對(duì)促進(jìn)黃土高塬溝壑區(qū)水土保持的高質(zhì)量發(fā)展具有重要意義。

土壤侵蝕監(jiān)測(cè)方法主要包括經(jīng)驗(yàn)估算法、淤地壩淤積反演法、徑流小區(qū)法、對(duì)比流域觀測(cè)法、土壤侵蝕模型模擬法以及水文站監(jiān)測(cè)資料推演法等[2-5]。經(jīng)驗(yàn)估算法是基于研究者的既有經(jīng)驗(yàn)知識(shí)并結(jié)合地形地貌及氣候特征對(duì)土壤侵蝕狀況進(jìn)行估算的方法,多用于觀測(cè)資料有限的區(qū)域。淤地壩淤積反演法耗時(shí)且成本較高,常用于小尺度流域的侵蝕反演[6]。徑流小區(qū)法研究坡面尺度土壤侵蝕特征,對(duì)小區(qū)均質(zhì)性要求較高,無(wú)法反映復(fù)雜地形條件下多種類(lèi)型侵蝕同時(shí)發(fā)生的區(qū)域性土壤侵蝕狀況[7-8]。受限于時(shí)間和空間上的低相似性和高異質(zhì)性,對(duì)比流域法多應(yīng)用于小流域,難以用于較大尺度區(qū)域[9]。土壤侵蝕模型模擬法及水文站監(jiān)測(cè)資料推演法是大尺度流域及區(qū)域尺度的常用方法,但模型模擬所需參數(shù)較多且部分參數(shù)難以獲取。相較之下,水沙資料監(jiān)測(cè)規(guī)范、時(shí)間序列長(zhǎng),故監(jiān)測(cè)資料推演法可謂大尺度區(qū)域侵蝕產(chǎn)沙研究較為有效的手段。此外,黃土高塬溝壑區(qū)塬面平坦開(kāi)闊,塬邊溝壁陡峭,雖然溝壑內(nèi)建淤地壩困難,但這種地貌可有效避免壩庫(kù)對(duì)流域把口站水沙監(jiān)測(cè)數(shù)據(jù)的影響,數(shù)據(jù)精度更高。鑒于此,本研究選取監(jiān)測(cè)資料推演法探究黃土高塬溝壑區(qū)水沙產(chǎn)輸特征。已有研究表明,在黃土高塬溝壑區(qū),水土保持能使小流域產(chǎn)洪次數(shù)減少,使地表徑流模數(shù)和徑流系數(shù)減小[10-11],是造成年尺度和次降雨尺度上水沙關(guān)系變化顯著的原因[12]。目前對(duì)黃土高塬溝壑區(qū)長(zhǎng)時(shí)間序列水沙關(guān)系演變特征和歸因的研究還較少。

董志塬位于甘肅省慶陽(yáng)市,是黃土高原最為典型、面積最大的黃土塬。硯瓦川流域位于董志塬的腹地,兼具黃土高塬溝壑區(qū)的塬、坡、溝主要地貌類(lèi)型。本研究選取硯瓦川流域?yàn)檠芯繀^(qū),獲取流域內(nèi)1976—2019年的降雨、徑流、泥沙數(shù)據(jù),探究區(qū)域產(chǎn)流產(chǎn)沙特征,旨在為區(qū)域侵蝕防治提供數(shù)據(jù)支撐。

1 資料與方法

1.1 硯瓦川流域概況

硯瓦川流域(35°31′N(xiāo)～35°44′N(xiāo),107°37′E～107°55′E)屬馬蓮河一級(jí)支流、涇河二級(jí)支流、渭河三級(jí)支流,發(fā)源于甘肅省慶陽(yáng)市西峰區(qū)董志塬腹地,在甘肅省寧縣匯入馬蓮河,干流全長(zhǎng)75 km,流域面積為376.69 km2,主要支流有慶豐溝、白玉溝、小崆峒溝等(圖1)。流域年均氣溫為8.1 ℃,年均降雨量為520 mm,植被屬暖溫帶森林草原地帶植被類(lèi)型,天然草原植被以鐵桿蒿(Artemisia gmelinii)以及白羊草(Artemisia gmelinii)和本氏針茅(Stipa capillata)混交群落為主。流域地質(zhì)構(gòu)造單一,除下游的溝谷底部砂巖外,其余被第四紀(jì)黃土所覆蓋。流域內(nèi)水土保持措施主要包括梯田、淤地壩和林草措施,2010年之前淤成壩地僅12.3 hm2,溝道工程甚少[13]。

圖1 硯瓦川流域地理位置及水文站、氣象站分布

1.2 資料來(lái)源

本研究采用西峰氣象站的日監(jiān)測(cè)資料(http://data.cma.cn/),同時(shí)搜集整理了硯瓦川流域內(nèi)硯瓦川、永豐、張鐵、元馬、朱寨、董志、胡同趙、吳家莊等雨量站的日降雨資料。徑流和泥沙數(shù)據(jù)來(lái)源于西峰水土保持科學(xué)試驗(yàn)站硯瓦川水文站,其控制面積為329 km2(2014年后修正為333.78 km2)。隨著甘肅省慶陽(yáng)市的城市快速發(fā)展,硯瓦川流域自2009年以來(lái)承納了城市排放的污水。鑒于慶陽(yáng)市污水處理廠于2008年12月正式運(yùn)行,根據(jù)其污水排放量,將硯瓦川流域2009—2019年的實(shí)測(cè)徑流量減去污水排放量,進(jìn)行了徑流還原處理。

1.3 研究方法

在無(wú)溝道壩庫(kù)攔蓄條件下,黃土高原中小流域侵蝕產(chǎn)沙量和河流輸沙量之比(即泥沙輸移比)接近于1[14-15],即河流水文站所測(cè)的輸沙量約等于或稍小于流域侵蝕產(chǎn)沙量。硯瓦川流域內(nèi)壩庫(kù)工程甚少,水文站實(shí)測(cè)河流輸沙量近似為流域產(chǎn)沙量。

研究中采用年產(chǎn)流模數(shù)及年產(chǎn)沙模數(shù)來(lái)表征硯瓦川流域的產(chǎn)水及產(chǎn)沙強(qiáng)度,以反映流域產(chǎn)水及侵蝕產(chǎn)沙程度。年產(chǎn)流模數(shù)指某年某流域水文站實(shí)測(cè)徑流量與該水文站所控制的流域面積之比。年產(chǎn)沙模數(shù)指某年某流域水文站實(shí)測(cè)輸沙量與該水文站所控制的流域面積之比。

采用年產(chǎn)流系數(shù)及產(chǎn)沙系數(shù)以綜合反映流域降雨的產(chǎn)流能力及侵蝕產(chǎn)沙能力。年產(chǎn)流系數(shù)指某年單位降雨條件下單位面積產(chǎn)流量,即年徑流深與年降雨量之比。產(chǎn)沙系數(shù)是指某年單位降雨量條件下單位面積土地上所產(chǎn)生的泥沙量,采用年輸沙模數(shù)與年降雨量之比求得。

雙累積曲線法由SEARCH J K和HARDISON C H[16]提出,現(xiàn)廣泛用于長(zhǎng)時(shí)間序列的水文氣象數(shù)據(jù)一致性和準(zhǔn)確性的檢驗(yàn),尤其在估算降雨和人類(lèi)活動(dòng)對(duì)徑流量和輸沙量的影響時(shí),是檢驗(yàn)兩個(gè)參數(shù)之間關(guān)系一致性和變差的常用方法。采用雙累積曲線法可識(shí)別流域徑流和輸沙在人類(lèi)活動(dòng)影響下發(fā)生突變的時(shí)間年份,同時(shí)據(jù)此來(lái)劃分研究時(shí)段并量化分析降雨和人類(lèi)活動(dòng)對(duì)徑流和輸沙變化影響的貢獻(xiàn)率[17-18]。對(duì)于給定的降雨時(shí)間序列Pi(i=1、2、…、N)和徑流量(輸沙量)時(shí)間序列Ri(i=1、2、…、N),時(shí)間序列的突變點(diǎn)為t(1

2 結(jié)果與分析

2.1 高塬溝壑區(qū)產(chǎn)流模數(shù)和產(chǎn)沙模數(shù)的階段性

黃土高塬溝壑區(qū)硯瓦川流域的年產(chǎn)流模數(shù)和年產(chǎn)沙模數(shù)的雙累積曲線檢驗(yàn)結(jié)果如圖2所示。

圖2 硯瓦川流域產(chǎn)流模數(shù)及產(chǎn)沙模數(shù)與年降雨量雙累積曲線

由圖2可知,黃土高塬溝壑區(qū)硯瓦川流域的產(chǎn)流模數(shù)及產(chǎn)沙模數(shù)均可劃分為兩個(gè)階段,均在1997年發(fā)生突變。由于硯瓦川流域在1994年黃土高原水土保持世界銀行貸款馬蓮河流域治理項(xiàng)目的實(shí)施下得到全面治理,植被面積迅速增加,梯田、谷坊、淤地壩等的數(shù)量增多,而水土保持措施的效益發(fā)揮具有一定的滯后性。因此,結(jié)合硯瓦川流域的治理過(guò)程,可將流域產(chǎn)流及產(chǎn)沙模數(shù)以1997年為界劃分為兩個(gè)階段,即1976—1996年的突變前時(shí)期及1997—2019年的突變后時(shí)期。

2.2 高塬溝壑區(qū)產(chǎn)流模數(shù)及產(chǎn)沙模數(shù)變化特征

黃土高塬溝壑區(qū)的硯瓦川流域不同年代產(chǎn)流模數(shù)及產(chǎn)沙模數(shù)的差異較大,詳見(jiàn)表1和表2。多年平均產(chǎn)流模數(shù)和產(chǎn)沙模數(shù)分別為22 991 m3/km2和1 889 t/km2。不同年代的產(chǎn)流模數(shù)及產(chǎn)沙模數(shù)均以1980—1989年最大、2000—2019年最小。21世紀(jì)黃土高塬溝壑區(qū)整體侵蝕強(qiáng)度處于中度以下水平。

表1 硯瓦川流域不同年代產(chǎn)流模數(shù)基本特征

表2 硯瓦川流域不同年代產(chǎn)沙模數(shù)基本特征

黃土高塬溝壑區(qū)硯瓦川流域年產(chǎn)沙模數(shù)與產(chǎn)流模數(shù)的年際變化過(guò)程和趨勢(shì)基本一致,如圖3所示。1976—2019年硯瓦川流域年產(chǎn)流模數(shù)和產(chǎn)沙模數(shù)的相關(guān)性呈極顯著水平(P<0.01),單相關(guān)系數(shù)達(dá)0.845,說(shuō)明高塬溝壑區(qū)侵蝕產(chǎn)沙強(qiáng)度隨產(chǎn)流強(qiáng)度的增減而基本同步增減。

圖3 硯瓦川流域產(chǎn)流模數(shù)和產(chǎn)沙模數(shù)的年際變化過(guò)程

硯瓦川流域年產(chǎn)流模數(shù)及年產(chǎn)沙模數(shù)的波動(dòng)都較為劇烈,尤以產(chǎn)沙模數(shù)更明顯,表明流域產(chǎn)沙較產(chǎn)流對(duì)降雨變化更為敏感(圖3及表2)。硯瓦川流域年產(chǎn)流模數(shù)最大值為47 933 m3/km2(1984年),最小值為11 510 m3/km2(2009年),最大值與最小值之比為4.16;年產(chǎn)沙模數(shù)最大值為12 030 t/km2(1984年),最小值為51 t/km2(2009年),最大值與最小值之比達(dá)234.96。產(chǎn)流模數(shù)及產(chǎn)沙模數(shù)的年際變差系數(shù)分別為0.36和1.51,表明產(chǎn)沙模數(shù)的年際變幅較產(chǎn)流模數(shù)的更為劇烈。

突變年后,年均產(chǎn)流模數(shù)及產(chǎn)沙模數(shù)顯著降低。在高塬溝壑區(qū)硯瓦川流域,突變年前的平均產(chǎn)流模數(shù)由25 646 m3/km2(1976—1996年)減少到突變后的20 566 m3/km2(1997—2019年),而年均產(chǎn)沙模數(shù)由3 051 t/km2減少到829 t/km2,減少幅度則分別達(dá)19.8%和72.8%,表明高塬溝壑區(qū)流域產(chǎn)流模數(shù)顯著減小,產(chǎn)沙模數(shù)更是大幅降低。

產(chǎn)沙模數(shù)與產(chǎn)流模數(shù)關(guān)系的變化一定程度上反映了流域徑流挾沙力的變化。高塬溝壑區(qū)的硯瓦川流域產(chǎn)沙模數(shù)與產(chǎn)流模數(shù)增減變化過(guò)程的一致性很高,進(jìn)一步分析得二者呈冪函數(shù)關(guān)系。不同階段硯瓦川流域產(chǎn)沙模數(shù)(Y)與產(chǎn)流模數(shù)(X)間的冪函數(shù)及其參數(shù)變化如圖4和表3所示。經(jīng)統(tǒng)計(jì)檢驗(yàn),擬合方程均達(dá)到0.001信度水平。

表3 不同內(nèi)生階段硯瓦川流域產(chǎn)沙模數(shù)與產(chǎn)流模數(shù)的冪函數(shù)參數(shù)

圖4 硯瓦川流域產(chǎn)沙模數(shù)與產(chǎn)流模數(shù)關(guān)系擬合曲線

采用表3的擬合公式,設(shè)定硯瓦川流域的年產(chǎn)流模數(shù)從10 000 m3/km2增至50 000 m3/km2時(shí),突變年前理論產(chǎn)沙模數(shù)從149 t/km2增至14 582 t/km2,而對(duì)應(yīng)突變年后的理論產(chǎn)沙模數(shù)從45 t/km2增加到4 840 t/km2,相應(yīng)減沙從70.0%變?yōu)?4.0%。即在相同的產(chǎn)流模數(shù)下,突變年份后相應(yīng)的產(chǎn)沙模數(shù)顯著減少,表明隨著硯瓦川流域水土保持高質(zhì)量規(guī)?；卫?其地表產(chǎn)流的挾沙力顯著降低。

2.3 高塬溝壑區(qū)產(chǎn)流系數(shù)和產(chǎn)沙系數(shù)變化

黃土高塬溝壑區(qū)的硯瓦川流域年產(chǎn)流系數(shù)及年產(chǎn)沙系數(shù)的年際變化過(guò)程如圖5所示。

圖5 硯瓦川流域年均產(chǎn)流系數(shù)及產(chǎn)沙系數(shù)年際變化過(guò)程

1976—2019年硯瓦川流域年均產(chǎn)流系數(shù)為0.044,最大為0.076(1984年),最小為0.024(2009年)。突變年之前年均產(chǎn)流系數(shù)為0.050(1976—1996年),之后為0.038(1997—2019年),表明突變年之后產(chǎn)流能力降低了24.0%。硯瓦川流域年均產(chǎn)沙系數(shù)為3.424 t/(km2·mm),最大為19.202 t/(km2·mm) (1984年),最小為0.106 t/(km2·mm)(2009年)。突變年后年均產(chǎn)沙系數(shù)為1.488 t/(km2·mm),較突變年前的5.544 t/(km2·mm)平均降低73.0%,流域產(chǎn)沙能力大幅度下降。產(chǎn)流系數(shù)與產(chǎn)沙系數(shù)的前后變化,表明硯瓦川流域水土保持措施的減沙作用遠(yuǎn)遠(yuǎn)大于減水作用。

2.4 高塬溝壑區(qū)河流含沙量年際變化

含沙量是河流徑流與輸沙量的綜合反映。硯瓦川流域的含沙量變化經(jīng)歷了3個(gè)時(shí)期:1998年之前,高含沙量,大振幅波動(dòng);1998—2010年,低值區(qū);2010年之后,有所反彈但增幅有限(圖6)。硯瓦川流域1976—2019年平均含沙量為64.5 kg/m3,變化于2009年的4.5 kg/m3和1984年的250.5 kg/m3之間。突變年之前平均含沙量為98.3 kg/m3(1976—1996年),突變年之后為33.6 kg/m3(1997—2019年)特別是2001—2009年(2006年除外)均不足10 kg/m3,突變年份前較突變年份后的平均含沙量降低了67%。從含沙量的年際過(guò)程看,硯瓦川流域近年受暴雨影響仍有反彈至均值的態(tài)勢(shì),某些年份如2006年則出現(xiàn)高含沙水流[19],說(shuō)明仍需推進(jìn)水土保持高質(zhì)量發(fā)展;區(qū)域高強(qiáng)度、規(guī)?；帘３种卫硪扬@著降低了區(qū)域土壤侵蝕能力;因?yàn)辄S土之質(zhì)地粗、松散等性質(zhì),以及黃土高原大量存在立地條件差的“溝谷坡”地,黃河及其支流的“河清”如何界定有待商榷。

圖6 硯瓦川流域含沙量年際變化過(guò)程

2.5 高塬溝壑區(qū)侵蝕產(chǎn)沙與降雨的關(guān)系及其演變

硯瓦川流域年降水量呈遞增趨勢(shì),如圖7所示,以2.44 mm/年的速率增加,流域多年平均降雨量為526.58 mm,最大值為817.3 mm(2003年),最小值為288.1 mm(1995年),極值比和變差系數(shù)Cv分別為1.84和0.20。

圖7 1976—2017年硯瓦川流域降水量年際變化過(guò)程

年內(nèi)不同日數(shù)的最大產(chǎn)沙量占年總產(chǎn)沙量的比例反映區(qū)域侵蝕產(chǎn)沙的季節(jié)分配特征。硯瓦川流域年內(nèi)最大1 d、最大3 d、最大7 d及最大10 d產(chǎn)沙量占年產(chǎn)沙量的比例變化如圖8及表4所示。平均最大3 d產(chǎn)沙量約占年產(chǎn)沙量的一半,且突變年份前后差異不大。而突變年前最大7 d及最大10 d的侵蝕產(chǎn)沙量占比較突變年后有增加趨勢(shì),表明年侵蝕產(chǎn)沙的集中度有所降低。就最大3 d產(chǎn)沙量而言,1976—2019年的44年中有44%的年份其最大3 d產(chǎn)沙量占年產(chǎn)沙量50%以上,24%的年份占比超過(guò)80%。最大3 d產(chǎn)沙量占年產(chǎn)沙量比例年際間差異顯著,最高的2001年為93%,最低的1984年為2%。高塬溝壑區(qū)是典型的水力侵蝕區(qū),年內(nèi)降雨集中度高,土壤侵蝕產(chǎn)沙主要因3～5場(chǎng)暴雨所引起,導(dǎo)致年內(nèi)侵蝕產(chǎn)沙的集中度極高。

表4 最大日產(chǎn)沙量占年產(chǎn)沙量的比例

圖8 不同最大日產(chǎn)沙量占年輸沙量的比例

硯瓦川流域產(chǎn)流模數(shù)隨年降雨量的增加而增加,呈顯著的線性關(guān)系。對(duì)1976—1996年及1997—2019年兩個(gè)時(shí)期內(nèi)產(chǎn)流模數(shù)與年降雨量以及日降雨量分別≥10、12、25、50 mm的年累計(jì)降雨量的相關(guān)關(guān)系進(jìn)行分析,結(jié)果表明:產(chǎn)流模數(shù)與年降雨量的直線回歸關(guān)系最優(yōu),兩個(gè)時(shí)段相關(guān)系數(shù)分別為0.671和0.734,達(dá)0.001信度的極顯著水平。線性回歸系數(shù)由54.07降至43.03,表明單位降雨量下的產(chǎn)流強(qiáng)度降低約21%。

硯瓦川流域產(chǎn)沙模數(shù)與降雨量的關(guān)系在不同時(shí)期表現(xiàn)不同。1976—1996年產(chǎn)沙模數(shù)與年降雨量的簡(jiǎn)單相關(guān)關(guān)系最好,相關(guān)系數(shù)為0.439,達(dá)顯著水平。1997—2019年產(chǎn)沙模數(shù)與日降雨量≥50 mm的年累計(jì)降雨量呈顯著指數(shù)函數(shù)關(guān)系,相關(guān)系數(shù)為0.676,亦達(dá)顯著水平。與各等級(jí)降雨量的關(guān)系變化表明,引起流域侵蝕產(chǎn)沙的臨界降雨強(qiáng)度顯著增大。高塬溝壑區(qū)是典型的水力侵蝕區(qū),引起土壤侵蝕產(chǎn)沙的因素可分為降雨(量、強(qiáng)度等)及地表覆蓋和微地形改變兩大類(lèi)。侵蝕產(chǎn)沙與降雨的關(guān)系變化表明,流域土地利用與覆被和微地形改變引起的侵蝕阻力增大,侵蝕產(chǎn)沙能力減弱。水沙變化及其與降雨關(guān)系的變化表明,黃土高原的土壤侵蝕防治應(yīng)從傳統(tǒng)的侵蝕性降雨災(zāi)害防治轉(zhuǎn)變?yōu)閷?duì)強(qiáng)暴雨侵蝕的防治。

3 結(jié)論

以高塬溝壑區(qū)典型支流硯瓦川為研究對(duì)象,通過(guò)水文站實(shí)測(cè)資料,采用產(chǎn)流(產(chǎn)沙)模數(shù)及產(chǎn)流(產(chǎn)沙)系數(shù)指標(biāo),剖析了黃土高塬溝壑區(qū)流域產(chǎn)流及侵蝕產(chǎn)沙變化過(guò)程及其特征,主要結(jié)論如下:

硯瓦川流域1976—2019年產(chǎn)流模數(shù)和產(chǎn)沙模數(shù)的年際變化波動(dòng)趨勢(shì)基本一致,侵蝕產(chǎn)沙強(qiáng)度隨產(chǎn)流強(qiáng)度的增減而同步變化,相關(guān)性高達(dá)0.845(P<0.01)。硯瓦川流域年產(chǎn)沙模數(shù)的變化幅度較產(chǎn)流模數(shù)更大,最大年產(chǎn)流模數(shù)約為最小值的4倍,最大年產(chǎn)沙模數(shù)約為最小值的235倍。

參照大規(guī)模水土保持措施實(shí)施年份,將整個(gè)研究時(shí)段劃分為治理前(1976—1996年)與治理后(1997—2019年)兩個(gè)時(shí)段。相比治理前,治理后硯瓦川流域平均年產(chǎn)流模數(shù)及產(chǎn)流強(qiáng)度分別減少23.9%、75.8%,年均產(chǎn)沙模數(shù)與產(chǎn)沙強(qiáng)度銳減,減小量高達(dá)72.8%,流域產(chǎn)流能力降低24.0%,產(chǎn)沙能力則大幅降低,降幅達(dá)73.0%。

降雨是影響硯瓦川流域產(chǎn)流產(chǎn)沙能力的重要因素。突變前后產(chǎn)流模數(shù)均與年降雨量的線性回歸關(guān)系最優(yōu);突變年之前產(chǎn)沙模數(shù)與年降雨量的相關(guān)關(guān)系最好,而突變年之后與日降雨量≥50 mm的年累計(jì)降雨量的相關(guān)關(guān)系最好。因此,在高塬溝壑區(qū)土壤侵蝕得到規(guī)?；卫砗?預(yù)防暴雨侵蝕災(zāi)害應(yīng)成為重點(diǎn)。