陳莉蘋(píng) 張舟

（1. 中節(jié)能?chē)?guó)禎環(huán)保科技股份有限公司，安徽合肥 230088；2. 珠江水資源保護(hù)科學(xué)研究所，廣東廣州 510611）

1 引言

科學(xué)合理地計(jì)算水環(huán)境容量及總量分配是流域水環(huán)境達(dá)標(biāo)的前提［1-3］，早期水環(huán)境容量的計(jì)算依靠傳統(tǒng)的試算法，即根據(jù)經(jīng)驗(yàn)調(diào)整污染負(fù)荷分布進(jìn)行反復(fù)試算，直到水域環(huán)境功能區(qū)達(dá)標(biāo)為止，該方法適用于簡(jiǎn)單排污情況，但在排污口眾多的時(shí)候，采用試算法工作量巨大且不是最優(yōu)解［4］。為解決此問(wèn)題，李適宇等［5］于1999 年提出分區(qū)達(dá)標(biāo)法，該方法通過(guò)引進(jìn)貢獻(xiàn)度系數(shù)的概念將二維容量問(wèn)題轉(zhuǎn)化為線性規(guī)劃問(wèn)題，避免了在海域排污口眾多時(shí)，傳統(tǒng)試算法帶來(lái)的巨大工作量，通過(guò)將水質(zhì)模型與線性規(guī)劃模型耦合，可以合理而精確地求取海域水環(huán)境容量。

2010 年鄭洪波等［6］將此方法應(yīng)用在大連金港區(qū)，給出了不同污染物（如COD、無(wú)機(jī)氮等）的海洋環(huán)境容量，建立了計(jì)算海洋環(huán)境容量的普遍方法；2015 年王宇駿等［7］綜合考慮效益、公平與水質(zhì)要求，改進(jìn)了分區(qū)達(dá)標(biāo)法，并將改進(jìn)的模型應(yīng)用于流溪河水環(huán)境容量的計(jì)算。2018 年王程等［8］選取多種因素如人口、GDP、水資源量、工業(yè)產(chǎn)值、工業(yè)廢水排放量、環(huán)保投資等建立優(yōu)化分配模型，提出削減方案，指導(dǎo)蘭江流域水質(zhì)管理。

以上主要采用分區(qū)達(dá)標(biāo)法來(lái)求解水環(huán)境容量，對(duì)于貢獻(xiàn)度系數(shù)的定義偏向于理論假設(shè)，未能同實(shí)際相結(jié)合。本文擬通過(guò)分區(qū)達(dá)標(biāo)法中對(duì)貢獻(xiàn)度系數(shù)的理論描述，運(yùn)用MIKE11 水環(huán)境模型計(jì)算出各污染源貢獻(xiàn)度系數(shù)，再采用MATLAB 擬合，得出貢獻(xiàn)度系數(shù)與污染源濃度、流量以及污染源與考核斷面距離之間的定量化關(guān)系，據(jù)此提出單位負(fù)荷貢獻(xiàn)度系數(shù)的概念及其表達(dá)式，并將其應(yīng)用于新安涌黑臭治理中，量化河涌現(xiàn)狀污染狀況。

該研究可為后續(xù)總量分配提供計(jì)算依據(jù)，便于更準(zhǔn)確貼合實(shí)際的總量分配，為后續(xù)水環(huán)境信息化管理作量化依據(jù)。

2 貢獻(xiàn)度系數(shù)研究

貢獻(xiàn)度系數(shù)是分區(qū)達(dá)標(biāo)法中的關(guān)鍵概念，在分區(qū)達(dá)標(biāo)法中通過(guò)引入貢獻(xiàn)度系數(shù)，可以將容量問(wèn)題轉(zhuǎn)化為線性規(guī)劃問(wèn)題，從而求取最優(yōu)的環(huán)境容量。

貢獻(xiàn)度系數(shù)定義可以表示為［5］：假設(shè)一條河流共有n 個(gè)排污口，m 個(gè)考核斷面，在1 個(gè)排污口給定1 個(gè)單位負(fù)荷量，即xj=1，其他排污口無(wú)負(fù)荷量排出，即xk=0（k=1，…，n；k≠j），通過(guò)模型計(jì)算出這種情況下的濃度分布，確定出m 個(gè)水質(zhì)控制斷面的濃度值，即為此排污口對(duì)考核斷面的貢獻(xiàn)度系數(shù)aij（i=1，…，m），改變排污口，重復(fù)以上步驟，就可以求出每個(gè)排污口的aij。在以上所有計(jì)算中，開(kāi)邊界上的濃度為0，即不考慮背景濃度。

根據(jù)貢獻(xiàn)度系數(shù)定性定義，可以看出影響貢獻(xiàn)度系數(shù)的主要因子為排污口負(fù)荷量以及排污口與水質(zhì)控制斷面的距離，具體來(lái)說(shuō)主要影響因子有排污口流量（Q）、濃度（C）以及排污口與水質(zhì)控制斷面距離（L）。

以下以實(shí)際河涌為研究對(duì)象，通過(guò)MATLAB 擬合出貢獻(xiàn)度系數(shù)與其主要影響因子之間的定量化關(guān)系。

2.1 研究概況

選取中山市新安涌作為研究對(duì)象，河涌長(zhǎng)度為1.1 km，其水文條件受感潮影響不大。根據(jù)污染源調(diào)查，該河涌旱季有水排口共15 個(gè)，其中，分流制污水排口2 個(gè)，分流制雨污混接雨水排口10 個(gè)；旱季無(wú)水排口28 個(gè)。本次研究旨在定量化貢獻(xiàn)度系數(shù)與流量、濃度以及距離之間的關(guān)系，因此選用分流制雨污混接雨水排口，排口分布及編號(hào)如圖1 所示。

圖1 新安涌有水排口分布示意

2.2 研究方法及原理

原理主要包括兩部分，MIKE11 模型計(jì)算原理和MATLAB 擬合方法。

2.2.1 MIKE11 模型計(jì)算原理

MIKE11 水動(dòng)力學(xué)模型主要基于HD 模塊，以圣維南（Saint－Venant）方程組為控制方程，利用Abbott六點(diǎn)隱式格式離散圣維南控制方程組，之后再用河道方程、汊點(diǎn)方程組、外邊界條件等解出隱式方程組，得到水動(dòng)力參數(shù)。MIKE11 水質(zhì)模型以一維對(duì)流擴(kuò)散方程作為控制方程，方程可以表示為：

式（1）中，C 為污染物濃度，mg/L；D 為擴(kuò)散系數(shù)，m2/s；A 為橫截?cái)嗝婷娣e，m2；K 為衰減系數(shù)，s-1；C2為源匯濃度，mg/L；q 為旁側(cè)入流，m3/s；x 為空間坐標(biāo)，m；t為時(shí)間坐標(biāo)，s。

采用時(shí)間和空間中心隱式差分格式對(duì)上述控制方程進(jìn)行離散，再利用水質(zhì)邊界條件、各種水質(zhì)過(guò)程的模擬及水質(zhì)過(guò)程與對(duì)流擴(kuò)散耦合計(jì)算來(lái)求解離散方程。

2.2.2 MATLAB 擬合方法

根據(jù)數(shù)據(jù)特點(diǎn)，本研究選取多項(xiàng)式擬合或者自定義方法，求出變量之間的定量化關(guān)系。

2.2.2.1 多項(xiàng)式擬合

將原始數(shù)據(jù)散點(diǎn)用多項(xiàng)式函數(shù)進(jìn)行擬合，即為多項(xiàng)式擬合。設(shè)由數(shù)據(jù)得到的函數(shù)y=f（x）在n 個(gè)點(diǎn)x1，x2，…，xn的值為y1，y2，…，yn，采用多項(xiàng)式?m（x）來(lái)擬合，設(shè)?m（x）=a0+a1x1+a2x2+…+amxm（m

對(duì)式（2）求解，使得Q（a0，a1，…，am）達(dá)到最小，由此確定多項(xiàng)式系數(shù)。

在MATLAB 中多項(xiàng)式擬合的實(shí)現(xiàn)方式為polyfit函數(shù)或者polyval 函數(shù)，可以對(duì)其直接調(diào)用。

2.2.2.2 自定義法

若數(shù)據(jù)特點(diǎn)不符合多項(xiàng)式規(guī)律，則需要采用自定義方法確定，根據(jù)數(shù)據(jù)特點(diǎn)，函數(shù)可能為指數(shù)函數(shù)，其一般表達(dá)式為f（x）=a×exp（-b×x）+c，通過(guò)調(diào)用MATLAB 中的fit 函數(shù)可以進(jìn)行指定函數(shù)的數(shù)據(jù)擬合，確定系數(shù)a，b，c，從而得出函數(shù)關(guān)系式。

2.3 研究結(jié)果

2.3.1 貢獻(xiàn)度系數(shù)與流量的關(guān)系

設(shè)定濃度和距離不變，選取氨氮為指標(biāo)，研究貢獻(xiàn)度系數(shù)與流量的關(guān)系。具體來(lái)說(shuō)，選定新安涌排口12，給定1 個(gè)單位濃度（1 mg/L），增大排口流量，分析流量對(duì)水質(zhì)及貢獻(xiàn)度系數(shù)的影響。分兩種情況討論。

2.3.1.1 排口氨氮濃度大于背景濃度

當(dāng)排口氨氮濃度大于河涌背景濃度時(shí)，排口會(huì)加重河涌污染。以排口12 為例，河涌背景的氨氮濃度設(shè)置為0 mg/L，計(jì)算出下游考核斷面濃度及貢獻(xiàn)度系數(shù)，以貢獻(xiàn)度系數(shù)為因變量y，流量為自變量x，MATLAB 擬合結(jié)果如圖2 所示。從圖2 可以看出，隨著流量增加，貢獻(xiàn)度系數(shù)逐步增大，增加至數(shù)值0.9 附近后趨于平緩。貢獻(xiàn)度系數(shù)與流量關(guān)系式為：y=-0.792e-0.085x+0.911 3，相關(guān)性系數(shù)R2=0.987 8，可以看出，當(dāng)排口氨氮濃度大于河涌背景濃度時(shí)，排口流量與貢獻(xiàn)度系數(shù)呈指數(shù)關(guān)系，流量增加，貢獻(xiàn)度系數(shù)增大。經(jīng)過(guò)試驗(yàn)研究，流量增加極限為115 m3/s，往后會(huì)影響模型的穩(wěn)定性。

圖2 流量與貢獻(xiàn)度系數(shù)關(guān)系（C 排口＞C 背景）

2.3.1.2 排口氨氮濃度小于背景濃度

當(dāng)排口氨氮濃度小于河涌背景濃度時(shí)，排口的排放會(huì)增加河涌稀釋流量。本次河涌背景的氨氮濃度設(shè)定為上游實(shí)際來(lái)水濃度（1.5 mg/L），排口12 的氨氮排放濃度給定為1 mg/L，流量值依次增大，流量與貢獻(xiàn)度系數(shù)關(guān)系見(jiàn)圖3。從圖3 可以看出，當(dāng)排口氨氮濃度大于河涌背景濃度時(shí)，隨著流量的增加，貢獻(xiàn)度系數(shù)逐步減小，其減小幅度隨流量持續(xù)增加而趨于平緩，最終其數(shù)值接近-0.47。以貢獻(xiàn)度系數(shù)為因變量y，流量為自變量x，通過(guò)MATLAB 擬合得出貢獻(xiàn)度系數(shù)與流量關(guān)系式為：y=3.399 5e-0.084x-0.457 6，相關(guān)性系數(shù)R2=0.989 5，排口流量與貢獻(xiàn)度系數(shù)呈負(fù)相關(guān)關(guān)系，流量增加，貢獻(xiàn)度系數(shù)減小。

圖3 流量與貢獻(xiàn)度系數(shù)關(guān)系（C 排口

2.3.2 貢獻(xiàn)度系數(shù)與濃度的關(guān)系

設(shè)定流量和距離不變，選取氨氮為指標(biāo)，研究貢獻(xiàn)度系數(shù)與濃度的關(guān)系。具體來(lái)說(shuō)，選定排口12，給定1 m3/s 的流量，氨氮濃度值依次從1 mg/L 增大到10 mg/L，計(jì)算出下游考核斷面濃度及貢獻(xiàn)度系數(shù)，擬合結(jié)果如圖4 所示。從圖4 可知，排口氨氮濃度與下游考核斷面濃度呈標(biāo)準(zhǔn)線性關(guān)系，擬合函數(shù)關(guān)系式y(tǒng)=0.125 8x-10-8，排口氨氮濃度越大，貢獻(xiàn)度系數(shù)越大。

圖4 排口氨氮濃度與貢獻(xiàn)度系數(shù)的關(guān)系

2.3.3 貢獻(xiàn)度系數(shù)與距離的關(guān)系

設(shè)定排口流量、濃度不變，取河涌12 個(gè)排口進(jìn)行計(jì)算，其里程由小到大。

選取氨氮為指標(biāo)，模型背景氨氮濃度全部設(shè)置為0 mg/L，12 個(gè)排口的濃度和流量都給定1 個(gè)單位（1 mg/L），各排口貢獻(xiàn)度系數(shù)如圖5 所示。由圖5 可見(jiàn)，排口與考核斷面距離越大，貢獻(xiàn)度系數(shù)越小，擬合的函數(shù)關(guān)系式y(tǒng)=-9.95×10-8x+0.125 9，表明距離和貢獻(xiàn)度系數(shù)呈負(fù)相關(guān)，線性斜率k=-9.95×10-8，其值非常小，因此距離對(duì)貢獻(xiàn)度系數(shù)影響很小，計(jì)算時(shí)不做定量考慮。

圖5 距離和貢獻(xiàn)度系數(shù)的關(guān)系

2.4 單位負(fù)荷貢獻(xiàn)度系數(shù)

根據(jù)上述擬合結(jié)果，貢獻(xiàn)度系數(shù)與距離和濃度均呈線性關(guān)系，與流量呈指數(shù)關(guān)系，其中影響貢獻(xiàn)度系數(shù)的主要因素為濃度和流量。

分區(qū)達(dá)標(biāo)法中貢獻(xiàn)度系數(shù)計(jì)算的設(shè)定不考慮背景濃度，偏向比較理想狀態(tài)，而在實(shí)際計(jì)算應(yīng)用中，需要加入背景濃度以更貼合實(shí)際情況，由此引入單位負(fù)荷貢獻(xiàn)度系數(shù)（aij），其表達(dá)式可表示為：

式（3）中，i 為考核斷面?zhèn)€數(shù)，共m 個(gè)；j 為排口個(gè)數(shù)，共n 個(gè)；Cij為第j 個(gè)排口單獨(dú)排放時(shí)目標(biāo)考核斷面濃度，mg/L；Cbj為目標(biāo)考核斷面背景濃度，mg/L；Cj為第j 個(gè)排口的排放濃度，mg/L；Qj為第j 個(gè)排口的排放流量，m3/s。

式（3）中，Cij主要通過(guò)MIKE11 水質(zhì)模型計(jì)算得出，aij的單位為d·L－1，從單位上看，單位負(fù)荷貢獻(xiàn)度系數(shù)主要與流量有關(guān)。

3 實(shí)際案例應(yīng)用

選擇新安涌作為實(shí)際案例應(yīng)用，采用MIKE11模型計(jì)算各排口單獨(dú)排放時(shí)下游考核斷面濃度，再根據(jù)式（3）計(jì)算各排口單位負(fù)荷貢獻(xiàn)度系數(shù)。

3.1 新安涌排口基本資料

根據(jù)污染源調(diào)查及實(shí)測(cè)資料，新安涌排口的基本情況（濃度、流量、與上游距離）如圖1 和表1 所示。從表1 可以看出，排口水質(zhì)部分相差很大（排口16，30，35），大多屬于分流制混接排口，有學(xué)校、政府等生活污水匯入。

表1 新安涌有水排口概況

3.2 MIKE11 水環(huán)境模型

根據(jù)污染源調(diào)查，新安涌主要污染物為氨氮，因此計(jì)算以氨氮為指標(biāo)，采用MIKE11 準(zhǔn)確計(jì)算新安涌排口單獨(dú)排放時(shí)下游考核斷面濃度。

模型設(shè)置：上游設(shè)置流量邊界，新安涌上游來(lái)水為新安村大小水坑和城桂路排洪渠，流量為7.27 m3/s；下游設(shè)置水位邊界，新安涌下游與白石涌相連，以河道下游常水位作為河道水動(dòng)力的下游邊界值，新安涌下游常水位取1.5 m。

新安涌治理目標(biāo)為消除黑臭，上游來(lái)水為黑臭水質(zhì)，則氨氮初始值為8 mg/L，氨氮的降解系數(shù)取0.12 d-1。

3.3 單位負(fù)荷貢獻(xiàn)度系數(shù)計(jì)算

各排口下游考核斷面濃度根據(jù)公式（3）計(jì)算的單位負(fù)荷貢獻(xiàn)度系數(shù)見(jiàn)表2。單位負(fù)荷貢獻(xiàn)度系數(shù)為負(fù)值的共有4 個(gè)排口，分別為排口12，14，24，26，這4 個(gè)排口對(duì)新安涌水質(zhì)變好有促進(jìn)作用；其他排口單位負(fù)荷貢獻(xiàn)度系數(shù)均為正值，這些排口會(huì)加重新安涌水質(zhì)污染。貢獻(xiàn)度系數(shù)值越大，對(duì)新安涌水質(zhì)污染越大，表2 中排口37 水質(zhì)污染最大。

表2 考核斷面濃度和單位負(fù)荷貢獻(xiàn)度系數(shù)

在河涌黑臭治理中，需要著重治理貢獻(xiàn)度系數(shù)大的排口，貢獻(xiàn)度系數(shù)為負(fù)值的排口，對(duì)河涌的水質(zhì)有一定的正向作用，后續(xù)可通過(guò)量化分配，將水環(huán)境容量分配至每個(gè)排口，實(shí)現(xiàn)優(yōu)化分配。

4 結(jié)論及展望

（1）對(duì)分區(qū)達(dá)標(biāo)中的貢獻(xiàn)度系數(shù)進(jìn)行了詳細(xì)研究，分析了貢獻(xiàn)度系數(shù)及其影響因子，MATLAB 擬合結(jié)果表明，貢獻(xiàn)度系數(shù)與污染源濃度、污染源與考核斷面距離之間呈線性關(guān)系，與流量呈指數(shù)關(guān)系，貢獻(xiàn)度系數(shù)主要影響因子為濃度和流量。

（2）根據(jù)貢獻(xiàn)度系數(shù)定義以及影響因子研究結(jié)果，提出了單位負(fù)荷貢獻(xiàn)度系數(shù)的概念并給出其計(jì)算表達(dá)式，該公式考慮了背景濃度，可以更貼合實(shí)際情況。

（3）通過(guò)新安涌實(shí)際案例計(jì)算，驗(yàn)證了單位負(fù)荷貢獻(xiàn)度系數(shù)表達(dá)式的實(shí)用性，通過(guò)計(jì)算單位負(fù)荷貢獻(xiàn)度系數(shù)可以對(duì)河涌排口治理有初步判斷，后續(xù)拓展至流域可以量化流域水質(zhì)現(xiàn)狀，對(duì)于流域水環(huán)境治理有實(shí)際的指導(dǎo)意義。

（4）貢獻(xiàn)度系數(shù)的研究可以為后續(xù)總量分配提供計(jì)算依據(jù)，便于將分區(qū)達(dá)標(biāo)法求取的最優(yōu)水環(huán)境容量分配至每個(gè)排口或污染源。該研究考慮了河涌背景濃度，有利于與運(yùn)營(yíng)實(shí)際相結(jié)合，可以更準(zhǔn)確貼合實(shí)際的總量分配，便于后續(xù)流域水環(huán)境信息化管理。