缺資料流域水污染物總量分配方法研究

張新華，李紅霞，肖玉成，陳 奕

（四川大學 水力學與山區(qū)河流開發(fā)保護國家重點實驗室，四川 成都 610065）

1 研究背景

流域水污染防控中，點源污染的防治得到了充分的重視，并取得了一定的治理效果。但是，非點源污染引起的水體污染等問題日益突出，如湖庫水體富營養(yǎng)化、水華和水藻爆發(fā)都與流域內(nèi)的非點源污染有密切關(guān)系。由于非點源污染物的成分復雜、類型多樣、排放分散、量大面廣，定量分析較為困難［1-2］。如何綜合考慮點源和非點源污染負荷進行總量控制與分配，是國內(nèi)外研究的熱點。

在污染總量控制與分配方面，負荷歷時曲線（LDC）由于其簡單易用，對資料要求較低，并且可以對點源和非點源污染同時分析，得到廣泛應用，尤其在資料缺乏地區(qū)。如程艷等［3］以洱海彌苴河流域作為典型研究區(qū)域利用LDC進行水質(zhì)分析，計算出超標負荷的總削減率，丁京濤［4］利用LDC進行了大寧河巫溪段水體總磷估算與分配研究。但LDC方法研究大都應用在整個流域，沒有細致計算流域內(nèi)部單元的污染物分配情況。

贛江袁河流域流量、水質(zhì)監(jiān)測資料都較少，進行污染模擬及總量控制與分配較為困難。本文以該流域為例，在進行流域控制單元劃分的基礎(chǔ)上，利用SWAT模型的模擬獲得各控制單元的流量，然后利用LDC方法進行資料缺乏流域的污染總量控制與分配方案研究。本文通過對污染負荷在不同水文條件下、不同控制單元和不同水質(zhì)目標下的點源和非點源之間的分配的研究，構(gòu)建袁河流域分區(qū)（控制單元）、分期（豐水期、平水期、枯水期）、分類（點源和非點源）和分級（不同水質(zhì)目標）的容量總量方案。

2 方法簡介

2.1 SWAT模型SWAT（Soil and Water Assessment Tool）模型，由美國農(nóng)業(yè)部農(nóng)業(yè)研究中心（US?DA-ARS）開發(fā)，是在SWRRB模型基礎(chǔ)上發(fā)展起來的一個大、中尺度的連續(xù)的流域環(huán)境模擬模型。該模型具有很強的水文物理機制，可以對不同的土壤條件、土地利用類型、氣候狀況和在人類活動的干擾下做出有效的徑流、水質(zhì)模擬預測［5］，SWAT模型在很多流域得到了成功應用［6-10］。

SWAT模型將研究流域劃分為若干個水文響應單元（Hydrologic Response Units，HRU），以HRU為最小的水文模擬單元進行模擬［11］。每一個HRU內(nèi)的水量平衡是基于降水、蒸散、滲透、地表徑流、壤中流、地下水回流和河道運移損失來計算。模型可以模擬流域內(nèi)部的多種水文循環(huán)物理過程：水的運動、泥沙的輸移、植物的生長以及營養(yǎng)物質(zhì)的遷移轉(zhuǎn)化等。模型的整個模擬過程分為兩大部分：子流域模塊（產(chǎn)流和坡面匯流部分）和匯流演算模塊（河道匯流部分）。前者控制著每個子流域內(nèi)主河道的水、沙、營養(yǎng)物質(zhì)和化學物質(zhì)的輸入量，后者決定了物質(zhì)從河網(wǎng)向流域出口的輸移運動及負荷的演算匯總過程。

2.2 LDC模型負荷歷時曲線法是Searcy在1959年首次提出，近年來在美國污染總量控制中得到廣泛應用［12-13］，污染負荷歷時曲線法只考慮水體的稀釋過程，無需水力學數(shù)據(jù)（流速、斷面形狀、坡度、糙率、離散系數(shù)等），而只需要水文數(shù)據(jù)（流量或水位）。對于單一控制斷面（有水文水質(zhì)同步信息的斷面）的容量分析，其容量分布特征可以直接從負荷歷時曲線考察得出。將流量系列乘以特定污染物濃度標準得到負荷通量作為縱坐標，建立的歷時曲線即為該污染物對應流量系列的LDC。其步驟如下。

（1）繪制流域內(nèi)控制站點的日流量歷時曲線。根據(jù)流量資料（無或缺資料地區(qū)可以采用水文學上的面積比或相似流域法獲得），將流量Q從大到小排序就可以得到流量歷時曲線。

（2）將流量曲線轉(zhuǎn)變?yōu)槿肇摵蓺v時曲線。確定控制站所在區(qū)域的水功能區(qū)劃及其對應的水質(zhì)要求，再將各種需要控制的污染物的容許濃度與流量歷時曲線各縱坐標相乘得到允許負荷量，將允許負荷量與相應的保證率建立一一對應的關(guān)系，繪制允許負荷頻率曲線圖，即得到負荷歷時曲線。

（3）總量控制與分配。將實測（或預測的規(guī)劃水平年）的水質(zhì)數(shù)據(jù)轉(zhuǎn)化為日負荷并點繪制到日負荷歷時曲線圖上，比較實測負荷量與允許負荷量，即可知不同季節(jié)情況下是否超標以及相應的削減量。

3 研究區(qū)域及數(shù)據(jù)

袁河屬贛江水系一級支流，發(fā)源于江西省萍鄉(xiāng)縣羅霄山脈武功山七盤嶺金頂峰北麓，自西南向東北流經(jīng)萍鄉(xiāng)、宜春、分宜、新余以及樟樹部分地區(qū)，至樟樹鎮(zhèn)金鳳洲注入贛江（圖1）。該河全長279km，河床比降0.585‰，集水面積6 262km2。多年平均降水量為1 602.6mm，多年平均氣溫為17.3℃，多年平均蒸發(fā)量為1 282.9mm。流域內(nèi)植被良好，土壤大多為紅壤。

所收集到的數(shù)據(jù)包括DEM、土壤類型、土地利用類型、降雨、氣象、徑流、泥沙、水質(zhì)等。本文中袁河DEM為30m×30m；土壤類型圖為1∶20萬；土地利用類型圖為1∶100萬；降雨數(shù)據(jù)（1999—2008年）來自袁河流域上16個雨量站（圖1），氣象數(shù)據(jù)（1999—2008年）來源于袁河流域上4個氣象站（圖1）。

徑流數(shù)據(jù)來自袁河上中游的蘆溪站和茅洲2個水文站，其中蘆溪徑流數(shù)據(jù)從1978—2008年，而茅洲數(shù)據(jù)為1978—1994年。考慮到數(shù)據(jù)長度的一致性，本文中模型率定階段為1999—2008年，因此需要將茅洲數(shù)據(jù)進行延長。本文通過建立茅洲月平均徑流和蘆溪月平均徑流的線性相關(guān)方程（R2=0.91），將茅洲徑流延長到2008年。

袁河上的非點源污染來源主要考慮農(nóng)業(yè)化肥、農(nóng)村生活污水以及禽畜養(yǎng)殖。點源污染中主要考慮城市生活污水和工業(yè)污水，袁河沿干支流分布有62個排污口。袁河水質(zhì)監(jiān)測結(jié)果表明，袁河主要以Ⅲ類水為主，污染較為嚴重的區(qū)域主要在新余工業(yè)區(qū)。

4 SWAT模擬結(jié)果分析

根據(jù)研究區(qū)自然情況以及研究需要，本文將袁河流域共劃分為38個子流域（見圖1）。采用試錯法對SWAT模型進行參數(shù)率定，即通過比較實測與模擬的徑流擬合情況來調(diào)整相關(guān)模型參數(shù)。由于袁河數(shù)據(jù)較少，本文將所獲得數(shù)據(jù)都用來進行參數(shù)率定。通過對參數(shù)的敏感性分析以及參照其他文獻確定需要率定的參數(shù)，主要徑流參數(shù)包括：CN2（SCS曲線數(shù)）、ESCO（土壤蒸發(fā)補償系數(shù)）、AWC（土壤可供水量）、SLOP（子流域坡度）等。首先利用蘆溪和茅洲（延長后的）的徑流數(shù)據(jù)進行徑流參數(shù)率定，然后輸出各控制單元的流量。

本文采用兩個評價指標：線性相關(guān)系數(shù)R2和納西效率系數(shù)NSE，對模擬結(jié)果進行分析。R2和NSE的計算公式如下：

式中：Xobs和Xsim分別為實測值和模擬值；分別為實測值和模擬值的算術(shù)平均值。

R2和NSE越接近于1表示實測與模擬值擬合越好，模擬精度越高。

圖2和圖3為蘆溪站和茅洲站的徑流的模擬結(jié)果。其中蘆溪徑流模擬的R2=0.92，NSE=0.78，茅洲徑流模擬的R2=0.89，NSE=0.73（表1）。一般認為，當R2＞0.6和NSE＞0.5時，模型模擬結(jié)果可以接受［12］。因此可以看出，SWAT對兩站的徑流模擬效果較好。此外，蘆溪站徑流模擬值比實測值略低，茅洲站徑流模擬值比實測值略高。但二者的模擬精度都遠遠超過了月徑流模擬的標準。

5 總量分配結(jié)果

5.1 控制單元劃分為了更有針對性地制定水環(huán)境保護方案，本文首先對研究流域進行分區(qū)，然后進行容量總量分配。分區(qū)原則主要是基于流域水環(huán)境生態(tài)系統(tǒng)的特征差異，即結(jié)合水文過程的集水范圍、行政分區(qū)、水功能分區(qū)以及劃分后單元便于管理統(tǒng)計等，體現(xiàn)流域水文特性、集水特性以及水功能保護。本文最終將袁河流域范圍劃分為4個一級控制單元（蘆溪、宜春、新余、樟樹），每個一級控制單元下又劃分為2～3個二級控制單元，即整個袁河流域共分為10個控制單元，從上游到下游依次為蘆溪保護區(qū)、蘆溪工業(yè)區(qū)、蘆溪保流區(qū)、宜春保留區(qū)1、宜春工業(yè)區(qū)、宜春保留區(qū)2、新余保留區(qū)1、新余工業(yè)區(qū)、新余保留區(qū)2和樟樹保留區(qū)，見圖4。

5.2 LDC總量分配結(jié)果根據(jù)水功能區(qū)劃可知各控制單元的水質(zhì)目標，其中，保護區(qū)水質(zhì)目標為Ⅱ類，保留區(qū)水質(zhì)目標為Ⅲ類，工業(yè)區(qū)水質(zhì)目標為Ⅳ類。通過比較污染負荷現(xiàn)狀值和水質(zhì)目標值即判斷污染物的超標情況。

流域內(nèi)控制單元內(nèi)均設(shè)有水質(zhì)監(jiān)測站，相應站點的流量由SWAT模型模擬獲得，由此可建立其相應的負荷歷時曲線。LDC計算結(jié)果如圖5所示（由于受篇幅限制只列出所有控制單元的氨氮及超標控制單元的總氮圖）。

從各控制單元的負荷歷時曲線上可以看出，明顯超標的有新余工業(yè)區(qū)和新余保留區(qū)2的總氮以及新余工業(yè)區(qū)的氨氮。為了判斷最關(guān)鍵水質(zhì)條件，采用上述代表值分別估算不同流量與實測水質(zhì)樣本的日負荷通量系列第90百分位數(shù)的濃度值相乘，與標準水質(zhì)對比得到削減率。

此外，LDC法可反映污染時的水文條件，其中，特豐條件：0～10%（歷時范圍）；豐水條件：10%～40%；平水條件：40%～60%；枯水條件：60%～90%；特枯條件：90%～100%。從圖中查看負荷超標發(fā)生的季節(jié)，即可判斷污染類型：若發(fā)生在豐水期，則表示主要為非點源污染；若水質(zhì)超標發(fā)生在枯水期，則表示主要為點源污染；若發(fā)生在平水期，則為點源與非點源的綜合。

由已計算出的LDC圖，按照上述算法算得各控制區(qū)各項指標削減率并分析其污染類型，見表1。從圖5和表1可以看出，新余工業(yè)區(qū)內(nèi)枯水期氨氮負荷超標率約為76%，豐水期超標率約為47%。由于枯水期基本是點源污染，豐水期以非點源為主，同時包括點源污染。因此說明新余工業(yè)區(qū)氨氮超標以點源為主；同樣，新余保留區(qū)枯水期總氮超標約68%，豐水期超標約53%，較新余工業(yè)區(qū)的非點源污染比重有所上升，但仍以點源為主。

在進行削減分配計算中，視枯水期污染超標全部為點源污染造成，豐水期污染為點源污染與非點源污染的綜合，即豐水期污染超標總量扣除枯水期污染超標量即為非點源污染分配量。同時應考慮到，工業(yè)區(qū)位于保留區(qū)上游，當工業(yè)區(qū)點源削減后，保留區(qū)點源削減量應當扣除工業(yè)區(qū)已經(jīng)削減的數(shù)值，若扣除后削減值小于等于零，則點源污染經(jīng)由上游削減后，下游可不用再削減；若扣除后大于零，則點源污染還需削減。此外，同一個控制分區(qū)內(nèi)，如有不同污染負荷削減，取其中最大值作為此控制分區(qū)的最終削減值，具體計算結(jié)果見表2。

從表2可以看出，對于新余工業(yè)區(qū)，豐水期和枯水期氨氮超標，其中點源削減75.66%，非點源削減0%；同時豐水期和枯水期TN超標，其中點源削減45.63%，非點源削減6.75%；取上述氨氮和TN削減比例的最大值，則新余工業(yè)區(qū)最終點源削減75.66%，非點源削減量為6.75%。對于新余保留區(qū)2，除去上游新余工業(yè)區(qū)的削減量后，新余保留區(qū)2的TN還應削減：點源削減35.22%，非點源削減19.46%?？傮w來看，袁河新余區(qū)點源超標嚴重，非點源略有超標。

表1 LDC結(jié)果分析

表2 點源和非點源分配結(jié)果

以上分析表明，袁河流域污染物超標河段主要在新余工業(yè)區(qū)及其下游的新余保留區(qū)，超標污染物主要是點源污染，因此袁河流域應嚴格控制流域內(nèi)點源的排放，同時對非點源污染進行一定的防治。

6 結(jié)論

本文利用分布式SWAT模型和污染負荷歷時曲線法，構(gòu)建了資料缺乏的贛江袁河流域非點源污染模擬與總量控制與分配方案，實現(xiàn)了對贛江袁河流域的主要污染物的模擬以及基于“分區(qū)、分期、分類、分級”思想的水污染物總量分配。研究結(jié)果表明該方法在無資料流域應用的可行性和適用性。本文所用技術(shù)和方法對其它資料缺乏流域的非點源污染模擬及總量分配也具有一定的借鑒意義和參考價值。

［1］郝芳華，程紅光，楊勝天.非點源污染模型——理論方法與應用［M］.北京：中國環(huán)境科學出版社，2006.

［2］李海杰.吉林省雙陽水庫匯水區(qū)農(nóng)業(yè)非點源污染研究［D］.長春：吉林大學，2007.

［3］程艷，李炳花，和壽芳，等.水污染LDC在洱海彌苴河流域水質(zhì)分析中的應用［J］.水電能源科學，2008，8：27-30.

［4］丁京濤.大寧河巫溪段水體總磷估算及分配研究［D］.北京：北京化工大學，2009.

［5］Arnold J G，Allen P M，Muttiah R，et al.Automated base-flow separation and recession analysis techniques［J］.Ground Water，1995，33（6）：1010-1018.

［6］Bosch N.The influence of impoundments on riverine nutrient transport：An evaluation using the Soil and Water Assessment Tool［J］.Journal of Hydrology，2008，355：131-147.

［7］Grizzetti B，Bouraoui F，Granlund K，et al.Modelling diffuse emission and retention of nutrients in the Vantaan?joki watershed（Finland）using the SWAT model［J］.Ecological Modelling，2003，169（1）：25-38.

［8］Santhi C，Arnold J G，William J R，et al.Validation of the SWAT model on a large river basin with point and nonpoint sources［J］.Journal of the American Water Resources Association，2001，37（5）：1169-1188.

［9］Saleh A，Arnold J G，Gassman P W，et al.Application of SWAT for the Upper North Bosque River Watershed［J］.Transactions of the ASAE，2000，43（5）：1077-1087.

［10］Neitsch S L，Arnold J G，Kiniry J R，et al.Soil and Water Assessment Tool Theoretical Documentation：Version 2000［Z］.Agricultural Research Service，Temple，Texas，2002.

［11］龐靖鵬.非點源污染分布式模擬--以密云水庫水源地保護為例［D］.北京：北京師范大學，2007.

［12］An Approach for Using Load Duration Curves in the Development of TMDLs［Z］.EPA 841-B-07-006，2007.

［13］Singh K P，Malik A，et al.Multivariate statistical techniques for the evaluation of spatial and temporal variations in water quality of Gomti River（India）——a case study［J］.Water Research，2004，38：3980-3992．