黃楨宇,沈良朵,陳浩民,黃楨蕾
(浙江海洋大學(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)行研究。
不考慮對流輸移擴(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)
根據(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)
深度平均的二維對流輸移擴(kuò)散方程為
(19)
式(19)中:H為波高;SW為污染物源匯項(xiàng)。
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)
在一個(gè)面積很大的海域,周邊邊界對污染物的影響非常小,可用以下邊界條件:
1) 當(dāng)流體流入計(jì)算區(qū)域時(shí),邊界處濃度可置為域外值,即C=0。
2) 當(dāng)流體流出計(jì)算區(qū)域時(shí),根據(jù)實(shí)際情況,可給定邊界條件為
(24)
式(24)中:v為邊界處流速矢量。
一般情況下可采用更簡潔的連續(xù)邊界條件,即
(25)
式(25)中:n為邊界外法線方向。
當(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)半徑處的污染物濃度下降。
模型計(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)定值。
本文結(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)位置處。
本文結(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)定值。