1981-2013年桃江流域徑流與泥沙模擬研究

李志強(qiáng), 齊述華, 劉旗福, 仝興慶, 劉貴花, 鄭海金

(1.江西師范大學(xué) 鄱陽湖濕地與流域研究教育部重點(diǎn)實(shí)驗(yàn)室/地理與環(huán)境學(xué)院, 江西 南昌 330022;2.贛州市水文局, 江西 贛州 341000; 3.江西省水土保持科學(xué)研究院, 江西 南昌 330029)

河道的徑流和泥沙是河流水文特征的最重要參數(shù)，同時(shí)受地形、氣候、地質(zhì)、土壤、植被等多種自然因素和人類活動(dòng)的影響[1]。開展人類活動(dòng)或自然條件變化對(duì)流域水文過程的影響評(píng)價(jià)，對(duì)區(qū)域的生態(tài)工程建設(shè)、維持區(qū)域可持續(xù)發(fā)展提供參考和依據(jù)。

目前，利用水文模型模擬人類活動(dòng)[2-3]、氣候變化[4-5]、植被覆蓋變化[6-7]等因素對(duì)水文過程的影響研究取得了顯著的進(jìn)展。按照水文模型的類型可以劃分為集總式模型和分布式模型[8]，集總式模型利用概念和經(jīng)驗(yàn)關(guān)系來描述水流運(yùn)動(dòng)狀態(tài)，具有所需參數(shù)較少、參數(shù)率定相對(duì)簡(jiǎn)單等優(yōu)勢(shì)，在諸多流域的徑流和泥沙模擬中取得了較好的效果[9-10]。雖然集總式水文模型應(yīng)用廣泛，但是模型的物理基礎(chǔ)較少，只能模擬出水文的宏觀表現(xiàn)，難以表達(dá)水文現(xiàn)象的本質(zhì)和物理機(jī)制。分布式水文模型因具有明確的物理意義，可以準(zhǔn)確的描述水文過程，且具有較強(qiáng)的適應(yīng)性，目前主要有AnnAGNPS[11]，VIC[12]，WATLAC[13]等模型應(yīng)用較為普遍，其中SWAT模型是由美國(guó)農(nóng)業(yè)部開發(fā)的流域分布式水文模型，其主要包括水文過程模型、土壤侵蝕模型和污染負(fù)荷模型3大子模型，目前均已得到較為廣泛的運(yùn)用[14]。

鄱陽湖流域地處東亞季風(fēng)氣候區(qū)，自1980s年代以來，鄱陽湖流域中主要河流泥沙含量呈下降趨勢(shì)[15]，這與1980s年代以來江西省發(fā)起以“植樹造林”為主的“山江湖”生態(tài)恢復(fù)工程顯著改變鄱陽湖流域植被覆蓋特征有關(guān)[16]，同時(shí)也與流域范圍內(nèi)修建水庫的攔沙效應(yīng)有關(guān)。鑒于此，本研究擬以贛江上游的桃江流域?yàn)檠芯繀^(qū)，運(yùn)用SWAT模型對(duì)該流域水文泥沙過程進(jìn)行模擬，分析不同時(shí)期土地利用情景對(duì)徑流與泥沙的影響，并根據(jù)泥沙變化情況評(píng)價(jià)居龍灘水利樞紐工程的攔沙效應(yīng)，以期界定桃江水文特征變化的主要原因。

1 研究區(qū)概況

桃江古稱彭水，又名信豐江，發(fā)源于贛粵交界的九連山脈，自西南向東北流經(jīng)江西省贛州市的全南縣、龍南縣和信豐縣，在贛縣注入貢水，經(jīng)貢水注入贛江，流域介于北緯24°30′—25°55′，東經(jīng)114°10′—115°20′，面積7 864 km2，是贛江的主要源區(qū)。桃江流域?qū)儆谥衼啛釒駶?rùn)季風(fēng)氣候，該區(qū)年均氣溫為19.94 ℃，年均降水量為1 497 mm，降水量年內(nèi)變化較大，一般4—6月為豐水期，占全年降水量的43.8%，11月至翌年1月為枯水期，占全年降水量的10.4%。流域年平均徑流量6.45×109m3,占贛江流域徑流量的9.4%，年平均輸沙量1.35×106t，占贛江流域輸沙量的13.8%。流域內(nèi)有水庫9個(gè)，其中居龍灘水利樞紐工程位于桃江下游，庫容量為7.36×107m3，為桃江流域內(nèi)的最大水庫，是一座以發(fā)電為主，兼有水庫養(yǎng)殖、改善航運(yùn)等綜合效益的中型水利工程，水庫正常蓄水位122.0 m，發(fā)電死水位117.0 m。出庫下行1.4 km為居龍灘水文站，是桃江流域的水文控制站。

桃江流域地勢(shì)南高北低，以山地為主,流域內(nèi)主要土地利用類型為：以馬尾松、杉木松、毛竹松為主的林地(77%)和耕地(17.5%)，特色經(jīng)濟(jì)作物為臍橙、柑橘和板栗等。流域內(nèi)主要土壤類型為紅壤(65.10%)，黃紅壤(16.83%)和潴育水稻土(13.41%)。

2 數(shù)據(jù)與方法

2.1 數(shù)據(jù)收集與預(yù)處理

SWAT模型的運(yùn)行需要針對(duì)研究區(qū)建立各要素?cái)?shù)據(jù)庫，包括數(shù)字高程模型(DEM)、土地利用和土壤等空間數(shù)據(jù)以及土壤屬性數(shù)據(jù)庫。

(1) DEM。由地理空間數(shù)據(jù)云(http:∥www.gscloud.cn/)下載的ASTER GDEM，該數(shù)據(jù)空間分辨率為30 m，利用ENVI軟件進(jìn)行校正并拼接生成研究區(qū)的DEM。

(2) 土地利用數(shù)據(jù)。土地利用矢量數(shù)據(jù)來自于江西省土地利用/覆蓋數(shù)據(jù)庫。該數(shù)據(jù)庫利用秋季獲取的無云Landsat衛(wèi)星遙感影像，參照《環(huán)境狀況評(píng)價(jià)技術(shù)規(guī)范》(HJ/T 192-2006)中土地利用分類方法建立的土地利用數(shù)據(jù)庫。該數(shù)據(jù)庫2005年以來的土地利用數(shù)據(jù)是經(jīng)過內(nèi)業(yè)校核和外業(yè)核查，Ⅰ級(jí)分類的總體精度達(dá)到90%以上，Ⅱ級(jí)分類精度大于85%，2005年之前的數(shù)據(jù)是按照前溯變化監(jiān)測(cè)方式完成土地利用制圖。依照SWAT模型運(yùn)行要求，將土地利用歸并為如表1 所示的耕地、林地、草地、水域、建設(shè)用地和未利用土地6類。

表1 桃江流域1980，2000，2015年各主要土地利用類型面積 km2

(3) 土壤數(shù)據(jù)。包括土壤類型數(shù)據(jù)和土壤屬性數(shù)據(jù)。土壤類型數(shù)據(jù)采用1∶50萬土壤圖，桃江流域主要分布有13類土壤。土壤屬性數(shù)據(jù)是記錄土壤層數(shù)、土壤粒徑組成、土壤濕密度等屬性要素。1∶50萬土壤圖采用國(guó)際制土壤質(zhì)地分類，與ArcSWAT采用的USDA簡(jiǎn)化的美國(guó)制土壤質(zhì)地不同，本文運(yùn)用3次樣條插值法將國(guó)際制土壤質(zhì)地標(biāo)準(zhǔn)轉(zhuǎn)換為適用于ArcSWAT模型支持的土壤質(zhì)地標(biāo)準(zhǔn)[17-18]，并建立土壤屬性數(shù)據(jù)庫。

(4) 氣象數(shù)據(jù)。包括1959—2013年溫度、降水量、相對(duì)濕度、風(fēng)速與日照時(shí)數(shù)等數(shù)據(jù)。ArcSWAT模型的運(yùn)行需要站點(diǎn)日降水量和日最高、最低溫度，其余數(shù)據(jù)利用桃江流域及其周邊共8個(gè)氣象站點(diǎn),由模型天氣發(fā)生器生成。

(5) 水文數(shù)據(jù)。收集了1980—2013年桃江下游的居龍灘水文站觀測(cè)的日徑流和泥沙含量資料。

2.2 方 法

(1) 子流域和水文響應(yīng)單元?jiǎng)澐?。ArcSWAT能夠依據(jù)DEM提取河網(wǎng)并劃分子流域，但由于模型通過默認(rèn)閾值提取的河網(wǎng)精度不高，文中采用Google Earth平臺(tái)矢量化的實(shí)際河網(wǎng)，然后根據(jù)提取的河網(wǎng)劃分子流域。

水文響應(yīng)單元(HRUs)為模型模擬的最小空間單元，依據(jù)土地利用、土壤和坡度等，通過設(shè)定閾值由ArcSWAT生成，將桃江流域劃分為43個(gè)子流域和185個(gè)水文響應(yīng)單元。

(2) 模型運(yùn)行與參數(shù)敏感性分析。由SUFI2，GLUE，ParaSol和MCMC等算法集成的SWAT-CUP是一個(gè)公共程序，用于ArcSWAT模型參數(shù)的敏感性分析，并對(duì)模擬結(jié)果進(jìn)行率定和驗(yàn)證。

按照先徑流、后泥沙的校準(zhǔn)順序調(diào)試運(yùn)行[19]，以ArcSWAT模型給定的參數(shù)范圍為基準(zhǔn)，利用SWAT-CUP的參數(shù)校準(zhǔn)工具和ArcSWAT中的sensitivity analysis工具確定各參數(shù)敏感性，得出表2中與徑流和泥沙有關(guān)的23個(gè)敏感參數(shù)。其中CN2.mgt，GW_REVAP.gw，GWQMN.gw和REVAPMN.gw等參數(shù)對(duì)徑流最為敏感，SPCON.bsn，SPEXP.bsn和USLE_P.mgt等參數(shù)對(duì)泥沙最為敏感。

表2 桃江流域徑流、泥沙敏感參數(shù)

2.3 精度評(píng)價(jià)

采用相對(duì)誤差(Re)，決定系數(shù)(R2)和Nash-Suttcliffe效率系數(shù)(Ens)分別評(píng)價(jià)模擬精度。其中Re評(píng)價(jià)模擬值與實(shí)際值的偏差，Re趨于0說明模擬結(jié)果越好；R2描述模擬值與觀測(cè)值的序列變化趨勢(shì)的一致性，R2趨于1，說明模擬精度高；Ens表示模擬值與實(shí)際值的接近程度：

Ens越接近于1，說明模擬值與實(shí)際值的偏離程度越小。一般來說，Ens>0.5時(shí)可以認(rèn)為模擬結(jié)果有一定的可信度，Ens>0.65時(shí)，模擬結(jié)果較好[20-21]。

3 結(jié)果與討論

3.1 徑流和泥沙模擬

將SWAT模型模擬時(shí)段劃分為預(yù)熱期，率定期和驗(yàn)證期，其中預(yù)熱期用于確定模型初始變量；率定期用來率定模型的敏感參數(shù)；驗(yàn)證期是對(duì)率定期確定的敏感參數(shù)和模擬結(jié)果的評(píng)價(jià)[22]。根據(jù)1980，1990，2000年3期土地利用數(shù)據(jù)，1980—1990年土地利用變化很小，將1980年設(shè)為模擬的預(yù)熱期，1981—1990年為模擬的率定期，1991—2000年為模擬的驗(yàn)證期。

由表3可以看出，經(jīng)過SWAT模型的預(yù)熱和率定，桃江流域在率定期和驗(yàn)證期的徑流和泥沙模擬效果良好。

表3 桃江流域率定期與驗(yàn)證期的徑流、泥沙模擬精度評(píng)價(jià)

圖1給出了1981—2013年的月徑流模擬值與觀測(cè)值變化趨勢(shì)，由表4可以看出模擬效果很好。從分時(shí)間段的模擬評(píng)價(jià)結(jié)果可以看出：2007—2013年的月徑流模擬Ens和R2都表明水電樞紐運(yùn)行沒有顯著改變桃江月徑流。這是由于居龍灘水電樞紐是以發(fā)電為主要功能的河道型水力工程，正常蓄水位122.0 m至起調(diào)水位118.0 m，在汛期遇有相應(yīng)量級(jí)的洪水時(shí)才進(jìn)行預(yù)泄，將庫水位降至118 m；在非汛期，為維持正常蓄水位，上游來水量基本等于水庫泄水量。這樣的水庫運(yùn)行方案，會(huì)影響調(diào)蓄洪水后幾天的徑流，但對(duì)月尺度的徑流影響較小。因此月徑流模擬的結(jié)果受居龍灘水電樞紐工程運(yùn)行的影響較小。

圖1 桃江流域月徑流模擬值與水文觀測(cè)值比較

從1981—2013年泥沙模擬效果的圖2和表4來看，模擬值與觀測(cè)值的變化趨勢(shì)基本一致，受2007年居龍灘水電站的建成運(yùn)行的攔沙影響，Ens只有0.64。以2007年為界，1981—2006年的泥沙模擬效果良好，泥沙模擬的Ens為0.81，模擬值與觀測(cè)值的決定系數(shù)R2為0.81，模擬值與觀測(cè)值達(dá)到顯著相關(guān)水平；而2007—2013年，盡管模擬值與觀測(cè)值變化趨勢(shì)一致決定系數(shù)達(dá)到極顯著，但模擬值明顯高于觀測(cè)值，模擬相對(duì)誤差顯著增加。

圖2 桃江流域月輸沙量模擬值與水文觀測(cè)值比較

評(píng)價(jià)目標(biāo)模擬期ReEnsR21981—20060．170．920．93流量2007—20139．320．910．921981—20131．920．920．921981—20061．460．810．81輸沙量2007—2013345．63-5．080．821981—201321．110．640．69

3.2 土地利用變化對(duì)徑流、泥沙的影響

分別以1980，2000，2015年3期的土地利用數(shù)據(jù)，模擬不同土地利用情景下的月尺度桃江流域的徑流和泥沙。由表5可以看出：不同土地利用情景下模擬的徑流和泥沙變化很小，這主要是由于按照SWAT模型運(yùn)行需要輸入的6大類土地利用類型中1980，2000，2015年這3個(gè)時(shí)期的土地利用變化很小，導(dǎo)致森林或草地的覆蓋度變化對(duì)徑流和泥沙的影響難以得到表達(dá)。盡管如此，仍取得較好的模擬結(jié)果，可以認(rèn)為：1980年以來，桃江流域的森林和草地的植被覆蓋度變化不是導(dǎo)致桃江徑流和泥沙變化的主要原因。

表5 不同土地利用情景下桃江流域年均徑流、泥沙比較

3.3 土壤侵蝕分析

根據(jù)SWAT模型的模擬結(jié)果，對(duì)桃江流域土壤侵蝕做進(jìn)行進(jìn)一步的分析：①參考土壤侵蝕強(qiáng)度分級(jí)標(biāo)準(zhǔn)表(SL190-96)，依據(jù)1980—2013年的模擬結(jié)果計(jì)算各子流域侵蝕模數(shù)，對(duì)各子流域土壤侵蝕強(qiáng)度進(jìn)行分級(jí)，得出桃江流域輕度侵蝕區(qū)域面積為3 877.06 km2，占桃江流域總面積的51.68%，中度侵蝕區(qū)域面積為3 624.90 km2，占桃江流域總面積的48.32%。②表6分析了桃江流域內(nèi)185個(gè)水文響應(yīng)單元不同土地利用方式對(duì)侵蝕產(chǎn)沙的影響，表明不同土地利用類型的土壤侵蝕強(qiáng)度順序?yàn)椋翰莸?耕地>建設(shè)用地>林地。

表6 桃江流域不同土地利用類型侵蝕模數(shù)

3.4 水利工程對(duì)河道輸沙的影響

居龍灘水電站于2007年3月第一臺(tái)機(jī)組正式運(yùn)行，2007年5月全部投入運(yùn)行[23]，導(dǎo)致桃江流域居龍灘水文站觀測(cè)的輸沙量大幅下降。從圖2可以看出，2007年以后，泥沙模擬值普遍大于泥沙觀測(cè)值，這正是由于難以獲得居龍灘水電站的調(diào)蓄水方案中的諸多參數(shù)所引起。結(jié)合水文站觀測(cè)泥沙含量序列和模擬結(jié)果，可以定量評(píng)價(jià)居龍灘水電站每年的攔沙量，結(jié)果如圖3所示：居龍灘水電站的建成與運(yùn)行使桃江輸入貢江的泥沙量顯著降低，2007—2013年期間居龍灘水電站的攔沙量約7.78×106t。

圖3 1981-2013年桃江累積輸沙量

4 結(jié)論與展望

(1) SWAT模型適用于桃江流域的徑流和泥沙模擬，率定期和驗(yàn)證期的模擬精度較高，達(dá)到顯著相關(guān)水平；

(2) 土地利用變化不是引起桃江徑流和泥沙變化的主要原因；

(3) 桃江流域輕度侵蝕區(qū)域面積為3 877.06 km2，占桃江流域總面積的51.68%，中度侵蝕區(qū)域面積為3 624.90 km2，占桃江流域總面積的48.32%；

(4) SWAT模型能夠較好地用于水利樞紐工程攔沙效應(yīng)的定量評(píng)價(jià)，桃江下游的居龍灘水利樞紐工程在2007—2013年總的攔沙量約為7.78×106t，占桃江流域輸沙總量的77.5%，顯著改變了桃江的輸沙特征。

由于1980,2000和2015年3個(gè)年份的森林和草地面積相近，土地利用變化對(duì)徑流和泥沙的情景模擬結(jié)果不能表達(dá)20世紀(jì)80年代以來植被覆蓋度變化對(duì)徑流和泥沙的影響，在未來的研究中，有必要進(jìn)一步強(qiáng)化SWAT模型針對(duì)不同覆蓋度的林、草地的水文特征參數(shù)數(shù)據(jù)庫的建設(shè)和完善；另外，在實(shí)現(xiàn)桃江流域徑流和輸沙模擬的基礎(chǔ)上，有必要進(jìn)一步對(duì)土地利用變化的環(huán)境效應(yīng)進(jìn)行模擬研究，特別是農(nóng)業(yè)生產(chǎn)和果園開發(fā)中大量使用化肥和農(nóng)藥對(duì)水質(zhì)的影響。

[1] 韋杰,賀秀斌.流域侵蝕產(chǎn)沙人類活動(dòng)影響指數(shù)研究:以長(zhǎng)江上游為例[J].地理研究,2012,31 (12):2259-2269.

[2] 呂樂婷,彭秋志,郭媛媛,等.基于SWAT模型的東江流域徑流模擬[J].自然資源學(xué)報(bào),2014,29(10):1746-1757.

[3] 鄭捷,李光永,韓振中,等.改進(jìn)的SWAT模型在平原灌區(qū)的應(yīng)用[J].水利學(xué)報(bào),2011,42(1):88-97.

[4] 郭軍庭,張志強(qiáng),王盛萍,等.應(yīng)用SWAT模型研究潮河流域土地利用和氣候變化對(duì)徑流的影響[J].生態(tài)學(xué)報(bào),2014,34(6):1559-1567.

[5] 丁相毅,周懷東,王宇暉,等.基于分布式水文模型的三峽庫區(qū)污染負(fù)荷對(duì)氣候變化的響應(yīng)研究[J].環(huán)境科學(xué)學(xué)報(bào),2012,32(8):1991-1998.

[6] 謝媛媛.基于SWAT模型的羅玉溝流域森林植被變化的水文生態(tài)響應(yīng)研究[J].水土保持研究,2009,16(4):197-202.

[7] 楊大文,雷慧閩,叢振濤.流域水文過程與植被相互作用研究現(xiàn)狀評(píng)述[J].水利學(xué)報(bào),2010,41(10):1142-1149.

[8] 徐宗學(xué).水文模型[M].北京:科學(xué)出版社,2009.

[9] 劉貴花,齊述華,朱婧瑄,等.氣候變化和人類活動(dòng)對(duì)鄱陽湖流域贛江徑流影響的定量分析[J].湖泊科學(xué),2016,28(3):682-690.

[10] 徐夏楠,高建華,賈建軍,等.氣候變化和人類活動(dòng)對(duì)鄱陽湖流域入湖輸沙量影響的定量估算[J].地理研究,2015,34 (5):838-850.

[11] 田耀武,黃志霖,肖文發(fā).基于AnnAGNPS模型的三峽庫區(qū)秭歸縣非點(diǎn)源污染輸出評(píng)價(jià)[J].生態(tài)學(xué)報(bào),2011,31(16):4568-4578.

[12] 馬俊超,杭慶豐,李瓊芳,等.基于VIC模型的灤河流域綜合干旱指數(shù)的構(gòu)建與應(yīng)用[J].水資源與水工程學(xué)報(bào),2015,26(2):79-84.

[13] 李云良,張奇,李相虎.鄱陽湖流域分布式水文模型的多目標(biāo)參數(shù)率定[J].長(zhǎng)江流域資源與環(huán)境,2013,22(5):565-572.

[14] 鄒松兵,陸志翔,龍愛華,等.ARCSWAT2009用戶指南[M].鄭州:黃河水利出版社,2012.

[15] 齊述華,熊夢(mèng)雅,廖富強(qiáng),等.人類活動(dòng)對(duì)鄱陽湖泥沙收支平衡的影響[J].地理科學(xué),2016,36(6):888-894.

[16] 劉沖,齊述華,湯林玲,等.植被恢復(fù)與氣候變化雙重影響下的鄱陽湖流域蒸散時(shí)空特征研究.地理研究,2016,35(12):2373-2383.

[17] 蔡永明,張科利,李雙才.不同粒徑制間土壤質(zhì)地資料的轉(zhuǎn)換問題研究[J].土壤學(xué)報(bào),2003,40(4):511-517.

[18] 陸穎, Sam BUCHANAN,何大明.流域水文模型中的土壤質(zhì)地轉(zhuǎn)換與飽和導(dǎo)水率Ks值確定[J].云南地理環(huán)境研究,2008,20(5):29-32.

[19] Neitsch S L, Arnold J G, Kiniry J R, et al. Soil and Water Assessment Tool User’s Manual(Version2000)[M]. Grassland, Soil and Water Research Laboratory, Agricultural Research Service, Temple, Texas, USA, 2002.

[20] Moriasi D N, Arnold J G, Van Liew M W, et al. Model evaluation guidelines for systematic quantification of accuracy in watershed simulations [J]. Transactions of the ASABE, 2007,50(3):885-900.

[21] 潘登,任理.分布式水文模型在徒駭馬頰河流域灌溉管理中的應(yīng)用Ⅰ.參數(shù)率定和模擬驗(yàn)證[J].中國(guó)農(nóng)業(yè)科學(xué),2012,45(3):471-479.

[22] 肖軍倉,周文斌,羅定貴,等.非點(diǎn)源污染模型:SWAT用戶應(yīng)用指南[M].北京:地質(zhì)出版社,2010.

[23] 周永忠,許榮貴.居龍灘水電站水力機(jī)械設(shè)計(jì)[J].小水電,2008,39(1):24-27.