殷 杰， 熊 庭

(武漢理工大學(xué) a. 交通學(xué)院； b． 能動學(xué)院， 武漢 430063)

固液兩相流管道水力輸送是現(xiàn)代疏浚工程中疏浚物料運輸最為常用的方法。在疏浚物泥漿輸送系統(tǒng)中，局部組件如閥門與彎頭等作為管路系統(tǒng)的重要組成部分，是整體管線阻力損失的重要組成部分。[1]N?RGAARD等[2]已通過試驗和理論計算的方法得出90°彎管、漸擴(kuò)管、漸縮管及分流三通內(nèi)顆粒流體流過時的具體局部阻力系數(shù)。但目前對于彎管流動形態(tài)及阻力特性的理論研究較少，試驗相對困難，理論模型并不全面。

隨著計算機技術(shù)的日趨成熟，計算流體動力學(xué)(Computational Fluid Dynamics, CFD)數(shù)值模擬方法可對局部組件內(nèi)的流動特性進(jìn)行較好地描述。Euler-Euler雙流體模型作為經(jīng)典模型被廣泛使用， CHEN等[3]利用基于CFD-DEM耦合的數(shù)值模擬方法對45°、60°和90°彎管內(nèi)液-固兩相流動的腐蝕磨損率、湍流強度和二次流速度矢量進(jìn)行預(yù)測。ZHANG等[4]使用離散單元法(Discrete Element Method, DEM)描述顆粒運動軌跡和顆粒間相互作用，并以顆粒與壁面間相互作用力大小描述壁面磨損情況。CAO等[5]使用歐拉-拉格朗日模型研究垂直彎管中的液固兩相流動流態(tài)、沖蝕磨損率等情況。眾多CFD數(shù)值仿真方法均可描述彎管內(nèi)流動進(jìn)行一定范圍內(nèi)描述，但鮮有研究具體表述各工況下泥漿彎管流動特性、阻力特性的具體分析。

基于兩相流宏觀連續(xù)介質(zhì)理論及顆粒動力學(xué)理論，本文采用無相變過程的歐拉雙流體模型，運用CFD 模擬方法研究垂直彎管內(nèi)泥漿流體的三維流動及阻力特性，并給出阻力特性評判方法。

1 數(shù)學(xué)描述

為精確地描述液固兩相流的流動特征，適應(yīng)較高濃度泥漿流的模擬，并且在準(zhǔn)確性和計算量之間達(dá)到平衡，采用基于顆粒流動力學(xué)理論的歐拉-歐拉雙流體模型進(jìn)行研究。

1.1 控制方程

歐拉雙流體模型將固相顆粒擬流體化，把顆粒相也看作連續(xù)流體相一樣的連續(xù)介質(zhì)。質(zhì)量守恒方程為

(1)

(2)

(3)

(4)

式(2)～式(4)中：p為流體靜壓力，Pa；g為重力加速度，m/s2；μf為液相黏度，Pa·s；μs為固相剪切黏度，Pa·s。

1.2 輸送方程

使用顆粒動力學(xué)理論[6](Kinetic Theory of Granular Flow，KTGF)模型來封閉守恒方程。顆粒流動力學(xué)理論引入“顆粒擬溫度”參量描述顆粒的相對無序化的運動過程。顆粒擬溫度輸運方程為

(5)

式(5)中：kΘs為顆粒能量擴(kuò)散系數(shù)；γΘs為碰撞能量耗散項；ρs為固相壓力。

1.3 湍流模型

使用基于雷諾時均方法的k-ε混合物湍流模型來求解液固兩相流動的湍流特征。湍動能k的輸運方程為

(6)

湍動能耗散率的ε輸運方程為

(7)

式(7)中：μt,m為混合物湍流黏度；Gk,m為湍動能的生成項。在k-ε混合物湍流模型中，各湍流常數(shù)的取值分別為：C1ε=1.44,C2ε=1.92,Cμ=0.09,σk=1.0,σε=1.3。

2 仿真方法

2.1 物理模型

本文所研究的排泥管線爬坡管段流體域物理模型見圖1。上部水平段為橡膠自浮管，自浮管與一段半自浮管(彎頭1)相連過度到水下，半自浮管與波紋管(彎頭1與彎頭2連接段及彎頭2均為波紋管)連接構(gòu)成爬坡管，爬坡管彎頭2與沉于海底的金屬沉管相連接，各部分具體參數(shù)見表1。當(dāng)爬坡高度h變化時，彎頭1和彎頭2發(fā)生彎曲，彎頭曲率半徑改變；連接段直管長度不變，傾斜角度隨之改變。

在同時兼顧計算結(jié)果精確性與計算資源經(jīng)濟(jì)性的前提下，通過對仿真結(jié)果的分析及對比，經(jīng)過網(wǎng)格獨立性驗證，在爬坡管數(shù)值模擬中所使用的計算網(wǎng)格及網(wǎng)格橫截面見圖2，管道內(nèi)流體域采用六面體網(wǎng)格，劃分網(wǎng)格單元數(shù)約為174萬。另外，考慮以下基本物質(zhì)的性質(zhì)：主相，海水，密度ρf=1 025 kg/m3，動力黏度μf=0.001 174 Pa·s；第二相，顆粒，粒徑0.5 mm，密度ρs=2 470 kg/m3。內(nèi)摩擦角選擇30°。

幾何參數(shù)符號參數(shù)值管道內(nèi)徑/mD0.9水平浮管長度/mL150水平沉管長度/mL250爬坡管傾角/(°)β45爬坡高度/mh5彎頭1曲率半徑/mRc12.7彎頭2曲率半徑/mRc22.7

2.2 邊界條件

速度入口條件應(yīng)用于入口。兩個相的速度和濃度是給定的值。在出口處施加大氣壓力的壓力出口條件。壁粗糙度設(shè)定為0.02 mm，液相使用無滑移邊界條件，粒子使用Johnson和Jackson部分滑移邊界條件，并將鏡面反射系數(shù)和顆粒-壁面恢復(fù)系數(shù)分別設(shè)置為0.05和0.97。[7]具體參數(shù)設(shè)置見表2。

2.3 求解過程和收斂方案

使用商用CFD軟件Fluent 16.1來求解上述連續(xù)方程和邊界條件。使用均方根殘差，收斂殘差設(shè)為10-4。選擇相位耦合SIMPLE算法以確保結(jié)果穩(wěn)定且準(zhǔn)確并且獲得收斂。采用QUICK法求解動量方程。壓力松弛因子和動量松弛因子分別為0.2，體積分?jǐn)?shù)為0.3，并使用其他因子的默認(rèn)值，時間步長為0.1 s，總計算時長100 s。

表2 計算工況參數(shù)匯總

3 結(jié)果與討論

3.1 阻力特性分析

本文選取局部阻力損失系數(shù)ζ為表征爬坡管段阻力損失特性的參數(shù)，其計算為

(8)

式(8)中:Δz為爬坡管段排高水頭損失；Δpb為彎頭1入口面平均壓力與彎頭2出口面平均壓力之差；ρm為泥漿密度；Vm為輸送速度。

此外使用狄恩數(shù)來分析爬坡管段的湍流流動變化。狄恩數(shù)為流體流經(jīng)彎管時的離心力與黏性力之比，其數(shù)值等于雷諾數(shù)與管道半徑、彎頭曲率半徑比的平方根之積[8]，其表達(dá)式為

(9)

式(9)中：Rc為彎頭的曲率半徑，m。泥漿懸液的黏度由費祥俊公式[9]計算得到，即

(10)

式(10)中：μm、μ0分別為泥漿混合物和載流體的黏度；Cv、Cv,max分別為顆粒體積分?jǐn)?shù)和顆粒體積比濃度極限。

通過對爬坡管的爬坡高度h= 3.6 ～ 5.0 m時不同狄恩數(shù)(即彎頭曲率半徑變化、混合物流速和泥漿體積分?jǐn)?shù)不同)下管道內(nèi)泥漿液固兩相流的數(shù)值模擬，得到爬坡管段的局部阻力損失系數(shù)隨狄恩數(shù)的變化曲線見圖3。由圖3可知：因海洋潮汐引起的爬坡高度變化對管道的阻力特性有顯著的影響，當(dāng)?shù)叶鲾?shù)在2.06×106～ 3.73×106范圍內(nèi)變化時，爬坡管的局部阻力損失系數(shù)在1.20 ～ 3.25范圍內(nèi)變化；對同一爬坡高度在不同工況下，其局部阻力損失系數(shù)ζ均隨狄恩數(shù)De增加而減?。划?dāng)混合物流速Vm不變(圖中虛線連接的工況點流速相同)，局部阻力損失系數(shù)ζ隨De增加而增加；泥漿體積分?jǐn)?shù)Cv=20%的局部阻力損失系數(shù)與Cv=30%對應(yīng)的工況相比更大。

不同工況下局部阻力系數(shù)隨爬坡高度的變化曲線見圖4。由圖4可知：當(dāng)混合物流速Vm和泥漿體積分?jǐn)?shù)Cv相同時，局部阻力損失系數(shù)ζ隨爬坡高度h增加而增加。因為狄恩數(shù)表征狄恩渦的強度[10]，垂直主流方向的渦流運動加強會使得泥漿的主流運動受阻。另外，爬坡高度h越大，泥漿流過彎管時其狄恩流越強，泥漿流動受到的切向擾流越劇烈，從而增大了爬坡管處的局部阻力損失。[11]

3.2 二次流現(xiàn)象及對流場的影響

爬坡高度h=5 m、混合物流速Vm=6 m/s、泥漿體積分?jǐn)?shù)Cv=20%工況下，爬坡管段部分截面處混合物切向速度云圖和矢量圖見圖5，其中圖5a)～圖5f)分別為彎頭1入口、彎頭1出口、彎頭2入口、彎頭2出口、X=5D和X=20D處的截面。速度云圖與矢量圖相對應(yīng)，云圖反映垂直主流動方向的分速度數(shù)值大小，矢量圖反映分速度的方向和相對大小，明顯出現(xiàn)疊加于主流之上的流動，這稱為二次流現(xiàn)象，在管道內(nèi)橫截面內(nèi)形成反向的渦流，渦核的位置在管道截面中對稱分布。[12]

a) 彎頭1入口b) 彎頭1出口c) 彎頭2入口

d) 彎頭2出口e) X=5D處的截面f) X=20D處的截面

圖5 爬坡管段部分截面處混合物切向速度云圖和矢量圖

在彎頭1入口處泥漿流體產(chǎn)生切向分速度，二次流開始發(fā)展；在彎頭1出口處可以觀察到完全發(fā)展的二次流；在彎頭2出口處分速度值最大，二次流強度最強；泥漿在離開彎頭部分以后，不再受到離心力的作用，混合相垂直分速度逐漸減小，但由圖5e)、圖5f)可看出X=5D處二次流仍有一定存留，在X=20D處二次流已基本消失。爬坡管內(nèi)泥漿所受到的離心力沿流動方向不斷變化，當(dāng)流經(jīng)彎頭2時與彎頭1中離心力方向相反，彎頭2入口面二次流強度較彎頭1出口明顯降低，在流過彎頭2后渦流方向改變。

隨著二次流的進(jìn)一步減弱，沉管內(nèi)的泥漿在湍流作用下流動形態(tài)逐漸恢復(fù)[13]，在X=20D處二次流已經(jīng)幾乎消失，管道內(nèi)流態(tài)恢復(fù)到直管中充分發(fā)展形態(tài)。[14]同樣，在其他工況下，各二次流渦流處也均可觀察到類似現(xiàn)象，這說明管道內(nèi)流動形態(tài)受到二次流的影響并與其有對應(yīng)關(guān)系，二次流現(xiàn)象是爬坡管中十分重要的現(xiàn)象。

3.3 體積分?jǐn)?shù)分布

粒徑dp= 0.5 mm時混合物流速和爬坡高度不同，泥漿流經(jīng)爬坡管段后截面垂直中心線上顆粒的體積分?jǐn)?shù)見圖6和圖7。

從圖6和圖7中可知：當(dāng)混合物流速Vm或者爬坡高度h變化時，管段后泥沙顆粒的體積分?jǐn)?shù)分布存在明顯的差異，其中混合物流速的影響較強，爬坡高度的影響較弱。當(dāng)混合物流速或者爬坡高度增加時，泥漿在垂直主流方向的二次流動加強泥漿的混合，故顆粒的體積分?jǐn)?shù)梯度降低。

由上述分析可知，在X= 20D處管道內(nèi)的二次流已經(jīng)基本消失，顆粒的體積分?jǐn)?shù)分布基本恢復(fù)到直管段充分發(fā)展態(tài)，對比圖6b)、圖6c)和圖7b)、圖7c)X= 5D和X= 20D處各曲線可看出，在X= 5D處的體積分?jǐn)?shù)曲線較X= 20D處仍有較明顯的變形。由于在X= 5D處泥漿不再受到離心力的作用，顆粒分布的改變主要由殘留的二次流渦流作用引起，故從顆粒體積分?jǐn)?shù)分布的梯度可反映出殘留二次流的強度。

4 結(jié)束語

本文以挖泥船排泥管道的爬坡管內(nèi)泥漿流動過程為研究對象，基于歐拉雙流體模型，利用CFD方法探究爬坡管內(nèi)泥漿流動的阻力特性和以及二次流現(xiàn)象及其對流動的影響，得到以下結(jié)論。

1) 當(dāng)爬坡高度不變時，爬坡管段的局部阻力損失系數(shù)均隨狄恩數(shù)增加而減?。划?dāng)泥漿成分和輸送速度不變時，爬坡管段的局部阻力損失系數(shù)隨爬坡高度近似為線性增長，且增長速率隨混合物流速或泥漿體積分?jǐn)?shù)的減小而增加。

2) 泥漿在爬坡管內(nèi)產(chǎn)生二次流，在爬坡管的彎頭1和彎頭2內(nèi)泥漿受方向相反的離心力作用，在兩彎頭的出口處二次流渦流方向相反。管道內(nèi)泥漿的流動形態(tài)受到二次流的影響，顆粒體積分?jǐn)?shù)分布變化規(guī)律與二次流渦流流向一致，其中在渦核處此變化趨勢尤為明顯。

3) 當(dāng)混合物流速或爬坡高度增加時爬坡管內(nèi)二次流強度變強，泥漿混合更充分，爬坡管出口處顆粒體積分?jǐn)?shù)分布的不對稱度降低。當(dāng)泥漿體積分?jǐn)?shù)增加時，雖然二次流強度進(jìn)一步增加，但重力的影響起到了主導(dǎo)作用，使得高體積分?jǐn)?shù)泥漿在沉管中較快恢復(fù)到直管中充分發(fā)展態(tài)。