交匯角對明渠交匯口污染物輸運特性影響的數值模擬分析

徐潔如,朱 海,毛頌平,武顯揚,王玲玲,余真真,王 寅

(1.河海大學水利水電學院,江蘇 南京 210098; 2.上?？睖y設計研究院有限公司,上海 200335;3.南京林業大學風景園林學院,江蘇 南京 210037;4.黃河水資源保護科學研究院,河南 鄭州 450004;5.江西省水利科學院防災減災與水工程安全研究所,江西 南昌 330029)

明渠交匯水流普遍存在于自然界和人工河道中。與順直明渠不同,明渠交匯口匯集不同方向的來流,水流頂托效應顯著、紊動摻混劇烈,復雜的水動力條件對水流中物質的輸運產生重要影響。因此,明渠交匯口水動力及物質輸移規律一直備受國內外學者的關注。

在明渠交匯口水動力特性方面,茅澤育等[1]將理論分析和物理模型試驗結合,研究發現明渠交匯角影響分離區的幾何尺寸,但對其平面幾何形狀的改變不明顯,與Best等[2]的水槽試驗及劉盛赟等[3]的水氣兩相流三維紊流k-ε數學模型研究結論一致;郭維東等[4]針對Y型對稱明渠交匯口建立了三維數學模型,得到了不同交匯角工況下分離區長度和寬度與流量比的線性擬合公式;Shakibainia等[5]針對干支流寬度比對交匯水流結構的影響開展數值模擬研究,結果表明,隨著支流相對寬度的減小,分離區變大,加速區變窄,水面變化加劇。雷諾應力模型(Reynolds stress model,RSM)在交匯水流水動力模擬方面具有廣泛的應用,周舟等[6-8]基于此類模型,結合流體體積(VOF)自由水面捕捉方法,對具有不同交匯角的明渠交匯口進行了數值模擬,研究了交匯角對交匯口壅水特性、流速分布特性、分離區幾何特征和位置以及二次流強度的影響,并繪制了分離區的三視圖。胥宗強等[9]采用大渦模擬方法建立了明渠交匯口的三維模型,定量研究了特定流量比條件下交匯角對交匯口下游流向渦以及分離區三維幾何特性的影響。已有的針對不同交匯角條件下明渠交匯水流水動力特性的研究表明,分離區的尺寸、形狀對交匯水流的頂托效應、下游有效過水寬度及流速分布均具有重要影響,但量化分離區幾何特征的參數尚不全面,有待進一步研究。在明渠交匯口物質輸運方面,魏娟等[10]采用RSM建立水氣兩相流數學模型,模擬了非等寬明渠交匯流的污染物分布并揭示了其三維特性,定性分析了污染帶形狀隨流量比、交匯角的變化規律;王平等[11-12]采用動床物理模型分析了清水和濁水在具有不同交匯角明渠中的沖淤特性,得出了減輕高含沙量河流交匯口淤堵程度的最優交匯角;Tang等[13]基于RSM建立了90°等寬明渠交匯的三維數學模型,研究發現流量比和床面形態是通過影響剪切層和螺旋流的強弱和位置,進而影響污染物混合的;陳凱霖等[14-15]建立了非等寬河流交匯物理模型,運用PIV技術分別研究了明渠交匯口流速分布和污染物擴散情況隨交匯角的變化;張濤[16]結合室內水槽試驗、數值模擬,研究了非等寬河流交匯區的水流結構與污染物輸運規律,得出了交匯角、動量比和寬深比與混合速率的函數關系,其結論在對汾河入黃口的實地觀測數據中得到了印證。雖然已有學者開展了明渠交匯口污染物輸運規律的相關研究,但對于污染物分布特征參數的研究尚不全面,與交匯角間的定量關系尚不明晰。

綜上所述,交匯角是影響交匯水流水動力、物質輸運的重要因素之一。本文針對明渠交匯流中的保守物質輸運問題,基于RSM建立明渠交匯流三維水動力-污染物耦合數學模型,旨在研究不同交匯角工況下明渠交匯口分離區的幾何特征和交匯口污染物分布,分析交匯角對污染物混合速率的影響,以期為交匯口污染物影響范圍劃定、水污染責任劃分等實際應用問題提供參考。

1 數學模型

1.1 控制方程

本文基于RSM建立三維明渠交匯流水動力-污染物耦合數學模型,模擬交匯流流速及污染物質量濃度分布。水流運動滿足質量及動量守恒方程,污染物質量濃度滿足物質輸運方程。RSM主要應用在雷諾應力具有各向異性的流動中,對復雜流動的模擬有較大優勢,已經在明渠交匯水流的數值模擬方面有所應用。該模型根據時均化法則,不采用渦粘性假設,直接建立并模化雷諾應力的輸運方程,并與耗散率ε的輸運方程聯立求解[17-18]。參照現有研究成果[6-8,18],模型參數取值分別為:Cs=0.11,C1=1.8,C2=0.60,Cε=0.15,Cε1=1.44,Cε2=1.92。

本文采用質量濃度通量的二階矩模式封閉物質輸運方程[17-18]。參照現有研究成果[18],模型參數取值分別為:Csc=0.11,C1c=3.0,C2c=0.33。

1.2 數值方法及邊界條件

控制方程的離散采用有限體積法,擴散項的離散采用中心差分格式,水動力模塊對流項的離散采用QUICK格式,污染物模塊對流項的離散采用二階迎風格式。采用PISO算法進行壓力場和速度場的解耦。

水氣交界面的模擬采用流體體積法,定義體積函數aw表示計算網格內水的體積占網格體積的比例,滿足如下方程:

(1)

式中:ui為流速;t為時間;xi為笛卡坐標,i=1,2,3分別代表x、y、z坐標方向。對于一個網格,若aw=1,則表示該網格被水充滿;若aw=0,則表示該網格內無水;若0

在明渠床面和側壁采用無滑移邊界條件,頂部采用對稱邊界。主渠出口的流速、雷諾應力和污染物質量濃度及其通量均采用零梯度邊界條件。主、支渠進口均采用速度邊界,主渠出口采用滿足靜水壓力分布的壓力邊界。根據經驗公式[10],確定入口邊界處的湍動能和耗散率值分別為

k=0.00375u2

(2)

ε=k1.5/0.4H

(3)

式中:u為斷面平均速度;H為斷面水深。

2 模型驗證

2.1 水動力模型驗證

本文采用Weber等[19-20]的90°等寬明渠交匯流物理模型試驗進行數學模型水動力模塊的驗證。計算區域(圖1)干支流渠道均為矩形斷面,主渠、支渠的長度分別為20W和10W(W為渠道寬度),渠道適當向上游方向延長以消除邊界條件的影響,并得到充分發展的紊流。采用結構化六面體網格剖分計算區域,并在交匯口、渠底、側壁和水氣交接界面附近進行了網格加密,總網格數為24.6萬個。

圖1 數學模型計算區域

在驗證工況下,流入交匯口的流量比q=Qu/Qd=0.75(Qu為主渠道的入口流量,Qd為下游出口流量,Qd=0.17m3/s),尾水水位H0=0.31m,弗勞德數Fr=0.37,雷諾數Re=186000。數學模型的Fr和Re與物理試驗保持一致,模擬結果以渠道寬度W和主渠出口的平均縱向流速u0進行無量綱化處理。經過計算,將縱向斷面y/W=0.167處的水面線、橫向斷面x/W=6處測點流速的物理模型測量值與數學模型計算值進行對比,結果分別見圖2、圖3。結果表明,水面線相對誤差為1.7%,3個測點的流速相對誤差分別為3.6%、1.4%和2.6%,說明本模型模擬獲得的交匯口水面線與流速分布與物理模型基本一致,能夠比較準確地反映明渠交匯口的水動力特性。

圖2 y/W=0.167縱向斷面水面線測量值與計算值對比

圖3 x/W=6處測點流速測量值與計算值對比

2.2 污染物輸運模型驗證

本文采用陳凱霖[21]的明渠交匯流物理模型試驗進行保守物質輸運模塊數學模型的驗證,不考慮污染物的降解作用。主渠、支渠寬度分別為W、W/3,其他幾何參數與水動力模塊驗證所用模型相同,并用19.2萬個結構化六面體網格劃分計算區域。在驗證工況下q=0.81,Fr=0.10,Re=6357.79。支流進口污染物質量濃度ρ0為2000μg/L,主流進口處無污染物,與物理模型試驗保持一致。

圖4為驗證工況下y/W=0.167和y/W=0.5縱斷面上z/W=0.427深度處污染物質量濃度的物理模型測量值與數學模型計算值對比,兩者分布趨勢基本一致,兩縱向斷面數學模型計算值和物理模型測量值平均誤差分別為6.8%、7.3%,數學模型可以比較準確地模擬明渠交匯口的污染物分布。

圖4 z/W=0.427處污染物質量濃度測量值與計算值對比

3 不同交匯角工況下明渠交匯水流水動力與污染物輸運特性

3.1 工況設置

為研究明渠交匯口污染物質量濃度三維分布規律及交匯角對其影響,在水動力模塊、污染物輸運模塊驗證成果的基礎上,針對45°、60°、75°和90°交匯角建立等寬矩形斷面明渠交匯口數學模型,各工況水動力參數均設置為Fr=0.37、Re=186000、Qd=0.17m3/s,交匯口流量比q=0.75。支渠入口污染物質量濃度ρ0均為1μg/L,主渠入口處無污染物進入計算區域。

3.2 縱向流速分布及分離區三維幾何特性

圖5 交匯口下游斷面縱向流速分布

圖6 交匯口下游斷面平均縱向流速梯度沿程分布

圖7為不同交匯角工況下近水面處的縱向流速分布云圖和平面流速矢量圖(為方便顯示,每4個網格顯示1個矢量箭頭)。不同交匯角工況下流場存在不同特征:在明渠交匯口處,交匯角決定了支渠水流的偏轉角度,交匯角越大,支渠水流流向偏折得越明顯、流速沿主渠展向(y方向)的流速分量越大,對主渠水流的擠壓效應越強;在交匯口下游,由于支渠水流對主渠的擠壓,形成了分離區、最大流速區、水流恢復區等流動特征區域,交匯角越大,支流側的分離區尺度越大,對側最大流速區的收窄效應越明顯、流速越大;交匯口下游較遠處為水流恢復區,流速的橫向梯度沿程減小,即流速在橫向上有恢復均勻分布的趨勢,交匯角越大,流速分布恢復均勻需要的距離越長。

圖7 交匯口近水面處(z/W=0.3)流速分布

在實際工程中,河道順主流流速對物質輸運有較大影響,通常采用順主流流速為0的等值線定義分離區的幾何邊界,這種方法稱為等值線法[22]。等值線法已在分離區的幾何特性研究中有所應用[1,8-9],為確保與已有試驗結果一致,本文也采用等值線法劃定分離區范圍,如圖8所示。等值線法劃定的分離區呈翼形,θ=45°工況未見明顯分離區。定義分離區邊界線在主渠流向的投影長度為分離區長度L,在主渠展向的投影長度為分離區寬度B,本文通過改進傳統的分離區幾何對稱系數R[9],使得在評估分離區的幾何對稱性的同時,也考慮了長度方向中心位置與最大寬度位置的相對關系,改進后的對稱系數Rm定義如下:

圖8 表層水體分離區形狀

Rm=exp{[(x1+x2)/2-xB]/L}

(4)

式中:x1、x2分別為分離區投影與主渠左岸側壁交點的坐標值;xB為分離區最大寬度位置的x軸坐標值。改進對稱系數的目的是定量表示分離區中心位置和分離區最大寬度位置xB的相對位置關系,若Rm=1,則分離區最大寬度位置與其長度方向中心位置重合,此時可認為其在x方向上為幾何對稱圖形;若Rm>1,則分離區最大寬度位置在其長度方向中心位置的下游側;若Rm<1,,則分離區最大寬度位置在其長度方向中心位置的上游側。

根據前文提到的分離區定義方法,從圖7、圖8可以得到定性結論:分離區的幾何尺寸隨著交匯角的增大而增大。在θ=45°工況下,分離區尺寸太小難以捕捉;在θ=60°工況下,分離區的絕大部分存在于中上部水體中,下部水體中的分離區同樣尺寸極小;在θ=75°和θ=90°工況下,均形成了貫穿水體各個深度的分離區,且具有一定的規模。圖9為θ=60°工況下的分離區幾何參數。由圖9(a)可知,分離區長度L在z/W≥0.05以上隨著高度的減小緩慢減小,靠近床面時(z/W≤0.05),近壁流體受床面摩擦影響導致L逐漸增大,并在壁面達到最大值。由圖9(b)可知,分離區寬度B總體呈現上大下小的特征,在水面和床面附近隨水深的變化較小,在水體中部變化較明顯。θ=60°時,L、B總體上均隨著水深的增加而減小,直至z/W≈0.13處減小到0。由圖9(c)可知,分離區寬長比B/L隨高度增大而增大,僅在近水面(z/W≥0.31)處有小幅度減小,即分離區在近床面處更狹長。除下部水體,交匯角越大,分離區寬長比也越大,可能與壁面摩擦有關。另外,陳凱霖[21]的研究中提到,支渠中靠近床面的水流與主渠流向夾角比表面小,因此下部水體分離區比中上部更為狹長。由圖9(d)可知,床面和自由水面附近分離區最大寬度位置偏向中心位置的上游側,中部水體分離區最大寬度位置偏向中心位置的下游側。同一水平位置上,交匯角的增大會造成分離區最大寬度位置向其中心位置的上游側移動,總體來說,θ=90°工況下各水平面的分離區對稱系數Rm更接近1,在x軸方向上的對稱性更強。

圖9 交匯口分離區幾何參數特征

以上分析與曾誠等[8-9]的研究結論一致,當交匯角越大,分離區尺寸越大,沿長度方向的中心位置越偏向下游。分離區的存在顯著影響了明渠的水流結構,如使局部流速減小、甚至出現回流,占據主渠交匯口下游的有效過水斷面、使主渠水流被迫收窄和加速等。不同交匯角工況下,分離區幾何尺寸和形狀上的差異對交匯口污染物輸運的影響不可忽視。

3.3 污染物質量濃度三維分布

在交匯水流中,污染物質量濃度分布受水動力影響呈現明顯的三維特征。圖10分別展示了交匯角為45°、60°、75°和90°時支流污染物進入主渠后下游各橫斷面的污染物質量濃度分布。支流高濃度污染物進入主渠后,部分聚集在支渠側,部分隨著支流的初始動能沖向主渠右岸,沿主渠橫向斷面形成明顯的質量濃度梯度;污染物質量濃度等值線在靠近交匯口處(x/W>1)總體保持垂直,垂向上變化不明顯,僅在底層略向右岸凸出;隨著干支流的進一步摻混,在遠離交匯口的下游區域,質量濃度等值線逐漸發生傾斜、扭曲,下部水體中污染物質量濃度等值線的凸出更加明顯,出現了垂向質量濃度梯度。交匯口下游下部水體污染物質量濃度等值線的凸出與斷面環流有密切聯系。

圖10 交匯口下游橫斷面污染物質量濃度分布

圖11為主渠x/W=8斷面流速矢量和污染物質量濃度分布。由圖11可見,在不同交匯角工況下主渠交匯口下游均形成了逆時針的斷面環流,靠近渠道左岸的污染物在二次流的作用下向底層水體輸運,下層污染物向右岸方向輸運,上層污染物向左岸方向輸運;質量濃度等值線在水體底層向右岸傾斜,在水體上部靠近左岸。在靠近渠道底面位置,受壁面摩擦阻力的影響污染物輸運受阻,橫向輸運距離較小。隨著交匯角變大,斷面環流的位置由斷面左側、下側向中部移動,使得質量濃度等值線傾斜、扭曲的程度更大,混合層的寬度越寬。

圖11 主渠x/W=8斷面流速矢量和污染物質量濃度分布

3.4 污染物質量濃度橫向輸運特性

不同交匯角工況下交匯區污染物橫向輸運特性不同。圖12為不同交匯角工況下交匯口表層水體的污染物質量濃度分布,定義5%ρ0質量濃度等值線為污染帶邊界[23],污染物質量濃度大于5%ρ0的區域稱為污染帶,其最大寬度為b。模擬結果表明,污染帶最大寬度隨著交匯角的增大逐漸變寬。符合陳凱霖[21]的研究中交匯角越小污染帶越狹長的結論,與前文中分離區的幾何形狀變化規律也一致,說明分離區具有引起污染物滯留的效應。當θ≥75°時,分離區回流現象明顯,挾帶部分污染物沿主渠流向的反向運動,使高污染物質量濃度區域(ρ/ρ0≥0.9)在分離區附近顯著加寬,隨后在其下游收窄,但污染帶最大寬度b在下游沿程基本不變,如圖12(c)(d)所示。在其他水平面上,不同交匯角對污染物橫向輸運距離的影響也表現出相似的規律。圖13為不同水深處污染帶寬度b/W與交匯角的關系曲線。計算結果表明,越靠近水面,污染帶寬度b/W越大,且在不同深度處b/W均隨著交匯角的增大而增大,增大的幅度與交匯角有密切關系。當45°≤θ≤60°和75°≤θ≤90°時,b/W隨交匯角增大的幅度較小;當60°<θ<75°時,b/W隨交匯角增大的幅度較大。這說明交匯角在60°到75°之間時,b/W對交匯角的變化較為敏感。污染物混合界面[23]定義為10%質量濃度等值線與90%質量濃度等值線包圍區域(見圖11)。計算結果表明,不同交匯角工況下無量綱化的污染物混合界面寬度bMI/W沿程變化曲線呈現相似的變化趨勢,即隨著干支流的混合,混合界面寬度沿程逐漸增大;交匯角越大,同一位置上的污染物混合界面寬度越大(見圖14)。這與張濤[16]、陳凱霖[21]的試驗研究結果基本一致。

圖12 交匯口表層水體污染物質量濃度分布

圖13 不同水深處污染帶寬度與交匯角的關系

圖14 近水面處混合界面寬度沿程變化

不同交匯角工況下污染物質量濃度橫向輸運特性具有顯著差異的主要原因是,在流量比相同時,交匯角差異使支渠來流動量在主渠橫向上的分量不同,主渠交匯口下游形成了尺寸不一的流動特征區域。其中分離區流速顯著下降,出現回流現象,使污染物容易在該區域滯留;分離區和流動加速區之間流速差大,存在剪切層,對污染物的橫向輸運起到抑制作用;流速恢復區流速分布逐漸趨于均勻,污染物質量濃度的橫向梯度逐漸減小。

3.5 污染物質量濃度縱向輸運特性

圖15 交匯口下游斷面污染物平均質量濃度沿程變化

本文采用完全混合偏差δ和斷面標準差σ結合的方式量化給定斷面的污染物混合程度[24]。完全混合偏差可定量表示網格質量濃度值ρi與完全混合后的理論質量濃度值ρp間的偏差[25]。理論質量濃度與完全混合偏差的計算公式如下:

(5)

(6)

式中:Qm、Qt分別為主流和支流中的流量;ρm、ρt分別為主流和支流中的污染物質量濃度。本文計算工況下ρp=0.25,δi=300%和δi=-100%分別表示該網格中的水為完全未混合的支渠污水和干流清水,δi=0%表示經完全混合的水。圖16為交匯口下游一系列斷面上完全混合偏差極值δmin和δmax的統計數據,結果表明,隨著污染物從支流進入干流,δmax沿程先升高后下降,δmin持續升高,兩者約在x/W=1斷面下游逐漸靠近并趨向于0;在計算區域出口斷面(x/W=10)上污染物質量濃度的差距已經很小,但δmax和δmin的絕對值均稍大于10%,尚未達到完全混合的標準[25]。

圖16 交匯口下游污染物質量濃度完全混合偏差極值沿程變化

污染物質量濃度標準差σ可定量表示給定斷面上污染物質量濃度差異。標準差越接近0,代表該斷面污染物質量濃度的差異越小,污染物分布越均勻、混合越充分。圖17為不同交匯角工況下交匯口下游斷面標準的沿程分布。由圖17可見,不同交匯角工況下,標準差的沿程變化呈現相似的變化規律:在交匯口上游頂點斷面(x/W=0),σ較小但不為零,且隨交匯角的增大而減小;隨著污染物進入干流,σ沿程增大,并在x/W=1斷面附近達到峰值,峰值隨交匯角的增大而增大,其位置隨交匯角增大有向上游移動的趨勢,這與坐標原點定義在交匯口上游頂點,故下游頂點的x坐標隨交匯角的增大而減小有關;最后σ沿程逐漸減小并趨于0。根據統計,在計算區域出口(x/W=10)處各工況的σ已小于0.03,此時,污染物已經進行了一定程度的混合,同一斷面污染物質量濃度差距較小。

圖17 交匯口下游污染物質量濃度標準差沿程變化

污染物的混合速率可用斷面標準差沿程梯度的絕對值|Δσ/Δx|定量表示[24]。支流污染物完全匯入主渠后,統計各工況沿程斷面的混合速率,如圖18所示。結果表明,混合速率沿程減小,數量級在10-2以下,與Ramn等[24]的研究結果一致。交匯角的大小僅能影響分離區及其上游有限范圍的混合速率,在0

圖18 交匯口下游污染物混合速率沿程變化

|Δσ/Δx|=0.039exp(-0.3995x/W)

(2.5

(7)

4 結 論

a.較小交匯角工況下分離區中上部更狹長,較大交匯角工況下分離區下部更狹長;由改進的對稱系數量化結果,較大交匯角工況下分離區幾何對稱性更強,且會引起分離區最寬位置向上游偏移。

b.交匯口下游斷面污染物質量濃度水平方向存在明顯梯度,在交匯口下游逐漸形成垂向質量濃度梯度,污染物質量濃度分布存在顯著三維特性;不同交匯角工況下交匯口下游產生斷面環流的位置不同,導致污染物空間分布存在顯著差異。

c.交匯角越大,污染帶最大寬度越大,混合界面寬度增大;分離區及其上游污染物混合速率隨交匯角的增大而增強,在流速恢復區及其下游,混合速率基本不受交匯角影響,其沿程變化規律符合指數函數關系。