王 月，徐 博,2，王艷麗

(1.北華大學(xué)機械工程學(xué)院，吉林吉林 132021； 2.燕山大學(xué)機械工程學(xué)院，秦皇島河北 066000；3. 吉林省工業(yè)技師學(xué)院，吉林吉林 132021)

引言

微細氣泡的等效直徑通常小于100 μm相對常規(guī)氣泡體積甚小，具有許多優(yōu)良特性，如氣泡比表面積較大、在兩相間傳質(zhì)效率高、上浮速率慢、表面自帶活性自由基等[1-2]。微細氣泡因具有如上所述特性而逐漸被廣泛應(yīng)用于臭氧投加處理污水、醫(yī)學(xué)檢測、美容保健、去果蔬農(nóng)殘等諸多領(lǐng)域，且在各領(lǐng)域中微細氣泡都表現(xiàn)出很好的潛在應(yīng)用及科研價值[3-5]。

近年來，國內(nèi)外有許多學(xué)者進行了大量的研究，并得到大量的研究成果。譚立新[6]研究了水利工程中水-氣二相流的一些基本理論，并提出了水-氣二相流的數(shù)值模擬方法與實測數(shù)值驗證的方式。趙鐸[7]通過在大型多相流試驗環(huán)道上進行氣-水二相流試驗及不同流量的數(shù)值模擬得出氣-液二相的分布及流速。柏超[8]通過理論分析方法解釋微細泡(液滴)生成過程的一些現(xiàn)象與流型轉(zhuǎn)換機理，并得出預(yù)測微細氣泡(液滴)大小的數(shù)量級關(guān)系式。袁希鋼等[9]針對T形微通道，運用數(shù)值模擬與試驗的方式研究了黏度和表面張力對氣-液二相流型分布的影響。楊麗等[10]針對流動聚焦模型，通過數(shù)值模擬研究了微管道內(nèi)液滴的形成機理及影響因素。馮俊杰[11]研究了2個共軸線上升氣泡之間的連接和聚并行為，以及氣泡在湍流場中的上升和破碎行為。

國外研究主要以剪切式研究為主，較少涉及同軸式的研制。CHEN W B等[12]運用界面描述法研究微細氣泡(液滴)形成時的表面形狀，仿真研究了在浸沒(同軸)式噴嘴處的兩相氣泡的形成和相變。整個氣泡由蒸氣和液相組成，假定內(nèi)部氣相被薄的冷凝層包圍，當兩相氣泡增長時蒸氣部分冷凝。界面描述法用于描述氣泡形成的動力學(xué)，氣泡界面處的氣液傳熱效率與不同流量下的氣液壓差有關(guān)。QIAN D等[14]在針對T形管的研究中得出，氣體段塞長度隨著表面氣體速度的提高而增加，并且表面液體速度降低。LI X等[15]通過數(shù)值模擬與試驗驗證研究T形管道中微液滴形成機制。WONG V L等[16]運用數(shù)值模擬方法研究T形管道微流體黏度對氣泡體積的影響，得出選擇聚合物濃度和分子量來控制剪切變形和彈性，可以調(diào)整液滴尺寸和形成機制。

本研究以同軸式微通道中的微細氣泡生成過程為研究對象，基于COMSOL開展微通道中的氣、液二相流相互作用過程數(shù)值模擬，獲得氣泡生成過程中不同時刻的輪廓圖像，分析不同氣體壓強對氣泡生成頻率、生成體積等生成特性的影響趨勢，并基于微細氣泡生成特性測試平臺，驗證微細氣泡生成過程數(shù)值模擬結(jié)果的準確性。

1 微細氣泡生成過程仿真建模

1.1 COMSOL簡介

COMSOL對模型定義是非常靈活的，可以用常數(shù)、通用函數(shù)、表達式以及插值函數(shù)等來表達材料屬性和邊界條件等。COMSOL具有多物理場應(yīng)用模式，可以解決氣液二相間多物理場耦合問題，用戶可以自主選擇需要的物理場進行數(shù)值計算。針對流體仿真，該軟件為用戶提供了許多優(yōu)化完善的控制方程，供用戶計算使用。

1.2 控制方程

由于微細氣泡流的雷諾數(shù)很小，因此將連續(xù)相設(shè)置為層流流動，在COMSOL Multiphysics中可用二相流水平集方法研究兩種流體之間的變化。數(shù)值模擬中涉及到的控制方程包括3個，即連續(xù)性方程、Navier-Stokes方程、水平集方程。

=▽·[-pI+μ(▽u+(▽u)T)+σkδn]

(1)

▽·u=0

(2)

(3)

式中，ρ—— 流體密度

u—— 流體速度向量

t—— 時間

p—— 總壓力

I—— 單位矩陣

μ—— 流體動力黏度

σ—— 表面張力系數(shù)

k—— 流體界面的曲率

δ—— 二相流界面函數(shù)

n—— 指向液滴的界面單位法線向量

σkδn—— 作用在二相流界面上的表面張力

φ—— 水平集函數(shù)

γ—— 重新初始化參數(shù)

ε—— 界面厚度控制參數(shù)

式(1)中的流體密度ρ和動力黏度μ可表示為:

ρ=ρl+(ρg-ρl)φ

(4)

μ=μl+(μg-μl)φ

(5)

式中，ρl—— 液體密度

ρg—— 氣體密度

μl—— 液體動力黏度

μg—— 氣體動力黏度

針對微通道內(nèi)二相流體的流動和相互作用，利用COMSOL的水平集模塊求解方程。在利用COMSOL進行數(shù)值模擬過程中，要求用戶針對潤濕壁邊界條件設(shè)定滑移長度，根據(jù)摩擦力定義，在此暫不考慮液體在潤濕壁面的相對滑動，遂將滑移長度初值取0。

1.3 幾何模型建立

在實際的試驗過程中，氣體入口是直徑60 μm的毛細金屬管，液體通道是邊長為10 mm的方形玻璃管。為減輕仿真計算量和仿真周期，先建立二維仿真模型，再在仿真中分別運用軸對稱旋轉(zhuǎn)和陣列得到三維立體仿真計算結(jié)果，其模型的每部分尺寸如圖1所示。

圖1 幾何模型尺寸

1.4 邊界條件設(shè)置

在同軸式微通道中，忽略能量交換，設(shè)定初始流場中A區(qū)域為空氣、B區(qū)域為液體。由于微細通道中流體雷諾數(shù)很小，Re<<1，因此流體流態(tài)為層流，針對同軸式微通道中氣泡生成的仿真模型，設(shè)置邊界條件如下：

1) 入口邊界條件

氣液二相流體在入口邊界的速度均勻，無梯度分布，氣體和液體均以各自恒定的流速進入微通道。

2) 出口邊界條件

氣液二相混合流體從固定出口流出時的壓力恒定，為1個大氣壓，pout=0 Pa或pout=1 Pa。

3) 對稱邊界條件

對于整體同軸微通道仿真模型而言，其中進行的流體流動數(shù)值模擬計算量較大，因此為減輕計算量，取整體模型的1/4部分作為仿真模型。另外，為方便控制方程識別，需要將仿真模型設(shè)置為軸對稱。

4) 潤濕壁邊界條件

將同軸微通道壁面添加潤濕性屬性，壁面潤濕程度可通過設(shè)置接觸角θ值確定。潤濕壁的滑移長度可設(shè)置為1 μm，如圖2所示。

圖2 同軸聚焦結(jié)構(gòu)邊界條件設(shè)置

1.5 初始條件設(shè)置

將毛細管道設(shè)置為氣體，外部管道設(shè)置為液體，圖2的初始界面為氣液分界面。根據(jù)實際試驗中獲得的液體流速和氣體流速設(shè)置仿真中的初始條件，液流流進管道的速度u1=0.049 m/s，氣流流入毛細管道的速度u2=2.97 m/s。設(shè)置接觸角為90°，重新初始化參數(shù)設(shè)置為0.5，界面厚度控制參數(shù)為ls.ep_default(表示細化網(wǎng)格后最小網(wǎng)格尺寸的1/2)，初始條件參數(shù)如表1所示。

表1 同軸聚焦模型初始條件

表1中取用的值為室溫(20 ℃)下的數(shù)值，接觸角選擇為氣泡最初形成時氣液二相分界面的切線與同道壁面的夾角。此外，仿真過程中的氣體流量由氣泡生成體積除以氣泡的生成時間求得。

1.6 網(wǎng)格單元劃分

由于連續(xù)相與分散相交界處網(wǎng)格單元尺寸較小，不能以較大尺寸的網(wǎng)格單元進行劃分，而全局網(wǎng)格若均設(shè)置成很小的網(wǎng)格單元，會增加仿真的計算量，致使仿真效率降低。綜合考慮以上條件，在此全局上采用物理場控制網(wǎng)格、單元格大小選擇及細化。局部分散相與連續(xù)相的交界處選擇自動細化網(wǎng)格，減少不必要的計算。在仿真之前，預(yù)先劃分一次網(wǎng)格，共59695個三角形網(wǎng)格單元，30311個網(wǎng)格頂點。在仿真過程中，系統(tǒng)自動細化網(wǎng)格，最終有116676個三角形網(wǎng)格單元和58944個網(wǎng)格頂點。圖3表示仿真前后模型的網(wǎng)格局部示意圖。

圖3 仿真前后網(wǎng)格對比圖

1.7 求解器設(shè)置

在同軸聚焦模型的求解中主要運用瞬態(tài)求解器，方法選擇向后迭代求解器，一致初始化部分選用向后歐拉法。在誤差估計結(jié)束后，完成數(shù)值模擬部分相關(guān)設(shè)置。

2 微細氣泡生成特性試驗

2.1 試驗設(shè)計及試驗平臺搭建

為探究同軸式微通道中微細氣泡的生成效果和生成特性，搭建了一套微細氣泡制備試驗系統(tǒng)。試驗中氣泡生成的載體為同軸機械式微流控芯片，普通氮氣瓶為芯片提供氣體，并利用YQD-4雙級氮氣減壓閥對氣體壓強進行調(diào)節(jié)；Harvard 704500微量注射泵作為芯片的液體供給單元，按需調(diào)整芯片液體入口流量。將Olympus CX21顯微鏡和Phantom v12.1高速攝像機組合使用，同時具備高速拍攝和顯微觀測的功能，用于觀察和采集氣泡生成過程中的圖像，試驗系統(tǒng)的原理示意圖如圖4所示。

圖4 微細氣泡生成特性測試平臺示意圖

由于數(shù)值模擬結(jié)果需要與試驗形成對比關(guān)系，因此試驗中需要收集到規(guī)律性較明顯的幾組參數(shù)，防止由于數(shù)值模擬算法的問題出現(xiàn)不收斂或無規(guī)律現(xiàn)象。試驗采取固定液體流量，改變氣體壓強的試驗條件，探究氣體壓強對微細氣泡脫離體積和生成頻率的影響規(guī)律。

2.2 試驗過程及方法

為探究不同液體流量和氣體壓強下的氣泡生成特性，利用所搭建的微細氣泡制備試驗系統(tǒng)開展氣泡生成特性研究試驗，所設(shè)定氣壓與液體流量記錄在表2中。

試驗過程中，液體流量由Harvard 704500微量注射泵提供，該注射泵可精確調(diào)節(jié)液體流量。氣體壓強可通過連接氮氣瓶的減壓閥和精密氣壓調(diào)節(jié)器精準調(diào)節(jié)。設(shè)置高速攝像機拍攝幀率為10000 fps，記錄每個微細氣泡的生成演變時間，試驗中的相關(guān)參數(shù)及數(shù)值，如表3所示。

表2 試驗變量參數(shù)及其數(shù)值

表3 試驗中的相關(guān)參數(shù)及其數(shù)值

2.3 試驗數(shù)據(jù)的選擇

為了可以和數(shù)值模擬結(jié)果進行清晰地對比，需要選擇規(guī)律性較為明確的試驗數(shù)據(jù)作為仿真的輸入控制參數(shù)。表4記錄了液體流量為150 mL/h時，不同氣體壓強pg下氣泡的脫離體積Vg和生成時間t。可以明顯看出，氣體壓強變化對氣泡的生成時間及體積影響較大，可以作為數(shù)值模擬輸入量仔細研究氣泡形成過程的變化。

表4 不同氣體壓強下氣泡生成時間及體積

3 仿真與試驗結(jié)果對比分析

3.1 微細氣泡形成過程

細致觀察氣泡的形成過程是研究微細氣泡生成的良好方法，而使用數(shù)值模擬方法可以更便捷的觀察微細氣泡的形成過程。圖5a為相同的控制量下通過數(shù)值模擬方法計算出的氣泡形成過程，圖5b為試驗中相同控制量下氣泡的形成過程，其中為了與試驗圖像形成對比選取了相同氣泡生長階段的圖形比較。

圖5 不同時刻下微細氣泡生成過程

將這些圖像進行對比分析，可以發(fā)現(xiàn)整個氣泡生成過程同樣可分為3個階段：膨脹階段、軸向拉伸階段和頸縮脫離階段。經(jīng)過對比容易發(fā)現(xiàn)，相同的氣泡生成階段，數(shù)值模擬中的形成時間較短，并且頸縮脫離階段的頸縮量相對較長。

3.2 氣泡生成頻率

由于數(shù)值模擬方法需要在實際工程中運用，因此必須找到數(shù)值模擬方法與實際氣泡生成頻率之間的區(qū)別。研究中選取了液體流量為140 mL/h，通過改變氣體壓強探究不同工況下氣泡生成頻率的變化規(guī)律。表5記錄了不同氣體壓強下氣泡的生成頻率，觀察表中氣泡生成頻率的試驗值fgt和仿真值fgs可以發(fā)現(xiàn)，隨著氣體壓強的增大，氣泡的生成頻率變快。

表5 不同氣體壓強下氣泡的生成頻率

為了更清晰的觀察實際試驗與數(shù)值模擬計算的氣泡生成頻率之間的區(qū)別，根據(jù)表5繪制實際氣泡生成頻率與數(shù)值模擬計算的氣泡生成頻率隨氣體壓強增加時的變化趨勢，如圖6所示?？梢悦黠@看出，經(jīng)過數(shù)值模擬計算出的氣泡生成頻率相對實際試驗氣泡的生成頻率較大，因此實際工程運用數(shù)值模擬方法分析氣泡生成頻率時應(yīng)考慮相應(yīng)的計算差值。

圖6 不同氣體壓強下的氣泡生成頻率

3.3 氣泡的體積

將液體流量設(shè)置為恒定，通過調(diào)節(jié)氣體壓強獲得固定流量不同氣壓下微細氣泡的脫離體積，并將不同工況下氣泡的體積記錄在表6中。

表6 恒定流量不同氣壓下氣泡的體積

表6中的數(shù)據(jù)記錄了不同氣體壓強下氣泡脫離體積的試驗值Vgt和仿真值Vgs，為了鮮明對比二者的變化趨勢，將表中數(shù)據(jù)進一步處理，得到不同氣體壓強下氣泡體積的變化曲線，如圖7所示。

圖7 不同氣體壓強下氣泡生成體積

從圖7中可以明顯看出，恒定液體流量不同氣體壓強下氣泡體積的試驗值和仿真值均隨氣體壓強的增大而減小，而且仿真結(jié)果普遍略小于實際試驗值。

4 結(jié)論

通過本研究可以得出：采用COMSOL Multiphysics軟件數(shù)值模擬機械式同軸型微管道內(nèi)氣泡的形成，與實際試驗規(guī)律明顯一致，但數(shù)值上具有一定差值。數(shù)值模擬在計算中的優(yōu)勢明確，可以顯著的觀察氣泡形成中的具體變化，經(jīng)過研究分析得到數(shù)值模擬方法在微細氣泡生成的研究分析上有著較大的運用前景。