豎直小通道內彈狀流氣彈長度的計算模型

閆超星，閻昌琪，孫立成，王 洋，張小寧

(哈爾濱工程大學 核安全與仿真技術國防重點學科實驗室，黑龍江 哈爾濱 150001)

鑒于緊湊式換熱器的發展需求，小通道內部的兩相流動特性得到了廣泛的關注。由于彈狀流內部氣彈和液彈互相尾隨交替出現，造成了通道內很大的密度差和流體的可壓縮性，易產生流動不穩定性，因此，針對彈狀流的實驗和理論研究備受國內外學者的關注[1-3]。

對以往的研究總結發現，學者們大都關注氣彈上升速度，且提出了較為準確的計算模型，但針對氣彈長度的研究僅給出了實驗現象和其影響因素，關于氣彈長度的計算模型還十分有限，且不完善。氣彈長度關系著氣液相的分布情況，進而影響傳熱和阻力特性，此外氣彈長度與兩相流型間的轉變密切相關。因此，有必要建立小通道內氣彈長度的計算模型，為今后的研究和工程應用提供技術支持。

1 實驗系統

1.1 實驗裝置及實驗內容

實驗以空氣和水為工質，在常溫、常壓條件下進行，實驗過程中空氣和水均保持在20～22 ℃。圖1為實驗回路示意圖。整個實驗裝置由供水系統、供氣系統、實驗段、數據采集系統和高速攝像系統5部分組成。實驗段為有機玻璃矩形通道，截面尺寸為3.25 mm×43 mm，總長2 m。實驗數據測量和采集系統詳見文獻[4]。高速攝像系統由高速攝像儀、光源和計算機組成。實驗中采用Photron公司的FASTCAM SA5型高速攝像儀垂直于實驗段寬邊進行拍攝，以漫射的背光為光源。圖像拍攝處位置距通道入口986 mm，流動已充分發展。

氣液流量均通過調節閥來調節，兩相混合物流經實驗段后，進入氣水分離器，在其內部依靠重力自然分離。實驗過程中，固定水流量，氣流量由小到大依次調節，待每個工況穩定后，記錄實驗數據，完成1個循環；然后再改變水流量，進行下組實驗。

圖1 實驗回路示意圖

1.2 圖像處理

圖像處理前首先確定標度因子γ。實驗段寬邊尺寸已知，可通過計算其圖像像素點的方式確定γ：

(1)

圖像處理過程示于圖2。圖2a、b和c分別代表n1幀、n2幀和n3幀對應的圖像。由圖2a、b得到氣彈上升速度，再由圖2a、c確定的時間間隔得到氣彈長度Lb：

Lb=γΔx(n3-n1)/(n2-n1)

(2)

圖2 實驗圖像處理過程

式中，Δx為n2幀與n1幀圖像氣彈頭部界面沿軸向的位移。需指出，窄矩形通道內彈狀流的氣彈長度具有很大隨機性，即使同一流動工況下，流過的氣彈長度也會有差異，這體現了兩相彈狀流的流動不穩定性。為減小主觀測量誤差，每個氣彈重復測量10次，對每個工況下15個氣彈進行處理后取平均值。

2 氣彈長度模型

豎直通道內彈狀流物理模型示于圖3。圖3中，h為距氣彈頭部的距離，α(h)和αb分別為距氣彈頭部h處和氣彈尾部的局部空泡份額，ub為氣彈上升速度，j為兩相折算速度，ufb為液膜末端的平均流速，wb和w分別為氣彈寬度(或直徑)和通道寬度(或直徑)。

圖3 豎直通道內彈狀流物理模型

Mishima等[5]針對圓形通道內彈狀流應用伯努利方程，得到距氣彈頭部h處的α(h)：

(3)

Xu等[6]和Hibiki等[7]將式(3)應用于矩形通道，其中分布參數C0和漂移速度Vgj分別采用下式計算：

(4)

(5)

式中：ρg、ρf和Δρ分別為氣、液相密度和兩相密度差；s和w分別為通道窄邊和寬邊尺寸。

Hibiki等[7]認為距氣彈頭部某一位置處，氣彈周圍液膜所受的重力與其受到的壁面剪切力相等，此時液膜的流動已充分發展。在此位置以下，向下流動的液膜與向上流動的液彈混合，液膜減速直至不再向下流動，液相攪混劇烈導致此處液膜不穩定，氣液界面的擾動將會打斷氣彈。矩形通道內，貼近寬邊處的液膜很薄，與窄邊處的液膜厚度相比可忽略，但寬邊處的液膜同樣受到壁面剪切力的作用；圓形通道內，氣彈周向液膜均勻分布。矩形通道和圓形通道內液膜的簡化模型示于圖4。假設同一橫截面上氣彈周圍液膜速度相等，本文以矩形通道為例進行推導。

圖4 彈狀流液膜簡化模型

對液膜受力分析可知[6]：

Δρgws(1-αb)dz

(6)

其中：dz為選取液膜微元段沿軸向的長度；f為壁面摩擦系數，其計算式如下：

(7)

其中：De為通道的當量直徑；υf為液相運動黏度；層流區n=1，湍流區n=0.25。Mishima等[5]提出，由于液膜內部存在小氣泡，且該處液膜緊靠氣彈尾部，液相攪混劇烈，因此該處液膜的流動用湍流模型更接近實際情況，本文暫取n=0.25。

Shah等[8]提出了計算矩形通道內單相層流區摩阻系數的關系式：

λ=CflRe-1=96(1-1.355 3ε+1.946 7ε2-

1.701 2ε3+0.956 4ε4-0.253 7ε5)

(8)

其中：λ為單相摩阻系數；ε為矩形通道橫截面的窄邊與寬邊之比，即s/w。

Sadatomi等[9]提出了考慮通道幾何形狀的湍流摩阻系數關系式：

Cft=Cft0[(0.015 4Cfl/Cfl0-0.012)1/3+0.85]

(9)

式中，Cfl0=16和Cft0=0.079分別為圓管內層流區和湍流區的摩阻系數。通過式(8)及(9)計算矩形通道內層流區和湍流區單相摩阻系數。將式(7)代入式(6)得到氣彈尾部液膜的平均速度：

(10)

此外，ufb還與αb存在以下關系[5]：

ufb=(αbub-j)/(1-αb)

(11)

又由漂移流模型可知：

ub=C0j+Vgj

(12)

結合式(10)～(12)可得：

(13)

觀察式(13)可發現，該式為計算αb的1個迭代關系式，并不能直接得到計算αb的表達式，因此，需對等式右邊1-αb的指數形式進行簡化處理，即：

≈β(1-αb)

(14)

圖5 簡化關系式驗證

當n=0.25、β分別取0.25和0.45時，式(14)兩側計算結果的逼近情況如圖5所示。實驗獲得的氣彈寬度與通道寬度的比值大多位于0.6～0.75范圍內，窄矩形通道內可近似認為αb等同于氣彈寬度與通道寬度的比值，因此，只要保證式(14)左右兩邊關系式在αb處于0.6～0.75范圍內具有較好的一致性即可。通過Mishima等[10]的實驗觀察和Wilmarth等[11]基于電導探針對液膜厚度的實驗研究，發現寬邊處的液膜與窄邊處相比可忽略，Hibiki等[7]得出當窄邊寬度小于2.4 mm時寬邊處液膜可忽略。隨矩形通道窄邊寬度的減小，氣彈受到壁面擠壓作用越明顯，氣彈寬度占通道寬度的比例增大。以窄邊尺寸2.4 mm為界，當s≥2.4 mm時，β=0.45；當s<2.4 mm時，β=0.25。

將式(14)代入式(13)可得：

αb=(j+β((De/υf)-nCfρf(w+s)/

(15)

根據式(3)可得，距氣彈頭部Lb處的局部空泡份額αb為：

(16)

結合式(15)和(16)可得：

(17)

從而得到氣彈長度的關系式為：

(18)

其中：n=0.25；s≥2.4 mm時，β=0.45；s<2.4 mm時，β=0.25。

圓形通道與矩形通道的機理相同，本文簡化推導過程，得到圓形通道計算關系式：

(19)

式(18)和(19)適用于豎直條件下常溫絕熱空氣-水兩相彈狀流動。

3 結果分析與討論

3.1 氣彈長度的影響因素

觀察式(18)和(19)發現，對于給定的流體物性參數，氣彈長度僅與通道尺寸和流動狀況(j、Cf)相關，實驗條件下Cf變化十分有限。因此，本文僅討論通道尺寸和兩相折算速度對氣彈長度的影響。以矩形通道為例，圖6示出不同高寬比ε時氣彈長度隨兩相折算速度j的變化趨勢。由圖6可知，氣彈長度隨高寬比及兩相折算速度的增加而增大。

3.2 計算模型的驗證

為驗證計算模型的實用性，將本文及Mishima等[10]和Cheng等[2]的實驗數據與計算值進行對比，結果示于圖7。所采用的實驗數據列于表1，實驗工質均為水和空氣，包括矩形和圓形通道。平均絕對誤差(MAE)定義為：

∑

(20)

其中：m為實驗數據點數；(Lb)pred和(Lb)exp分別為氣彈長度預測值和實驗值。

對于當前實驗數據，計算模型的MAE為26.8%，而對Mishima等和Cheng等全部實驗數據的MAE為34.9%，說明計算模型具有較好的適用性。模型對Mishima等的40 mm×1.07 mm通道的實驗數據預測效果稍差，可能是模型中認為氣彈周向液膜的流速相等引起，因窄邊僅為1.07 mm的矩形通道內寬邊液膜流動受到抑制，加劇了寬邊與窄邊液膜速度的不均勻程度，與模型假設存在一定偏差。模型的定性尺寸為表1列出的通道尺寸，需指出，小通道氣彈長度的模型十分有限，對其他通道尺寸的適用性有待進一步驗證，本文只是提供分析氣彈長度的一種方法。

圖6 高寬比ε對氣彈長度的影響

a——本文實驗數據；b——公開發表文獻的實驗數據

表1 氣彈長度實驗數據

4 結論

通過對彈狀流液膜的受力分析，結合漂移流模型建立了小通道內氣彈長度的計算模型，并最終給出了氣彈長度的計算關系式。通過可視化實驗數據及公開發表文獻中數據對計算方法進行了評價。計算模型對本實驗段的平均絕對誤差為26.8%，對Mishima等和Cheng等實驗數據的平均絕對誤差為34.9%；對于尺度為1 mm左右的小通道，模型的預測精度低于2 mm以上的通道。

參考文獻：

[1] 夏國棟，周芳德，胡明勝. 傾斜上升彈狀流中Taylor氣泡運動速度研究[J]. 化學工程，1997，25(5)：36-41.

XIA Guodong, ZHOU Fangde, HU Mingsheng. Study on the velocity of Taylor bubble in inclined upward slug flow[J]. Chemical Engineering, 1997, 25(5): 36-41(in Chinese).

[2] CHENG T W, LIN T L. Characteristics of gas-liquid two-phase flow in small diameter inclined tubes[J]. Chemical Engineering Science, 2001, 56(21-22): 6 393-6 398.

[3] BHUSAN S, GHOSH S, DAS G, et al. Rise of Taylor bubbles through narrow rectangular channels[J]. Chem Eng J, 2009, 155(1-2): 326-332.

[4] 王廣飛，閻昌琪，孫立成，等. 窄矩形通道內兩相流動壓降特性研究[J]. 原子能科學技術，2011,45(12)：1 331-1 336.

WANG Guangfei, YAN Changqi, SUN Licheng, et al. Resistance characteristics of two-phase flow through narrow rectangular duct[J]. Atomic Energy Science and Technology, 2011, 45(12): 1 331-1 336(in Chinese).

[5] MISHIMA K, ISHII M. Flow regime transition criteria for upward two-phase flow in vertical tubes[J]. Int J Heat Mass Transfer, 1984, 27(5): 727-737.

[6] XU J L, CHENG P, ZHAO T S. Gas-liquid two-phase flow regimes in rectangular channels with mini/micro gaps[J]. Int J Multiphase Flow, 1999, 25(3): 411-432.

[7] HIBIKI T, MISHIMA K. Flow regime transition criteria for upward two-phase flow in vertical narrow rectangular channels[J]. Nuclear Engineering and Design, 2001, 203(2-3): 117-131.

[8] SHAH R K, LONDON A L. Laminar flow forced convection in ducts: A source book for compact heat exchanger analytical data[M]. New York: Academic Press, 1978.

[9] SADATOMI M, SATO Y, SARUWATARI S. Two-phase flow in vertical noncircular channels[J]. Int J Multiphase Flow, 1982, 8: 641-655.

[10] MISHIMA K, HIBIKI T, NISHIHARA H. Some characteristics of gas-liquid flow in narrow rectangular ducts[J]. Int J Multiphase Flow, 1993, 19(1): 115-124.

[11] WILMARTH T, ISHII M. Interfacial area concentration and void fraction of two-phase flow in narrow rectangular vertical channels[J]. J Fluids Eng, 1997, 119(4): 916-922.