魏文禮，趙小軍，劉玉玲

(西安理工大學(xué) 陜西省西北旱區(qū)生態(tài)水利工程重點(diǎn)實驗室，陜西 西安 710048)

氣泡浮力羽流廣泛應(yīng)用于工業(yè)領(lǐng)域，如采用深層曝氣方式的污水處理，以及石油工程中的氣液分層等。在很多實際工程問題中，都需要對氣泡浮力羽流的流動特性進(jìn)行研究，討論其運(yùn)動特性以及對實際工程的影響。氣泡浮力羽流屬于典型的氣液兩相流動，且氣相和液相之間相互混摻劇烈。目前對氣泡羽流的研究主要采用試驗觀測的手段，在早期的研究中，主要是對氣泡羽流流動參數(shù)進(jìn)行測量，如速度、含氣率分布等[1-2]。近年來,隨著試驗條件以及測量技術(shù)的發(fā)展，對氣泡羽流流動特性的觀測和流動細(xì)節(jié)的捕捉也越來越精確，可以更準(zhǔn)確地得到不同試驗參數(shù)下氣泡羽流流型的變化和氣液混合效果的影響等[3]。試驗布置的不同(如羽流發(fā)生器傾斜角度不同)，也會使羽流產(chǎn)生不同程度的附壁效應(yīng),并會對羽流尾部的傾斜角度和范圍產(chǎn)生影響[4]。隨著計算流體動力學(xué)的發(fā)展，采用數(shù)值模擬的手段對氣泡浮力羽流進(jìn)行研究也被越來越多的學(xué)者所采用，模擬的結(jié)果也越來越精確。目前，模擬采用的方法主要有歐拉方法和拉格朗日方法，而以歐拉方法為基礎(chǔ)的雙流體模型得到了廣泛的運(yùn)用。根據(jù)相關(guān)學(xué)者的對比研究，認(rèn)為拉格朗日模型對低氣流率下的氣泡羽流模擬結(jié)果較好，而歐拉方法對不同氣流率下的氣泡羽流均能得到較好的模擬結(jié)果[5]。有學(xué)者通過對氣泡羽流模型進(jìn)行三維和二維簡化模擬，發(fā)現(xiàn)不同的簡化方式適應(yīng)的紊流模型也有差別[6]。Kerdouss等[7]對1個雙葉輪反應(yīng)器進(jìn)行了數(shù)值模擬，主要研究了反應(yīng)器內(nèi)部氣泡尺寸及氣泡的破滅與聚合過程，并與相關(guān)試驗進(jìn)行了比較，取得了良好的效果。苗一等[8]對渦輪攪拌器內(nèi)氣液混合過程進(jìn)行了大渦模擬，得到了混合時間的示蹤劑響應(yīng)曲線。李志鵬等[9]對攪拌槽內(nèi)單循環(huán)流動特性進(jìn)行了大渦模擬，并與相關(guān)試驗和模擬數(shù)據(jù)進(jìn)行了詳細(xì)比較，取得了良好的模擬結(jié)果。Kerdouss等[10]采用數(shù)值模擬的方法對攪拌反應(yīng)器內(nèi)的氣液傳質(zhì)系數(shù)進(jìn)行了預(yù)測，并與試驗值進(jìn)行了比較，表明試驗數(shù)據(jù)和數(shù)值模擬結(jié)果具有良好的一致性。Ahmed等[11]對多葉輪生物反應(yīng)器進(jìn)行了試驗和數(shù)值模擬研究，對反應(yīng)器內(nèi)流型變化進(jìn)行了觀測，并且測得了反應(yīng)器內(nèi)氣泡停留與氣液混合時間以及功率消耗等參數(shù)，最后與模擬結(jié)果進(jìn)行了比較，取得了很好的研究成果。李偉等[12]對Rushton攪拌釜進(jìn)行了數(shù)值模擬，分析了攪拌釜內(nèi)的氣液流動特性以及氣含率的分布，并與試驗數(shù)據(jù)進(jìn)行對比，表明試驗值和模擬值吻合較好。

綜觀現(xiàn)有研究成果，并未見有關(guān)利用歐拉雙流體模型(Euler)和混合模型(Mixture)模擬曝氣池氣泡浮力羽流動力特性的差異的討論。為此，本研究以文獻(xiàn)[1]的曝氣試驗為基礎(chǔ)，分別采用歐拉雙流體模型和混合模型結(jié)合標(biāo)準(zhǔn)κ-ε紊流模型對該試驗進(jìn)行模擬，通過模擬結(jié)果與相關(guān)試驗數(shù)據(jù)的對比與分析，選擇模擬氣泡浮力羽流的適宜模型， 計算分析不同氣流率下氣泡羽流的流動特性，以期為提高曝氣池中深層曝氣方式的污水處理效率提供依據(jù)。

1 數(shù)學(xué)模型

1.1 歐拉雙流體模型

歐拉雙流體模型的守恒方程由質(zhì)量守恒方程和動量守恒方程組成。

質(zhì)量守恒方程：

(1)

(2)

(3)

動量守恒方程：

(4)

(5)

q相的剪切黏度μq的計算公式為：

(6)

式中：Cμ為常數(shù)，其值為0.09；kq為q相的紊動動能；εq為q相的紊動動能的耗散率。

q相的紊動動能kq和耗散率εq可以通過k-ε紊流模型獲得，其方程如下:

(7)

(8)

式中：Gk,q為q相紊動動能的產(chǎn)生速率；σk、σε、C1ε、C2ε為常數(shù)，其值分別為1.0，1.3，1.44，1.92[13]。

1.2 混合模型

混合模型利用小空間尺度上的局部平衡來求解混合相的動量方程、連續(xù)性方程和第2相的體積分率。

混合模型連續(xù)性方程為：

(9)

(10)

(11)

混合模型的動量方程可以通過對所有相各自的動量方程求和獲得，其形式如下：

(12)

(13)

混合模型第2相(用p表示)的體積分?jǐn)?shù)方程的表達(dá)式為：

(14)

混合模型中相對速度(滑移速度)被定義為第2相(用p表示)相對于第1相(用q表示)的速度，即有：

(15)

對于任意一相，如k相的質(zhì)量分?jǐn)?shù)可定義為：

(16)

(17)

紊動黏性項采用κ-ε紊流模型求得。

混合相紊動黏度μm的計算公式為：

(18)

式中：k為混合流體的紊動動能；ε為混合流體紊動動能耗散率；Cμ為常數(shù)，其值為0.09。

混合流體的紊動動能k和耗散率ε的表達(dá)方程如下：

(19)

(20)

式中：Gk為混合流體紊動動能的產(chǎn)生速率。

以上方程組構(gòu)成了求解流場分布規(guī)律的封閉方程組，根據(jù)實際工況施加相應(yīng)的邊界條件后，就構(gòu)成該方程組的定解問題。

2 物理模型與邊界條件

2.1 物理模型區(qū)域與計算網(wǎng)格

模擬采用模型的尺寸為0.5 m×0.08 m×2 m(圖1)，矩形通氣孔位于底部距離左壁面0.15 m處，其邊長為0.04 m。初始液面高度為1.5 m。通氣量為1.6 L/min，通氣速度為0.017 m/s。文獻(xiàn)[14-15]討論了初始液面高度與頂面空氣層高度之間比例關(guān)系對模擬的影響，根據(jù)其結(jié)論本次模擬采用的空氣層厚度對模擬結(jié)果的干擾可以忽略不計。根據(jù)氣泡羽流的流動特性以及多相流模型的適應(yīng)性，適合的多相流模型有混合模型和歐拉模型。本研究將采用這2種模型分別進(jìn)行模擬。模擬采用的模型區(qū)域圖和網(wǎng)格圖如圖1和圖2所示。

2.2 邊界條件

由于氣相在水中屬于分散相，而溢出水面以后，在空氣層中又屬于連續(xù)相。因此對頂面出口邊界的設(shè)置有幾種選擇，包括自由出流、壓力出口以及速度進(jìn)口。其中自由出流和壓力出口邊界都比較好理解，而在頂面采用速度進(jìn)口的情況卻比較少。在底部通氣速度很大，空氣從頂面溢出以后造成水面的劇烈波動，且底部空氣進(jìn)入空氣層區(qū)域以后變成連續(xù)相。在這種情況下，為避免模擬結(jié)果產(chǎn)生較大的誤差，可在頂部設(shè)置速度進(jìn)口，給一個較小的速度，以減小模擬誤差。就本研究而言，由于通氣速度較小，空氣溢出水面以后對水面的影響也很小。因此，結(jié)合試驗的布置方式和過程，采用速度進(jìn)口對通氣口進(jìn)行描述，頂面采用壓力出口邊界，模型底面和四周區(qū)域采用無滑移固壁邊界。

速度與壓力耦合求解使用了壓力隱式算子分裂PISO算法。

圖1 模擬模型區(qū)域圖

3 結(jié)果與分析

3.1 模型比較

本研究分別采用混合模型和歐拉模型進(jìn)行模擬。為了討論2種模型的精確性，就2種模型對氣泡羽流中氣相和液相之間滑移速度的捕捉效果進(jìn)行判別。由于氣泡羽流是基于底部通入空氣帶動區(qū)域里面的水進(jìn)行流動而形成的，因此氣相和液相之間有速度差異。圖3和圖4分別為采用混合模型和歐拉模型計算得到的z/H=0.35(H為物理模型的水體深度)時氣相的滑移速度。

圖3 采用混合模型和歐拉模型模擬的氣相滑移速度的分布

對比圖3中2種模型的模擬結(jié)果可以看出，用混合模型時，計算得到的滑移速度比較小，比采用歐拉模型時得到的滑移速度值小1個數(shù)量級，幾乎可以忽略不計?？梢?，歐拉模型較混合模型可以更好地捕捉氣液兩相之間的速度差異。

3.2 流動特性分析

圖4為20,60和100 s時計算得到的液相速度矢量圖。從圖4可以發(fā)現(xiàn)：速度矢量沿水深方向的增加趨勢并不穩(wěn)定，特別是在上部區(qū)域，氣泡羽流出現(xiàn)了擺動；此外，在靠下部區(qū)域，氣泡羽流的擺動不明顯，且這一部分的速度相對較大，說明該通氣速度對底部區(qū)域羽流形狀的影響不大。結(jié)合60和100 s時的速度矢量圖可以看出，擺動只存在于上半部區(qū)域，且羽流的尾部形狀時刻變化，這一規(guī)律與文獻(xiàn)[16]的結(jié)果一致。

圖5為20,60和100 s時的氣相體積分?jǐn)?shù)等值面分布圖。從圖5可以看出，隨著時間的推移，氣泡羽流上部的擺動對氣相的分布有較大影響?？諝怏w積分?jǐn)?shù)隨著羽流的擺動時刻改變，而在羽流的下部區(qū)域空氣體積分?jǐn)?shù)的分布比較穩(wěn)定。

圖4 不同時刻液相瞬時速度矢量圖

圖6為不同時刻液相軸向速度的分布圖。從圖6可以看出，在z/H=0.625的平面內(nèi)，不同時刻液相速度存在波動性，但是波動的峰值基本一致，且沿水平方向左右移動，這是由于羽流尾部的左右擺動造成的。

圖7表示在z/H=0.22和z/H=0.35 2個平面內(nèi)氣相速度的模擬值與試驗值的對比圖。從圖7-(a),(b)可以看出，氣相速度模擬值與試驗值吻合較好，從而驗證了歐拉雙流體模型結(jié)合標(biāo)準(zhǔn)κ-ε紊流模型的精確性。

圖6 不同時刻液相軸向速度的分布

圖7 不同z/H時氣相速度模擬值與試驗值的比較

4 結(jié) 論

本研究分別采用歐拉雙流體模型和混合模型結(jié)合標(biāo)準(zhǔn)κ-ε紊流模型對曝氣池中氣泡浮力羽流的動力特性進(jìn)行數(shù)值模擬，得到了氣泡羽流的流動特性為：速度矢量沿水深方向的增加是不穩(wěn)定的，特別是在上部區(qū)域氣泡羽流出現(xiàn)了擺動；在靠下部區(qū)域，氣泡羽流的擺動不明顯，但這一部分的速度相對較大。通過對模擬結(jié)果的比較與分析，認(rèn)為歐拉模型較混合模型能夠更加精確地捕捉氣液兩相之間的速度差異。通過模擬結(jié)果與試驗數(shù)據(jù)的比較，進(jìn)一步驗證了歐拉模型的精確性和適用性。

[參考文獻(xiàn)]

[1] Becker S,Sokolichin A,Eigenberger G.Gas-liquid flow in bubble columns and loop reactors:Ⅱ.Comparison of detailed experiments and flow simulations [J].Chemical Engineering Science，1994,49:5747-5762.

[2] Pfleger D,Gomes S,Gilbert N,et al.Hydrodynamic simulations of laboratory scale bubble columns fundamental studies of the Eulerian-Eulerian modelling approach [J].Chemical Engineering Science，1999,54(21):5091-5099.

[3] Amiri T Y,Moghaddas J S,Moghaddas Y.A jet mixing study in two phase gas-liquid systems [J].Chemical Engineering Research and Design，2011,89(3):352-366.

[4] 肖柏青，張法星，劉春艷，等.曝氣池內(nèi)氣泡羽流附壁效應(yīng)的試驗研究 [J].水力發(fā)電學(xué)報，2012(4):104-107.

Xiao B Q,Zhang F X,Liu C Y,et al.Coanda effect of bubble plume in aeration tanks [J].Journal of Hydroelectric Engineering,2012(4):104-107.(in Chinese)

[5] Ashraf Ali B,Pushpavanam S.Analysis of unsteady gas-liquid flows in a rectangular tank: Comparison of Euler-Eulerian and Euler-Lagrangian simulations [J].International Journal of Mu-ltiphase Flow，2011,37(3):268-277.

[6] 徐禮嘉，陳彩霞，夏梓洪，等.鼓泡塔內(nèi)氣泡羽流周期性擺動的數(shù)值模擬 [J].化學(xué)工程，2012(9):48-51.

Xu L J,Chen C X,Xia Z H,et al.Numerical simulation of bubble plume periodic oscillation in bubble column [J].Chemical Engineering,2012(9):48-51.(in Chinese)

[7] Kerdouss F,Bannari A,Proulx P.CFD modeling of gas dispersion and bubble size in a double turbine stirred tank [J].Chemical Engineering Science，2006,61:3313-3322.

[8] 苗 一，潘家禎，閔 健，等.渦輪槳攪拌槽內(nèi)混合過程的大渦模擬 [J].華東理工大學(xué)學(xué)報:自然科學(xué)版,2006,32(5):624-628.

Miao Y,Pan J Z,Min J，et al.Large eddy simulation of mixing process in stirred tank with rushton turbine [J].Journal of East China University of Science and Technology:Natural Science Edition,2006,32(5):624-628.(in Chinese)

[9] 李志鵬，高正明.渦輪槳攪拌槽內(nèi)單循環(huán)流動特性的大渦模擬 [J].過程工程學(xué)報,2007,7(5):900-904.

Li Z P,Gao Z M.Large eddy simulation of single loop flow field in a rushton impeller stirred tank [J].The Chinese Journal of Process Engineering,2007,7(5):900-904.(in Chinese)

[10] Kerdouss F,Bannari A,Proulx P,et al.Two-phase mass transfer coefficient prediction in stirred vessel with a CFD model [J].Computers and Chemical Engineering,2008,32:1943-1955.

[11] Ahmed S U,Ranganathan P,Pandey A,et al.Computational fluid dynamics modeling of gas dispersion in multi impeller bioreactor [J].Journal of Bioscience and Bioengineering,2010,109(6):588-597.

[12] 李 偉，劉躍進(jìn)，孔麗娜，等.Rushton攪拌釜內(nèi)的氣含率分布及其流動特性的模擬 [J].化工學(xué)報,2011,62(10):2691-2698.

Li W,Liu Y J,Kong L N,et al.Simulation for gas hold-up distribution and flow behavior in Rushton stirred tank [J].CIESC Journal,2011,62(10):2691-2698.(in Chinese)

[13] 魏文禮，李 超.錐閥附近水力特性的三維數(shù)值模擬研究 [J].計算力學(xué)學(xué)報，2011,28(1):108-112.

Wei W L,Li C.3D numerical simulation of hydraulic characteristics around fixed-cone valve [J].Chinese Journal of Computational Mechanics,2011，28(1):108-112.(in Chinese)

[14] Hu G,Celik I.Eulerian-Lagrangian based large-eddy simulation of a partially aerated flat bubble column [J].Chemical Engineering Science,2008,63(1):253-271.

[15] Bravo H R,Gulliver J S,Hondzo M.Development of a commercial code-based two-fluid model for bubble plumes [J].Environmental Modelling and Software,2007,22(4):536-547.

[16] Simiano M.Experimental investigation of large-scale three dimensional bubble plume dynamics [D].Zurich，Swiss:Swiss Federal Institute of Technology,2005.