靜水條件下垂直圓管負(fù)浮力射流的數(shù)值模擬

胡玉珊，張仁貢，孫 燦，鄭福祥，王 盛，尹則高

(1.浙江同濟(jì)科技職業(yè)學(xué)院，浙江 杭州 311231；2.浙江禹貢信息科技有限公司，浙江 杭州 310000；3.中國海洋大學(xué) 工程學(xué)院，山東 青島 266100；4.中交疏浚技術(shù)裝備國家工程研究中心有限公司，上海 201208)

1 研究背景

負(fù)浮力射流是自然界和水利工程中常見的復(fù)雜水動力現(xiàn)象，例如污水在近海海域的排放。其主要運動特點是射流流體同時受到浮力和初始動量的雙重影響，且二者作用方向相反[1-2]。近年來，國內(nèi)外諸多學(xué)者對射流的機(jī)理和行為進(jìn)行了大量分析和研究。理論分析方面，Priestley[3]、Hirst[4]、槐文信等[5]假設(shè)高斯流速分布和線性卷吸率，結(jié)合整個流場的對稱性，研究了靜水環(huán)境下各種射流的時均特征；Lee等[6]建立了熱平流下單孔射流的高斯-渦積分模型，模型預(yù)測結(jié)果與試驗結(jié)果吻合較好。物理試驗方面，李晶金等[7]利用粒子圖像測速技術(shù)(particle image velocimetry,PIV)得到了負(fù)浮力射流在均勻環(huán)境水體下的最大上升高度，并給出了無量綱最大上升高度的表達(dá)式；李瓊[8]研究了均勻和分層環(huán)境流體中負(fù)浮力射流的運動規(guī)律，揭示了負(fù)浮力射流流動行為規(guī)律和卷吸特性；Milton-Mcgurk等[9]利用PIV技術(shù)和平面激光誘導(dǎo)熒光技術(shù)(planner laser induced fluorescence,PLIF)分別觀測了射流流體的速度和濃度的沿程衰減過程，發(fā)現(xiàn)負(fù)浮力射流的卷吸系數(shù)小于純水射流的卷吸系數(shù)；Burridge等[10]發(fā)現(xiàn)負(fù)浮力射流的穩(wěn)態(tài)高度在大尺度渦旋的影響下會發(fā)生波動；Abessi等[11]給出了不同傾角下負(fù)浮力射流的稀釋度公式；張玉玲等[12]采用PIV技術(shù)和PLIF技術(shù)對波流環(huán)境下多孔垂向射流近區(qū)的速度場和濃度場進(jìn)行了測量，對比分析了不同波浪條件下射流運動軌跡、速度場、濃度場及渦量場特性，重點探討了波流和純流環(huán)境下射流運動特性的差異以及波浪作用對多孔射流發(fā)展的影響；丁宏偉等[13]利用PIV技術(shù)對靜水環(huán)境下雙孔浮射流的速度場進(jìn)行了測量，得到了軸線垂向速度的沿程分布。數(shù)值模擬方面，李蒙等[14]采用大渦數(shù)學(xué)模型對橫流環(huán)境下正弦振蕩射流進(jìn)行了三維模擬，分析了非振蕩射流與振蕩射流運動特性和稀釋規(guī)律的差異，并探討了振蕩參數(shù)對射流的影響；袁麗蓉[15]建立了靜水環(huán)境中垂向圓形射流的三維各向異性紊流模型，分析了射流軸線速度和濃度衰減的規(guī)律；尹則高等[16]分析了垂直圓管射流流體擴(kuò)散特征，研究了規(guī)則波對負(fù)浮力射流三維時均速度以及瞬時濃度的影響，并得到了不同波參數(shù)下的中軸線射流垂向速度的沿程衰減規(guī)律；張松磊等[17]運用Fluent軟件對非靜態(tài)環(huán)境中的負(fù)浮力射流進(jìn)行了模擬和研究，得到了不同時刻的水流流態(tài)圖以及整個流動過程的流線分布圖；Ardalan等[18]開展了靜水環(huán)境下傾斜射流的數(shù)值模擬，發(fā)現(xiàn)稀釋度與密度弗勞德數(shù)的比值可表達(dá)為射流初始浮力的二次函數(shù)；方神光等[19]對靜水中正方形孔口浮力射流特性進(jìn)行了研究，給出了處于自相似區(qū)內(nèi)射流中軸線上的濃度分布和溫度分布；諶霞等[20]建立了三維泥沙和鹽度輸移擴(kuò)散的RNGk-ε紊流模型，該模型考慮了泥沙和鹽度引起的浮力射流影響，利用該模型對低鹽度含沙水流排入某海域后的水力特性以及鹽度和含沙量擴(kuò)散對環(huán)境水體的影響進(jìn)行了模擬預(yù)測。

本文建立了靜水條件下垂直圓管負(fù)浮力射流的三維水動力數(shù)學(xué)模型和室內(nèi)物理試驗?zāi)Ｐ停梦锢砟Ｐ驮囼炗^測的垂向速度對數(shù)學(xué)模型可靠性進(jìn)行了驗證，進(jìn)而探討了射流初始速度、初始密度、噴口直徑和密度弗勞德數(shù)對負(fù)浮力射流最大上升高度、穩(wěn)態(tài)高度和中軸線垂向速度的影響，可以為工程應(yīng)用提供參考依據(jù)。

2 三維水動力數(shù)學(xué)模型

2.1 控制方程

本文基于FLOW-3D軟件，利用考慮密度變化的連續(xù)性方程、Reynolds時均N-S方程、RNGk-ε湍流模型建立了靜水條件下垂直圓管負(fù)浮力射流三維水動力數(shù)學(xué)模型，并利用VOF(volume of fluid)方法進(jìn)行自由液面的位置捕捉，具體的方程形式可參見文獻(xiàn)[16]。

2.2 初始和邊界條件

數(shù)值模擬在靜水水槽內(nèi)展開，計算開始前設(shè)定靜水部分水體體積分?jǐn)?shù)為1；水槽底板和邊壁均為墻邊界條件；基于質(zhì)量動量法生成負(fù)浮力射流源項，具體步驟可參見文獻(xiàn)[16]，本文不再贅述。

3 物理模型試驗及數(shù)學(xué)模型驗證

3.1 物理模型試驗

在中國海洋大學(xué)長25.0 m、寬1.0 m、高1.2 m的水槽中進(jìn)行室內(nèi)物理模型試驗，設(shè)水槽長為x方向，寬為y方向，高為z方向。試驗靜水深d=0.6 m，射流噴口距離水槽造波端和底板分別為10.0和0.1 m，噴口內(nèi)徑D=0.015 m，定義射流噴口圓形中心為坐標(biāo)原點。利用進(jìn)水管安裝的CN-8300型水泵將水箱中密度為ρj=1050 kg/m3的鹽水垂直射入密度為ρa=1000 kg/m3的清水中。在進(jìn)水管上安裝流量計進(jìn)行射流鹽水流量監(jiān)測，同時采用ADV(acoustic Doppler velocimetry)測量典型測點流速。物理模型試驗布置如圖1所示。

圖1 物理模型試驗布置示意圖(單位：m)

3.2 數(shù)學(xué)模型驗證

圖2(a)給出了射流初始速度wj=1.57 m/s，D=0.015 m，ρj=1050 kg/m3工況的無量綱垂向速度w0/wj的驗證，其中w0為中軸線上射流垂向速度。圖2(b)為距離噴口不同位置剖面水平線上w/w0分布，其中w為流場中某一點的平均垂向速度。由圖2可知，數(shù)學(xué)模型計算結(jié)果與物理模型試驗數(shù)據(jù)基本吻合。

4 結(jié)果分析與討論

4.1 數(shù)值模擬工況

利用驗證后的數(shù)學(xué)模型，對不同射流初始速度、初始密度和噴口直徑下負(fù)浮力射流水動力特性進(jìn)行探究，具體工況見表1，其中靜水深d=0.6 m，靜止清水密度ρa=1000 kg/m3。

4.2 射流場ρ/ρa值分布

圖3給出了射流初始速度wj=1.40 m/s、初始射流密度ρj=1050 kg/m3、射流噴口直徑D=0.015 m工況下，不同時刻xoz剖面上流體密度與清水密度的比值ρ/ρa分布。由圖3(a)～3(e)可知，射流主體垂直進(jìn)入環(huán)境水體后，不斷減速上升到最大上升高度Zm，射流流體在中軸線兩側(cè)呈對稱分布，以倒錐形向外擴(kuò)展；整個上升過程中，射流主體邊緣和頂部的密度比值不斷減小，這是由于射流主體的運動會帶動周圍水體的運動，促進(jìn)周圍水體與射流流體摻混。到達(dá)最大上升高度后，射流頂端流體在負(fù)浮力作用下開始下降，但下降動量被黏性效應(yīng)抑制，在射流上部形成了彌散的云狀結(jié)構(gòu)。由圖3(f)～3(i)可以發(fā)現(xiàn)，重力作用下彌散的云狀結(jié)構(gòu)緩慢下降，并與周圍水體混合，云狀結(jié)構(gòu)逐漸消失。隨后，在射流主體兩側(cè)形成對稱的下降流，射流主體坍塌，射流高度降低并保持到一個穩(wěn)態(tài)高度Zss。其他工況的xoz剖面上射流流體的分布也都有類似的規(guī)律，只是在射流初始速度和初始密度不同的工況下，射流主體的最大上升高度和穩(wěn)態(tài)高度不同。

圖2 垂直圓管負(fù)浮力射流的無量綱垂向速度驗證

表1 數(shù)值模擬工況表

4.3 射流初始速度wj對w0的影響

圖4(a)和4(b)分別給出了不同wj條件下最大上升高度時刻和穩(wěn)態(tài)高度時刻w0的沿程變化，其中ρj=1050 kg/m3，D=0.015 m。由圖4(a)可知，當(dāng)z<0.075 m(5.0D)時，不同工況中軸線上的w0衰減不大，這是由于此區(qū)域受射流初始動量的影響較大，浮力影響相對較小，同時射流流體與環(huán)境水體之間的黏性效應(yīng)對w0的衰減影響較??；當(dāng)z>0.075 m(5.0D)時，射流初始動量的影響減弱，浮力的影響逐漸增強(qiáng)，中軸線上w0衰減趨勢加快。由于各工況射流初始密度相同，浮力的影響差別不大，w0衰減趨勢也基本一致。隨著wj從1.0 m/s增加到1.5 m/s，射流初始動量不斷增加，射流最大上升高度由0.31 m增加到0.50 m。圖4(b)給出的穩(wěn)態(tài)高度時刻w0沿射流中軸線變化規(guī)律與圖4(a)相似；噴口斷面上的w0數(shù)值與wj基本相等；隨著遠(yuǎn)離噴口，穩(wěn)態(tài)高度時刻的w0值顯著降低；隨著wj由1.0 m/s增加到1.5 m/s，射流的穩(wěn)態(tài)高度由0.19 m增加到0.35 m。

圖4(c)和4(d)分別給出了最大上升高度時刻和穩(wěn)態(tài)高度時刻射流中軸線上w0/wj與z/D關(guān)系。當(dāng)z<5.0D時，射流受初始動量的影響較大，數(shù)據(jù)相對分散，規(guī)律性較差；這里只給出了z>5.0D區(qū)域內(nèi)，浮力占主導(dǎo)作用區(qū)域內(nèi)的數(shù)據(jù)。w0/wj沿z/D的衰減可用下列經(jīng)驗公式來描述：

(1)

式中：B和n為衰減系數(shù)和衰減指數(shù)。

由圖4(c)可以看出，各擬合曲線衰減系數(shù)B值介于8.35～9.83之間，n值介于-1.42～-1.25之間，變化幅度相對較小，說明射流初始速度wj對最大上升高度時刻的垂向速度衰減早晚和衰減率的影響有限。與文獻(xiàn)[16]得到的B值和平均n值相比，本文最大上升高度時刻的B值較小，n值較大，說明靜水條件下垂直圓管負(fù)浮力射流中軸線上垂向速度衰減較早，衰減率較小。由圖4(d)可以看出，隨著wj從1.0 m/s增加到1.5 m/s，擬合曲線的B值從25.86減小至7.95，n值則從-2.15增大至-1.23，說明穩(wěn)態(tài)高度時刻射流中軸線上垂向速度衰減較早，衰減率減小。圖4(c)和4(d)給出的n值均小于文獻(xiàn)[21]給出規(guī)則波工況下純水射流的n值-1.07，說明穩(wěn)態(tài)高度時刻射流中軸上垂向速度的衰減率大于規(guī)則波作用下純水射流的對應(yīng)值。

圖3 xoz剖面上不同時刻流體密度與清水密度比值ρ/ρa分布(ρj=1050 kg/m3，wj=1.40 m/s，D=0.015 m)

4.4 射流初始密度ρj對w0的影響

圖5(a)和5(b)分別給出了不同ρj條件下最大上升高度時刻和穩(wěn)態(tài)高度時刻w0的沿程變化，其中wj=1.30 m/s，D=0.015 m。由圖5(a)可知，射流到達(dá)最大上升高度時刻，w0沿程分布可近似分為3段：z<0.075 m(5.0D)時，射流初始動量作用較大，浮力作用相對較小，不同工況中軸線上w0衰減趨勢均較緩；0.075 m(5.0D)0.180 m(12.0D)時，射流初始動量的作用進(jìn)一步減弱，浮力對射流的作用開始占據(jù)主導(dǎo)地位，導(dǎo)致w0分布出現(xiàn)明顯差異；并且隨著遠(yuǎn)離噴口，不同工況速度分布差異增加。隨著射流初始密度由1 030 kg/m3增加到1 070 kg/m3，射流最大上升高度由0.50 m降低到0.33 m。圖5(b)給出穩(wěn)態(tài)高度時刻的w0沿射流中軸線變化規(guī)律與圖5(a)相似；噴口位置上的w0數(shù)值與wj基本相等；隨著遠(yuǎn)離噴口，穩(wěn)態(tài)高度時刻的w0值顯著減??；隨著ρj從1 030 kg/m3增加到1 070 kg/m3，射流穩(wěn)態(tài)高度由0.38 m降低到0.28 m。

圖5(c)和5(d)分別給出了z>5.0D范圍內(nèi)最大上升高度時刻和穩(wěn)態(tài)高度時刻射流中軸線上w0/wj與z/D關(guān)系。由圖5(c)可以看出，各擬合曲線衰減系數(shù)B值介于7.57～8.88之間，n值介于-1.29～-1.19之間，變化幅度相對較小，說明射流初始密度ρj對最大上升高度時刻的垂向速度衰減早晚和衰減率的影響不大。由圖5(d)可以看出，隨著ρj從1 030 kg/m3增加到1 070 kg/m3，穩(wěn)態(tài)高度時刻擬合曲線的B值從7.79增加到16.84，n值則從-1.21減小到-1.77，說明此時射流中軸線上垂向速度衰減較晚，衰減率增加。與文獻(xiàn)[16]得到的平均n值相比，本文最大上升高度時刻和穩(wěn)態(tài)高度時刻的n值均較大，說明靜水條件下垂直圓管負(fù)浮力射流中軸線上垂向速度衰減率較小。

圖5 不同ρj條件下射流中軸線上垂向速度w0沿程變化(wj=1.30 m/s，D=0.015 m)

4.5 射流噴口直徑D對w0的影響

圖6(a)和6(b)分別給出了不同噴口直徑D條件下最大上升高度時刻和穩(wěn)態(tài)高度時刻w0的沿程變化，其中ρj=1050 kg/m3，且射流流量保持恒定，均為0.277 L/s。由圖6(a)和6(b)可知，隨著噴口直徑D的增加，射流中軸線上w0相應(yīng)減小。隨著噴口直徑D由0.014 m增加到0.016 m，射流最大上升高度由0.50 m降低至0.42 m，射流穩(wěn)態(tài)高度由0.50 m降低至0.25 m。

圖6(c)和6(d)分別給出了z>5.0D范圍內(nèi)最大上升高度時刻和穩(wěn)態(tài)高度時刻射流中軸線上w0/wj與z/D關(guān)系。由圖6(c)和6(d)可知，隨著D從0.016 m減小至0.014 m，最大上升高度時刻擬合曲線的B值從8.76降低到6.27，n值則從-1.29增加到-1.07；穩(wěn)態(tài)高度時刻擬合曲線的B值從14.36降低到7.32，n值則從-1.61增加到-1.12，即衰減較早，衰減率減小。與文獻(xiàn)[16]得到的平均n值相比，本文最大上升高度時刻和穩(wěn)態(tài)高度時刻的n值均較大，說明靜水條件下垂直圓管負(fù)浮力射流中軸線上垂向速度衰減率較小。

圖6 不同D條件下射流中軸線上w0沿程變化(ρj=1050 kg/m3)

4.6 密度弗勞德數(shù)的影響

Zm∝M03/4/B01/2∝D1/2wj/(g′)1/2=D·Fr

(2)

圖7給出了射流無量綱最大上升高度Zm/D和無量綱穩(wěn)態(tài)高度Zss/D與密度弗勞德數(shù)Fr的關(guān)系以及擬合直線，其中Zm和Zss分別為射流最大上升高度和穩(wěn)態(tài)高度。由圖7可以看出，在本文密度弗勞德數(shù)研究范圍內(nèi)(10

圖8為射流最大上升高度與穩(wěn)態(tài)高度的比值λ(λ=Zm/Zss)與密度弗勞德數(shù)Fr的關(guān)系。由圖8可知，λ與Fr基本無關(guān)。圖8中方形點為Burridge等[10]給出0.4≤Fr≤45范圍內(nèi)的試驗數(shù)據(jù)，其擬合直線表達(dá)式為λ=1.43，與密度弗勞德數(shù)Fr無關(guān)。圓形點為10≤Fr≤23范圍內(nèi)本研究的數(shù)值計算值，對其進(jìn)行擬合得到λ=1.52，略大于Burridge等[10]給出的λ=1.43。原因是由于Fr<6.5范圍內(nèi)，Burridge等[10]給出的λ值很大一部分小于1而引起的。

圖7 無量綱最大上升高度Zm/D及穩(wěn)態(tài)高度Zss/D與密度弗勞德數(shù)Fr的關(guān)系

圖8 最大上升高度與穩(wěn)態(tài)高度比值λ與密度弗勞德數(shù)Fr關(guān)系

5 結(jié) 論

基于數(shù)值模擬等方法，對靜水條件下垂直圓管負(fù)浮力射流特性進(jìn)行研究，得到以下結(jié)論：

(1)隨著射流初始速度的增加，射流初始密度的降低和噴口直徑的減小，射流最大上升高度和穩(wěn)態(tài)高度均增加。

(2)靜水條件下，射流初始速度和射流初始密度對最大上升高度時刻的中軸線垂向速度衰減早晚和衰減率的影響不大；隨著射流初始密度的增加，穩(wěn)態(tài)高度時刻的射流中軸線上垂向速度衰減較晚，衰減率增加；垂直圓管負(fù)浮力射流量一定時，隨著噴口直徑的增加，最大上升高度時刻和穩(wěn)態(tài)高度時刻射流中軸線上垂向流速衰減較早，衰減率減小。與規(guī)則波作用下的垂直圓管負(fù)浮力射流相比，本研究中最大上升高度時刻射流中軸線上垂向速度衰減較早，衰減率也較??；與規(guī)則波工況下的純水射流相比，本研究中穩(wěn)態(tài)高度時刻射流中軸上垂向速度衰減率相對較大。

(3)建立了密度弗勞德數(shù)表達(dá)的無量綱最大上升高度和穩(wěn)態(tài)高度線性關(guān)系式；發(fā)現(xiàn)本研究中射流最大上升高度和穩(wěn)態(tài)高度比值趨近于常數(shù)1.52，與密度弗勞德數(shù)基本無關(guān)。需要指出的是，為研究方便，本文靜水深均取0.6 m，導(dǎo)致噴口直徑較小和射流初始速度工況下的最大上升高度和穩(wěn)態(tài)高度基本相等，未出現(xiàn)較大靜水深條件下前者大于后者的現(xiàn)象。因此，下一步需要深入研究不同水深對負(fù)浮力射流上升高度的復(fù)雜影響，并闡明其影響規(guī)律。