基于格子Boltzmann方法的微通道內(nèi)氣液兩相流流型和壓力降特性研究*

高宇航 馮 凱 張會臣

(大連海事大學(xué)船舶與海洋工程學(xué)院 遼寧大連 116026)

微流控技術(shù)以其樣品消耗量小、高通量、易控制、傳熱傳質(zhì)效率高等優(yōu)點引起了人們的廣泛關(guān)注[1]。微通道內(nèi)氣液兩相流技術(shù)應(yīng)用廣泛，包括傳熱傳質(zhì)[2]、生物醫(yī)學(xué)[3]、化工生產(chǎn)[4]和食品加工[5]等諸多領(lǐng)域。隨著微流控器件的發(fā)展，對微通道內(nèi)氣液兩相流的流動特性進行研究尤為重要。

T型微通道是最常見的微通道構(gòu)型，是設(shè)計復(fù)雜結(jié)構(gòu)的基礎(chǔ)，形成氣-液、液-液兩相流或多相流，用于傳熱傳質(zhì)和混合強化等。微通道內(nèi)氣液兩相流流型和流場決定微器件的設(shè)計和功能實現(xiàn)。氣泡的頭部和尾部的流場、形狀和液膜厚度是影響傳熱傳質(zhì)的主要因素[6]。壓力降是反映流體流動性的重要參數(shù)之一。氣泡的生成和運動，與微通道內(nèi)的界面力和液體的表面張力密切相關(guān)，微通道內(nèi)壓力變化復(fù)雜，在氣相入口和液相入口，甚至于整個通道內(nèi)，壓力呈周期性的變化[7]。

格子Boltzmann方法(LBM)是一種介觀尺度的數(shù)值模擬方法，能夠很好地應(yīng)用于單相流及多相流的模擬。鄭馨瑤等[8]通過將非牛頓流體冪律模型引入牛頓流體格子Boltzmann模型的方法，數(shù)值模擬了矩形截面微通道內(nèi)剪切稀化流體的流動特性。白冰等人[9]利用單組分汽液兩相LBM模型與Laplace定律進行研究，結(jié)果表明液滴表面張力隨溫度升高而下降，與實際流體表面張力與溫度的關(guān)系定性上一致。婁欽等人[10]給出了在不同毛細數(shù)條件下，分離的牛頓流體液滴尺寸，并與數(shù)值結(jié)果進行對比，得到了一致的結(jié)果。FALLAH和RAHNI[11]采用LBM對普通和改進T型結(jié)構(gòu)微通道中液滴的形成過程進行了數(shù)值模擬，發(fā)現(xiàn)改進結(jié)構(gòu)中的擴展壁導(dǎo)致小回流區(qū)的產(chǎn)生和連續(xù)相方向的改變，增強了剪切過程，在相同條件下產(chǎn)生較小的液滴和較短的液滴間距。BA等[12]采用LBM模擬研究了對稱和非對稱T型結(jié)構(gòu)中液滴的形成。研究發(fā)現(xiàn)，在連續(xù)相中，由于分散相的阻塞，上游壓力大于下游壓力，壓差是液滴破碎的主要驅(qū)動力，主要由連續(xù)相產(chǎn)生的擠壓壓力所致。

然而，當(dāng)兩相流密度比較大時，模擬變得更為復(fù)雜。LEE和LIN[13]提出了一種非理想氣體的穩(wěn)定離散化方法，用于模擬具有高密度和高黏度比的不可壓縮兩相流，結(jié)果表明，在密度比為1 000的條件下與解析解吻合得很好，液滴碰撞后發(fā)射的液膜擴散系數(shù)也遵循已知的擴散冪律。SARITHA和BANERJEE[14]提出了一種基于偽勢的密度比達到1 000單組分多相LBM流動求解器，研究了雷諾數(shù)為1和10定常和非定常泊肅葉流動中單個空泡的流動特性。章詩婷等[15]提出一種相場簡化多相流格子Boltzmann方法，有效地模擬密度比為1 200和黏度比為500的復(fù)雜界面算例，如拉普拉斯定律等。綜上可以看出，目前對氣-液兩相流的數(shù)值模擬研究較為常見，但采用格子Boltzmann方法模擬具有大密度比的氣-液兩相流的研究較少，并且所建模型普適性差，模擬問題單一且不易收斂，許多問題亟待深入探討。

本文作者基于格子Boltzmann方法建立了密度比達到500的氣-液兩相流模型，研究了二維T型微通道內(nèi)氣液流速、毛細數(shù)、壁面潤濕性等因素對兩相流形成、壓力場和流速的影響，更加深入地解釋了微通道內(nèi)氣液兩相流流型和壓力降特性變化的規(guī)律和影響機制。

1 兩相流格子Boltzmann模型

格子Boltzmann方法是一種介于微觀分子動力學(xué)和宏觀流體動力學(xué)之間的一種介觀尺度的數(shù)值模擬方法，它可以看作是一個空間、時間和粒子速度都是離散的簡化虛擬分子動力學(xué)模型[16]。LBM回避了N-S方程中壓力項的隱式，減少了計算所需的時間。碰撞和遷移過程都是局部的，易于編程，并具有良好的并行性。在微觀層面引入分子間相互作用力，能夠模擬傳統(tǒng)CFD方法難以處理的多相流和組分流中的相分離和界面動力學(xué)。

基于ZHENG等[17]提出的大密度比格子Boltzmann模型，文中對T型微通道內(nèi)氣液兩相流運動特性進行模擬。通過分布函數(shù)gi和fi來求解流體界面和速度場，其對應(yīng)的演化方程為

(1)

τn、τφ是量綱一松弛時間，τn取決于流體運動黏度：

式中：υ是運動黏度；q是一個常數(shù)。

μφ為化學(xué)勢，可通過公式(4)計算：

μφ=A(4φ3-4φ*2φ)-κ?2φ

(4)

系數(shù)A和κ可通過公式(5)計算：

式中：σ為表面張力系數(shù)；W為界面層厚度；φ*為序參數(shù)的平衡狀態(tài)，可由初始流體密度確定，φ*=(ρl-ρg)/2。

式中：Γ為控制界面移動速率的遷移率系數(shù)；n=(ρl+ρg)/2。

宏觀參數(shù)，如序參數(shù)φ、密度n以及速度u可通過公式(10)—(12)計算。

模擬中使用的參數(shù)數(shù)值如表1所示。

表1 模擬參數(shù)

2 模擬結(jié)果分析

2.1 單液滴模擬及模型驗證

將液滴滴到水平壁面上，根據(jù)壁面的親水性或者不同液滴的性質(zhì)形成不同的接觸角。文中采用建立的大密度比LBM模型模擬了液滴的接觸角，模型如圖1所示。

圖1 液滴接觸角模擬示意

數(shù)值模擬的計算域為200×100，液氣兩相的密度比為500∶1，初始圓形液滴的直徑為20，與下壁面相切。通過調(diào)整序參數(shù)比形成不同的接觸角大小。當(dāng)達到穩(wěn)定狀態(tài)時，呈現(xiàn)出不同的接觸角形狀，如圖2所示。

圖2 靜止液滴接觸角模擬結(jié)果

當(dāng)序參數(shù)比為負時，液滴接觸角較大，隨著序參數(shù)比增大，液滴接觸角逐漸減小。采用GHASSEMI和PAK[18]的方法計算接觸角，如圖3所示。

圖3 接觸角的定義

接觸角計算公式[19]為

序參數(shù)比與接觸角之間的理論關(guān)系[19]為

通過測量模擬結(jié)果中的液滴與固體壁接觸的長度L和液滴的高度H，根據(jù)公式(13)計算的接觸角模擬值與根據(jù)公式(14)計算的接觸角理論值之間的關(guān)系如圖4所示。

圖4 不同序參數(shù)比條件下接觸角計算值和模擬值比較

通過比較，發(fā)現(xiàn)模擬值與理論值之間具有較好的一致性，可以證明此模型能夠準(zhǔn)確地模擬液滴與壁面之間的接觸角。

2.2 T型微通道內(nèi)氣液兩相流

T型微通道的模型如圖5所示，計算域的格子尺寸為100×600，微通道格子寬度為20。左側(cè)為液體入口，下面是氣體入口，兩相流在T型微通道內(nèi)混合生成氣液兩相流，通過右側(cè)出口流出微通道。兩相密度比為500∶1，與接觸角的模擬相同，壁面處施加反彈邊界，入口處施加速度邊界，出口施加充分發(fā)展邊界。

圖5 T型微通道的建模

2.2.1 接觸角對兩相流流型的影響

在入口速度和毛細數(shù)相同時，對5種不同接觸角(接觸角為57.5°、72.9°、90°、108.9°、129.5°)條件下氣泡的形成過程進行模擬。在模擬時間為2×105時，氣泡呈現(xiàn)出不同的形態(tài)，如圖6所示。

圖6 不同接觸角條件下的兩相流流型

從圖6(a)可以看出，當(dāng)接觸角為57.5°時，還未生成氣泡；而隨著接觸角增大，通道內(nèi)已有氣泡生成，并產(chǎn)生變形，如圖6(b)—(e)所示。這表明液體與壁面的接觸角是影響氣泡生成的重要因素，壁面接觸角越大越有利于氣泡的生成。這是由于壁面接觸角的增大，導(dǎo)致壁面附近的液體出現(xiàn)滑移速度，增大了液體的速度，進而截斷氣體形成氣泡。

通過圖6中液相的形態(tài)，可以觀察到，接觸角較小時，液相(黃色)尾部內(nèi)凹，頭部平坦，如圖6(b)所示；隨著接觸角不斷增大，液相尾部凹陷的程度逐漸減小，如圖6(d)—(e)所示，尾部變平，甚至有向外凸的趨勢。液相的頭部逐漸由略微凹陷到向外凸起，這與模擬單液滴的接觸角度變化相一致。因此，由單靜止液滴與壁面的接觸角可以推測兩相流動的氣泡/液滴的形態(tài)。

2.2.2 毛細數(shù)對兩相流流型的影響

毛細數(shù)也是影響兩相流流型的重要因素之一。為了快速生成液滴，接觸角選擇129.5°，氣液兩相入口量綱一速度均為0.002。將毛細數(shù)分別設(shè)置為0.000 25、0.000 19、0.000 17、0.000 16進行模擬，在模擬時間為2×105時，流型圖如圖7所示。

圖7 不同毛細數(shù)條件下的兩相流流型

從圖7(a)中可以看出，毛細數(shù)越大，氣泡形成的分層流段越長，還未形成氣泡，需要繼續(xù)向前運動才會生成氣泡。由圖7(b)—(d)可以觀察到，在相同的時間步數(shù)條件下，隨著毛細數(shù)的減小，氣泡脫離之前的聯(lián)結(jié)部位(圖中標(biāo)識處)逐漸變薄。因此，毛細數(shù)的減小有利于氣泡的生成，該結(jié)果與SHI等[20]的結(jié)論相一致。該現(xiàn)象是由于毛細數(shù)的減小使表面張力的影響增大，加快了氣泡的生成。

2.2.3 氣液流速比對兩相流流型和壓力降的影響

實際微通道的壁面接觸角一般小于90°，故下面分析模擬接觸角選擇為57.5°。為使氣泡快速生成，毛細數(shù)為0.000 16，界面層厚度為9。文中選擇了6種氣液流速比進行模擬，分別為1∶1(vG∶vL=0.001∶0.001和0.002∶0.002)、1∶2(vG∶vL=0.001∶0.002)、1∶3(vG∶vL=0.001∶0.003)、2∶1(vG∶vL=0.002∶0.001)和3∶1(vG∶vL=0.003∶0.001)。模擬時間為1×105時的流型如圖8所示。

圖8 不同氣液流速條件下兩相流流型

從圖8(a)—(c)中可知，在液相速度一定時，氣相速度越大，生成氣泡的位置距離T型微通道交叉處越遠，甚至不會再形成氣泡，即當(dāng)氣液流速比達到3∶1時，生成的液體部分不能完全充滿整個微通道，其寬度約為通道寬度的1/2，而氣相一直是連續(xù)的。從圖8(c)—(e)可知，當(dāng)氣相速度一定時，隨著液相速度的增加，生成的液滴越來越長，相應(yīng)的氣泡長度略微減小，由1∶1時的25格子數(shù)減小到1∶3時的16格子數(shù)，減小了36%。從圖8(c)、(f)可知，即使在相同的氣液流速比的情況下，由于速度的不同，流型也存在差別，如氣液流速比1∶1時，氣液流速的增加，使氣泡的長度有所減小，由氣液入口速度0.001時的25個格子數(shù)減小到入口速度0.002時的18個格子數(shù)，減小了28%，但相比不同氣液流速比條件下氣泡長度變化較小。同時還發(fā)現(xiàn)，當(dāng)氣相流速較大，生成液滴/氣泡的位置遠離T型微通道交叉處，分層流的長度增加。

對不同氣液流速條件下的兩相流進行模擬，模擬時間為1×105，探究不同氣液流速條件的壓力及流速變化情況。微通道中心線上壓力的變化如圖9(a)—(c)所示。圖9(d)所示為氣液流速比為1∶3情況下壓力、流型和流速的對照。

圖9 不同氣液流速條件下通道中心位置壓力變化

通過圖9可以發(fā)現(xiàn)，由單相區(qū)域到兩相區(qū)域時，壓力首先降低，然后經(jīng)歷兩相交匯區(qū)域，壓力存在波動，再到兩相區(qū)域，壓力繼續(xù)降低。與兩相交匯區(qū)域不同，兩相區(qū)域是與流型有關(guān)的有規(guī)律的壓力波動，氣泡前后端的壓力幾乎相同，并高于液相壓力，而液相壓力下降，液相的壓力降在兩相壓力降中占主導(dǎo)地位。液相壓力較之前單相區(qū)域下降得更為明顯。從兩相和單相流的壓力變化斜率可以發(fā)現(xiàn)，兩相流的壓力變化要大于單相流的壓力變化。

微通道方向中心線速度變化曲線如圖10所示。沿微通道方向中心線速度的峰值存在“滯后”現(xiàn)象，分別出現(xiàn)在入口后的一段距離的位置(距離入口7～9個格子數(shù)范圍內(nèi))和兩相交匯之后的一段距離的位置(距離兩相交匯位置3～12個格子數(shù)范圍內(nèi))，并不是入口位置速度或者兩相交匯之處速度出現(xiàn)峰值。之后，速度逐漸下降。速度波動幅度隨氣液流速比增大而增大。從速度峰值出現(xiàn)的范圍也可以發(fā)現(xiàn)，氣相的加入使得峰值出現(xiàn)的范圍增大，也證明了氣相的加入使兩相區(qū)域的速度波動增加。所以氣相速度的大小，對流動的穩(wěn)定性具有較大的影響。

圖10 微通道中心線的速度分布

Fig.10 Velocity distribution along the center of microchannel

vG：vL為0.002∶0.002時，微通道各個軸向位置橫截面上的速度分布如圖11所示。

圖11 微通道各個軸向位置橫截面速度分 布(vG∶vL=0.002∶0.002)

在單相流區(qū)域(Xaxis=11)流體速度為一條最大速度約為1.5倍入口速度的拋物線。在氣相入口處的兩側(cè)(Xaxis=41和Xaxis=61)觀察到靠近氣體入口一側(cè)出現(xiàn)較大的速度，并在下游(Xaxis=61)的整體速度最大。在兩相流區(qū)域，速度逐漸減小，速度受氣相的影響出現(xiàn)非拋物線分布的情況。

3 結(jié)論

(1)隨著毛細數(shù)的增大，分散相脫離點遠離兩相入口，形成更長的分層流。增大接觸角，減小毛細數(shù)有利于氣泡的生成。當(dāng)氣相流速較大，生成液滴/氣泡的位置遠離T型微通道交叉處，分層流的長度增加。

(2)沿微通道方向壓力逐漸減小，兩相區(qū)域內(nèi)液相的壓力降占主導(dǎo)地位。兩相流的壓力降大于單相流的壓力降，兩相流在氣液兩相交匯區(qū)域和兩相流區(qū)域壓力都存在波動。

(3)流速的峰值出現(xiàn)“滯后”現(xiàn)象，氣相的加入使滯后出現(xiàn)的范圍有所增加，氣液流速比越大，速度波動越大。在氣相入口處下方氣液兩相運動速度較大，氣液兩相整體運動速度在下游最大。