黃代忠，方神光，張義樂

(珠江水利科學(xué)研究院，水利部珠江河口海岸工程技術(shù)研究中心，廣東廣州510611)

水是生命之源，但隨著工業(yè)和第三產(chǎn)業(yè)的粗獷式發(fā)展，水污染事故正成為危害人民生存健康的重大環(huán)境問題。由于突發(fā)性水污染事故具有發(fā)生、發(fā)展、危害的不確定性，流域性、影響的長期性和應(yīng)急主體不明確等顯著特征[1-3]，一旦發(fā)生，不僅對人群安全與生存帶來極大危害，還會對生態(tài)環(huán)境造成極大影響[4]。據(jù)統(tǒng)計從2000—2011年期間就發(fā)生1 176起突發(fā)性水污染事故。如2004年2月沱江支流由于約2 000 t純氨氮排入，造成江水嚴重污染；2005年11月苯類污染物污染松花江重大水污染事故和2005年12月粵北北江流域嚴重鎘超標水污染事故等。這一系列水污染事故的發(fā)生，給人民的身體健康和財產(chǎn)造成了巨大損失，因此突發(fā)性水體污染的應(yīng)對和解決方法已成為當前水環(huán)境研究的熱點問題。

在修建有梯級水庫的中上游河道，通過梯級水庫的聯(lián)合調(diào)度來應(yīng)對突發(fā)性水體污染事故的方法已經(jīng)在實際情況中得到初步使用并取得了不錯的結(jié)果。如2005年北江鎘污染事件中采用了加大上游水庫排量稀釋水體污染物、利用人工小洪峰加快污染物運移到下游處置區(qū)的應(yīng)急措施以控制污染事故的惡化；2012年廣西龍江發(fā)生的錫污染事件中，通過調(diào)度柳江上游大埔、麻石等水庫水量來進行稀釋，保證下游柳州市的飲水安全。由于突發(fā)性水污染事件頻繁發(fā)生，對于該領(lǐng)域的研究也得到快速發(fā)展。方神光等[5]對西江中游突發(fā)性污水團沿河運動特性進行了分析，提出了其影響變化范圍主要分為隨流剪切擴散階段和稀釋影響衰減階段。辛小康等[6]建立MIKE21[7]模型，以宜昌江段為研究對象，計算選取了三峽水庫對三種排放類型不同的調(diào)水方式對污染物輸移傳播的影響；于磊等[8]對大寧水庫發(fā)生突發(fā)性水污染事故時，綜合考慮風(fēng)、水庫調(diào)度的因素，設(shè)定了5種場景來模擬污染物擴散的規(guī)律。為了解西江中游水庫群對污染物輸移擴散的影響，尋找合適的梯級水庫調(diào)度方案以應(yīng)對突發(fā)性水污染事故，本文建立了西江中游河網(wǎng)水動力及污染物輸移擴散數(shù)學(xué)模型，通過不同調(diào)度方式下污染物擴散結(jié)果的對比分析，探討梯級水庫在突發(fā)性水體污染下的合理調(diào)度模式。

1 數(shù)學(xué)模型

1.1 河道梯級水庫水動力及污染物輸移擴散數(shù)學(xué)模型

水動力及污染物輸移擴散數(shù)學(xué)模型如下：

(1)

(2)

(3)

式中s——距離，m；Q——流量，m3/ s；Z——水位，m；g——重力加速度，m/s2；t——時間，s；A——過水斷面面積，m2；R——水力半徑，m；n——河道糙率系數(shù)；D——擴散系數(shù)，m2/s；C——某斷面t時刻的濃度，mg/m3；q——網(wǎng)格單位長度上的源或匯，m3/(s·m)；cl——源項或匯項的污染物濃度。

為保障污染物輸移擴散模型的迎風(fēng)特性，在對式(3)中對流項進行離散和計算過程中，需要判斷每個網(wǎng)格節(jié)點位置的流向，同時河道水系中的污染物輸移擴散主要以對流剪切為主，紊動擴散影響相對微小，因此在計算過程中進行了忽略處理。一維水動力數(shù)學(xué)模型采用Preissmann四點隱式差分格式進行離散，通過追趕法對離散得到的大型對角稀疏方程組進行求解。

1.2 模型驗證

通過水庫調(diào)度來改變污染物輸移轉(zhuǎn)化的特性，能在很大程度上減輕污染物的危害程度和危害范圍。但由于突發(fā)性水污染事故具有發(fā)生、發(fā)展等特性，真正通過大型水庫及梯級水庫群調(diào)水的方式來處理水污染的案例和研究很少。為驗證此污染物輸移擴散數(shù)學(xué)模型的準確性和可靠性，本文通過對相同工況下數(shù)值模擬和解析公式的結(jié)果分析，對模型進行驗證。

設(shè)置一個長10 km，底寬67.5 m，邊坡1∶2.5，底坡為0.15‰的棱柱形梯形河道，河道糙率系數(shù)取0.027，上邊界條件為保持恒定流量2 000 m3/s，下邊界根據(jù)明渠均勻流公式推算到的水深確定水位并保持恒定不變。假定在河道上游某斷面突然泄露了1 t的可溶性難降解污染物，試估算下游5、10 km斷面處污染物濃度的變化過程以及0.5、1 h后河道的污染物濃度分布[9]。此類問題可給出輸移擴散式(3)的解析解如下：

(4)

式中u——斷面平均流速，m/s；x——下游距離污染源的距離，m；M0——進入水體污染初始面源強度，g/m2；k——污染物衰減系數(shù)，s-1。

忽略污染物在垂向和橫向的擴散，推算斷面摩阻流速0.11 m/s；初始面源強度934.929 g/m，縱向擴散系數(shù)D=6.01hu*=7.4 m2/s。圖1給出了下游不同控制斷面位置污染物濃度隨時間的變化以及不同時刻污染物濃度沿程分布的模型解和解析解比較，驗證了建立的污染物輸移擴散數(shù)學(xué)模型的準確性和可靠性。

a) 斷面濃度變化過程線

b) 不同時間點污染物濃度沿程分布圖1 模擬解與解析解的對比

1.3 西江中游河網(wǎng)水動力及污染物輸移數(shù)學(xué)模型建立

西江為典型的樹狀河網(wǎng)，主干河流有紅水河、黔江、潯江；柳江、郁江為主要支流。由于西江中游梯級水庫的建設(shè)，河道上下游水流失去了水力連續(xù)性，圣維南方程組在此處不能適用，為達到同步計算的目的，并融入水庫的調(diào)度方案，此處引入子河網(wǎng)的概念，將水利樞紐作為河網(wǎng)細分為子河網(wǎng)的特殊節(jié)點，以便計算過程中考慮水利樞紐的調(diào)度作用，從而將各個河段連接成一個計算整體，詳細可參見相關(guān)文獻[10]。

2 紅水河突發(fā)性水污染事故下的水庫調(diào)度分析

2.1 突發(fā)性水體污染地點及調(diào)度水庫的選定

大化瑤族自治縣位于廣西壯族自治區(qū)中部偏西北的紅水河中游，隸屬于河池市。大化鎮(zhèn)是該縣的政治、經(jīng)濟、文化、交通中心，是集建材、冶煉、建筑、加工等產(chǎn)業(yè)為一體的工業(yè)城鎮(zhèn)。大化工業(yè)集中區(qū)是由原來的大化縣城南工業(yè)園區(qū)和巖灘工業(yè)區(qū)合并而成，在園內(nèi)集中著以金屬硅、硅錳合金冶煉，電解錳等為主的企業(yè)，給當?shù)氐纳鷳B(tài)環(huán)境帶來巨大的隱患。

本次計算分析假定在大化瑤族自治縣大化鎮(zhèn)發(fā)生突發(fā)性水污染事故，事故地點上距天峨水文站147 km，下距巖灘水電站2.2 km，具體位置見圖2。假設(shè)污染物從2月1日8時開始進入河道，在10時流量達到50 m3/s的最大流量，一直以該流量維持到18時，然后再20時降為0。污水所含的污染物濃度值初始值設(shè)為1，計算范圍內(nèi)污染物濃度的本底值是0，設(shè)河道內(nèi)的污染物濃度值超過0.000 1時，視為超標。

選取參與調(diào)度的水庫有巖灘、紅花、大藤峽，這些水庫在枯水期調(diào)度的規(guī)則見表1，在不影響水庫運行的前提條件下調(diào)度一定的水量去稀釋及加快污染物的輸移，讓其危害程度降到最低。

表1 突發(fā)性污染物事故應(yīng)對方案

2.2 選取的水文條件及還原

枯季天然來流小，斷面流速小，污染物停留在河道中時間變長，影響距離增大。因此，選擇了2014年2月1日至3月13日的流量、水位實測資料作為水文邊界條件(圖3)。

區(qū)間入流采用分段還原，將西江中游河網(wǎng)分為天峨—橋鞏段，橋鞏、柳州—武宣段，武宣、南寧—大湟江口、大湟江口—梧州段，梧州—高要段共計5部分。結(jié)合上下游洪量差值還原區(qū)間入流，區(qū)間入流分別設(shè)于刁江、洛清江、蒙江、桂江和賀江入流口。天峨—橋鞏分段還原時，上、下邊界條件分別由天峨站、橋鞏電站實測資料確定，巖灘、大化、百龍灘、樂灘壩前水庫都采用實測資料值，區(qū)間入流設(shè)在刁江口。橋鞏、柳州—武宣分段還原時，上、下游邊界條件分別由遷江、柳州和武宣站實測資料確定，紅花水電站壩前水位為正常蓄水位，區(qū)間入流設(shè)在洛清江口。武宣、南寧—梧州分段還原時，上、下游邊界條件分別由武宣、南寧和梧州站實測資料確定，西津、貴港、桂平長洲樞紐壩前水位由實測資料確定，區(qū)間入流分別設(shè)在蒙江口、桂江口。梧州—高要分段還原時，上、下游邊界條件根據(jù)梧州和高要站實測資料確定，區(qū)間入流設(shè)在賀江口。各河段還原獲得的水位、流量計算值與實測值在允許誤差范圍內(nèi)，其中天峨—橋鞏段還原與驗證結(jié)果見圖3。

圖2 突發(fā)性水污染發(fā)生地點示意

a) 上下游邊界條件

b) 區(qū)間入流過程還原

c) 橋鞏流量驗證圖3 天峨—橋鞏段還原與驗證

d) 天峨水位驗證續(xù)圖3 天峨—橋鞏段還原與驗證

2.3 無應(yīng)對措施下污染物輸移特性

假設(shè)污染物進入河道后，下游各梯級水庫未采取調(diào)度措施進行應(yīng)對，事故發(fā)生地點下游巖灘電站、大化電站、樂灘電站、橋鞏電站、石龍三江口、桂平三江口、梧州水文站點位置污染物濃度變化曲線見圖4。表2給出了各特征斷面污染物濃度變化特征值。

圖4 無應(yīng)對措施下各控制站點污染物濃度變化過程線

典型斷面開始超標時刻開始達標時刻出現(xiàn)峰值時刻峰值濃度/(mg·L-1)超標時長/h巖灘電站2月1日15時2月3日9時2月2日4時0.033 644大化電站2月5日17時2月8號9時2月6日17時0.023 866樂灘電站2月10號10時2月12號18時2月11號12時0.019 159橋鞏電站2月13號9時2月15日14時2月14號10時0.017 262石龍三江口2月17號4時2月19號21時2月18號2時0.011 167武宣水文站2月21日3時2月24日6時2月22日13時0.011 076桂平三江口2月27號5時3月2日1時2月28日11時0.004 970梧州水文站3月5號9時3月8日23時3月7號3時0.006 176

巖灘電站從2月1日15時開始超標，在2月2日4點達到最大值0.033 6，至2月3日9時消失，總共超標歷時44 h；大化電站壩址污染物濃度從2月5日17時開始超標，在2月6日17時達到最大值0.023 8，至2月8號9時開始達標，總共超標歷時66 h；樂灘電站壩址污染物濃度從2月10號10時開始超標，在2月11號12點達到最大值0.019 1，至2月12號18時消失，總共超標歷時59 h；橋鞏電站壩址處濃度在2月13號9時開始超標，在2月14號10時達到最大值0.017 2，至2月15日14時消失，總共超標歷時62 h；石龍三江口污染物濃度從2月17號4時開始超標，在2月18號2時達到最大值0.011 1，至2月19號21時消失，總共超標歷時67 h；武宣水文站位置從2月21日3時開始超標，在2月22日13點達到最大值0.011，至2月24日6點消失，總共超標歷時76 h；桂平三江口從2月27號5時開始超標，在2月28日11時達到最大值0.004 9，至3月2日1時消失，總共超標歷時70 h；梧州水位站位置污染物濃度從3月5號9時開始超標，在3月7號3時達到最大值0.006 1，至3月8日23時消失，總共超標歷時76 h。

從圖4中看出各典型斷面位置污染物濃度變化呈拋物線狀，拋物線從上游向下游呈現(xiàn)由尖瘦型向矮胖型過渡，顯示各控制站點污染物濃度都呈現(xiàn)某個時間段逐漸增大，達到最大值后逐漸減少直至消失，即污染團通過了該站點。各斷面污染物濃度的超標時間從上游的巖灘電站44 h增大到下游梧州水文站76 h，峰值濃度從0.033 6降低到0.006 1，在樂灘電站總超標歷時減少是由于刁江支流的匯入加快了污染物的輸移，桂平三江口位置總超標歷時減少是郁江水流的匯入加快了污染物的輸移傳播。

2.4 水庫調(diào)度方案設(shè)置

2.4.1梯級水庫基本特性

事故發(fā)生地下游沿線各梯級水庫情況見表3，由表中可知，事故發(fā)生地下游各梯級水庫興利庫容較大的有巖灘、紅花和大藤峽，由于本次模擬事故發(fā)生時間為枯水期，考慮水庫沖污、興利等綜合因素，本次調(diào)度方案主要對巖灘、紅花和大藤峽3個梯級水庫進行調(diào)度。

2.4.2水庫應(yīng)對調(diào)度方案設(shè)置

針對該突發(fā)性水污染事故，設(shè)置了4種水庫調(diào)度方案，見表4。

a) 由于水污染事故發(fā)生的突然性，現(xiàn)實發(fā)生水污染事件時采取的應(yīng)對手段多采用單庫調(diào)度進行應(yīng)對，因此水庫調(diào)度方案1為當巖灘電站壩址處的污染物濃度從2月1日15時開始超標時，巖灘水庫開始調(diào)度。

表3 事故地點下游各梯級水庫情況

b) 水庫調(diào)度方案2在方案1調(diào)度巖灘水庫的基礎(chǔ)上，加入紅花水庫一起調(diào)度，為并聯(lián)式雙庫調(diào)度，考慮到紅花調(diào)電站下泄流量與巖灘流量在石龍三江口匯合以增強調(diào)度效果的需要，紅花水庫的起調(diào)時間定為2月12日9時。

c) 水庫調(diào)度方案3在方案2同時調(diào)度巖灘和紅花水庫的基礎(chǔ)上，加入大藤峽水庫，大藤峽水庫的起調(diào)時間為當壩址處污染物濃度開始超標時，即從2月22日2時。

d) 大藤峽水庫興利庫容較大，如有需要可進一步增加下泄流量促進污染物輸移擴散，因此水庫調(diào)度方案4在方案3調(diào)度巖灘、紅花和大藤峽水庫的基礎(chǔ)上，將大藤峽水庫調(diào)度下泄流量從方案3的5 000 m3/s增大為6 000 m3/s。

表4 水庫調(diào)度方案設(shè)置

2.5 調(diào)度結(jié)果分析

2.5.1水庫調(diào)度方案1結(jié)果分析

水庫調(diào)度方案1實施后的計算結(jié)果見表5。方案1單獨調(diào)度巖灘水庫，選擇以增大下泄流量至4 000 m3/s，直至庫內(nèi)水位降至219 m的最低限制水位，調(diào)度時長共計26 h。圖5、6給出了典型斷面在調(diào)度方案1實施前后污染物濃度峰值削減比例和超標時長減少值(正值表示超標時長減少，負值表示增加)。結(jié)果顯示：巖灘電站壩址處水質(zhì)達標的時間從2月3日9時提前到2月2日的1時，污染物濃度的峰值降低了22.3%，污染物超標時間縮短32 h。武宣水文站處的污染物濃度超標時刻由2月21日3時提前到2月16日9時，達標時刻從2月22日13時提前到2月20日8時，污染物濃度的峰值降低了41.8%，但污染物超標時長增加了20 h；梧州水文站處污染物濃度超標時刻由3月5號9時提前到3月2日3時，污染物濃度達標時間從3月8日23時提前到3月6日12時，污染物濃度峰值降低了49.1%，但污染物超標時長增加了32 h。分析調(diào)度方案1得到以下結(jié)論。

表5 典型斷面污染物濃度特征值(水庫調(diào)度方案1)

a) 在突發(fā)性水體污染發(fā)生在巖灘庫區(qū)時，通過調(diào)度巖灘水庫，加大下泄流量，可以顯著降低庫區(qū)及下游典型斷面位置處的污染物濃度峰值，起到稀釋污染物的作用，典型污染物峰值削減幅度總體呈現(xiàn)向下游增大趨勢。

b) 由于加大巖灘出庫下泄流量，相應(yīng)增加了巖灘壩址下游河道的流速，下泄水流攜帶污染物以更快的速度向下游運動，巖灘水庫下游典型斷面位置污染物濃度超標時刻及達標時刻都顯著提前。

c) 巖灘水庫調(diào)度，可以顯著削減污染物在巖灘庫區(qū)和大化庫區(qū)滯留的時長，但大化壩址以下河道典型斷面的污染物超標時長呈現(xiàn)顯著增加趨勢。

圖5 調(diào)度方案1下的典型斷面污染物濃度峰值削減比例

圖6 調(diào)度方案1下的典型斷面污染物濃度超標時長削減值

2.5.2其他水庫調(diào)度結(jié)果分析

表6給出了水庫調(diào)度方案2實施后的主要斷面濃度特征值結(jié)果。由于調(diào)度方案2是在方案1調(diào)度巖灘水庫調(diào)度的基礎(chǔ)上，加入了紅花電站，紅花電站控制下泄流量為2 000 m3/s，于石龍三江口與巖灘下泄流量匯合，起到降低污染物濃度并快速向下游推移的作用。因此，從武宣水文站斷面污染物濃度變化過程來看，污水團抵達該斷面的時間為2月15日18時，較方案1提前了15 h，污染物峰值濃度也較方案1下降了21.9%，但該斷面污染物濃度超標時長則較方案1增加了11 h。污水團抵達梧州水文站斷面的時間同樣較方案1提前了18 h，污染物峰值濃度較方案1下降了19.4%，但污染物超標時長則增加了7 h。

表7給出了水庫調(diào)度方案3實施后的主要斷面污染物濃度特征值。該方案是在方案2巖灘和紅花電站聯(lián)合調(diào)度的基礎(chǔ)上，加入即將建成的大藤峽水庫，大藤峽水庫調(diào)度開始時間為污水團抵達壩址斷面的時刻，大藤峽水庫調(diào)度期間控制下泄流量5 000 m3/s。在大藤峽水庫參與調(diào)度的情形下，可以看到，污水團抵達大藤峽斷面的時間與方案2一致，但由于大藤峽水庫通過加大下泄流量的方式加入調(diào)度，加快了庫區(qū)水體的流動，因此污水團峰值較方案2提前25 h抵達壩址斷面，污染物峰值濃度與方案2相差不大，但該斷面污染物濃度超標時長較方案2顯著減少了44 h；污水團抵達梧州水文站斷面較方案2提前了26 h，污染物濃度峰值較方案2提前18 h抵達，峰值濃度相差不大，但超標時長較方案2顯著減少17 h，在峰值濃度差別不大的情況下，顯然梧州斷面污染物濃度超標時長減少是由于大藤峽水庫加大下泄流量從而對污染物起到極大稀釋作用形成的。

表8給出了調(diào)度方案4下的大藤峽壩址斷面和梧州水文站斷面的污染物特征值。該方案與調(diào)度方案3基本一致，只是大藤峽電站調(diào)度期間下泄流量比方案3多1 000 m3/s?？傮w來看，該方案實施后的變化趨勢與方案3一致，污水團及濃度峰值更快的抵達大藤峽壩址斷面和梧州水文站斷面，峰值濃度基本一樣，兩斷面污染物濃度超標時長都進一步縮短。

表6 典型斷面污染物濃度特征值(方案2)

表7 典型斷面污染物濃度特征值(方案3)

表8 典型斷面污染物濃度特征值(方案4)

3 結(jié)論

本文在假定大化瑤族自治縣大化鎮(zhèn)發(fā)生突發(fā)性水污染事故的前提下，設(shè)置不同的調(diào)度方案，聯(lián)合調(diào)度西江流域上一座或某幾座水庫來應(yīng)對該事故，得到的結(jié)論如下。

a) 單獨調(diào)度巖灘水庫，加大下泄流量?？杉铀傥廴疚锵乱扑俣龋@著降低庫區(qū)及下游典型斷面位置處的污染物濃度峰值，削減污染物在巖灘庫區(qū)和大化庫區(qū)滯留時長；同時使巖灘水庫下游典型斷面位置污染物濃度超標時刻及達標時刻都顯著提前；但大化壩址以下河道典型斷面的污染物超標時長呈現(xiàn)顯著增加趨勢。當增大巖灘下泄流量至4 000 m3/s時，巖灘電站壩前水質(zhì)污染物超標時間縮短32 h，污染物濃度的峰值降低了22.3%；武宣水文站斷面污染物濃度的峰值降低41.8%，超標時長增加了20 h；梧州水文站斷面污染物濃度峰值降低49.1%，污染物超標時長增加32 h。

b) 同時調(diào)度巖灘水庫和紅花電站?？蛇M一步降低石龍三江口下游典型斷面位置處的污染物濃度峰值；促進下游典型斷面位置污染物濃度超標時刻及達標時刻的提前；也會進一步延長下游河道典型斷面的污染物超標時長。當增加紅花電站控制下泄流量為2 000 m3/s時，武宣水文站斷面污染物超標時間較方案1提前15 h，超標時長增加11 h，污染物峰值濃度下降21.9%；梧州水文站斷面污染物超標時間提前18 h，超標時長增加7 h，污染物峰值濃度下降19.4%。

c) 同時調(diào)度巖灘水庫、紅花電站和大藤峽水庫。大藤峽上游污染物移動未發(fā)生改變，但是大藤峽壩址處污染物峰值到達時間提前，持續(xù)時間降低；大藤峽水庫下游各典型斷面污染物到達時間提前，持續(xù)時間降低，峰值濃度不變。當增加大藤峽水庫控制下泄流量5 000 m3/s時，大藤峽壩前水質(zhì)超標時間較方案2提前25 h，超標時長減少了44 h，污染物峰值濃度基本不變；梧州水文站斷面水質(zhì)超標時間較方案2提前26 h，超標時長減少17 h，峰值濃度相差不大。

d) 同時調(diào)度巖灘水庫、紅花電站和大藤峽水庫，增加大藤峽水庫下泄流量。大藤峽壩址及下游各典型斷面污染物移動速度進一步加快，但是峰值濃度不變。當大藤峽水庫控制下泄流量6 000 m3/s時，大藤峽壩前水質(zhì)超標時間較方案3提前11 h，超標時長減少了8 h，污染物峰值濃度基本不變；梧州水文站斷面水質(zhì)超標時間較方案3提前22 h，超標時長減少6 h，峰值濃度相差不大。

e) 對比以上4種調(diào)度方式中污染物運動情況可知，當西江中游發(fā)生突發(fā)性水體污染時：①對于干流沿線的水庫，充分利用各梯級水庫的可利用庫容，根據(jù)污染物移動的具體情況，上游水庫向下游逐級補水；②對于支流上的水庫，通過調(diào)控使支流水庫的調(diào)度流量過程傳遞到三江匯合口與污染物輸移傳播的過程保持同步。通過聯(lián)合調(diào)度西江中游各梯級水庫，可以取得比單庫調(diào)度更好的效果。