劉 超，李建柱，徐彬鑫

(1.河北省邢臺水文勘測研究中心，河北 邢臺 054000；2.天津大學(xué)水利工程仿真與安全國家重點實驗室，天津 300350；3.水利部海河水利委員會，天津 300170)

水資源是人類生活、生產(chǎn)和生態(tài)環(huán)境的重要保障，中國的水資源僅為世界平均水平的1/4，且降水的多年變化和年內(nèi)集中程度均很大[1]，徑流多以洪水的形式發(fā)生，較大一部分徑流不能被利用。隨著社會經(jīng)濟的發(fā)展，點源污染和非點源污染導(dǎo)致的水環(huán)境問題也逐漸凸顯，進一步加劇了我國水資源的供需矛盾。由于治理能力的不斷提升與管理措施的不斷嚴格，點源污染基本得到控制，非點源污染已成為造成水體污染的主要原因。

流域作為自然形成的水資源動態(tài)系統(tǒng)，是水循環(huán)過程及水行政管理中的重要單元。流域問題的研究，是水資源開發(fā)利用，合理配置優(yōu)化，綜合決策管理等宏觀問題研究的落腳點[2]。首要任務(wù)是獲取流域詳實的水文、氣象信息，但受限于地理條件與歷史因素，同時考慮到經(jīng)濟成本，觀測數(shù)據(jù)往往集中在流域的中下游區(qū)域[3]，且在時間上也不具有連續(xù)性，客觀上限制了流域決策管理的精細化。因此，流域水文、水動力、水質(zhì)模型成為補充缺測數(shù)據(jù)的重要手段，以及流域水資源綜合管理的有效工具。

潘家口水庫是天津市、唐山市重要的水源地，承擔(dān)著天津、唐山工業(yè)、農(nóng)業(yè)、城市生活和環(huán)境用水的重任[4]。截至2018年底，累計向天津市供水187億 m3，為天津市經(jīng)濟可持續(xù)發(fā)展做出了重要貢獻。2000年之前，潘家口水庫各斷面水質(zhì)均為III類水標(biāo)準，近年來，隨著社會經(jīng)濟的快速發(fā)展，潘家口水庫上游工業(yè)和城市生活污水、農(nóng)業(yè)生產(chǎn)中大量農(nóng)藥化肥的施用、水庫網(wǎng)箱養(yǎng)魚及水庫周邊旅游設(shè)施造成的污染日益嚴重。2016-2017年，潘家口水庫由于水質(zhì)較差，天津中斷了引灤供水。因此，構(gòu)建潘家口水庫控制流域污染物負荷模型，模擬分析潘家口水庫的污染物負荷，對于精準分析污染物來源，并提出針對性措施和對策可提供科學(xué)依據(jù)。

1 研究區(qū)概況

灤河發(fā)源于河北省承德市境內(nèi)的巴彥古爾圖山北麓，流經(jīng)河北省、內(nèi)蒙古自治區(qū)及遼寧省等地區(qū)，最終于河北省樂亭縣匯入渤海[5]。灤河流域位于39°10′N～42°30′N，115°30′E～119°15′E，控制面積為44 750 km2。潘家口水庫是灤河上的大型水庫，是引灤入津的水源地，目前承擔(dān)著保障天津市河道生態(tài)需水的重要任務(wù)，流域位置及水系情況如圖1。

根據(jù)2018年國控斷面相關(guān)河湖水質(zhì)監(jiān)測資料，灤河流域重要斷面水質(zhì)為Ⅲ～劣Ⅴ類。主要是一些月份的水質(zhì)情況不滿足各自水環(huán)境控制單元與水功能區(qū)的控制目標(biāo)。其中，瀑河大桑園控制斷面的水質(zhì)最差時為劣Ⅴ類，化學(xué)需氧量(COD)、五日生化需氧量(BOD)等指標(biāo)明顯超標(biāo)；灤河干流的上板城大橋控制斷面水質(zhì)為Ⅴ類，高錳酸鹽指數(shù)、化學(xué)需氧量等指標(biāo)超標(biāo)。研究流域內(nèi)水質(zhì)控制斷面2018年水質(zhì)情況見表1。

圖1 研究區(qū)域地理位置及水系圖

2 流域污染物負荷模型構(gòu)建

根據(jù)流域水文氣象、土地利用及污染物來源等數(shù)據(jù)，構(gòu)建SWAT水文水質(zhì)模型對流域徑流過程和污染物負荷進行模擬分析。

2.1 SWAT模型原理

SWAT模型(Soil and Water Assessment Tool)是美國農(nóng)業(yè)部開發(fā)的流域尺度分布式水文模型，模型可用于模擬單一流域或多個具有水文聯(lián)系流域的水文過程[6]。SWAT模型內(nèi)包含有土地利用、地下水、流域管理、河道水質(zhì)情況等多個模塊，能夠預(yù)測長時間序列下流域的水文、泥沙與營養(yǎng)物質(zhì)的變化情況[7]。

表1 研究流域重要水域控制斷面水質(zhì)

SWAT模擬的水文循環(huán)基于水量平衡方程

(1)

式中：SWt為最終土壤水含量，mm；SW0為第i天初始土壤含水量，mm；t為天數(shù)，d；Rday為第i天降水量，mm；Qsurf為第i天地表徑流，mm；Ea為第i天蒸發(fā)量，mm；Wseep為第i天土壤底層滲透量，mm；Qgw為第i天地下水出流，mm。

SWAT在降雨徑流計算中則采用SCS曲線法，該方程自上世紀50年代起開始得到廣泛應(yīng)用，其方程式為：

(2)

式中：Qsurf為地表徑流量，mm；Rday為降雨量，mm；Ia為初損量，mm；S為滯留系數(shù)，mm；表達式如下：

(3)

式中：CN為模擬日的曲線數(shù)，初損Ia通常估計為0.2S。

表2 SWAT模型建模所需基礎(chǔ)數(shù)據(jù)

2.2 SWAT模型構(gòu)建

SWAT模型在構(gòu)建過程中需要輸入數(shù)字高程數(shù)據(jù)、土壤數(shù)據(jù)、土地利用數(shù)據(jù)與氣象數(shù)據(jù)(降水、氣溫、太陽輻射、風(fēng)速及相對濕度等)，具體輸入數(shù)據(jù)見表2。在模型建立后，還需要利用實測徑流數(shù)據(jù)對模型進行率定。具體建模過程如下：首先通過八流向法根據(jù)DEM與設(shè)定的集水面積閾值生成流域內(nèi)河網(wǎng)，并劃定出若干子流域；其次，定義流域內(nèi)土地利用情況、土壤類型與坡度分級，通過設(shè)定各自對應(yīng)的閾值，進而形成閾值組合以在各子流域內(nèi)生成水文響應(yīng)單元(HRUs)；之后在模型中加載氣象數(shù)據(jù)并運行，SWAT以HRU為最小計算單元各自獨立進行水文循環(huán)過程模擬，之后通過SCS曲線法計算得到地表徑流，河道水流演算則采用馬斯京根法。

2.2.1 數(shù)字高程(DEM)數(shù)據(jù)處理

潘家口水庫控制流域DEM由NASA提供的SRTM的DEM數(shù)據(jù)裁剪得到，分辨率為30 m×30 m(如圖2)，將研究區(qū)域劃分了34個子流域。

圖2 研究區(qū)DEM及子流域劃分情況

2.2.2 土地利用及土壤數(shù)據(jù)

選用2005年土地利用數(shù)據(jù)，并將土地利用數(shù)據(jù)進行重分類，分為林地、草地、耕地、水域、建設(shè)用地及未利用土地六大類，

土壤數(shù)據(jù)資料來源于HWSD土壤數(shù)據(jù)庫(http://www.fao.org)，該數(shù)據(jù)庫由聯(lián)合國糧農(nóng)組織(FAO)、中國科學(xué)院南京土壤研究所(ISSCAS)、歐盟委員會聯(lián)合研究中心(JRC)等國內(nèi)外多家機構(gòu)聯(lián)合構(gòu)建。研究區(qū)域土地利用情況與土壤類型見圖3。

圖3 研究區(qū)域土壤類型及土地利用情況

2.2.3 水文氣象數(shù)據(jù)

降水?dāng)?shù)據(jù)來自三道河子、寬城、下河南等10個雨量站1998-2018年逐日數(shù)據(jù)，氣象數(shù)據(jù)選用多倫、圍場、承德及青龍4個氣象站同時段逐日數(shù)據(jù)，數(shù)據(jù)來源為中國氣象科學(xué)數(shù)據(jù)共享服務(wù)網(wǎng)(http://cdc.cma.gov.cn)。

2.2.4 污染物數(shù)據(jù)概化

農(nóng)村生活污染與畜禽污染年內(nèi)污染負荷產(chǎn)生量較為固定，因此在SWAT模型中將兩種污染分別根據(jù)其各自TN、TP、NH3-N的比例，定義新的化肥類型，以連續(xù)施肥的形式進行概化。施肥持續(xù)時間設(shè)置為1年，施肥頻次設(shè)置為每日一次，同時按照各子流域人口與面積加權(quán)，確定子流域污染排放量。

對于化肥污染，承德市主要種植作物為玉米、大豆、馬鈴薯，生長季節(jié)約在每年3-10月份，因此設(shè)置施肥時間為每年3-9月，施肥頻次設(shè)置為每月一次。

點源污染數(shù)據(jù)采用各排污口排放數(shù)據(jù)輸入模型，由于收集到資料為年排放總量，因此將其概化為每月輸出量相等，即年內(nèi)12個月均勻排放至流域。

2.3 SWAT模型適用性評價

采用相關(guān)系數(shù)R2及納什效率系數(shù)(NSE)兩種指標(biāo)來評價模擬結(jié)果，其計算公式分別見式(4)、(5)。

(4)

(5)

2.3.1 徑流量率定驗證結(jié)果分析

采用SWAT-CUP對模型參數(shù)進行率定，經(jīng)過2 000次迭代得到最終模擬結(jié)果，計算各水文站率定期與驗證期的相關(guān)系數(shù)與納什效率系數(shù)，見表3，潘家口水庫入庫徑流量模擬結(jié)果見圖4。可見流域內(nèi)各控制站R2、NSE均高于0.50，說明模型能較好地反映流域內(nèi)水文變化過程。尤其是在潘家口水庫入庫徑流模擬時，率定期和驗證期NSE分別為0.80和0.79，表明構(gòu)建的SWAT模型對潘家口水庫水量模擬結(jié)果很好。

表3 各水文站模擬效果評價

圖4 潘家口水庫徑流模擬與實測對比

2.3.2 水質(zhì)率定驗證結(jié)果

選擇2018年對模型進行率定，2017年對模型進行驗證。對各站TN、TP模擬值與實測值的相關(guān)系數(shù)進行計算，結(jié)果如表4、表5所示，R2均在0.57以上。從圖5可以看出，率定后的模型能較好地反映研究流域營養(yǎng)物遷移轉(zhuǎn)化情況，但對于汛期的極值點模擬效果并不理想，主要原因為汛期氮磷負荷主要受泥沙影響，而本研究中由于缺乏泥沙監(jiān)測數(shù)據(jù)，未對其進行率定驗證，因而造成了誤差。

表4 各監(jiān)測斷面TN模擬效果評價

表5 各監(jiān)測斷面TP模擬效果評價

圖5 各控制站TN、TP模擬值與實測值對比

3 結(jié)語

基于DEM、土地利用、水文氣象與污染物等數(shù)據(jù)構(gòu)建了潘家口水庫控制流域的SWAT水量水質(zhì)模型，并針對水量、水質(zhì)等指標(biāo)對模型進行了率定與驗證，采用R2與NSE評價了模型在該地區(qū)的適用性，主要結(jié)論如下：

(1)基于水文氣象與水質(zhì)資料，選取潘家口水庫控制流域，劃分了34個子流域，1 473個水文響應(yīng)單元，構(gòu)建了水量與污染物負荷的SWAT模型，

(2)采用2000-2010年實測徑流數(shù)據(jù)對水量模型進行率定，2011-2015年數(shù)據(jù)進行驗證，各水文站徑流模擬的相關(guān)系數(shù)與NSE均較高，說明模型能較好地反映流域內(nèi)的水量變化過程。

(3)選取總氮與總磷兩個水質(zhì)指標(biāo)，以2017、2018年實測斷面數(shù)據(jù)進行了污染物負荷模型的率定和驗證，結(jié)果表明模型能較好地反映污染物的遷移轉(zhuǎn)化過程。