洋河水庫(kù)流域非點(diǎn)源污染遷移特性及流域數(shù)學(xué)模型的構(gòu)建

李大鳴,栗琪程，卜世龍，張弘強(qiáng)，李彥卿,顧利軍，姚志帆，傅長(zhǎng)鋒

(1 天津大學(xué) 水利工程仿真與安全國(guó)家重點(diǎn)實(shí)驗(yàn)室，天津 300350；2 中國(guó)市政工程西北設(shè)計(jì)研究院有限公司，甘肅 蘭州 730000；3中興長(zhǎng)天信息技術(shù)有限公司，江西 南昌330000；4河北省水利水電勘測(cè)設(shè)計(jì)研究院，天津 300250)

據(jù)報(bào)道，非點(diǎn)源污染目前已成為水環(huán)境的主要污染源，其中有60%左右的河流和50%左右的湖泊污染與非點(diǎn)源污染有關(guān)[1]。我國(guó)密云水庫(kù)、洱海等水域非點(diǎn)源污染比例已超過(guò)點(diǎn)源污染，每年流入洱海的非點(diǎn)源污染部分TN、TP污染比例分別占了97.1%和92.5%[2-3]。目前，中國(guó)的農(nóng)業(yè)非點(diǎn)源污染形勢(shì)十分嚴(yán)峻，全國(guó)化肥投入總量從1990年的2 590.3 萬(wàn)t增加到2008年的5 239.0萬(wàn)t，化肥投入量也從1990年的270.75 kg/hm2增加到2008年的430.43 kg/hm2，遠(yuǎn)遠(yuǎn)超過(guò)了發(fā)達(dá)國(guó)家為防止過(guò)度使用化肥對(duì)水體造成污染所設(shè)定的安全上限[4]。

SWAT模型在非點(diǎn)源污染控制評(píng)價(jià)方面的有效性，在國(guó)外的諸多流域[5-9]已經(jīng)得到了驗(yàn)證。而在國(guó)內(nèi)，學(xué)者們利用長(zhǎng)期實(shí)測(cè)資料對(duì)SWAT模型進(jìn)行率定與驗(yàn)證，并在密云水庫(kù)、大寧河流域、沙河水庫(kù)、南四湖流域、阿什河流域[10-15]等非點(diǎn)源污染的研究中得到應(yīng)用，且已取得了良好的效果，證明了SWAT模型在我國(guó)的適用性。非點(diǎn)源污染的遷移規(guī)律受到地區(qū)、降雨、氣候、植被等多種因素的影響，張展羽等[16]利用SWAT模型研究了不同灌溉方式下小流域非點(diǎn)源污染的氮磷遷移規(guī)律。Qi等[17]開(kāi)發(fā)了一種基于物理特性的土壤溫度模塊，用來(lái)替代 SWAT中使用的經(jīng)驗(yàn)土壤溫度模塊，進(jìn)一步增強(qiáng)了SWAT模型的適用性；郭藝[18]研究了流域非點(diǎn)源磷污染輸出系數(shù)的時(shí)空變異性，證明了輸出系數(shù)模型對(duì)于非點(diǎn)源污染研究的有效性。

洋河水庫(kù)作為河北秦皇島市主要的城市供水源之一，2002-2013年其水質(zhì)整體標(biāo)準(zhǔn)為Ⅳ-Ⅴ類[19-21]。其中，秦皇島水務(wù)局2013年9月的水質(zhì)監(jiān)測(cè)結(jié)果顯示，根據(jù)Ⅲ類水質(zhì)標(biāo)準(zhǔn)要求，洋河水庫(kù)流域TN整體超標(biāo)2.5～22.5倍，TP整體超標(biāo)1.6～2.9倍，非點(diǎn)源氮、磷為流域內(nèi)主要污染物，為研究洋河水庫(kù)流域中非點(diǎn)源氮、磷的污染特性，本研究利用SWAT模型對(duì)洋河水庫(kù)流域非點(diǎn)源污染進(jìn)行模擬計(jì)算，從空間角度分析該地區(qū)非點(diǎn)源氮、磷的分布特點(diǎn)及遷移特性，同時(shí)為簡(jiǎn)化模型計(jì)算，本研究借助SWAT模型中部分模塊功能，利用Fortran語(yǔ)言建立流域的污染物遷移擴(kuò)散數(shù)學(xué)模型，模擬并驗(yàn)證2015年流域水質(zhì)狀況，以期快速高效地完成庫(kù)區(qū)污染狀況的模擬，為洋河水庫(kù)非點(diǎn)源污染入庫(kù)前的治理提供參考。

1 研究區(qū)域概況

洋河水庫(kù)于1959年10月修建于秦皇島市撫寧縣大灣子村北，控制流域面積755 km2，流域共有4條主要水系，分別為東洋河水系、迷霧河水系、麻姑營(yíng)河水系以及西洋河水系。多年平均降雨量750 mm。流域降水量年內(nèi)分配不均，夏季降水量占全年降水量的75%以上，并多以暴雨形式出現(xiàn)，容易造成洪澇災(zāi)害。

1.1 非點(diǎn)源污染負(fù)荷統(tǒng)計(jì)與分析

將污染源按形成特點(diǎn)劃分為生活源污染、農(nóng)業(yè)污染、分散式禽畜養(yǎng)殖污染、城鎮(zhèn)地表徑流攜帶物污染、水土流失污染5類，各類污染負(fù)荷估算系數(shù)根據(jù)文獻(xiàn)[22]及《河北省水資源質(zhì)量評(píng)價(jià)》(2004年)確定，各污染源年均負(fù)荷產(chǎn)生量的統(tǒng)計(jì)結(jié)果如表1所示。

1.2 土地利用類型及土壤數(shù)據(jù)庫(kù)的建立

SWAT模型中所用到的土地利用類型數(shù)據(jù)來(lái)源于河北省秦皇島市1∶5萬(wàn)空間矢量數(shù)據(jù)土壤覆蓋圖層，土地利用類型包括耕地、林地、草地、水域、居民用地和未利用土地6個(gè)一級(jí)類型以及25個(gè)二級(jí)類型。通過(guò)ArcGIS分析模塊，結(jié)合SWAT所需土地利用類型數(shù)據(jù)結(jié)構(gòu)，對(duì)洋河水庫(kù)流域土地利用類型重新分類，共分為6類，分別為草地、園地、林地、農(nóng)村居民點(diǎn)、耕地、水域，具體土地利用類型對(duì)應(yīng)編碼及其空間分布如圖1所示。

表1 洋河水庫(kù)流域各污染源年均負(fù)荷產(chǎn)生量

土地利用類型數(shù)據(jù)存儲(chǔ)格式為shp格式。從圖1可以看出，研究區(qū)域主要土地利用類型為耕地，占總面積的44.05%。

土壤是降雨轉(zhuǎn)化為徑流的介質(zhì)基礎(chǔ)，通過(guò)土壤，降雨進(jìn)行下滲、側(cè)向流動(dòng)、表面蒸發(fā)等多次分配過(guò)程。本研究的土壤數(shù)據(jù)來(lái)源于FAO構(gòu)建的世界和諧土壤數(shù)據(jù)庫(kù)，土壤類型分布如圖2所示。

由圖2可知，洋河水庫(kù)流域主要土壤類型為雛形土和高活性淋溶土。土壤屬性數(shù)據(jù)庫(kù)中部分物理參數(shù)可直接從FAO-90分類中查得，水分參數(shù)可借助SPAW軟件的Soil Water Characteristics 模塊計(jì)算得到，化學(xué)屬性參數(shù)則采用模型默認(rèn)值[23]。

1.3 氣象數(shù)據(jù)的來(lái)源

水文過(guò)程是確定非點(diǎn)源污染數(shù)學(xué)模型的基礎(chǔ)。本研究采用“中國(guó)地面國(guó)際交換站氣候資料日值數(shù)據(jù)集”省、市站點(diǎn)數(shù)據(jù)，站點(diǎn)分別為河北省懷來(lái)、承德、樂(lè)亭3個(gè)基本站點(diǎn)和天津基本站點(diǎn)。降雨及相關(guān)水文數(shù)據(jù)來(lái)源于洋河水庫(kù)上游流域12個(gè)雨量站2006-2014年的日降雨監(jiān)測(cè)結(jié)果以及1個(gè)同期水文站(洋河水庫(kù)站)實(shí)測(cè)結(jié)果，水文站和雨量站的分布如圖3所示。

2 SWAT模型構(gòu)建及模擬結(jié)果分析

2.1 模型構(gòu)建

2.1.1 子流域劃分及水文響應(yīng)單元生成 基于DEM數(shù)據(jù)，運(yùn)用軟件ArcSWAT自動(dòng)劃分子流域。在劃分子流域時(shí)，事先導(dǎo)入了矢量格式的數(shù)字河網(wǎng)文件，反復(fù)嘗試不同子流域面積閾值，以保證基于DEM生成的河網(wǎng)與真實(shí)河網(wǎng)保持一致。在生成河網(wǎng)后，根據(jù)實(shí)際情況與研究需要對(duì)河網(wǎng)節(jié)點(diǎn)采取人工添加、刪除、修改等操作。最終劃分成35個(gè)子流域，其中 15#、18#、28#子流域直接入庫(kù)(圖4)。

圖3 洋河水庫(kù)流域水文站、雨量站分布

圖4 洋河水庫(kù)流域子流域的劃分圖

在確定水文響應(yīng)單元(hydrologic response unit,HRU)時(shí)，為確保精度，本研究采用每個(gè)子流域設(shè)置多個(gè)HRU。由于研究區(qū)域內(nèi)的主要土地類型是耕地，占總面積達(dá)44.05%；其次是林地和草地，分別占28.19%，18.43%，為清除次要土地利用類型，最終確定土地利用閾值為20%；同理可得土壤閾值為10%，坡度閾值為20%，生成HRU，共生成240個(gè)HRU。

2.1.2 模型參數(shù)敏感性分析 敏感性分析及參數(shù)率定采用SWAT-CUP軟件中SUFI-2算法。該算法中敏感分析方法采用全局敏感性分析。在全局敏感性分析中，參數(shù)敏感性取決于多元回歸系統(tǒng)，其對(duì)拉丁超立方生成參數(shù)與目標(biāo)函數(shù)值進(jìn)行回歸，有：

(1)

式中：g為預(yù)測(cè)值，α為不確定性隨機(jī)殘差，m為迭代次數(shù)，βi為回歸系數(shù)，bi為參數(shù)變量。

根據(jù)所要校核的目標(biāo)對(duì)象，選取適當(dāng)?shù)拿舾行苑治鰠?shù)，模型參數(shù)與敏感性分析的參數(shù)及含義如表2所示。

表2 敏感性分析參數(shù)及其含義

其中，對(duì)于徑流模擬而言，CN2、SOL_AWC、ESCO 3個(gè)參數(shù)最為敏感；而對(duì)于N、P等污染物的模擬，CN2、ESCO、ALPHA_BF、SURLAG、NPERCO、BIOMIX、RSDCO等參數(shù)最為敏感。因此,分別以28#和24#子流域?yàn)槔?，選取以上8個(gè)參數(shù)作為敏感性參數(shù)對(duì)徑流和總氮(TN)進(jìn)行分析驗(yàn)證，敏感性分析結(jié)果如表3,4所示 。

表3 洋河水庫(kù)流域28#子流域徑流敏感性分析結(jié)果

表4 洋河水庫(kù)流域24#子流域總氮(TN)敏感性分析結(jié)果

表3和表4中，t值表示參數(shù)相對(duì)顯著性，t的檢驗(yàn)概率值P體現(xiàn)t統(tǒng)計(jì)量的顯著性，其中t的絕對(duì)值越大，敏感性越高；P值越接近0，顯著性越大。由表3和表4可知，對(duì)于徑流而言，基流分割系數(shù)ALPHA_BF以及正常濕潤(rùn)情況下徑流曲線值CN2最為敏感；對(duì)TN而言，CN2最為敏感。

2.1.3 模型率定及驗(yàn)證 SUFI-2算法是通過(guò)拉丁超立方體隨機(jī)采樣法隨機(jī)生成一組參數(shù)代入SWAT中進(jìn)行目標(biāo)函數(shù)的計(jì)算方法，驗(yàn)證模型的可靠性。評(píng)價(jià)SWAT模型模擬適用性時(shí)，一般采用目標(biāo)函數(shù)R2與Nash-Sutcliffe效率系數(shù)(Ens)，且判別依據(jù)是Ens≥0.50，R2≥0.6視為滿意。為使包括土壤含水率等初始值在內(nèi)的一系列參數(shù)的初始值不為0，模型設(shè)置2007年為模型預(yù)熱期，模型校準(zhǔn)期為2008-2012年，率定時(shí)間步長(zhǎng)為月，洋河水庫(kù)流域流量校準(zhǔn)結(jié)果如圖5所示。

圖5 洋河水庫(kù)流域流量校準(zhǔn)過(guò)程(2008-2012年)

利用校準(zhǔn)好的參數(shù)修改模型初始值，采用實(shí)測(cè)的2013-2014年各月流量統(tǒng)計(jì)值進(jìn)行驗(yàn)證，驗(yàn)證結(jié)果如圖6所示。

模型校準(zhǔn)結(jié)果顯示，R2=0.96，Ens=0.95；模型驗(yàn)證結(jié)果顯示,R2=0.98，Ens=0.98。可知模擬值與實(shí)測(cè)值擬合效果較好，表明所建模型在一定程度上能很好地模擬洋河水庫(kù)流域水文過(guò)程，且模型模擬水文過(guò)程具有較高的可靠性。洋河水庫(kù)上游沒(méi)有較大規(guī)模的集中排污點(diǎn)，主要污染源為非點(diǎn)源污染，故模型構(gòu)建時(shí)未設(shè)置點(diǎn)源排放，水質(zhì)模擬結(jié)果全部來(lái)源于非點(diǎn)源污染。受水質(zhì)資料限制，本研究根據(jù)東洋河入口處2013年TN入庫(kù)量的月監(jiān)測(cè)數(shù)據(jù)進(jìn)行率定及驗(yàn)證，結(jié)果如圖7所示。

模型率定結(jié)果顯示，R2=0.80，Ens=0.79，滿足R2≥0.50，Ens≥0.60的適用性模擬要求，污染負(fù)荷模擬總體趨勢(shì)基本一致，SWAT模型模擬效果較為滿意。

圖6 洋河水庫(kù)流域流量驗(yàn)證過(guò)程(2013-2014年)

圖7 2013年?yáng)|洋河入口處TN入庫(kù)量的率定及驗(yàn)證結(jié)果

2.2 模擬結(jié)果分析

通過(guò)對(duì)2008-2014年各年的洋河水庫(kù)流域徑流、TN和TP負(fù)荷的模擬計(jì)算，得出2008-2014年洋河水庫(kù)流域35個(gè)子流域的年均非點(diǎn)源負(fù)荷總量：TN為845.96 t，流失率10.6%，其中有機(jī)氮772.5 t，對(duì)TN貢獻(xiàn)率91%；TP為240.56 t，流失率6.4%，其中有機(jī)磷95.8 t，對(duì)TP貢獻(xiàn)40%，可見(jiàn)有機(jī)態(tài)污染物為流域內(nèi)主要非點(diǎn)源污染物；流域年均產(chǎn)沙量為865 357.40 t。由圖8和圖9可以看出，有機(jī)氮與有機(jī)磷負(fù)荷量在空間分布上與產(chǎn)沙量基本保持一致，且西洋河水系有機(jī)氮與有機(jī)磷負(fù)荷量明顯大于東洋河水系。其中33#、16#子流域分別為有機(jī)氮、有機(jī)磷負(fù)荷量最大單元，達(dá)到66.87和8.27 t；3#子流域有機(jī)氮、有機(jī)磷負(fù)荷量均最小，分別為0.126和0.015 t。這主要是因?yàn)槲綉B(tài)污染負(fù)荷的運(yùn)移主要依靠土壤流失進(jìn)入河道，西洋河水系人口密度遠(yuǎn)大于東洋河水系，相應(yīng)的生活污水、畜禽養(yǎng)殖、淀粉加工業(yè)等污染源大于東洋河水系。從土壤植被覆蓋情況來(lái)看，西洋河水系植被類型多為耕地或裸露地，化肥、農(nóng)藥施用量較大，土壤流失較為嚴(yán)重；而東洋河水系植被覆蓋較好，多為林地、草地等，是良好的水土保持天然屏障。

圖8 2008-2014年洋河水庫(kù)流域各子流域非點(diǎn)源有機(jī)污染物負(fù)荷量多年平均分布情況

圖9 2008-2014年洋河水庫(kù)流域各子流域產(chǎn)沙量、降雨量多年平均分布情況

此外，流域入庫(kù)口處為污染負(fù)荷高風(fēng)險(xiǎn)區(qū)，這主要與非點(diǎn)源污染的富集特性有關(guān)，有機(jī)氮、有機(jī)磷等污染物受降雨徑流過(guò)程影響，污染物逐漸由上游向下游匯聚，最終造成庫(kù)區(qū)水體的嚴(yán)重污染。

3 污染物遷移擴(kuò)散數(shù)學(xué)模型的構(gòu)建

農(nóng)業(yè)非點(diǎn)源污染是造成我國(guó)流域污染的主要原因之一，而數(shù)學(xué)模型是流域污染研究、水質(zhì)分析與預(yù)測(cè)的主要技術(shù)手段[24]。SWAT模型結(jié)構(gòu)復(fù)雜，模型自身包含的數(shù)據(jù)庫(kù)是以北美洲水文環(huán)境為基礎(chǔ)進(jìn)行構(gòu)建的，在國(guó)內(nèi)應(yīng)用時(shí)，各種水文水質(zhì)資料要求很高，且模型構(gòu)造較為繁瑣，不易達(dá)到模擬精度要求[25]，因此本研究借助SWAT模型中部分模塊功能，利用Fortran語(yǔ)言構(gòu)建洋河水庫(kù)流域的污染物遷移擴(kuò)散數(shù)學(xué)模型。模型主要結(jié)構(gòu)如圖10所示。

圖10 洋河水庫(kù)流域污染物遷移擴(kuò)散數(shù)學(xué)模型結(jié)構(gòu)示意圖

3.1 流域計(jì)算單元?jiǎng)澐帜K的構(gòu)建

模型首先借助SWAT模型中子流域自動(dòng)劃分功能對(duì)洋河區(qū)域進(jìn)行計(jì)算單元?jiǎng)澐?，由于模型中缺乏?duì)水文響應(yīng)單元的相關(guān)定義，為滿足模型的精度要求，本研究規(guī)定單元?jiǎng)澐衷瓌t為每個(gè)單元包含唯一一條主要分支河道，河道分流、匯流點(diǎn)位于單元邊界處，單元內(nèi)匯流與單元間匯流由地形條件及上下游條件確定。基于此，本研究將單元?jiǎng)澐譃閮深?，一類是源流單元，一類為匯流單元，共劃分單元237個(gè)，單元?jiǎng)澐旨跋鄳?yīng)匯流關(guān)系見(jiàn)圖11。

圖11 洋河水庫(kù)流域單元?jiǎng)澐旨皡R流級(jí)別

借助Fortran語(yǔ)言自編程序?qū)τ?jì)算單元進(jìn)行等級(jí)劃分，將無(wú)匯入計(jì)算單元作為1級(jí)單元，規(guī)定N個(gè)1級(jí)單元匯入的計(jì)算單元定義為2級(jí)單元，以此類推，最終全流域得到的最高計(jì)算單元等級(jí)為20級(jí)。

3.2 污染物分布模塊的構(gòu)建

模型首先統(tǒng)計(jì)了洋河水庫(kù)流域各村莊污染物的產(chǎn)生量，其中對(duì)于生活污水以及禽畜排泄物等排放點(diǎn)相對(duì)穩(wěn)定的污染源，采用直接疊加法將污染源分配到計(jì)算單元；而對(duì)于化肥農(nóng)藥及固體廢棄物等污染源需要根據(jù)其影響面積加以相應(yīng)分配。由于缺乏具體的村莊范圍數(shù)據(jù)，本研究利用ArcGIS對(duì)導(dǎo)入的各個(gè)村莊制作泰森多邊形，作為村莊的模擬范圍，根據(jù)各計(jì)算單元中村莊模擬范圍所占面積的比值分配相應(yīng)的污染物，進(jìn)而完成污染物分布由村莊向計(jì)算單元的轉(zhuǎn)化，分配后的各水系污染物負(fù)荷量面密度的計(jì)算結(jié)果如表5所示。

表5 洋河水庫(kù)流域各水系污染物負(fù)荷量面密度統(tǒng)計(jì)結(jié)果

表5(續(xù)) Continued table 5

3.3 降雨匯流模塊的構(gòu)建

降雨匯流的計(jì)算模式稱為擬序匯流計(jì)算方法，規(guī)定由源流計(jì)算單元向匯流計(jì)算單元中匯集，由低級(jí)別計(jì)算單元向高級(jí)別計(jì)算單元匯集，匯流需滿足水量平衡條件，計(jì)算公式為：

(2)

(3)

單元產(chǎn)流系數(shù)估算公式為：

(4)

3.4 污染物濃度遷移擴(kuò)散模塊的構(gòu)建

根據(jù)污染物的平衡，可以得到污染物遷移擴(kuò)散公式為：

(5)

(6)

污染物推移系數(shù)的計(jì)算公式為：

(7)

3.5 產(chǎn)流系數(shù)及污染物推移系數(shù)的率定

根據(jù)公式(4)、單元植被占比、阻水系數(shù)以及中心坡度的數(shù)據(jù)，首先暫定各計(jì)算單元產(chǎn)流系數(shù)的率定系數(shù)為1，根據(jù)洋河水庫(kù)流域2013年雨量站數(shù)據(jù)，可以得到洋河水庫(kù)流域4大水系算數(shù)平均面雨量，其中東洋河水系為662.40 mm，西洋河水系為700.46 mm，迷霧河水系為681.30 mm，麻姑營(yíng)河水系為774.73 mm。運(yùn)用公式(2)、(3)可以計(jì)算得到各計(jì)算單元的匯流水量及出流水量，其中4個(gè)入庫(kù)計(jì)算單元的出流水量即為對(duì)應(yīng)的4個(gè)水系的入庫(kù)流量，依據(jù)洋河水庫(kù)2013年入庫(kù)流量實(shí)測(cè)數(shù)據(jù)，對(duì)比4個(gè)入庫(kù)計(jì)算單元的出流水量及入庫(kù)流量的實(shí)測(cè)值，可得到4個(gè)水系的率定系數(shù)Xn(n為水系編碼)，根據(jù)每個(gè)計(jì)算單元所在的水系，便可得到各計(jì)算單元的率定系數(shù)Xim，并以此重新率定產(chǎn)流系數(shù) 。

將污染物推移系數(shù)初值及推移系數(shù)的率定系數(shù)都取為1，根據(jù)公式(7)可以初步計(jì)算得到各計(jì)算單元的污染物推移系數(shù)，再根據(jù)污染分布模塊中得到的2013年各計(jì)算單元污染物元素入河量以及公式(5)和(6)，在不考慮河道工程的情況下，可以計(jì)算得到各計(jì)算單元的污染物元素入流量及污染物元素出流量，其中4個(gè)入庫(kù)計(jì)算單元的污染物出流量，即為對(duì)應(yīng)的4大水系的污染物入庫(kù)量，進(jìn)而得到各污染物的質(zhì)量濃度，對(duì)比4個(gè)入庫(kù)計(jì)算單元的污染物元素出流質(zhì)量濃度及4大水系污染物元素入庫(kù)濃度的實(shí)測(cè)值，便可得到4個(gè)水系的污染物率定系數(shù)Kn，根據(jù)每個(gè)計(jì)算單元所在的水系便可得到各計(jì)算單元的污染物推移系數(shù)的率定系數(shù)Kim，并以此重新率定污染物推移系數(shù)Yim。

3.6 模型計(jì)算結(jié)果及驗(yàn)證

3.6.1 2015年污染物出流量 利用降雨匯流模塊及污染物濃度遷移擴(kuò)散模塊計(jì)算得出2015年各水系入庫(kù)計(jì)算單元出流量，根據(jù)入庫(kù)口計(jì)算單元的出流水量和污染物TN、TP出流量，可計(jì)算得出東洋河水系污染元素入庫(kù)的質(zhì)量濃度，模型計(jì)算得出東洋河水系各單元出口污染物出流量，結(jié)果如圖12和13所示。

3.6.2 模擬結(jié)果驗(yàn)證 經(jīng)污染物遷移擴(kuò)散數(shù)學(xué)模型計(jì)算得出，東洋河入庫(kù)口處TN年污染物入庫(kù)量為478.778 t，TP年污染物入庫(kù)量為1.302 t。將2015年?yáng)|洋河入庫(kù)口的污染物質(zhì)量濃度計(jì)算值與實(shí)測(cè)值進(jìn)行對(duì)比，結(jié)果如表6所示。

表6 2015年?yáng)|洋河入庫(kù)口污染物質(zhì)量濃度實(shí)測(cè)值與計(jì)算值的對(duì)比

由表6可知，模型計(jì)算精度最高的是TN質(zhì)量濃度，相對(duì)誤差低于10%；TP質(zhì)量濃度的計(jì)算精度較低，相對(duì)誤差在20%以上，主要是由于同種土地利用類型磷的輸出系數(shù)具有較大的時(shí)空變異性，對(duì)于輸出系數(shù)模型而言，其不僅與徑流、施肥量有關(guān)，同時(shí)還與離河距離遠(yuǎn)近、土壤類型等因素息息相關(guān)。TP質(zhì)量濃度的相對(duì)誤差為33.30%，而絕對(duì)誤差僅為0.008 mg/L，總體上看，模型模擬結(jié)果較為可靠。

4 結(jié) 論

以洋河水庫(kù)流域?yàn)檠芯繀^(qū)域，通過(guò)建立SWAT模型對(duì)水庫(kù)2008-2014年的徑流和N、P負(fù)荷進(jìn)行初步模擬計(jì)算，然后以Fortran語(yǔ)言為工具，利用SWAT模型中的部分模塊建立洋河水庫(kù)流域污染物遷移擴(kuò)散數(shù)學(xué)模型。得到以下研究結(jié)論：

1)有機(jī)態(tài)污染源為洋河水庫(kù)流域主要污染流失源；在空間分布上，洋河水庫(kù)流域有機(jī)態(tài)污染負(fù)荷呈現(xiàn)“西高東低”的特征，其中有機(jī)氮、有機(jī)磷年均負(fù)荷量最大子流域?yàn)槲餮蠛铀祪?nèi)的33#、16#子流域，分別為66.87和8.27 t；有機(jī)氮、有機(jī)磷年均負(fù)荷量最小子流域均為東洋河水系內(nèi)3#子流域，分別為0.126和0.015 t。受降雨和徑流的影響，入庫(kù)流域區(qū)為污染高風(fēng)險(xiǎn)區(qū)，此外有機(jī)態(tài)污染物負(fù)荷量與泥沙流失量具有空間一致性，可見(jiàn)有機(jī)氮和有機(jī)磷的遷移主要吸附于土壤顆粒，隨水土流失過(guò)程進(jìn)入河道。

2)洋河水庫(kù)流域污染物遷移擴(kuò)散數(shù)學(xué)模型采用擬序算法，根據(jù)水量平衡條件及污染物平衡規(guī)定模擬了洋河水庫(kù)流域徑流及非點(diǎn)源污染遷移狀況，極大地減少了水文水質(zhì)資料的需求，同時(shí)大大提高了模型的運(yùn)算效率。從模擬結(jié)果來(lái)看，TN質(zhì)量濃度的相對(duì)誤差為7.59%，模擬效果較好；TP質(zhì)量濃度的絕對(duì)誤差為0.008 mg/L?？傮w來(lái)看，經(jīng)驗(yàn)?zāi)Ｐ湍M效果符合精度要求。實(shí)際工作中，可以以各計(jì)算單元出口作為流域污染的控制點(diǎn)，通過(guò)輸入不同時(shí)段下的降水條件及相應(yīng)的污染物產(chǎn)生狀況，可以快速得到流域中各個(gè)控制點(diǎn)的污染現(xiàn)狀，這為流域污染治理奠定了基礎(chǔ)。