瞬時(shí)點(diǎn)源和連續(xù)源隨流輸移擴(kuò)散

黃楨宇，沈良朵，陳浩民，黃楨蕾

(浙江海洋大學(xué) 港航與交通運(yùn)輸工程學(xué)院,浙江 舟山 316022)

近年來，隨著航運(yùn)業(yè)的不斷發(fā)展，由船舶引發(fā)的燃油泄漏事故時(shí)有發(fā)生，對海洋環(huán)境有十分嚴(yán)重的影響。如何處理這類隨流擴(kuò)散污染物已成為人們關(guān)注的熱點(diǎn)，其中污染物輸移擴(kuò)散的燃油泄漏規(guī)律研究尤為受到重視。近岸環(huán)流主要出現(xiàn)在緩坡海域波浪破碎之后形成的波浪破波區(qū)，近岸環(huán)流速度沿水深的鉛直分布是不均勻的。由于近岸緩坡海域波浪破碎區(qū)水深較淺，近岸環(huán)境速度的這種不均勻性可忽略，認(rèn)為流速沿水深混合充分、均勻，把運(yùn)動(dòng)簡化為平面二維運(yùn)動(dòng)。[1]以往的研究普遍認(rèn)為潮流是污染物輸移擴(kuò)散的主要驅(qū)動(dòng)力[2]，這種分析方式適用于在遠(yuǎn)海和擁有巖石陡坡地形的近海環(huán)境。然而，在近岸緩坡地形條件下不能忽略波浪對物質(zhì)輸移的影響。目前，國內(nèi)對近岸污染物擴(kuò)散問題只進(jìn)行了試探性研究。例如：孫濤等[3]和陶建華等[4]主要研究基于數(shù)學(xué)模型的波浪作用對近岸污染物擴(kuò)散的影響；丁雷[5]主要在解析法上對污染物紊動(dòng)擴(kuò)散系數(shù)進(jìn)行推導(dǎo)研究；彭云[6]采用物理模型研究船舶燃油泄漏之后水體污染物的擴(kuò)散特征；武周虎等[7]研究變擴(kuò)散系數(shù)下三維濃度的擴(kuò)散問題；江春波等[8]主要針對二維物質(zhì)輸移擴(kuò)散問題提出有限體積算法。國外學(xué)者對比也做了相應(yīng)的研究。[9-12]

基于上述研究成果，本文以瞬時(shí)點(diǎn)源和連續(xù)源污染物為研究對象，給出無限區(qū)域點(diǎn)源和連續(xù)源污染物濃度解析解、數(shù)值解和墨水?dāng)U散試驗(yàn)。在此基礎(chǔ)上，對簡單波流場下點(diǎn)源和連續(xù)源污染物濃度輸移擴(kuò)散特征進(jìn)行研究。

1 瞬時(shí)點(diǎn)源解析解

不考慮對流輸移擴(kuò)散情況下的二維擴(kuò)散方程為

(1)

式(1)中:C為深度平均的污染物濃度；Dcwx和Dcwy分別為在波流作用下污染物沿x和y方向的擴(kuò)散系數(shù)。

由于Oxy平面上的某點(diǎn)濃度在x方向和y方向相互獨(dú)立，故某點(diǎn)的濃度C(x,y,t)可由2部分乘積構(gòu)成，即

C(x，y，t)=C1(x，t)×C2(y，t),

C1≥0;C2≥0

(2)

式(2)中:C1和C2分別為x方向和y方向深度平均的污染物濃度。

瞬時(shí)點(diǎn)源質(zhì)量M為

(3)

S={(x,y)|x2+y2≤R2}

(4)

式(3)和式(4)中：M1和M2分別為瞬時(shí)點(diǎn)源在深度平均的x和y方向的分質(zhì)量。

將式(2)代入式(1)可得

(5)

由濃度分布的對稱性，有

(6)

式(6)中的瞬時(shí)點(diǎn)源解析式分別為

(7)

由式(2)可得二維擴(kuò)散方程的瞬時(shí)點(diǎn)源解為

(8)

各向同性情況下有

Dcwx=Dcwy=D

(9)

將式(9)代入式(8)可進(jìn)一步化簡為

(10)

考慮點(diǎn)源位置(x0,y0)隨流輸移速U和v可得濃度解析解為

(11)

2 連續(xù)源解析解

根據(jù)瞬時(shí)點(diǎn)源解析解推導(dǎo)過程，對于無限區(qū)域，將式(8)進(jìn)行疊加可得到連續(xù)源解析解，經(jīng)過時(shí)間t后連續(xù)源微元濃度的解析解為

(12)

對等式兩邊積分，有

(13)

由于積分結(jié)果無法用基本初等函數(shù)表示，因此需借助三維空間，采用由三維到二維的方法將該問題設(shè)想為有一系列厚度為dx的薄片(見圖1),假設(shè)波流場以一個(gè)和沿x正方向的速度u經(jīng)過與坐標(biāo)原點(diǎn)重合的源點(diǎn)O，此時(shí)每一薄片接受的污染物質(zhì)量為Mdτ(dτ為每一薄片通過源點(diǎn)的歷時(shí)dτ=dx/u)，然后在平面上做橫向擴(kuò)散(如圖1所示)，根據(jù)二維擴(kuò)散的瞬時(shí)點(diǎn)源解，薄片中單位面積上的質(zhì)量為

(14)

式(14)中:Dcwz和Dcwy分別為波流作用下污染物沿z和y方向的擴(kuò)散系數(shù)。

圖1 連續(xù)源三維推導(dǎo)圖

薄片的位置由x=ut給出，且三維的濃度是薄片單位上的質(zhì)量除以厚度，因此其解為

(15)

由于是在簡單波流暢情況下,因此有Dcwz=Dcwy=D，令r2=z2+y2，并將其代入式(14)，化簡得

(16)

根據(jù)式(16)，令z為0，可將其轉(zhuǎn)化為二維解析式

(17)

令無量綱濃度為

原解析式即為

(18)

3 數(shù)值求解

3.1 基本控制方程

深度平均的二維對流輸移擴(kuò)散方程為

(19)

式(19)中:H為波高；SW為污染物源匯項(xiàng)。

3.2 擴(kuò)散系數(shù)

3.2.1純水流作用下的擴(kuò)散系數(shù)

純水流作用下深度平均的污染物擴(kuò)散系數(shù)為

(20)

(21)

式(20)和式(21)中：Dcx和Dcy為流作用下污染物沿x和y方向的擴(kuò)散系數(shù)；K為謝才系數(shù);k1和k2分別取5.93和0.15。

3.2.2純波浪作用下的擴(kuò)散系數(shù)

純波浪作用下污染物擴(kuò)散系數(shù)為

(22)

式(22)中:Dwx和Dwy分別為波流作用下污染物沿x和y方向的擴(kuò)散系數(shù);h為當(dāng)?shù)厮?；H=h+η;η為波高與當(dāng)?shù)厮畹牟钪?；T為波浪周期;破碎區(qū)α1=5.5η/h，非破碎區(qū)α1=1.0。

3.2.3波流共同作用下的擴(kuò)散系數(shù)

在波浪和潮流作用下污染物擴(kuò)散系數(shù)可采用比較簡單的線性迭代法得到，即

Dcw=Dc+Dw

(23)

3.3 邊界條件

在一個(gè)面積很大的海域，周邊邊界對污染物的影響非常小，可用以下邊界條件：

1) 當(dāng)流體流入計(jì)算區(qū)域時(shí)，邊界處濃度可置為域外值，即C=0。

2) 當(dāng)流體流出計(jì)算區(qū)域時(shí)，根據(jù)實(shí)際情況，可給定邊界條件為

(24)

式(24)中:v為邊界處流速矢量。

一般情況下可采用更簡潔的連續(xù)邊界條件，即

(25)

式(25)中:n為邊界外法線方向。

4 結(jié)果比較

4.1 點(diǎn)源擴(kuò)散

當(dāng)計(jì)算區(qū)域足夠大時(shí)，可忽略邊界影響，與理論推導(dǎo)的無限區(qū)域情況相近。擴(kuò)散范圍300 000 m×300 000 m；初始點(diǎn)源放在區(qū)域中心位置；為簡化計(jì)算，取各向同性，水平速度為U=v=0.5 m/s；D=Dcwx=Dcwy=10 000 m2/s，M=10 000 000個(gè)單位質(zhì)量；為使計(jì)算穩(wěn)定并快速得到計(jì)算結(jié)果，空間步長取Δx=Δy=3 000 m,時(shí)間步長取Δt=10 s。

為便于比較分析初始階段和穩(wěn)定階段點(diǎn)源解析解和數(shù)值解，當(dāng)點(diǎn)源解析解和數(shù)值解為1 000 s和5 000 s時(shí)，相應(yīng)濃度分布三維圖分別見圖2和圖3。當(dāng)點(diǎn)源解析解和數(shù)值解分別為1 000 s和5 000 s時(shí)，相應(yīng)濃度分布等值線比較見圖4和圖5。

a) 解析解

a) 解析解

圖4 點(diǎn)源解析解數(shù)值解1 000 s等值線對比圖 圖5 點(diǎn)源解析解數(shù)值解5 000 s等值線對比圖

由圖2～圖5可知:在簡單各向同性波流場下，點(diǎn)源隨流擴(kuò)散數(shù)值解與解析解吻合良好，呈現(xiàn)出圓環(huán)狀分布；濃度值由中心向外部遞減，隨著時(shí)間的推移，相同圓環(huán)半徑處的污染物濃度下降。

4.2 連續(xù)源擴(kuò)散

模型計(jì)算區(qū)域和參數(shù)選取與點(diǎn)源情況一致。當(dāng)連續(xù)源解析解和數(shù)值解為1 000 s時(shí)，相應(yīng)濃度分布三維圖見圖6。連續(xù)源解析解和數(shù)值解無量綱化等值線對比見圖7；連續(xù)源污染物擴(kuò)散中某點(diǎn)濃度隨時(shí)間變化圖見圖8。

a) 解析解

圖7 連續(xù)源解析解和數(shù)值解無量綱化等值線對比

圖8 連續(xù)源某點(diǎn)濃度隨時(shí)間的變化曲線

由圖6和圖7可知:在各向同性簡單波流場下，連續(xù)源隨流輸移擴(kuò)散數(shù)值解與解析解吻合良好，濃度隨空間呈橢圓分布，距離源點(diǎn)越遠(yuǎn)，相對濃度越小。由圖8可知:隨著時(shí)間的推移，距離投放源點(diǎn)固定位置處的污染物濃度趨于一個(gè)穩(wěn)定值。

5 試驗(yàn)數(shù)值結(jié)果比較

本文結(jié)合文獻(xiàn)[1]中的點(diǎn)源和連續(xù)源墨水?dāng)U散試驗(yàn)進(jìn)行比較分析，其坐標(biāo)示意見圖9。取工況3和工況6，坡度為1∶40，波浪入射位置的靜水深h=0.45 m，波浪入射角θ=30°，入射波高H0=0.05 m，入射波周期T=2 s，瞬時(shí)點(diǎn)源取源項(xiàng)=C(3，4)=1，連續(xù)源取源項(xiàng)Sm(3，4)=0.1(坐標(biāo)表示投放點(diǎn)距離靜水線3 m，距離沿岸流上游4 m)；計(jì)算參數(shù)參考文獻(xiàn)[1]中的數(shù)據(jù)(紊動(dòng)混摻系數(shù)N=0.001 2，波浪破碎指標(biāo)γ=0.7，波流場底摩擦因數(shù)cf=0.006，空間步長根據(jù)要求為Δx=Δy=0.2 m，時(shí)間步長Δt=0.1 s)，點(diǎn)源和連續(xù)源試驗(yàn)及其計(jì)算結(jié)果分別見圖10和圖11，其中濃度等值線從內(nèi)到外依次為100%、30%和5%。

圖9 物質(zhì)擴(kuò)散模型坐標(biāo)示意

a) t=22 s

c) t=32 s

a) t=32 s

c) t=52 s

由圖10和圖11可知:同一時(shí)刻點(diǎn)源計(jì)算結(jié)果和試驗(yàn)結(jié)果較為吻合，進(jìn)一步表明該數(shù)值模型能較好地模擬點(diǎn)源的實(shí)際擴(kuò)散效應(yīng)。與點(diǎn)源橢圓分布不同的是，連續(xù)源濃度最大值點(diǎn)并不在中心，而是在投放源點(diǎn)位置處。

6 結(jié)束語

本文結(jié)合二維對流輸移擴(kuò)散瞬間點(diǎn)源和連續(xù)源的解析解、數(shù)值模型以及墨水?dāng)U散試驗(yàn)得出瞬時(shí)點(diǎn)源和連續(xù)源在各向同性的簡單波流場環(huán)境下的濃度輸移擴(kuò)散分布特征，即：瞬時(shí)點(diǎn)源呈圓環(huán)分布濃度值由中心向外部遞減，隨著時(shí)間推移，相同圓環(huán)半徑處污染物濃度降低；連續(xù)源濃度隨空間呈橢圓分布，距離源點(diǎn)越遠(yuǎn)，相對濃度越小，隨著時(shí)間的推移，距離投放源點(diǎn)固定位置處污染物濃度將趨于一個(gè)穩(wěn)定值。