郝桂珍，宋鳳芝，徐 利，朱貴有

(1.河北建筑工程學(xué)院，河北 張家口 075000；2.張家口市水務(wù)局，河北 張家口 075000)

近年來(lái)，隨著流域污染防治工作的不斷深入，農(nóng)業(yè)非點(diǎn)源已成為國(guó)內(nèi)外水環(huán)境的首要污染源。從全球范圍來(lái)看，由污染導(dǎo)致退化的12 億hm2耕地中，農(nóng)業(yè)非點(diǎn)源污染導(dǎo)致的耕地退化約占12%[1]；我國(guó)《第一次全國(guó)污染源普查公報(bào)》顯示，農(nóng)業(yè)源排放的TN、TP分別占其排放總量的57.19%和67.27%。國(guó)內(nèi)外的實(shí)踐研究表明，模型能夠模擬流域復(fù)雜的污染過(guò)程、分析非點(diǎn)源污染負(fù)荷的時(shí)空分布特征，是量化和管控非點(diǎn)源污染最直接有效的方法[2]，因此自國(guó)外引進(jìn)以來(lái)便在國(guó)內(nèi)迅速得到推行。相關(guān)研究表明[3,10]，SWAT模型在我國(guó)北方和南方流域都具有良好的適用性，例如王中根等[3]對(duì)海河流域的研究是SWAT模型在我國(guó)大型流域的成功應(yīng)用；李爽等[4]利用SWAT模型對(duì)南四湖流域的非點(diǎn)源氮、磷污染進(jìn)行了定量模擬和時(shí)空分布特征分析；秦云[5]利用SWAT模型模擬及預(yù)測(cè)了梁子湖流域的非點(diǎn)源污染，并基于模型分析提出了非點(diǎn)源污染控制措施。

本文以清水河流域?yàn)檠芯繉?duì)象，構(gòu)建SWAT模型對(duì)流域內(nèi)的非點(diǎn)源污染進(jìn)行模擬?；谀Ｐ湍M結(jié)果統(tǒng)計(jì)流域的非點(diǎn)源總氮(TN)、總磷(TP)污染負(fù)荷，并從時(shí)間、空間和降水等方面分析非點(diǎn)源污染特征，以掌握清水河流域的非點(diǎn)源污染狀況，為流域的污染控制和管理提供科學(xué)依據(jù)。

1 研究區(qū)概況

1.1 清水河流域概況

清水河流域位于冬奧舉辦城市張家口。2017年1月，習(xí)總書(shū)記視察河北省張家口市時(shí)做出重要指示，強(qiáng)調(diào)張家口要加強(qiáng)生態(tài)建設(shè)，樹(shù)立生態(tài)優(yōu)先意識(shí)，建成首都水源涵養(yǎng)功能區(qū)和生態(tài)環(huán)境支撐區(qū)，這是深入實(shí)施京津冀協(xié)同發(fā)展戰(zhàn)略的重大舉措。冬奧會(huì)的舉辦，更是代表中國(guó)向世界展示山清水秀、生機(jī)盎然、欣欣向榮的美麗中國(guó)形象，因此該流域的環(huán)境形象是向世界展示中國(guó)走綠色可持續(xù)發(fā)展道路的樣板。

該流域內(nèi)的清水河是海河流域永定河水系洋河的主要支流。如圖1所示，清水河發(fā)源于張家口市崇禮區(qū)樺皮嶺南麓，流經(jīng)清三營(yíng)、獅子溝、白旗、西灣子、高家營(yíng)等五個(gè)鄉(xiāng)鎮(zhèn)，至張家口市匯入洋河，全長(zhǎng)109 km，上游分為東溝、正溝和西溝三大支流。

1.2 清水河水質(zhì)現(xiàn)狀

為保證冬奧期間的用水安全，張家口市已對(duì)污染企業(yè)進(jìn)行關(guān)閉或停產(chǎn)整治[7]，點(diǎn)源污染得到有效控制。因此，為進(jìn)一步提高水環(huán)境質(zhì)量，有必要對(duì)流域內(nèi)的非點(diǎn)源污染進(jìn)行分析研究。

2 流域非點(diǎn)源污染模型的應(yīng)用

2.1 SWAT模型概述

SWAT(Soil and Water Assessment Tool)模型是20世紀(jì)90年代中后期Jeff Arnold博士為美國(guó)農(nóng)業(yè)部農(nóng)業(yè)研究局(United States Department of Agriculture- Agricultural Research Service，USDA-ARS)開(kāi)發(fā)的流域尺度模型[10]。該模型是一種基于地理信息技術(shù)(Geographic Information System，GIS)的分布式流域水文模型，由水文、土壤侵蝕與泥沙運(yùn)輸、營(yíng)養(yǎng)物遷移轉(zhuǎn)化三個(gè)模塊構(gòu)成，不僅可以進(jìn)行長(zhǎng)期連續(xù)的模擬，而且計(jì)算效率較高，被廣泛應(yīng)用于非點(diǎn)源污染的預(yù)防和治理。因此，本次應(yīng)用該模型對(duì)清水河流域農(nóng)業(yè)非點(diǎn)源污染進(jìn)行模擬。

2.2 SWAT模型構(gòu)建

借助Arc GIS平臺(tái)構(gòu)建清水河流域SWAT非點(diǎn)源污染模型。SWAT模型作為一個(gè)分布式水文模型，在應(yīng)用于流域非點(diǎn)源污染模擬及控制時(shí)，所需數(shù)據(jù)包括數(shù)字高程模型(Digital Elevation Model，DEM)、土地利用數(shù)據(jù)、土壤數(shù)據(jù)、水文氣象數(shù)據(jù)和農(nóng)業(yè)管理數(shù)據(jù)等，將這些基礎(chǔ)數(shù)據(jù)轉(zhuǎn)換為模型所需的格式，建立SWAT模型基礎(chǔ)數(shù)據(jù)庫(kù)，數(shù)據(jù)來(lái)源及說(shuō)明見(jiàn)表1。

表1 數(shù)據(jù)來(lái)源及說(shuō)明Tab.1 Data source and description

2.2.1 劃分子流域

高精度水系流域的生成有助于提高模型的模擬精度，DEM數(shù)據(jù)是SWAT模型用來(lái)自動(dòng)提取流域邊界基礎(chǔ)，子流域劃分是模型運(yùn)行的重要組成部分?？紤]到清水河流域的實(shí)際情況，以張家口市水文站為流域總出口，將河道閾值設(shè)定為4 000 hm2，通過(guò)模型調(diào)用Arc GIS的水文模塊完成子流域劃分，共生成集水面積1 987.03 km2，并劃分為23個(gè)子流域，見(jiàn)圖3(a)。

2.2.2 水文響應(yīng)單元的劃分

子流域劃分完成后，需進(jìn)一步在每個(gè)子流域內(nèi)劃分水文響應(yīng)單元(Hydrologic Response Unit，HRU)，包括土地利用、土壤、地形數(shù)據(jù)加載及重分類(lèi)和HRU的定義。

例如：教師在開(kāi)展朗讀教學(xué)前，可以先將文章朗讀一遍，然后標(biāo)準(zhǔn)每一個(gè)段落需要采取的情感、語(yǔ)氣等，同時(shí)找出學(xué)生可能出現(xiàn)的問(wèn)題，對(duì)其進(jìn)行標(biāo)記，在開(kāi)展課堂教學(xué)過(guò)程中，也能夠更好地指正。

如圖3(b)所示，SWAT模型將清水河流域的土地利用類(lèi)型重分為AGRL(耕地)、FRST(林地)、PAST(草地)、WATR(灘地)、URMD(城鎮(zhèn)用地)和URML(城鄉(xiāng)居民點(diǎn))6類(lèi)，其中草地、林地和耕地分別占流域總面積的39.50%、31.60%和27.88%，是流域內(nèi)主要的土地利用類(lèi)型。

流域內(nèi)土壤被重分為13種類(lèi)型，其中滯水高活性淋溶土(LVj)和不飽和雛形土(CMd)是流域內(nèi)的主要土壤類(lèi)型，占流域總面積的62.52%。疊加子流域、土地利用和土壤數(shù)據(jù)，并按模型操作手冊(cè)的建議，將土壤最小閾值設(shè)置為5%定義HRU，最終將清水河流域劃分為398個(gè)水文響應(yīng)單元。

2.2.3 氣象和土地管理數(shù)據(jù)庫(kù)

利用天氣發(fā)生器將崇禮氣象站2007-2016年監(jiān)測(cè)的日氣象數(shù)據(jù)轉(zhuǎn)換為SWAT模型計(jì)算所需的.dbf格式，構(gòu)建SWAT模型氣象數(shù)據(jù)庫(kù)；而后依據(jù)實(shí)地調(diào)研和統(tǒng)計(jì)年鑒，構(gòu)建農(nóng)業(yè)管理數(shù)據(jù)庫(kù)。至此，SWAT模型基礎(chǔ)數(shù)據(jù)庫(kù)構(gòu)建完成。

2.3 模型率定與驗(yàn)證

2.3.1 敏感性參數(shù)選擇

利用SWAT-CUP軟件中的全局敏感性分析(Global Sensitivity)方法，篩選出10個(gè)對(duì)模型模擬結(jié)果影響最大的敏感性參數(shù)，以便對(duì)模型進(jìn)行率定及驗(yàn)證。參數(shù)的物理意義和取值范圍如表2所示。

表2 敏感性參數(shù)的物理意義和取值范圍Tab.2 The physical significance and value range of sensitivity parameters

2.3.2 模型的率定與驗(yàn)證

采用張家口水文站的實(shí)測(cè)徑流數(shù)據(jù)，利用SWAT-CUP軟件中的SUFI2算法對(duì)模型進(jìn)行率定和驗(yàn)證。以2007-2009年為預(yù)熱期、2010-2013年為率定期、2014-2016年為驗(yàn)證期，根據(jù)納什效率系數(shù)Ens和線(xiàn)性回歸系數(shù)R2的結(jié)果對(duì)模型進(jìn)行適用性評(píng)價(jià)，監(jiān)測(cè)站點(diǎn)實(shí)測(cè)值與模擬值擬合結(jié)果和校驗(yàn)結(jié)果見(jiàn)圖4。結(jié)果表明，模型率定期的Ens和R2分別為0.505、0.660，驗(yàn)證期的Ens和R2分別為0.674、0.738，均滿(mǎn)足Ens≥0.5、R2≥0.6的精度要求，表明SWAT模型在清水河流域的適用性較好。

3 非點(diǎn)源污染特征分析

3.1 時(shí)間分布特征

3.1.1 年際變化特征

對(duì)SWAT模型進(jìn)行率定和驗(yàn)證后，通過(guò)子流域輸出文件(output.sub)中的ORGN、NSURQ、LAT_Q_NO3和GWNO3累加得到非點(diǎn)源TN污染負(fù)荷，通過(guò)ORGP、SOLP和SEDP累加得到非點(diǎn)源TP污染負(fù)荷。如表3所示，2011、2013、2015和2016年的年均非點(diǎn)源污染輸出負(fù)荷表現(xiàn)為T(mén)N>TP，氮、磷輸出負(fù)荷表現(xiàn)為2013年最高、2015年最低，TN、TP的最高輸出負(fù)荷量分別為626.32和160.80 t/a，研究區(qū)域內(nèi)TN造成的非點(diǎn)源污染較為嚴(yán)重。

表3 清水河流域的年均模擬輸出負(fù)荷 t

3.1.2 季節(jié)變化特征

如圖5所示，清水河流域非點(diǎn)源污染輸出負(fù)荷變化具有明顯的季節(jié)性：夏季(6-8月)TN和TP的輸出負(fù)荷占全年總負(fù)荷的比例均大于60%，是非點(diǎn)源污染產(chǎn)生的高峰期；而冬季是非點(diǎn)源污染產(chǎn)生的低谷期，TN、TP的輸出負(fù)荷占比分別為2.31%和0.74%；受降雨量和徑流量的影響，春秋兩季的輸出負(fù)荷量有明顯的波動(dòng)變化，TN的平均輸出負(fù)荷分別為7.00%、15.23%，TP的平均輸出負(fù)荷分別為1.52%、16.80%。

3.1.3 各子流域月變化特征

研究區(qū)域內(nèi)23個(gè)子流域在2011、2013、2015和2016年不同月份TN和TP輸出負(fù)荷累計(jì)值如圖6所示?？梢钥闯?，各子流域非點(diǎn)源TN和TP的輸出負(fù)荷變化特征非常相似，5、6、16、20號(hào)子流域累計(jì)輸出負(fù)荷量較大，而23、17號(hào)子流域輸出量較小?？傮w看來(lái)，汛期的輸出負(fù)荷明顯大于非汛期，所有子流域在汛期(6-8月)的輸出負(fù)荷占比超過(guò)總負(fù)荷的50%，可見(jiàn)，降雨量對(duì)非點(diǎn)源污染負(fù)荷有較大影響。

3.2 空間分布特征

根據(jù)子流域輸出文件，分別作研究區(qū)內(nèi)非點(diǎn)源TN、TP輸出負(fù)荷強(qiáng)度(kg/hm2)的空間分布圖。如圖7和圖8所示，各子流域TN、TP的輸出負(fù)荷強(qiáng)度變化范圍分別為0.038～0.469和0.002～0.136 kg/hm2，其年均值分別為0.193、0.043 kg/hm2。2011年非點(diǎn)源TN、TP輸出負(fù)荷主要集中在5、6、10、12、14號(hào)子流域，普遍達(dá)到0.246和0.053 kg/hm2；2013年非點(diǎn)源TN和TP 輸出負(fù)荷主要集中在2、3、5、6、10、12、14號(hào)子流域，普遍達(dá)到0.393和0.106 kg/hm2；2015年非點(diǎn)源TN和TP 輸出負(fù)荷主要集中在1、2、3、5、6、14號(hào)子流域，普遍達(dá)到0.274和0.054 kg/hm2；2016年非點(diǎn)源TN和TP 輸出負(fù)荷主要集中在2、3、5、6、7、10、12、14號(hào)子流域，普遍達(dá)到0.314和0.068 kg/hm2?？梢?jiàn)，2013年的輸出負(fù)荷強(qiáng)度最高。

清水河流域地形東北高、西南低，河流順地勢(shì)流淌，地形差異使得流域內(nèi)河網(wǎng)密度分布不均，在東溝和西溝區(qū)域中河網(wǎng)較密集，耕地和草地也集中分布在河流兩岸坡度較小的平原和丘陵區(qū)；植被量由東溝至西溝逐漸減少，水土流失現(xiàn)象逐漸加重。非點(diǎn)源污染負(fù)荷強(qiáng)度較高的2、3、5、6、7、10、12、14號(hào)子流域約占流域總面積的30.24%，而TN、TP輸出負(fù)荷量分別約占總負(fù)荷量的50.41%和51.40%，是研究區(qū)域非點(diǎn)源污染負(fù)荷流失的關(guān)鍵源區(qū)。關(guān)鍵源區(qū)集中分布在東溝和西溝水系兩側(cè)，主要土地利用類(lèi)型為耕地和草地，由于水土流失和翻耕、施肥等農(nóng)業(yè)管理措施的影響，導(dǎo)致了較高強(qiáng)度的污染流失；其次，在正溝水系周?chē)簿哂休^高的輸出強(qiáng)度，該區(qū)域林地的植被覆蓋度較高，植物根系對(duì)營(yíng)養(yǎng)物質(zhì)的固定作用較好，有效緩了降雨徑流的沖刷和污染物的流失。

3.3 污染負(fù)荷與降水的關(guān)系

3.3.1 時(shí)間變化關(guān)系

為進(jìn)一步探究清水河流域非點(diǎn)源污染輸出負(fù)荷變化與降水的關(guān)系，如圖9所示，統(tǒng)計(jì)流域2011、2013、2015和2016年的月均降雨量及非點(diǎn)源污染逐月輸出負(fù)荷量。由圖9可見(jiàn)，研究區(qū)域非點(diǎn)源TN、TP負(fù)荷的變化趨勢(shì)基本一致，年際變化較大，并呈現(xiàn)出中間高、兩頭低的形狀變化，最高值均出現(xiàn)在6-8月。降雨量的形狀變化特征與其相似，峰值出現(xiàn)在6-8月，是研究區(qū)域的汛期。污染物年內(nèi)總變化統(tǒng)計(jì)如表4所示，污染物在汛期(6-8月)的輸出量占比較高，TN、TP污染負(fù)荷分別占全年總量的75.07%和79.56%，降雨量和徑流量在汛期占全年總量的比重也超過(guò)50%。在農(nóng)業(yè)管理措施中，6-8月是流域內(nèi)作物種植、施肥的重要時(shí)期，化肥積累量大；而汛期降雨量豐富，徑流量也較大，降雨徑流的沖刷和運(yùn)移對(duì)流域內(nèi)非點(diǎn)源污染物影響較大。

3.3.2 相關(guān)性分析

為進(jìn)一步探究清水河流域內(nèi)非點(diǎn)源輸出負(fù)荷與降雨量的關(guān)系，利用SPSS 19.0軟件對(duì)研究區(qū)域的降雨量、徑流量和非點(diǎn)源TN、TP污染負(fù)荷進(jìn)行相關(guān)性分析，如表5所示，徑流量的形成與降雨量呈正相關(guān)關(guān)系，其Pearson相關(guān)性為0.650；非點(diǎn)源污染負(fù)荷與降雨量、徑流量的相關(guān)性均>0.7，在0.01水平上顯著相關(guān)?？梢?jiàn)，降雨量和徑流量對(duì)非點(diǎn)源污染的產(chǎn)生和遷移至關(guān)重要。

表4 清水河流域汛期指標(biāo)及占全年比例Tab.4 The index of Qingshui river basin in flood season and its proportion for the whole year

表5 各值相關(guān)性分析Tab.5 Correlation analysis of each value

研究區(qū)域各子流域在模擬期內(nèi)的非點(diǎn)源污染輸出負(fù)荷與降雨量的相關(guān)系數(shù)空間分布如圖10所示，總體上看，降水作為污染物流失的驅(qū)動(dòng)因子，非點(diǎn)源污染輸出負(fù)荷均與降水呈正相關(guān)。非點(diǎn)源TN輸出負(fù)荷與降水相關(guān)性較強(qiáng)的區(qū)域主要集中分布在正溝、西溝兩側(cè)和東溝上游，其中相關(guān)性最高的區(qū)域?yàn)?、4、5、14號(hào)子流域，相關(guān)系數(shù)為0.769～0.857，因此降水是導(dǎo)致該區(qū)域內(nèi)非點(diǎn)源污染的流失主導(dǎo)因素。各子流域的非點(diǎn)源TP輸出負(fù)荷與降水相關(guān)性的空間分布與TN有部分差別，如3、6、15、16號(hào)子流域的TP輸出負(fù)荷與降水的相關(guān)性明顯弱于TN。耕地是這些區(qū)域中的主要土地利用類(lèi)型，并且多為坡度較小的平原和丘陵區(qū)，水土流失現(xiàn)象較嚴(yán)重；農(nóng)業(yè)生產(chǎn)施用的磷肥進(jìn)入土壤后較難溶解，而溶解的部分容易被土壤膠粒吸附形成較穩(wěn)定的物質(zhì)，不易被釋放和移動(dòng)，而氮素易被淋洗而被徑流攜帶流失，因此使得該區(qū)域的氮磷輸出負(fù)荷與降水和相關(guān)系數(shù)存在一定的差異。

相關(guān)系數(shù)高的區(qū)域，非點(diǎn)源污染物的流失與降水密切相關(guān)，高降雨量及強(qiáng)降水都會(huì)引起污染物的大量流失。此外，結(jié)合非點(diǎn)源污染負(fù)荷的空間分布可知，該區(qū)域自身非點(diǎn)源污染負(fù)荷強(qiáng)度較高，會(huì)加重非點(diǎn)源污染。因此，非點(diǎn)源污染物與降水相關(guān)性較高的區(qū)域一定要加強(qiáng)管理，特別是要嚴(yán)格控制農(nóng)業(yè)種植區(qū)的化肥施用量和施用強(qiáng)度，合理施肥；草地及林地區(qū)應(yīng)合理建造植被過(guò)濾帶，以加強(qiáng)水土保持效果，削減徑流量，減少非點(diǎn)源污染物的流失。

4 結(jié) 論

利用SWAT模型對(duì)清水河流域進(jìn)行非點(diǎn)源污染模擬及特征分析，主要結(jié)論如下：

(1)研究區(qū)域農(nóng)業(yè)非點(diǎn)源污染的時(shí)間分布特征表明：在年際變化上，2011、2013、2015和2016年的年均非點(diǎn)源污染輸出負(fù)荷均表現(xiàn)為T(mén)N>TP，TN、TP的最高輸出負(fù)荷量分別為626.32和160.80 t/a，TN造成的非點(diǎn)源污染較為嚴(yán)重；在季節(jié)變化上，非點(diǎn)源污染輸出負(fù)荷具有明顯的季節(jié)性，夏季和冬季分別是其產(chǎn)生的高峰期和低谷期；此外，各子流域的年累積負(fù)荷量在6-8月份占比超過(guò)總輸出負(fù)荷的50%，呈現(xiàn)出汛期明顯大于非汛期的特征。

(2)研究區(qū)域農(nóng)業(yè)非點(diǎn)源污染的空間分布特征表明：2、3、5、6、7、10、12、14號(hào)子流域TN、TP輸出負(fù)荷量分別約占總負(fù)荷量的50.41%和51.40%，是非點(diǎn)源污染負(fù)荷流失的關(guān)鍵源區(qū)，這些區(qū)域集中分布在東溝和西溝水系兩側(cè)，主要土地利用類(lèi)型為耕地和草地。

(3)根據(jù)非點(diǎn)源污染負(fù)荷的時(shí)空分布特征，進(jìn)一步探究污染輸出負(fù)荷與降雨量的關(guān)系，結(jié)果表明：在時(shí)間上，流域內(nèi)的非點(diǎn)源TN、TP的逐月輸出負(fù)荷與月均降雨量均呈現(xiàn)出中間高、兩頭低的形狀變化，汛期(6-8月)TN、TP污染負(fù)荷分別占全年總量的75.07%和79.56%，降雨量和徑流量的比重也超過(guò)50%；相關(guān)性分析表明，各子流域的污染輸出負(fù)荷與降雨量呈正相關(guān)關(guān)系，輸出負(fù)荷與降水相關(guān)性較強(qiáng)的區(qū)域集中分布在正溝、西溝兩側(cè)和東溝上游，該區(qū)域內(nèi)非點(diǎn)源污染物的流失與降水密切相關(guān)，高降雨量及強(qiáng)降水都會(huì)引起污染物的大量流失。

□