基于HSPF模型的東江流域氮磷污染研究

孫滔滔，趙 鑫，尹魁浩，彭盛華，林 青，吳 靜

(1.深圳市環(huán)境科學(xué)研究院國家環(huán)境保護(hù)飲用水水源地管理技術(shù)重點實驗室，廣東 深圳 518001；2.廣東粵港供水有限公司，廣東 深圳 518021；3.San Francisco Estuary Institute，CA，USA)

水污染問題目前已成為制約我國社會經(jīng)濟(jì)可持續(xù)發(fā)展的問題之一。近30年來，由于人口聚集、城市更新，氮、磷等污染物過量排放，很多重要水域的水體污染，尤其是氮、磷富營養(yǎng)化問題日益突出，已經(jīng)對很多地區(qū)人民生活用水構(gòu)成威脅，對生產(chǎn)和經(jīng)濟(jì)發(fā)展造成巨大損失。東江肩負(fù)著深港兩地的跨流域調(diào)水供應(yīng)，是廣東省重要的水源地，其水質(zhì)狀況直接影響著沿岸以及深港地區(qū)的供水安全。隨著東江流域經(jīng)濟(jì)的快速發(fā)展、城鎮(zhèn)化進(jìn)程的加快，東江水質(zhì)部分河段已趨于富營養(yǎng)化狀態(tài)[1]。如何確定氮、磷的確切來源以及各種污染因子的時空分布特征和遷移演化規(guī)律，是有效地提升污染防治水平，達(dá)成東江流域污染控制目標(biāo)和保持水質(zhì)的關(guān)鍵所在。

近年來，國內(nèi)學(xué)者對東江的污染問題從不同角度開展了一系列的研究。江濤等[2]根據(jù)東江干流博羅站、河源站和龍川站的監(jiān)測資料，對CODMn、BOD5和NH3-N等水質(zhì)參數(shù)的濃度變化趨勢進(jìn)行分析。王博等[1]結(jié)合東江干流區(qū)域差異性對干流的營養(yǎng)狀態(tài)進(jìn)行了分段研究。李星等[3]全面研究了東江干流水體氮的時空變化特征，并利用附生藻的δ15N值示蹤了其來源。溫美麗等[4]通過采樣分析氮磷質(zhì)量分?jǐn)?shù)，研究了新豐江水庫上游氮磷污染的時空變化。曾金鳳等[5]基于東江源區(qū)1975-2016年共42年的歷史數(shù)據(jù)，解析了源區(qū)水體的水質(zhì)及特征污染物氨氮的變化情況。這些研究都是基于對實測數(shù)據(jù)的分析。

流域模型是研究復(fù)雜的流域水文水質(zhì)輸移過程及各種污染因子的時空分布特征和遷移演化規(guī)律的有效途徑和方法，因此被廣泛用于點源及非點源氮、磷污染研究。近年來，流域模型也被用于研究東江流域降雨徑流及污染負(fù)荷特征。劉潔等從土地利用[6]和降雨[7]角度，用SWAT(Soil & Water Assessment Tool)模型對東江流域徑流變化進(jìn)行研究。呂樂婷等[8]用SWAT模擬了東江流域1960-2010年51年間月徑流與降水關(guān)系。董延軍[9]和石赟赟[10]等用HSPF(Hydrological Simulation Program FORTRAN)對東江流域建模發(fā)現(xiàn)，城鎮(zhèn)化會導(dǎo)致流域調(diào)蓄能力下降，進(jìn)而產(chǎn)流能力增強(qiáng)，尤其是汛期這種表現(xiàn)更為顯著。這些研究表明流域模型是研究東江這樣的大型流域上復(fù)雜的水文水質(zhì)問題的一個有效實用的工具。

本文研究的主要目的是用HSPF構(gòu)建一個高精度的東江流域模型，來研究流域內(nèi)降雨徑流及污染負(fù)荷特征，從而為東江流域的污染治理提供決策基礎(chǔ)。HSPF模型作為分布式綜合性流域模型的代表,是研究大型、混合型土地利用流域污染問題的最佳模型之一。目前，HSPF模型因為對輸入資料要求更高，建模專業(yè)門檻更高等原因，在國內(nèi)相對應(yīng)用較少。但是，近年來基于HSPF建模的美國最大的海灣Chesapeake Bay流域[11]、滇池流域[12]、苕溪流域[13]等研究區(qū)的案例表明，HSPF模型對流域長期連續(xù)的徑流量模擬具有較好的適用性，而且能夠完整地模擬氮磷遷移、轉(zhuǎn)化過程。因此我們采用HSPF來研究東江流域氮磷污染。

1 研究區(qū)域概況

東江流域是廣東、香港的水資源供應(yīng)地，流域內(nèi)多年平均降雨量約為1 800 mm，年徑流總量約為330 億m3。東江發(fā)源于江西省尋烏縣，干流全長約562 km，包括16條一級支流，流域總面積約3.5 萬km2，其中廣東省境內(nèi)面積約3.2 萬km2，占比約90%(見圖1)。

圖1 東江流域及主要水系[14]Fig.1 Dongjiang river basin and major tributary[14]

2 研究方法

2.1 HSPF模型綜述

HSPF是由美國環(huán)保署開發(fā)維護(hù)的一個半分布式綜合性流域模型。HSPF模型將常見的污染物和毒性有機(jī)物模擬納入到模型中，能夠?qū)崿F(xiàn)多種污染物地表、壤中流過程及蓄積、遷移、轉(zhuǎn)化的綜合模擬[15]。HSPF能夠?qū)崿F(xiàn)不同時空尺度的降雨、下滲等過程的動態(tài)和連續(xù)模擬?？臻g尺度方面，研究者可以將研究區(qū)劃分為具有承接關(guān)系的子流域，子流域可以進(jìn)一步劃分為響應(yīng)單元，子流域和響應(yīng)單元的大小可以根據(jù)不同需求進(jìn)行調(diào)整。時間尺度方面，HSPF一般用于模擬時間尺度為小時的產(chǎn)匯流過程，根據(jù)需要也可以模擬年、月、日、分鐘等時間尺度。過去三十多年，HSPF模型已被廣泛應(yīng)用于流域水文、水質(zhì)模擬研究，包括氣候及土地利用變化對流域產(chǎn)流的影響，流域點源或非點源污染負(fù)荷確定，泥沙、營養(yǎng)物質(zhì)、殺蟲劑傳輸模擬以及各種流域管理措施對河流水質(zhì)的影響等方面的研究。

2.2 流域空間的劃分

HSPF模型構(gòu)建的第一步是通過氣象分布以及流域物理特性的差異，把流域劃分成多個具有不同類型的水文響應(yīng)的子流域。根據(jù)GIS(Geographic Information System)地形地貌以及土壤數(shù)據(jù)，東江流域在模型中被劃分為407個子流域，面積從0.008 km2到412 km2。這是目前東江流域模型研究中空間離散精度最高的。這樣高度細(xì)分的子流域能充分反應(yīng)流域空間地形地貌的差異，為準(zhǔn)確模擬流域各個地區(qū)的水文水質(zhì)響應(yīng)和分布的提供了基礎(chǔ)。

2.3 模型輸入數(shù)據(jù)

HSPF模型需要大量的輸入數(shù)據(jù)，包括氣象、土地利用及污染源數(shù)據(jù)。氣象數(shù)據(jù)包括每小時降水量、蒸發(fā)量、日平均氣溫、平均風(fēng)速及太陽輻射等。東江流域模型的氣象數(shù)據(jù)來源于14個國家氣象站(見圖2)，2009到2015年7年的氣象數(shù)據(jù)用來驅(qū)動模型。土地利用數(shù)據(jù)采用2015年100 m空間分辨率的GIS數(shù)據(jù)，在模型中分為六類：農(nóng)業(yè)、森林、草地、城市、水體及其他。東江流域的點源污染包括生活污水和工業(yè)污水(見圖3)，面源污染包括禽畜養(yǎng)殖業(yè)和農(nóng)業(yè)肥料。生活污水和工業(yè)污水以及禽畜養(yǎng)殖業(yè)數(shù)據(jù)來自于廣東省2015年環(huán)境統(tǒng)計數(shù)據(jù)。農(nóng)業(yè)施肥根據(jù)廣東省農(nóng)業(yè)手冊的年施肥量以及農(nóng)作物生長季節(jié)，估算春秋兩季每畝的施肥量。

圖2 東江流域水文站和氣象站Fig.2 Hydrological station and weather station in Dongjiang river basin

2.4 模型校核

模型校核是一個通過調(diào)整模型參數(shù)來吻合模型結(jié)果和監(jiān)測數(shù)據(jù)的過程。東江流域設(shè)有14個長期的水文站(見圖2)，其中12個站2011至2014年的日流量數(shù)據(jù)用來水文校核(見表1，楓樹壩水庫壩下站和新豐江水庫站用實際泄流數(shù)據(jù)，所以不需要校核)。HSPF模型的水文校核通常用圖形及統(tǒng)計指標(biāo)來綜合判斷[16]，主要包括時間序列圖、確定性系數(shù)R2、總徑流量誤差以及Nash-Sutcliffe 效率系數(shù)(NSE)。確定性系數(shù)R2可用來評價實測值與模擬值變化趨勢的一致性，R2越接近于1，表示模擬值和實測值變化趨勢越吻合。總徑流量誤差和模型系數(shù)NSE表示了實測值與模擬值的偏離程度，可以用來判斷模型結(jié)果和監(jiān)測數(shù)據(jù)的總體吻合度。NSE計算公式如下：

(1)

式中：Yobs為觀測徑流；Ysim為模擬徑流；Ymean為觀測徑流的平均值。

一般R2大于0.6，NSE≥0.5時，模型的模擬結(jié)果被認(rèn)為是可信的。對HSPF模型，當(dāng)總徑流誤差絕對值<10%，NSE≥0.65時，模擬結(jié)果被認(rèn)為非常好[16]。由于數(shù)據(jù)限制，水質(zhì)校核采用了合理保證的方法。模型中主要的水質(zhì)參數(shù)選用了HSPF水質(zhì)校核指導(dǎo)文獻(xiàn)[17]中推薦的典型值，然后根據(jù)水質(zhì)模擬結(jié)果與水質(zhì)分類指標(biāo)和相關(guān)數(shù)據(jù)的對比進(jìn)行了微調(diào)。校核后的模型然后用來確定氮、磷污染的時空分布特征和遷移演化規(guī)律。

3 結(jié)果與討論

3.1 模型校核結(jié)果

對于每個校核站點，時間序列圖可以用來判斷模擬的流量在數(shù)量級以及時間上的分布是否和監(jiān)測數(shù)據(jù)吻合。圖4和圖5分別是東江干流河源和博羅站2011-2014年的日流量模擬結(jié)果與實測結(jié)果的對比圖。這兩個圖表明HSPF模擬的流量和實測數(shù)據(jù)十分吻合。再結(jié)合相應(yīng)的統(tǒng)計指標(biāo)(見表1)，結(jié)果表明模型在這兩個站的模擬效果很好。

圖5 東江博羅站日流量模擬結(jié)果與實測結(jié)果的對比圖Fig.5 Comparison of daily flow simulation results and measured results of Boluo Station

表1列出了東江模型在所有12個水文站的統(tǒng)計指標(biāo)。其中9個站的確定性系數(shù)R2>0.6；10個站的總徑流誤差達(dá)到了HSPF模型的流量校核標(biāo)準(zhǔn)(絕對值<10%)，只有九洲和平山兩個站略微超標(biāo)；7個站的模型系數(shù)NSE大于0.5，但上游兩個支流站(勝前和岳城)的模型系數(shù)小于0.3。對這兩個站的降雨徑流分析表明這個結(jié)果主要是受降雨分布影響，模型所用的氣象數(shù)據(jù)可能不能完全代表水文站當(dāng)?shù)氐慕涤昵闆r，所以導(dǎo)致了結(jié)果的偏差。但總體來講，綜合各個指標(biāo)，東江模型在大部分站都達(dá)到了良好標(biāo)準(zhǔn)，尤其是干流上的主要站點。

3.2 東江水質(zhì)的時空變化

校核后的HSPF模型可以輸出任意地點的水量和多種污染物指標(biāo)，包括水溫，溶解氧，BOD，泥沙，以及氮、磷濃度，因此可以用來分析東江水量和水質(zhì)在時間和空間上的變化。本文選取東江干流上龍川，河源，嶺下，博羅(見圖2)4個監(jiān)測站2010-2015年每月平均氨氮、總氮、總磷的濃度來探討水質(zhì)的變化。

表1 東江模型水文校核統(tǒng)計結(jié)果Tab.1 Dongjiang river model hydrological calibration result

東江干流上游龍川站氨氮濃度最低，河源站最高(見圖6)。龍川站月平均氨氮濃度為0.04 mg/L到0.16 mg/L,而河源站的氨氮濃度則為0.17 mg/L到1.17 mg/L，是龍川站的5倍多。這主要是由于河源站附近有多個污染源造成的(見圖2)。從河源到嶺下再到博羅，氨氮濃度依次降低，主要是沿途污染源負(fù)荷較低，加上大量的來水稀釋了上游來水的濃度。從時間上來說，由于污染源數(shù)據(jù)缺乏變化，每個站的氨氮濃度年際間變化不大，但是月與月之間濃度變化較大。河源站和嶺下站汛期(4-6月)濃度最高，表明這兩個站受面源的影響很大，但另外兩個站則是枯水期(11-1月)濃度最高，表明它們受徑流的影響更大。

圖6 東江干流2010-2015每月平均氨氮濃度Fig.6 Monthly average ammonia nitrogen concentration of Dongjiang main stream(2010-2015)

從空間分布上來講，東江干流的總氮濃度上下游變化不大(見圖7)，盡管上游龍川站的濃度低于其他3個站。龍川站總氮濃度的范圍是0.5到3.2 mg/L，河源站的范圍是0.8到5.6 mg/L，相差不大。而且和氨氮不同的是，總氮濃度從河源到嶺下再到博羅幾乎沒有變化，表明上游的氨氮部分轉(zhuǎn)化為硝氮。從時間分布上，各個站的總氮濃度雖然年際間變化不大，但是總體呈現(xiàn)上升的趨勢(見圖7)，而且汛期(4-9月)濃度高于非汛期，表明總氮受降雨徑流及面源污染的影響很大。

圖7 東江干流2010-2015每月平均總氮濃度Fig.7 Monthly average total nitrogen concentration of Dongjiang main stream(2010-2015)

東江干流總磷濃度的空間分布呈現(xiàn)出和氨氮相似的趨勢(見圖8)，上游龍川站濃度最低，河源站最高，從河源到嶺下再到博羅總磷濃度遞次降低?？偭诐舛仍邶埓ㄕ局挥?.05到0.15 mg/L，但河源站的濃度是龍川的10倍之多，從0.21 mg/L到1.93 mg/L。這主要是由于河源站附近有多個污染源以及汛期上游的來水?dāng)y帶的污染物造成的(見圖2)。從時間分布上，各個站的總磷濃度年際間變化不大，但是總體上有上升的趨勢(見圖8)，而且汛期(4-9月)濃度高于非汛期。這和總氮的時間分布一致，表明水質(zhì)受汛期來水和面源污染的影響很大。

圖8 東江干流2010-2015每月平均總磷濃度Fig.8 Monthly average total phosphorus concentration of Dongjiang main stream(2010-2015)

東江流域的HSPF模型還可以用來估算各地區(qū)的來水和污染物負(fù)荷。表2列出了模型測算的東江一級支流年平均水量和氮磷負(fù)荷及其百分比。東江支流年平均來水量為241 億m3，其中增江河、西枝江、新豐江水庫及楓樹壩水庫四條支流的來水占61.4%，其余38.6%來自于其他19條支流。氨氮年平均負(fù)荷約9 000 t，總氮約6 萬t，總磷約3 000 t。增江河、公莊河、新豐江水庫和西枝江四條支流總共貢獻(xiàn)了71.9%的氨氮，61.6%的總氮和73.8%的總磷。在這四條支流中，公莊河的水量占6.1%，但負(fù)荷卻很高(22.4%的氨氮，12.7%的總氮，15.9%的總磷)，分析是由于這個支流流域內(nèi)有多個禽畜養(yǎng)殖場排放氮磷(見圖3)，造成污染物濃度偏高。其他3個支流的高負(fù)荷則主要是由于來水量大造成的。東江上游以農(nóng)業(yè)經(jīng)濟(jì)為主，禽畜業(yè)和農(nóng)業(yè)面源污染是其主要污染源，應(yīng)該要著重控制。中下游則是要集中處理城市化帶來的生活和工業(yè)污水。

表2 東江流域主要支流年平均水量和氮磷負(fù)荷Tab.2 Annual average water volume and nitrogen and phosphorus load of major tributaries

4 結(jié) 論

本文構(gòu)建了一個基于HSPF的分布式模型來研究東江流域徑流和氮磷的時空分布。東江流域HSPF模型包括407個子流域，在12個水文站上進(jìn)行了校核。綜合各個指標(biāo)來看，東江模型在大部分站，尤其是干流上的主要站點，都達(dá)到了良好標(biāo)準(zhǔn)。校核后的模型可以用來確定氮、磷污染的時空分布特征和遷移演化規(guī)律，模型結(jié)果表明：

(1)東江干流上游龍川站氨氮濃度最低，河源站最高。從龍川到博羅氨氮濃度年際間變化不大，但是月與月之間濃度變化較大。

(2)東江干流的總氮濃度上下游變化不大，表明上游的氨氮部分轉(zhuǎn)化為硝氮。各個站的總氮濃度年際間變化不大，主要受降雨徑流及面源污染的影響。

(3)東江干流總磷濃度上游龍川站濃度最低，河源站最高，從河源到嶺下再到博羅總磷濃度遞次降低。各個站的總磷濃度年際間變化不大，汛期(4-9月)濃度高于非汛期，表明水質(zhì)受汛期來水和面源污染的影響很大。

(4)東江氨氮年平均負(fù)荷約9 000 t，總氮約6 萬t，總磷約3 000 t。增江河、公莊河、新豐江水庫和西枝江四條支流總共貢獻(xiàn)了71.9%的氨氮，61.6%的總氮和73.8%的總磷。因此東江負(fù)荷控制的重點應(yīng)該放在這幾條大的支流。

(5)分布式流域模型為研究東江流域復(fù)雜的水文水質(zhì)及各種污染因子的時空分布特征提供了一個有效的方法?？煽康臄?shù)據(jù)是構(gòu)建模型的關(guān)鍵，因此應(yīng)該把監(jiān)測和收集流域內(nèi)的水文水質(zhì)數(shù)據(jù)作為東江流域未來污染研究的一個重點。

□