周 磊，閆 曉，黃善仿，黃彥平，肖澤軍

（1.中核集團(tuán) 核反應(yīng)堆熱工水力技術(shù)重點(diǎn)實(shí)驗(yàn)室，四川 成都 610041；2.清華大學(xué) 工程物理系，北京 100084）

臨界熱流密度（CHF）發(fā)生時(shí)，冷卻劑對(duì)流換熱系數(shù)減小、燃料元件壁溫升高，嚴(yán)重的可使燃料元件燒毀從而引起放射性泄漏。因此，準(zhǔn)確預(yù)測(cè)CHF對(duì)于反應(yīng)堆設(shè)計(jì)和安全分析具有重要的意義。通常認(rèn)為在環(huán)狀流區(qū)，CHF是由液膜的干涸引起，據(jù)此建立的CHF預(yù)測(cè)模型常稱為液膜干涸（DRYOUT）模型。由于環(huán)狀流型所對(duì)應(yīng)的含汽率范圍廣泛（約0.1～1.0），DRYOUT不僅對(duì)于沸水堆，而且對(duì)于壓水堆的安全分析及瞬態(tài)過程均具有很重要的意義。近年來，DRYOUT模型倍受重視［1－2］。

板型燃料元件結(jié)構(gòu)緊湊，有利于提高堆芯功率密度、減小堆芯體積，在一些研究堆、新型核動(dòng)力堆和國(guó)外生產(chǎn)堆中已得到應(yīng)用。其冷卻劑流道是窄縫隙的矩形通道。矩形窄通道在一些電子設(shè)備的冷卻中也有重要應(yīng)用，因此相關(guān)熱工水力性質(zhì)已成為研究熱點(diǎn)［3］。

文獻(xiàn)中已有的DRYOUT模型主要是面向圓管常規(guī)通道開發(fā)的，對(duì)于矩形窄通道是否可用、性質(zhì)如何，缺乏足夠的實(shí)驗(yàn)數(shù)據(jù)驗(yàn)證。本文通過優(yōu)化DRYOUT模型中重要本構(gòu)關(guān)系式的組合、適當(dāng)選取模型中的沉積率和夾帶率等公式，以得到適用于矩形窄通道的DRYOUT模型。

1 流動(dòng)結(jié)構(gòu)和控制方程

圖1為矩形窄通道中DRYOUT模型的示意圖。由于熱量的加入，沿流動(dòng)方向含汽率逐步升高，當(dāng)達(dá)到某一轉(zhuǎn)變值xt時(shí)流型進(jìn)入環(huán)狀流；在環(huán)狀流區(qū)，蒸汽在通道中間流動(dòng)形成汽芯，液體則一部分沿壁面形成液膜，一部分以?shī)A帶液滴的形式存在于汽芯中。由于蒸發(fā)、液滴沉積及夾帶等過程的共同作用，液膜流量沿流向逐漸變小，液膜厚度逐漸變薄。對(duì)于通道均勻加熱的情況，若在通道出口處液膜厚度足夠小，此時(shí)的熱流密度即為臨界熱流密度。

圖1 環(huán)狀流區(qū)干涸模型示意圖Fig.1 Sketch of DRYOUT model in annular flow regime

取通道上長(zhǎng)度為Δz的一段，液膜的連續(xù)性方程為：

其中：Wlf為液膜質(zhì)量流量，kg／s；me、md、mv分別為夾帶率、沉積率和蒸發(fā)率，kg／（m2·s）；ξi和ξh分別為汽芯和液膜的交界面周長(zhǎng)和熱周，m；z為流動(dòng)方向的坐標(biāo)，m。

對(duì)于圓管通道，熱周和界面周長(zhǎng)相等，并可取為通道周長(zhǎng)。而矩形通道熱周取決于具體加熱情況，而界面長(zhǎng)度可近似取為通道周長(zhǎng)。令Δz→0，即得Wlf滿足如下控制方程：

在DRYOUT模型中通常假設(shè)液膜的溫度近似為對(duì)應(yīng)壓力下的飽和溫度，加入的熱量全部用于液體蒸發(fā)，因此有：

其中：qw為壁面熱流密度，W／m2；hfg為汽化潛熱，J／kg。

因此，只要給定環(huán)狀流起始點(diǎn)和沉積率、夾帶率公式及相應(yīng)的本構(gòu)方程，通過式（2）即可求得任意位置處的液膜質(zhì)量流速，進(jìn)而通過迭代過程可求得在出口發(fā)生干涸時(shí)對(duì)應(yīng)的CHF。

2 本構(gòu)方程

2.1 環(huán)狀流起始點(diǎn)條件

對(duì)于環(huán)狀流起始點(diǎn)，文獻(xiàn)中發(fā)表了大量的判定關(guān)系式。敏感性計(jì)算表明，在DRYOUT模型中，當(dāng)加熱板長(zhǎng)度較小或臨界含汽率較小時(shí)，環(huán)狀流起始點(diǎn)對(duì)CHF有顯著影響，因而環(huán)狀流起始點(diǎn)判據(jù)須慎重選取。本文選取著名的Mishima ＆ Hibiki關(guān)系式［4］作為判據(jù)：

其中：Jv為汽相表觀速度；σ為表面張力，N／m；ρv為汽相密度，kg／m3；Δρ為液相和汽相密度差，kg／m3；g為重力加速度，即9.81m／s2。

上式是專門針對(duì)矩形窄通道建立的。

無量綱粘度數(shù)Nμf定義如下：

其中：μf為液相粘度，Pa·s；ρf為液相密度，kg／m3。

王俊峰［5］通過流型可視化實(shí)驗(yàn)證實(shí)矩形窄縫通道中環(huán)狀流存在且可用式（4）進(jìn)行較好預(yù)測(cè)。

2.2 沉積率和夾帶率關(guān)系式

敏感性計(jì)算表明，沉積率和夾帶率對(duì)液膜質(zhì)量流量或液膜厚度有著全程性的影響，因而對(duì)計(jì)算的CHF有顯著影響。但文獻(xiàn)中公開發(fā)表的沉積率、夾帶率關(guān)系式數(shù)量眾多。本文僅選取5組應(yīng)用廣泛的關(guān)系式進(jìn)行簡(jiǎn)要介紹和比較，這5組關(guān)系式均是在大量的液滴沉積實(shí)驗(yàn)數(shù)據(jù)的基礎(chǔ)上擬合得到的，有的還被用于已有的反應(yīng)堆熱工水力分析程序中，具有較好的通道適用性。

1）Wurtz－Sugawara關(guān)系式（case 1）

由Sugawara對(duì)Wurtz關(guān)系式進(jìn)行改進(jìn)得到，以便更好地體現(xiàn)壓力的影響，改進(jìn)后的公式在COBAR－TF子通道分析程序中得到應(yīng)用［6］：

其中：ug為汽相平均速度，m／s；τi為汽液界面摩擦應(yīng)力，N／m2。

其中：uf為液膜平均速度，m／s。系數(shù)k為：

其中：Reg為汽相雷諾數(shù)。

ks為等效沙粒粗糙度（無量綱），由實(shí)驗(yàn)數(shù)據(jù)擬合得到：

式中：δ為液膜厚度，m。

液滴沉積率和液滴濃度成正比，即：

其中：C為液滴濃度，kg／m3；kd為沉積率系數(shù)，m／s，由下式計(jì)算：

其中：Sc為無量綱斯密特?cái)?shù)。

汽芯中液滴濃度C定義如下：

式中：αd、αg分別為液滴和汽相的體積份額。

2）Govan關(guān)系式（case 2）

Govan等［7］直接對(duì)沉積和夾帶實(shí)驗(yàn)數(shù)據(jù)進(jìn)

行擬合，得到：

其中：min表示取小值函數(shù)；D為通道等效直徑，m。

夾帶率為：

其中：Gg、Glf分別為汽相和液膜質(zhì)量流速，kg／（m2·s）。

臨界液膜質(zhì)量流速Glfc由下式給出：

3）Okawa ＆ Tsuyoshi關(guān)系式（case 3）

他們假設(shè)夾帶率和界面剪切應(yīng)力成正比，而和表面張力成反比，首先定義無量綱數(shù)π［8］：e

式中：fi為剪切應(yīng)力系數(shù)。

則夾帶率可表示為：

并且只有當(dāng)液膜雷諾數(shù)Relf足夠大時(shí)，才發(fā)生液滴的夾帶，即：

對(duì)于沉積率，Okawa和Tsuyoshi等的研究結(jié)論是高壓下Govan公式較好，而低壓下Sugawara公式更好，因此他們的沉積率關(guān)系式是二者的線性組合（對(duì)于Govan公式有少量修改以消除原公式的不連續(xù)性）：

權(quán)重w的大小取決于壓力以及含汽率的情況。

4）Kataoka ＆Ishii關(guān)系式（case 4）

Ishii ＆ Kataoka給出的沉積率關(guān)系式［9］如下：

其中：E為汽芯中夾帶液滴在全部液相質(zhì)量中的比重。

Ishii ＆Kataoka給出的夾帶率公式考慮了入口效應(yīng)等，因此結(jié)構(gòu)上較為復(fù)雜，即：

E∞為平衡夾帶份額，按下式計(jì)算：

而夾帶Webber數(shù)定義為：

當(dāng)液膜雷諾數(shù)小于某臨界值時(shí)，夾帶不會(huì)發(fā)生，臨界液膜雷諾數(shù)為：

5）Okawa ＆ Kataoka關(guān)系式（case 5）

Okawa和 Kataoka［10］對(duì)沉積實(shí)驗(yàn)數(shù)據(jù)分析后認(rèn)為，實(shí)驗(yàn)數(shù)據(jù)難以使用統(tǒng)一的函數(shù)形式表達(dá)。當(dāng)液滴濃度較高時(shí)，沉積率系數(shù)可視為液滴濃度的函數(shù)，即：

式中：C＊為無量綱液滴濃度；k＊d為無量綱沉積率系數(shù)。

當(dāng)液滴濃度很低時(shí)，沉積系數(shù)是摩擦速度的函數(shù)：

uτ為摩擦速度，m／s，定義為：

其中：τw為壁面摩擦應(yīng)力，N／m2。夾帶率公式的擬合思路和Okawa ＆Tsuyoshi等相同，只是函數(shù)形式略有改變：

無量綱量πe1為：

臨界汽相表觀速度Jgc為：

若汽相表觀速度低于臨界值，則夾帶率為0。

2.3 沸騰夾帶率

上述沉積率和夾帶率公式均是根據(jù)絕熱等溫實(shí)驗(yàn)數(shù)據(jù)擬合得到的，并未考慮汽化過程的影響。實(shí)際上對(duì)于有加熱的情況，汽泡不斷產(chǎn)生并進(jìn)入汽芯，也會(huì)形成液滴夾帶，因此對(duì)于有加熱的情況，總的夾帶率應(yīng)為：

meb為沸騰夾帶率，kg／（m2·s）。本文采用 Ueda公式［11］計(jì)算：

其中，液膜厚度可根據(jù)液膜質(zhì)量流量求得。

2.4 統(tǒng)計(jì)參量和評(píng)價(jià)實(shí)驗(yàn)數(shù)據(jù)庫(kù)

由于缺乏矩形窄通道沉積率和夾帶率實(shí)驗(yàn)數(shù)據(jù)，不能對(duì)上述關(guān)系式進(jìn)行直接評(píng)價(jià)。但可將這些關(guān)系式引入DRYOUT模型中，通過和CHF實(shí)驗(yàn)數(shù)據(jù)的對(duì)比，得到最優(yōu)的DRYOUT模型。為此，收集準(zhǔn)備了矩形窄通道CHF實(shí)驗(yàn)數(shù)據(jù)庫(kù)，參數(shù)范圍列于表1。表中DRYOUT數(shù)據(jù)是根據(jù)Mishima ＆Hibiki關(guān)系式［4］選出的。

按照習(xí)慣，定義如下統(tǒng)計(jì)參數(shù)。模型計(jì)算偏差CHFR，表示CHF模型計(jì)算值與實(shí)驗(yàn)值之比：

模型平均偏差，即CHFR的平均值為：

標(biāo)準(zhǔn)偏差僅體現(xiàn)CHFR的分散程度，若CHFR分布散度很小但總體傾斜很大，結(jié)果仍不好，因此引入CHFR的分布帶寬BW：

其中，N表示樣本總數(shù)。

標(biāo)準(zhǔn)偏差為：

表1 矩形窄通道CHF實(shí)驗(yàn)數(shù)據(jù)庫(kù)Table 1 Experimental database for CHF in narrow rectangular channel

3 模型開發(fā)結(jié)果

3.1 模型優(yōu)化過程

圖2、3示出了由5組沉積率和夾帶率關(guān)系式計(jì)算得到的CHFR的分布情況，本文僅給出CHFR隨質(zhì)量流速和含汽率的分布。表2列出了采用不同沉積率和夾帶率關(guān)系式計(jì)算結(jié)果的統(tǒng)計(jì)特征。

從圖2、3和表2可看出，case 1、4的計(jì)算結(jié)果總體上小于實(shí)驗(yàn)值，模型平均偏差僅為0.863 和 0.895，CHFR 的 分 散 度 也 很 大。case 2的結(jié)果絕大部分分布在±20%之內(nèi)，標(biāo)準(zhǔn)偏差也較小，但從CHFR－xC的分布圖上看出，當(dāng)含汽率大于0.25時(shí)，模型計(jì)算結(jié)果偏大。case 5的結(jié)果分布趨勢(shì)和case 2相似，但分散度更大，且有更多點(diǎn)超出了20%偏差線。只有case 3，全部結(jié)果的精度在±20%以內(nèi)，且CHFR的分布帶寬BW小于0.3，標(biāo)準(zhǔn)偏差為0.043，這 是 因 為 其 沉 積 率 公 式 是 Wurtz－Sugawara關(guān)系式和Govan關(guān)系式的線性組合，因此基本消除了case 2中CHF在高含汽率時(shí)偏大的現(xiàn)象，結(jié)果在5組關(guān)系式中最好。因此，Okawa ＆ Tsuyoshi關(guān)系式被用于本文的DRYOUT模型中。表3列出了本文DRYOUT模型最終采用的主要關(guān)系式。

3.2 與已有DRYOUT模型的比較

對(duì)于本文收集的實(shí)驗(yàn)數(shù)據(jù)，同時(shí)采用Celata模型和Kataoka模型進(jìn)行了計(jì)算。計(jì)算結(jié)果的對(duì)比見表4和圖4?？煽闯?，Celata模型計(jì)算結(jié)果全部大于實(shí)驗(yàn)值，總體偏大約30%。Kataoka模型雖引入了臨界液膜的概念，但計(jì)算結(jié)果總體上仍偏大，且含汽率較低時(shí)計(jì)算誤差更大，這是因?yàn)樵撃Ｐ筒捎昧瞬煌沫h(huán)狀流起始點(diǎn)判據(jù)，含汽率低于0.2時(shí)基本超出了模型的推薦使用范圍。而本文模型對(duì)于全部實(shí)驗(yàn)數(shù)據(jù)計(jì)算精度達(dá)到20%，且標(biāo)準(zhǔn)偏差最小，性能明顯優(yōu)于已有模型。

圖2 Wurtz－Sugawara（a、d）、Govan（b、e）和 Okawa ＆ Tsuyoshi（c、f）關(guān)系式計(jì)算結(jié)果Fig.2 Caculation results of Wurtz－Sugawara（a，d），Govan（b，e）and Okawa ＆ Tsuyoshi（c，f）correlations

圖3 Kataoka ＆Ishii（a、c）和 Okawa ＆ Kataoka（b、d）關(guān)系式計(jì)算結(jié)果Fig.3 Caculation results of Kataoka ＆Ishii（a，c）and Okawa ＆ Kataoka（b，d）correlations

表2 不同關(guān)系式得到的DRYOUT模型統(tǒng)計(jì)量Table 2 Statistics of DRYOUT models by different correlations

表3 本文DRYOUT模型采用的主要關(guān)系式Table 3 Correlations adopted in optimized DRYOUT model

表4 本模型和已有模型結(jié)果統(tǒng)計(jì)量對(duì)比Table 4 Comparison of statistics of results by this model and existing models

圖4 各模型計(jì)算值對(duì)比Fig.4 Comparison of calculated results by different models

4 結(jié)論

矩形窄通道以其特有的熱工性能而獲得越來越多的關(guān)注。本文對(duì)DRYOUT型臨界熱流密度模型進(jìn)行了研究，主要結(jié)論如下。

1）已有的DRYOUT模型主要是面向圓管常規(guī)通道開發(fā)的，能否適用于矩形窄通道缺乏足夠的實(shí)驗(yàn)支持；對(duì)于本文的實(shí)驗(yàn)數(shù)據(jù)，已有模型計(jì)算結(jié)果總體偏高，從安全的角度考慮，已有DRYOUT模型不宜直接用于矩形窄通道CHF的預(yù)測(cè)計(jì)算。

2）通過匹配DRYOUT模型中重要中間量的計(jì)算和合理選取沉積率、夾帶率關(guān)系式，得到了對(duì)于矩形窄通道適用的優(yōu)化模型。模型計(jì)算結(jié)果明顯優(yōu)于已有模型，對(duì)本文全部實(shí)驗(yàn)數(shù)據(jù)的計(jì)算精度在20%以內(nèi)，且主要系統(tǒng)參數(shù)對(duì)CHFR的分布無系統(tǒng)性影響，已用于相關(guān)的分析程序中。

［1］CELATA G P，MISHIMA K，ZUMMO G.Critical heat flux prediction for saturated flow boiling of water in vertical tubes［J］.Int J Heat Mass Transfer，2001，44（22）：4 323－4 331.

［2］OKAWA T，KOTANI A，KATAOKA I，et al.Prediction of the critical heat flux in annular regime in various vertical channels［J］.Nucl Eng Des，2004，229（2）：223－236.

［3］HUANG J，HUANG Y P，WANG Q W，et al.Numerical study on effect of gap width of narrow rectangular channel on critical heat flux enhancement［J］.Nucl Eng Des，2009，239（2）：320－326.

［4］HIBIKI T，MISHIMA K.Flow regime transition criteria for upward two－phase flow in vertical narrow rectangular channels［J］.Nucl Eng Des，2001，203（2－3）：117－131.

［5］王俊峰.矩形窄縫通道內(nèi)水流動(dòng)沸騰流型及其演變機(jī)理研究［D］.成都：中國(guó)核動(dòng)力研究設(shè)計(jì)院，2010.

［6］SUGAWARA S.Droplet deposition and entrainment modeling based on the three－fluid model［J］.Nucl Eng Des，1990，22（1－3）：67－84.

［7］GOVAN A H，HEWITT G F.An improved CHF modeling code［C］∥Proceedings of the sedond UK national heat transfer conference.London：Oxford University Press，1981：33－48.

［8］OKAWA T，KITAHARA T，KENJI Y，et al.New entrainment rate correlation in annular twophase flow applicable to wide range of flow condition［J］.Int J Heat Mass Transfer，2002，45（1）：87－98.

［9］KATAOKA I，ISHII M，NAKAYAMA A.Entrainment and deposition rates of droplets in annular two－phase flow［J］.Int J Heat Mass Transfer，2000，43（9）：1 573－1 589.

［10］OKAWA T，KATAOKA I.Correlations for the mass transfer rate of droplets in vertical upward annular flow［J］.Int J Heat Mass Transfer，2005，48（23－24）：4 766－4 778.