王青青，吳 昊，王臘春，季曉敏

(南京大學(xué)地理與海洋科學(xué)學(xué)院，江蘇 南京 210023)

近年來，我國的突發(fā)水污染事故進入了密集高發(fā)時段[1-2]。長江作為黃金水道，水上運輸發(fā)達，且沿江多工業(yè)企業(yè)，突發(fā)污染事故風(fēng)險高。據(jù)統(tǒng)計，自2000年以來長江流域較大的水污染事件就發(fā)生了50多起[3]。不同于常規(guī)污染，突發(fā)污染事故具有不確定性與隨機性[4]，污染物在瞬間或短期內(nèi)大量進入水體，在破壞水域生態(tài)環(huán)境的同時也對水源地的供水安全造成了嚴(yán)重威脅，影響社會經(jīng)濟穩(wěn)定。

為了降低突發(fā)水污染事故的危害，對事故的模擬研究也逐漸成為熱點，目前主要以水動力、水質(zhì)數(shù)學(xué)模型為基礎(chǔ)[5-9]，結(jié)合情景設(shè)定[10]來研究污染物的遷移擴散過程，確定污染物的運動規(guī)律和時空分布特征，進一步為突發(fā)污染事故預(yù)警與應(yīng)急提供依據(jù)[11]。不少學(xué)者對長江突發(fā)污染事故進行了模擬研究，例如袁玥等[12]應(yīng)用MIKE21FM二維模型對長江蘄春段排污事故進行模擬，量化污染事故的影響范圍與強度；張楚天等[13]采用平面二維水質(zhì)模型與GIS集成，實現(xiàn)對長江武漢段突發(fā)污染事故的動態(tài)可視化模擬。不同于長江中上游河段，長江下游河段多汊道地形[14]且受到潮汐作用的影響[15]，水動力條件復(fù)雜，目前對長江下游感潮河段這一復(fù)雜水文條件下的突發(fā)污染事故模擬較少。且由于地形資料、水文資料等的限制，突發(fā)污染事故的水動力基礎(chǔ)多選取事故發(fā)生的特定河段[16]，模擬范圍較小，構(gòu)建大范圍尺度的水動力模塊作為基礎(chǔ)來模擬突發(fā)污染事故的研究鮮有報道。

筆者構(gòu)建了大通到徐六涇段的二維水動力模型，在水動力模塊的基礎(chǔ)上，選取事故風(fēng)險高且水文條件復(fù)雜的南京段進行突發(fā)污染事故模擬。基于長江南京段潛在的事故風(fēng)險類型，同時設(shè)定工業(yè)和船舶2種突發(fā)污染事故情景，模擬污染團的遷移擴散過程以及到達各個水廠取水口的時間及濃度，研究結(jié)果可對長江水源地取水用水進行安全預(yù)警，同時為水源地應(yīng)急調(diào)水提供依據(jù)。

1 數(shù)據(jù)來源與研究方法

1.1 數(shù)據(jù)來源

研究區(qū)的地形文件由南京市航道局2009年編制的長江下游航行參考圖數(shù)字化生成，總長度約531 km，共用航道圖20張，比例尺皆為1∶40 000。流量數(shù)據(jù)采用2014年大通水文站實測日流量數(shù)據(jù)。潮位數(shù)據(jù)采用2014年南京站、鎮(zhèn)江(二)站、江陰站、徐六涇站實測日潮位數(shù)據(jù)，每日有4個潮位數(shù)據(jù)，包括2個高潮位和2個低潮位。

1.2 研究方法

1.2.1模型選取

長江下游干流的感潮河段范圍廣，且江面寬度遠大于深度，因此采用二維水動力模型作為突發(fā)污染事故的研究基礎(chǔ)。選取MIKE21 Flow Model FM模型進行模擬，該模型采用的是FM網(wǎng)格結(jié)構(gòu)，可對網(wǎng)格信息進行編輯，在研究區(qū)適用性更高。

二維水動力模塊建立于Navier-Stokes方程上，將三維流動控制方程簡化為二維淺水流動質(zhì)量和動量守恒控制方程(公式1～3)。水質(zhì)模塊建立于二維空間下的對流擴散方程(公式4)。

連續(xù)性方程：

(1)

x方向動量方程：

(2)

y方向動量方程：

(3)

對流擴散方程：

(4)

式(1)～(4)中,t為時間，s；(x，y)為直角或斜角坐標(biāo)系下的坐標(biāo)，m；(u，v)為(x，y)方向水流流速，m·s-1；h為總水深，m；η為潮位，m；S為源匯項，kg·m-3·s-1；(us，vs)為源匯項在(x，y)方向水流速度，m·s-1；g為重力加速度，m·s-2；ρ和ρ0分別為水體和參考水密度，kg·m-3；Pa為水表面大氣壓，Pa；τsx、τbx、τsy、τby為底部剪切應(yīng)力，N；Sxx、Sxy、Syx、Syy為輻射應(yīng)力分量，N；f為科氏力系數(shù)，m-1；Txx、Txy、Tyx、Tyy為粘滯切應(yīng)力分量，N；C為物質(zhì)質(zhì)量濃度，kg·m-3；(Dx，Dy)為(x，y)方向上的擴散系數(shù)，m2·s-1。字母上方有一橫杠的表示平均值。

1.2.2水動力模塊構(gòu)建與驗證

采用非結(jié)構(gòu)三角網(wǎng)格對河道地形進行處理，三角網(wǎng)格總數(shù)為79 685個，最小網(wǎng)格分辨率為100 m。如圖1所示，上邊界設(shè)為大通水文站，下邊界設(shè)為徐六涇潮位站，模型計算的時間步長設(shè)置為10 min。

圖1 模型地形文件和網(wǎng)格文件

長江中下游的糙率系數(shù)(n)一般取值在 0.014～0.030之間[17]。該研究的二維水動力模型中采用曼寧系數(shù)(M)，M與n之間存在倒數(shù)關(guān)系，即M=1/n。由于選取的河道較長，且上下游河道的地形差距較大，因此將河底地形考慮在內(nèi)，使用水深數(shù)據(jù)對公式進行優(yōu)化。改進后的曼寧系數(shù)公式為

(5)

式(5)中，M為曼寧系數(shù)；n0和n1為底部的粗糙程度，n0的取值范圍一般為0.01～0.02，n1的取值范圍一般為0.005～0.01；d為水深，m。在取值范圍內(nèi)選取不同的n0和n1值對水位進行模擬，經(jīng)過率定得出n0=0.014，n1=0.01。

選取2014年2月南京站、鎮(zhèn)江(二)站、江陰站潮位數(shù)據(jù)進行驗證，將這3個站點的實測潮位數(shù)據(jù)與模型的模擬水位數(shù)據(jù)進行對比，結(jié)果如圖2所示。模擬值與實測值吻合性良好，可作為水質(zhì)模塊的水動力基礎(chǔ)。

圖2 南京站、鎮(zhèn)江(二)站、江陰站潮位實測數(shù)據(jù)與模擬數(shù)據(jù)對比

1.2.3突發(fā)污染事故情景設(shè)計與參數(shù)設(shè)置

1989年揚子石化公司人員誤將設(shè)備剩余物料倒入馬汊河，造成附近水廠供水無法飲用；2014年，揚子石化煉油廠酸罐爆炸，導(dǎo)致5個儲罐泄漏著火，污染物通過排水口進入馬汊河；2013年“鑫川8號”船舶碰擦南京長江大橋后，行至八卦洲洲頭沉沒，導(dǎo)致沉船內(nèi)80多t燃料重油溢出。綜合來看，長江南京段水源安全威脅主要來自于2個方面：一是由于企業(yè)或其他污染源超常規(guī)排放；二是由于車輛、船舶運輸事故造成的污染物泄漏。突發(fā)污染物類型集中在油類、化學(xué)品(硫酸、苯酚和苯類等)和重金屬(以鎘為代表)上。結(jié)合長江南京段水源地的分布，設(shè)定2種突發(fā)污染事故情景(圖3)：一為馬汊河工業(yè)事故，二為南京長江三橋船舶事故，同時選取對長江南京段水源威脅大的苯酚作為特征污染物。

圖3 突發(fā)污染事故點位

考慮水動力條件最弱情況下污染物擴散對水源地產(chǎn)生的影響最大[18-19]，選取枯水期(2014年2月)作為水動力基礎(chǔ)，假定100 t苯酚在2 h內(nèi)進入水體?？紤]事故發(fā)生點到水源地的距離(表1)，在馬汊河工業(yè)事故中設(shè)定落潮和漲潮2個事故發(fā)生的時間點，分別為2月5日06：00與2月5日12：00；在三橋船舶事故中，設(shè)定事故發(fā)生時間點為落潮段2月5日06：00。

表1 突發(fā)污染事故點距水源地距離

Table 1 Distance between the accident point and water resources

事故地點水源地距離1)/km夾江江浦—浦口燕子磯上壩龍?zhí)恶R汊河-22.6-22.1-13.7-4.519.3三橋 5.913.123.130.743.0/54.82)

1)負值表示該水源地位于突發(fā)污染事故點上游；2)三橋事故點至龍?zhí)端吹亟?jīng)過八卦洲汊道，南汊距離43.0 km，北汊距離54.8 km。

在泄露污染初期污染源的擴散情況尚不清晰，且長江水量豐沛，故僅考慮苯酚的橫向擴散，以污染最大化來建立預(yù)測[20-21]。另由于苯酚屬于非常規(guī)污染物類型，無常規(guī)水質(zhì)監(jiān)測數(shù)據(jù)資料，故結(jié)合其他參考文獻設(shè)定苯酚的白天降解系數(shù)為0.20 d-1[22]，夜晚降解系數(shù)為0.176 d-1[23]，擴散系數(shù)為10 m2·s-1[24]。

2 結(jié)果與分析

根據(jù)GB 3838—2002《地表水環(huán)境質(zhì)量標(biāo)準(zhǔn)》，在水質(zhì)檢測中，Ⅲ類水允許的苯酚極限值為0.005 mg·L-1。因此在輸出濃度結(jié)果時，將濃度值超出0.005 mg·L-1的區(qū)域定義為受影響區(qū)域。

2.1 情景一：馬汊河工業(yè)事故

如圖4所示，在落潮和漲潮的水動力條件下上壩水源地均能檢測到污染物，但污染物濃度量級僅為10-11mg·L-1，因此可認為在潮汐作用中污染團會隨著潮汐上溯，但對上壩水源地基本無影響；龍?zhí)端吹卦诼涑焙蜐q潮下均受到影響，表現(xiàn)為落潮情況下污染團先到達龍?zhí)端吹?，濃度峰值低于漲潮，2種情況下污染物的整體運動特征一致。

圖4 漲落潮情況下上壩水源地和龍?zhí)端吹匚廴疚锏臐舛茸兓?/p>

以落潮情況為例進行分析(表2)：龍?zhí)端吹卦谑鹿拾l(fā)生后24～54.5 h內(nèi)受到影響，總時長約為30.5 h。在此過程中，污染物濃度峰值發(fā)生了5次連續(xù)波動?？紤]到水流的往復(fù)運動，可概化成2個連續(xù)的峰值過程：事故發(fā)生后30.6～35.3 h污染物質(zhì)量濃度達到首個峰值，為0.196 mg·L-1，波峰持續(xù)時間為4.7 h；事故發(fā)生后41.5～44.3 h污染物質(zhì)量濃度到達第2個較小的峰值，為0.094 mg·L-1，波峰持續(xù)時間為2.8 h。

表2 龍?zhí)端吹匚廴疚餄舛确逯底兓?/p>

Table 2 Peak concentration changes of pollutants in Longtan water source

波動時間(月-日-時:分)峰值質(zhì)量濃度/(mg·L-1)第1次02-06-10:400.180第2次02-06-12:100.193第3次02-06-17:200.196第4次02-06-23:30 0.087第5次02-07-02:200.094

第1次峰值過程污染團受徑流作用，從上游向下游運動〔圖5(a)〕；第2個峰值過程污染團受到漲潮流控制，從下游向上游回溯[25]〔圖5(b)〕。第2個峰值濃度比第1個小，且對龍?zhí)端吹氐挠绊憰r間更短，這是由于在污染團的來回運動中污染物被稀釋和降解，使得回溯過程對水源地的影響較小。

圖5 龍?zhí)端吹匚廴疚镞w移擴散時空分布

由此可見，由于潮流的控制，污染團會發(fā)生往復(fù)運動，水源地多次出現(xiàn)污染物濃度峰值，受影響的時間也隨之加長。在馬汊河工業(yè)事故情景下，污染團僅對龍?zhí)端吹禺a(chǎn)生影響，并未上溯到馬汊河上游的上壩水源地和燕子磯水源地。因此在龍?zhí)端吹厥艿接绊懙?0.5 h內(nèi)可就近使用上壩和燕子磯水源地作為龍?zhí)端吹氐膽?yīng)急水源。

2.2 情景二： 南京長江三橋處船舶事故

2.2.1江浦—浦口、燕子磯、上壩、龍?zhí)端吹?/p>

如表3所示，情景2南京長江三橋事故對下游5個水源地均有影響，除夾江水源地外，三橋北側(cè)、南側(cè)事故對其他4個水源地的影響基本相同。以北側(cè)事故為例，江浦—浦口水源地最先受到影響，持續(xù)時間為18.5 h；燕子磯為第2受影響水源地，持續(xù)時間為24.5 h；第3受影響的是上壩水源地，持續(xù)時間為27.8 h；最后受影響的是龍?zhí)端吹兀掷m(xù)時間為145.8 h。

表3 各水源地受到污染的時間統(tǒng)計表

Table 3 Pollution affected time of different water sources

水源地水廠三橋北側(cè)事故三橋南側(cè)事故開始受影響的時間(月-日-時:分)結(jié)束影響的時間(月-日-時:分)受影響總時長/h開始受影響的時間(月-日-時:分)結(jié)束影響的時間(月-日-時:分)受影響總時長/h夾江北河口02-05-23:5002-06-15:0015.202-05-23:5002-12-07:40151.8城南02-05-18:2002-06-16:2022.002-05-14:5002-08-05:5063.0江浦—浦口江浦02-05-10:1002-05-21:2011.202-05-10:3002-05-21:3011.0浦口02-05-15:2002-06-04:4012.302-05-15:4002-06-04:5013.2燕子磯上元門02-05-19:2002-06-18:4023.302-05-19:3002-06-19:0023.5城北02-05-20:4002-06-20:1023.502-05-20:5002-06-20:3023.7上壩遠古02-06-05:3002-07-09:2027.802-06-05:5002-07-09:3027.7龍?zhí)洱執(zhí)?2-06-18:1002-12-20:00145.802-06-18:2002-12-20:00145.7

各水廠在受影響的過程中污染物濃度均出現(xiàn)多次峰值。以北側(cè)事故為例(圖6)，其中北河口、城南、江浦、浦口、上元門水廠出現(xiàn)2次峰值，城北、遠古水廠出現(xiàn)3次峰值，龍?zhí)端畯S出現(xiàn)6次以上峰值。綜合水源地來看，江浦—浦口水源地在事故發(fā)生后11.5 h出現(xiàn)濃度峰值，為0.667 mg·L-1；燕子磯水源地在事故發(fā)生后16.8 h出現(xiàn)濃度峰值，為0.440 mg·L-1；遠古水源地在事故發(fā)生后29.3 h出現(xiàn)濃度峰值，為0.283 mg·L-1；龍?zhí)端吹卦谑鹿拾l(fā)生后47.8 h出現(xiàn)濃度峰值，為0.193 mg·L-1。

圖6 三橋北側(cè)事故中各水廠取水口污染物濃度變化

同情景1，江浦—浦口、燕子磯、上壩水源地水廠取水口出現(xiàn)2～3次峰值，為徑流、潮流共同作用的結(jié)果，峰值的出現(xiàn)是連續(xù)的。但值得注意的是，龍?zhí)端吹胤逯档某霈F(xiàn)不是完全連續(xù)的，這是由于污染團在經(jīng)過八卦洲汊道時分為2團分別沿著長江干流、八卦洲左汊進行運動(圖7)。長江干流攜帶的污染團先到達龍?zhí)端吹?，在徑流、潮流的共同作用下來回運動，使得龍?zhí)端吹爻霈F(xiàn)前5個峰值，峰值質(zhì)量濃度為0.045～0.189 mg·L-1；八卦洲左汊的污染團后到達龍?zhí)端吹?，使龍?zhí)端吹爻霈F(xiàn)后4個峰值，質(zhì)量濃度低于0.019 mg·L-1。

圖7 龍?zhí)端吹匚廴疚镞w移擴散時空分布

綜上可知，水源地距離事故點的位置越遠，受事故的影響越晚，但受影響的時間會加長；水源地出現(xiàn)的峰值個數(shù)隨距離增加而增加，但峰值濃度值減小。這是由于在污染團的運動過程中一方面污染帶長度隨著遷移擴散距離的增加而增加，在徑流、潮流共同作用下往復(fù)過程多，污染物濃度不斷下降；另一方面，長江南京段汊道較多，污染團通過汊道時分散，導(dǎo)致位于汊道下游的水源地經(jīng)過多個污染團，下游水源地受影響的時間加長、污染物濃度峰值個數(shù)增加。

2.2.2夾江水源地

如圖8所示，夾江水源地在2次事故中受到的影響差距很大。在北側(cè)事故中，夾江水源地受影響的時間為2月5日18:20—2月6日15:00，持續(xù)時間為20.7 h；位于水源地內(nèi)的北河口水廠出現(xiàn)2個質(zhì)量濃度峰值，最大值為0.032 mg·L-1，城南水廠也出現(xiàn)2個峰值，最大值為0.006 mg·L-1。在南側(cè)事故中受影響的時間為2月5日14:50—2月12日07:40，持續(xù)時間為160.8 h；北河口水廠出現(xiàn)6個以上質(zhì)量濃度峰值，最大值為0.055 mg·L-1，城南水廠出現(xiàn)6個峰值，最大值為0.199 mg·L-1。

圖8 城南和北河口水廠取水口污染物濃度變化

分析二維水質(zhì)模型，在三橋北側(cè)事故中污染團少量進入夾江水源地，僅對上游城南水廠產(chǎn)生影響，未到達下游北河口水廠，北河口水廠的污染團來自于長江主江段潮流回溯過程；在三橋南側(cè)事故中污染團進入夾江水源地，對城南水廠、北河口水廠均產(chǎn)生影響，并在夾江水源地內(nèi)來回運動。

夾江水源地地形狹窄，水流量小，一旦污染團進入，會在水源地內(nèi)來回運動，持續(xù)影響水源地水質(zhì)。因此，對于夾江水源地而言，要著重防止污染團進入夾江段，同時要防止由于潮流作用污染團上溯進入北河口水廠。

3 結(jié)論

基于構(gòu)建的大通到徐六涇二維水動力模型，結(jié)合長江南京段實際污染來源，設(shè)定馬汊河工業(yè)污染和長江三橋船舶污染2種突發(fā)污染事故情景，模擬污染物苯酚的遷移擴散過程。

(1)該研究構(gòu)建的大通到徐六涇二維水動力模型中，南京站、鎮(zhèn)江站和江陰站的模擬潮位與實測數(shù)據(jù)吻合度良好，為感潮河段下突發(fā)污染事故污染物遷移擴散的模擬奠定了可靠的基礎(chǔ)。

(2)在設(shè)定的2種典型突發(fā)污染事故情景中，馬汊河工業(yè)事故僅對下游龍?zhí)端吹赜?0.5 h的影響；三橋船舶事故中，江浦—浦口、燕子磯、上壩、龍?zhí)端吹厥苡绊懙臅r間分別為18、24、28和146 h，夾江水源地在三橋北側(cè)事故中受影響時間為18 h，在三橋南側(cè)事故中為154 h。

(3)長江南京段受到徑流和潮汐作用的共同影響，污染物在河道內(nèi)往復(fù)運動，隨著污染物的遷移擴散，污染帶長度增加，下游水源地受影響時間加長。長江南京段具有獨特的汊道地形，污染團進入不同汊道，對下游龍?zhí)端吹卦斐啥啻挝廴尽?/p>

該研究構(gòu)建的大通到徐六涇的二維水動力模型為長江下游突發(fā)污染事故的模擬提供了一定的研究基礎(chǔ)。同時模擬過程中發(fā)現(xiàn)潮汐作用對污染物的遷移擴散產(chǎn)生一定影響，但均未上溯到事故發(fā)生點之上的水源地，后續(xù)可進一步進行不同潮汐類型以及漲潮落潮情況下的研究，為突發(fā)污染事故的應(yīng)急、水源地選址等提供依據(jù)。