窄矩形通道內攪混流和環狀流相界面參數的計算方法

程林海，谷海峰,*，劉安泰，虞 想，張君毅

(1.哈爾濱工程大學 核安全與仿真技術國防重點學科實驗室，黑龍江 哈爾濱 150001;2.中國核動力研究設計院 中核核反應堆熱工水力技術重點實驗室，四川 成都 610213)

兩流體模型是目前兩相系統流體傳熱學中最詳細、精確的宏觀模型[1]。準確獲取窄矩形通道內相界面參數是構建和驗證窄矩形通道中兩流體模型的必要條件[2-3]，對板狀燃料的核反應堆堆芯傳熱分析具有重要意義?？张莘蓊~與界面濃度是至關重要的相界面參數[4]，同時由于窄矩形通道橫截面具有較大的長寬比，導致通道寬邊方向上局部參數具有一定的分布規律[5-7]，局部參數反映局部位置的傳熱傳質能力，局部位置發生傳熱惡化往往會威脅整個設備的安全運行。因此獲取窄矩形通道內空泡份額和界面濃度，以及兩參數在其寬邊方向上的局部分布，對設備設計、安全評估、構建和驗證兩流體模型十分重要。

目前在窄矩形通道內相界面參數的實驗研究中，因受通道尺寸限制，測量方法以非接觸式的高速攝像法為主，且研究對象主要是空泡份額較低的泡狀流和彈狀流工況。張偲琪[7]對橫截面為40 mm×4 mm的窄矩形通道進行拍攝，在研究脈動條件下泡狀流的相界面參數過程中提出了相界面參數計算方法。Chalgeri等[8]利用高速攝像儀對橫截面為66.6 mm×2.35 mm的窄矩形通道進行拍攝，利用Fij軟件計算了泡狀流下的相界面參數。Shen等[9]針對40 mm×0.993 mm的窄矩形通道內泡狀流和帽狀流工況，提出了空泡份額和界面濃度計算方法。然而在窄矩形通道內高含氣率下的攪混流和環狀流流型中，氣液相界面較為復雜，且氣相空間缺少封閉且規則的幾何形狀，難以通過可視化實驗測量相界面參數。這導致該流型下所對應的相界面參數數據較少[10]，給兩流體模型的驗證帶來困難。

本文通過在橫截面為65 mm×2 mm的豎直窄矩形通道內進行空氣-水向上兩相流動實驗，基于高速攝像技術拍攝攪混流和環狀流條件下兩相流動圖像，提出一種獲取該流型下通道內空泡份額和界面濃度以及兩參數在其通道寬邊方向上局部分布的計算方法。

1 實驗裝置及工況

實驗回路示意圖如圖1所示，窄矩形通道試驗段材質為高透光率的聚甲基丙烯酸甲酯(PMMA)，中部流道截面為65 mm×2 mm，可視化有效長度為2 m。

通過調節氣、水回路的閥門，使氣、水流量穩定在目標工況，通道內出現相應的流型。調節高速攝像儀，使其拍攝到清晰且完整的兩相流視頻圖像。拍攝速率為6 000幀/s，分辨率為1 280×960 pixel。實驗中對試驗段上、中、下游分別進行拍攝，拍攝中央分別距混合腔出口600、1 200、1 800 mm。

圖1 實驗回路示意圖Fig.1 Sketch of experimental loop

為使氣相和液相在入口處充分混合，試驗段入口設有多孔混合腔室，如圖2所示?；旌锨皇覂染鶆蛟O有直徑為1 mm、間距為1 mm的孔洞。

圖2 兩相混合腔結構Fig.2 Two-phase mixing cavity structure

所開展兩相流動實驗中，氣相折算速度jg=1～9 m/s，液相折算速度jf=0.1～1.5 m/s，涵蓋攪混流、環狀流流型區間。具體實驗工況如圖3所示。

2 相界面參數獲取方法

2.1 相界面特征

實驗拍攝圖像如圖4所示。從圖4可看出，高空泡份額下氣相空間內具有以下特征：1) 氣相空間具有連續性，無論是攪混流(圖4a)還是環狀流(圖4c)，圖像中央位置均存在較長的氣相空間；2) 氣相空間邊界較復雜，尤其是攪混流中某些局部位置(圖4b)。

圖3 實驗工況流型圖Fig.3 Flow pattern map of experimental case

圖4 高速攝像儀拍攝的圖像Fig.4 Image of high speed camera

氣相空間的連續性導致有限的拍攝區域內存在不完整的氣泡和氣相空間，邊界復雜性導致氣相空間沒有規則的幾何形狀，難以計算其表面積和體積，因此需探索合適的計算方法。同時復雜的相邊界會對光產生反射與折射[11]，導致數據處理過程中無法準確獲取相界面輪廓，這也給相界面邊界識別帶來了難度。

2.2 圖像處理

根據已有圖像處理的研究[12-15]，針對攝像法拍攝的兩相流圖像進行氣相空間識別，識別過程如下：增強對比、閾值分割、孔洞填充和分水嶺分割。

攪混流圖像處理過程如圖5所示，通過增強對比和閾值分割，可初步識別相邊界(圖5c)。閾值分割后的圖像，由于邊界復雜性導致某些位置相邊界不連續，進而導致孔洞填充后存在不能識別的氣相空間(圖5d)。邊界不連續的位置可稱為“斷點”，本文采用手動補全的方法進行優化調整，直到完全識別所有氣相空間(圖5e)。

手動補全過程如下：以“斷點”前后相界面邊界像素點寬度為補全線條的寬度，對照拍攝原圖邊界走向對“斷點”進行補全，補全后與原圖進行反復對比，以盡量減小誤差。

采用分水嶺分割的方法將相互粘連的氣泡分隔開，由于相界面的復雜性，存在將不規則氣相空間錯誤分割的現象(圖5f)，本文通過手動調整錯誤分割的方法進行調整(圖5g)。由于環狀流工況圖像同樣存在“斷點”的現象，使用手動補全的方法對閾值分割后的圖像進行處理，圖像處理過程如圖6所示。

圖5 攪混流圖像處理過程Fig.5 Image processing of churn flow

圖6 環狀流圖像處理過程Fig.6 Image processing of annular flow

2.3 空泡份額與界面濃度

通過圖像處理可得到清晰、準確的相界面二值化圖像。從圖5、6可看出，拍攝區域內同時存在較小的氣泡和較大的氣相空間。通過程序識別方法可得到每個氣泡、氣相空間的投影面積和投影周長。根據投影面積可計算投影面積直徑ds(mm)：

(1)

式中，S為識別出的每個氣泡、氣相空間的投影面積，mm2。

在窄矩形通道內的兩相流動中，根據氣泡的大小可將氣泡分為兩類[16]：1) 氣泡投影面積直徑ds小于或等于通道寬度s，即ds≤s(本文中s為2 mm)時，認為是球形氣泡；2) 氣泡或氣相空間較大(ds>s)時，將會受到通道寬邊的擠壓，變為壓扁的氣泡或氣相空間。由于二者形態不同，因此具有不同的體積和表面積計算方法。

1) 球形氣泡的體積和表面積(ds≤s)

ds≤s時，球形氣泡的體積Vb(mm3)和表面積Ab(mm2)由下式計算：

(2)

2) 壓扁氣泡和氣相空間的體積和表面積(ds>s)

Shen等[16]在窄矩形通道內兩相流動研究中，假設壓扁氣泡邊緣為半圓形。同樣地，本文假設受擠壓的氣泡和氣相邊緣為半圓形，如圖7所示。

圖7 壓扁氣泡(a)和氣相空間(b)示意圖Fig.7 Sketch map of compressed bubble (a) and gas-phase space (b)

圖7中，長方形區域為高速攝像拍攝區域，通過程序識別出的投影面積S是圖7a中藍色實線圍成的面積。圖7b中的投影面積S是藍色實線與上下拍攝邊界圍成的面積。程序識別出的投影周長L是圖7a、b中為藍色實線的周長。圖中藍色實線是壓扁氣相最外圈的邊界，虛線是氣相空間與窄矩形通道寬邊內壁面直接接觸的邊界。從圖7a、b可看出，假設壓扁氣相邊緣為半圓形的區域為虛線與實線圍成的部分。

從圖5～7可看出，氣相空間較復雜，且有不同的長徑比。擁有相同投影面積而不同長徑比的氣相，其周長不同，直接影響氣相體積和表面積。長徑比對周長L的影響示于圖8。由圖8可計算得，在投影面積相同的情況下，長徑比為4的規則橢圓形壓扁氣泡的周長，較投影面積為正圓形(長徑比=1)壓扁氣泡的周長大46%。因此，在相界面復雜的工況中，若將拍攝圖像內所有受擠壓的氣相都等效為投影面積為圓形的氣泡進行計算，表面積和體積的計算結果將產生較大的誤差。

為準確求得上述氣泡、氣相空間的體積和表面積，提出將其拆分求解的思想。將任意受擠壓的氣泡和氣相拆分為中央部分和外圈部分，如圖9所示，中央部分為直接與窄通道壁面接觸的氣相部分，外圈部分為受到擠壓、橫截面為半圓形的氣相部分。

圖8 相同投影面積下不同長徑比氣泡的周長對比Fig.8 Perimeter comparison of bubbles with different aspect ratio in same area

圖9 壓扁氣泡/氣相空間分解圖Fig.9 Decomposition diagram of compressed bubble/gas-phase space

對拆分后的中央部分采用面積等效法(即等效前后投影面積相等)，等效為圓柱體進行參數計算。對拆分后的外圈部分采用周長等效法(即等效前后周長相等)，等效為圓環進行參數計算。對等效后的圓柱體和圓環分別計算體積和表面積，再將結果相加，即可得到拆分前氣泡或氣相空間的體積和表面積。

這樣既保留了原氣相中央部分的實際體積與表面積，又保留了因氣相復雜形態對氣泡體積和表面積產生的影響。

經過拆分等效后的圓柱體和圓環的體積和表面積計算方法如下。

1) 圓柱體的體積和表面積

等效后圓柱體的體積和表面積與拆分前原氣相中央部分的體積和表面積相等，計算公式如下：

(3)

(4)

式中：V柱體為圓柱體的體積，mm3；A柱體為圓柱體的前后表面積之和，mm2；r為等效后圓柱體的半徑，mm。本實驗中窄矩形通道寬度s為2 mm。

2) 圓環的體積和表面積

等效后圓環的周長與原氣相周長相等，因此求得的圓環體積和表面積即為等效前氣相外圈部分的體積和表面積。圓環的體積和表面積通過旋轉體積分公式求解：

V圓環=π2RLs2/4+πs2/6

(5)

A圓環=π2RLs+πs2

(6)

式中：V圓環為圓環的體積，mm3；A圓環為圓環的外側表面積，mm2；RL為圓環的外圈半徑，mm。

(7)

式中，L為氣泡投影周長，mm。

因此，結合式(3)～(6)，則壓扁氣泡的總體積與表面積為：

(8)

通過圖像識別程序，可得到每幅圖像內每個氣泡和氣相的周長和投影面積。以投影面積直徑2 mm為分界線，將氣泡進行分類，分別代入式(2)和(8)，可求得每個氣泡的體積和表面積。

對某張圖像內所有氣泡體積求和，然后除以通道的體積即為此張圖像的空泡份額，對所有氣泡表面積求和后除以通道體積即為此張圖像的界面濃度。選取某工況下某位置處某段拍攝時間內的多張照片，對每張圖片的空泡份額和界面濃度求平均值，即為此位置在此工況下的時均空泡份額和界面濃度。

(9)

(10)

式中：n為每張圖像內氣泡的數目；N為某段時間內圖像的數目；Vb,i為第i個氣泡的體積，mm3；Ab,i為第i個氣泡的表面積，mm2；Vchannel為拍攝區域內通道的體積，mm3。

2.4 局部參數模型

由于窄矩形通道橫截面具有較大的長寬比，因此在通道寬邊方向上，空泡份額和界面濃度具有一定分布特性。為研究其分布特性，需將拍攝圖像沿通道寬邊方向劃分為多個區域，計算每個區域內的局部參數，如圖10所示，圖中a、b區域分別為攪混流和環狀流工況下的一處劃分區域。

圖10 劃分區域示意圖Fig.10 Schematic diagram of area division

通過圖像識別程序，可獲得局部區域內每個氣相空間的面積Sloc和周長Lloc。圖10中黑色區域為液相，白色區域為氣相。局部區域內每個白色區域面積為Sloc，每個白色區域與黑色區域接觸界面的周長為Lloc。

從圖10可看出，劃分區域內的氣泡往往不是完整的，當投影面積直徑小于等于窄矩形通道寬度(ds≤s)時，假設氣泡為球形，此時氣泡體積和面積按照球形公式(式(2))計算。此假設可能將大氣泡的一角當作是球形氣泡進行計算，會導致結果存在一定的偏差。

當局部氣相的投影面積直徑大于窄矩形窄邊寬度(ds>s)時，其體積和表面積計算方法如下：

Vloc=sSloc-(s(s/2)-π(s/2)2/2)Lloc

(11)

Aloc=2Sloc+(πs/2-s)Lloc

(12)

式中：Vloc為局部氣相的體積，mm3；Aloc為局部氣相的表面積，mm2。

由于假設中氣泡邊緣為半圓形，局部氣相的體積若按照柱體體積來計算，會忽略受擠壓的氣相邊緣為半圓形，將部分液相區域算作氣相，導致計算結果偏大。因此需減去多余的部分體積，如式(11)所示。局部氣相的表面積若按照2倍投影面積來計算，同樣忽略了氣相邊緣為半圓形的部分，導致計算結果偏小，因此需加上部分表面積，如式(12)所示。

體重增加的速度和多少是影響肚中黑線產生的重要因素，所以預防肚中黑線的有效方法就是不要快速發胖，孕期體重應該漸進式地增加。建議整個懷孕過程中的體重增加控制在11～14千克。

當分區內無氣泡邊界，即Lloc=0時，如圖10中環狀流條件下的b區域，式(11)、(12)變為：

Vloc=sSloc

(13)

Aloc=2Sloc

(14)

因此窄矩形通道寬邊方向上局部空泡份額αloc與局部界面濃度ai,loc(m-1)的計算式為：

(15)

(16)

式中：Vloc,i為分區內第i個氣相空間的體積，mm3；Aloc,i為分區內第i個氣相空間的表面積，mm2；Vchannel,loc為局部分區內的通道體積，mm3。

3 結果驗證

對每個實驗工況下的兩相流動視頻，以0.1 s為時間間隔進行圖像提取和數據處理，每個工況處理數據時長為60 s。對得到的空泡份額、界面濃度和局部參數進行驗證。

1) 空泡份額

Ishii[17]提出的適用于窄距形通道的空泡份額經驗關系式如下：

jg/α=C0j+Vgj

(17)

(18)

式中，Vgj為漂移速度，m/s。

Ishii給出的適用于圓管的漂移速度經驗關系式如下。

攪混流：

(19)

環狀流：

Vgj=(1-α)·

(20)

Woldesemayat等[18]證明了此經驗關系式適用于窄矩形通道。將實驗得到的上、中、下游空泡份額與經驗關系式進行對比，如圖11所示。可看出，實驗得到的空泡份額與經驗公式的相對偏差大部分在10%以內。

圖11 空泡份額實驗值與經驗公式計算值的對比Fig.11 Experiment value vs calculated valueof empirical formular for void fraction

2) 界面濃度

由于暫無窄通道內攪混流和環狀流界面濃度經驗關系式，自定義系數K來驗證界面濃度計算方法的合理性。K表示界面濃度與空泡份額的相關性，定義如下：

(21)

從圖5、6可看出，在攪混流和環狀流區間內，空泡份額越大，相界面邊界越簡單，中央氣空間的幾何形狀越接近于圓柱體。當拍攝空間內充滿氣體時，界面濃度和空泡份額可近似按照柱體的計算方法計算。則K的表達式如下：

(22)

當空泡份額較小時，攪混流相界面較復雜，此時相界面面積中側面積部分占比增加，因此K應大于1，且攪混程度越深，K值偏離1越大。

將不同工況下的K值作圖，如圖12所示。

圖12 空泡份額與界面濃度的相關性Fig.12 Correlation between void fraction and interfacial area concentration

從圖12可看出，隨著空泡份額的增加，K>1逐漸向K=1靠近，且總是大于1，符合上述理論分析，因此認為界面濃度計算方法是合理的。

3) 局部參數

為證明分區求解局部參數計算方法的合理性，將拍攝區域沿寬邊方向等分為10份。對每份區域內計算得到的空泡份額和界面濃度求平均值，與整體圖像所得結果進行對比，如圖13所示。

空泡份額和界面濃度分區計算結果與未分區計算結果整體相對偏差在5%以內。存在偏差的主要原因在于分區參數計算過程中，將所有投影面積直徑小于窄矩形窄邊寬度的氣相(ds≤s)，均假設為球形來計算。此假設下，可能將大氣泡的一角當作球形氣泡，偏離了實際情況。

圖13 整體圖像計算結果與分區計算結果的對比Fig.13 Whole image calculation compared with partition calculation results

4 結論

本文通過搭建窄矩形通道內空氣-水向上流動實驗回路，獲取了攪混流、環狀流流型下相界面圖像，提出了適用于復雜相界面工況下的空泡份額和界面濃度的計算方法，得到如下主要結論。

1) 基于氣相拆分的思想，實現了窄矩形通道內攪混流和環狀流下空泡份額的計算。得到的空泡份額與經驗公式計算值的相對偏差在10%以內，證明了空泡份額計算方法的合理性。

2) 基于氣相拆分的思想，實現了窄矩形通道內攪混流和環狀流下界面濃度的計算。自定義界面濃度與空泡份額的相關系數K，實驗中K值隨著空泡份額的增加逐漸趨近于1，符合理論分析情況。

3) 基于窄矩形通道的特點和分區計算的思想，提出了適用于計算分區通道內，局部位置空泡份額和界面濃度的計算方法。將拍攝區域沿通道寬邊等分為10份，求得分區內空泡份額與界面濃度的平均值，與未分區計算得到的參數進行了對比，相對偏差整體在5%以內。