桂民洋，田文喜，吳 迪，陳榮華，張 魁，蘇光輝，秋穗正

(西安交通大學(xué) 核科學(xué)與技術(shù)學(xué)院，陜西 西安 710049)

在核反應(yīng)堆的設(shè)計(jì)與安全分析中，臨界熱流密度是重要的限制性熱工參數(shù)，它直接影響核反應(yīng)堆的安全性和經(jīng)濟(jì)性。臨界熱流密度是指沸騰傳熱機(jī)理發(fā)生變化而使核燃料元件表面發(fā)生傳熱惡化時(shí)發(fā)熱元件表面單位面積產(chǎn)生的熱量，如果反應(yīng)堆燃料元件表面發(fā)生沸騰臨界，將會(huì)導(dǎo)致燃料元件表面溫度過高從而造成加熱壁面燒毀，放射性物質(zhì)泄漏，進(jìn)而造成嚴(yán)重的反應(yīng)堆運(yùn)行事故[1]。因此，準(zhǔn)確預(yù)測(cè)臨界熱流密度對(duì)反應(yīng)堆熱工水力設(shè)計(jì)和安全分析有著十分重要的意義。臨界熱流密度模型一般分為兩類：偏離泡核沸騰(DNB)型和干涸型。常用的預(yù)測(cè)方法主要有經(jīng)驗(yàn)關(guān)系式和臨界熱流密度查詢表(CHF look-up table)[2]。

但受限于特定流體種類、熱工邊界、幾何尺寸等實(shí)驗(yàn)條件，大多數(shù)實(shí)驗(yàn)關(guān)系式的適用范圍不能輕易外推，或外推后預(yù)測(cè)精度下降。因此，一些學(xué)者轉(zhuǎn)而從臨界發(fā)生時(shí)的微觀機(jī)理入手，力圖獲得適用范圍更廣的臨界熱流密度機(jī)理模型。如Weisman等[3]、Ying等[4]提出的近壁面氣泡壅塞模型，他們認(rèn)為，發(fā)生臨界時(shí)壁面附近大量汽泡聚集形成一個(gè)汽泡層，當(dāng)汽泡層外沿處的液相流動(dòng)減弱到不能穿透汽泡層到達(dá)加熱壁面時(shí)，汽泡層的汽泡會(huì)結(jié)合形成汽泡膜阻隔傳熱從而發(fā)生DNB；而Lee等[5]、Liu等[6]針對(duì)DNB提出了微液層蒸干模型，認(rèn)為壁面附近會(huì)有由多個(gè)小汽泡結(jié)合而形成一個(gè)拉長(zhǎng)的大汽塊，汽塊下存在一層很薄的微液層，當(dāng)微液層蒸干時(shí)即發(fā)生臨界。對(duì)于干涸，其主要發(fā)生在環(huán)狀流末期，此時(shí)液滴的作用顯著，因而大多數(shù)學(xué)者主要關(guān)注環(huán)狀流流型區(qū)液滴的夾帶沉積現(xiàn)象[7-12]。

作為反應(yīng)堆堆芯熱工水力分析的重要工具，子通道程序中臨界熱流密度基本通過實(shí)驗(yàn)關(guān)系式求得，而本文則基于臨界熱流密度機(jī)理模型，通過與子通道分析方法耦合，實(shí)現(xiàn)對(duì)棒束通道內(nèi)臨界熱流密度的預(yù)測(cè)，其中子通道程序用于計(jì)算通道中的焓分布和流量分布，為臨界熱流密度機(jī)理模型提供邊界條件，進(jìn)而判斷臨界的發(fā)生。

1 數(shù)學(xué)物理模型

本研究中考慮3類模型，即子通道模型、DNB型和干涸型臨界熱流密度機(jī)理模型。其中，子通道模型基于SACOS子通道程序[13]。

1.1 子通道模型

子通道內(nèi)流體的流動(dòng)包含軸向和橫向兩個(gè)方向，其控制方程組(連續(xù)方程、能量方程、軸向動(dòng)量方程和橫向動(dòng)量方程)如下。

連續(xù)方程：

(1)

式中：i、j為相鄰?fù)ǖ谰幪?hào)；Ai為通道i流通截面積，m2；ρi為通道i流體密度，kg/m3；mi為通道i軸向質(zhì)量流量，kg/s。

能量方程：

(2)

軸向動(dòng)量方程：

(3)

式中：pi為通道i的壓力，Pa；g為重力加速度，m/s2；θ為通道傾斜角度，(°)；fT為橫向動(dòng)量因子；ui為通道i的軸向流速，m/s；uj為通道j的軸向速度，m/s；u*為通道間橫向攪混流體流速，m/s；fi為通道i的摩擦系數(shù)；Di為通道i的等效水力直徑，m；Ks為局部阻力系數(shù)；Δz為軸向控制體高度，m。

橫向動(dòng)量方程：

(4)

式中：sij為通道i與j間隙長(zhǎng)度，m；lij為通道i與j間攪混長(zhǎng)度，m；KG為通道i與j間橫流阻力系數(shù)；ρ*為通道間橫流冷卻劑密度，kg/m3，為來流通道的冷卻劑參數(shù)。

1.2 DNB型臨界熱流密度機(jī)理模型

DNB型臨界熱流密度機(jī)理模型主要基于微液層蒸干模型，假設(shè)加熱壁面附近產(chǎn)生的小汽泡結(jié)合形成大汽塊，在汽塊下存在非常薄的液相層，稱為微液層，當(dāng)微液層蒸干即發(fā)生臨界，如圖1所示。

圖1 微液層蒸干模型示意圖

微液層蒸干所對(duì)應(yīng)的熱流密度可表示為：

qDNB=ρfδhfgUB/LB

(5)

式中：δ為汽塊下微液層厚度，m；hfg為汽化潛熱，kJ/kg；UB為汽塊移動(dòng)速度，m/s；LB為汽塊長(zhǎng)度，m。

在微液層蒸干模型中，汽塊下微液層厚度δ、汽塊移動(dòng)速度UB與汽塊長(zhǎng)度LB是模型求解的3個(gè)關(guān)鍵參數(shù)。其中，汽塊移動(dòng)速度可通過軸向方向施加在汽塊上的浮力FB和拖拽力FD的平衡得到，即：

FB=FD

(6)

最終可得到汽塊速度的表達(dá)式為：

(7)

式中：UBL為汽塊徑向位置處的主流速度，m/s；CD為拖拽系數(shù)。

汽塊的長(zhǎng)度LB假設(shè)為Helmholtz臨界波長(zhǎng)，可得：

(8)

微液層的厚度δ由徑向方向施加在汽塊上的力的平衡來確定，主要考慮蒸發(fā)力FE、側(cè)面提升力FL[5]、壁面潤(rùn)滑力FWL和Marangoni力FM[6]。最終可得其表達(dá)式為：

(9)

式中：DB為汽塊當(dāng)量直徑，m；q為加熱壁面熱流密度，W/m2；CWL為壁面潤(rùn)滑系數(shù)；CL為提升力系數(shù)；G為通道質(zhì)量流密度，kg/(m2·s)。

圖2為微液層蒸干模型計(jì)算流程圖。

圖2 微液層蒸干模型計(jì)算流程圖

1.3 干涸型臨界熱流密度機(jī)理模型

干涸型臨界熱流密度機(jī)理模型主要基于環(huán)狀流液膜蒸干模型。加熱通道內(nèi)飽和沸騰區(qū)兩相流流型如圖3所示，主要經(jīng)歷了泡狀流、攪混流，最后形成環(huán)狀流直到液膜逐漸燒干。通道內(nèi)環(huán)狀流的特點(diǎn)是液膜在沿著通道壁面流動(dòng)，同時(shí)夾帶有液滴的汽芯在通道中央流動(dòng)。

圖3 環(huán)狀流液膜蒸干模型示意圖

在環(huán)狀流區(qū)域存在液膜蒸發(fā)、液滴夾帶和液滴沉積的復(fù)雜現(xiàn)象，干涸的發(fā)生可認(rèn)為是三者之間相互競(jìng)爭(zhēng)的結(jié)果。考慮沿通道軸向液膜質(zhì)量流量的變化，可得到液膜質(zhì)量流量的控制方程：

(10)

式中：Wf為液膜質(zhì)量流量，kg/s；md為液滴沉積率，kg/(m2·s)；me為液滴夾帶率，kg/(m2·s)；mv為液相蒸發(fā)率，kg/(m2·s)；P為通道濕周(加熱周長(zhǎng))，m。

一旦確定了環(huán)狀流起始點(diǎn)，通過引入液滴夾帶沉積的本構(gòu)關(guān)系式，由式(10)可得到沿軸向壁面液膜質(zhì)量流量的變化趨勢(shì)，當(dāng)液膜質(zhì)量流量為0時(shí)，即認(rèn)為發(fā)生干涸。

1) 環(huán)狀流起始點(diǎn)

環(huán)狀流起始點(diǎn)采用Wallis[14]的理論，其推薦如下由攪混流到環(huán)狀流轉(zhuǎn)變的關(guān)系式：

(11)

(12)

即可通過通道內(nèi)含汽率確定環(huán)狀流起始點(diǎn)。

2) 液滴的夾帶和沉積率

液滴的夾帶率和沉積率關(guān)系式選用Okawa模型[9]，它是基于大量的液滴實(shí)驗(yàn)數(shù)據(jù)擬合而來，對(duì)于沉積率md，其與氣相中液滴濃度C有關(guān)：

md=kdC

(13)

式中，kd為液滴沉積系數(shù)，表達(dá)式為：

(14)

(15)

針對(duì)液滴夾帶率，從夾帶產(chǎn)生的機(jī)理出發(fā)，主要考慮兩類效應(yīng)：中心氣相對(duì)液相的剪切夾帶(液相雷諾數(shù)Ref>320)和氣液界面處氣泡的破裂造成的夾帶，即：

me=me1+me2

(16)

對(duì)于前者，Okawa提出了如下關(guān)系式：

(17)

式中：fi為相間剪切應(yīng)力系數(shù)；Jg為氣相表觀速度，m/s；δ為液膜厚度，m；σ為表面張力，N/m。

對(duì)于后者，采用了Ueda等[15]建立的夾帶關(guān)系式：

(18)

3) 液相的蒸發(fā)率

模型中假設(shè)環(huán)狀流區(qū)域液膜處于飽和狀態(tài)，則加熱壁面的熱量全部用于液相的蒸發(fā)，則液相的蒸發(fā)率為：

mv=q/hfg

(19)

2 子通道模型與機(jī)理模型的耦合

子通道程序已被證明可較好地計(jì)算棒束通道中的流量和焓，得到通道中的流量和溫度分布；而臨界熱流密度機(jī)理模型從微觀機(jī)理現(xiàn)象入手，可克服經(jīng)驗(yàn)關(guān)系式對(duì)幾何條件、熱工邊界條件的限制要求，通過兩者的耦合計(jì)算，可充分發(fā)揮各自的作用。但一般的機(jī)理模型是基于單通道提出的，對(duì)于子通道還需考慮通道間的質(zhì)量、能量攪混，因而對(duì)于兩者耦合還需做特殊的處理。

2.1 子通道模型與DNB型機(jī)理模型的耦合

一般情況下，DNB主要發(fā)生在過冷沸騰區(qū)，若通道軸向均勻加熱，則DNB最先發(fā)生在出口；但通道軸向功率非均勻分布時(shí)，DNB在通道出口之前某處就有可能發(fā)生，因此，為充分考慮非均勻加熱，軸向通道將被分為若干個(gè)控制體，由子通道程序計(jì)算得到通道的流量、焓參數(shù)傳遞給機(jī)理模型，每個(gè)控制體作為單獨(dú)通道均調(diào)用機(jī)理模型進(jìn)行計(jì)算?？紤]到徑向功率的不均勻性(如1個(gè)子通道為4根加熱棒中間圍成的區(qū)域)，保守來看應(yīng)取徑向熱流密度最大值與式(5)所求的qDNB進(jìn)行比較，從而判斷該位置處是否發(fā)生DNB，耦合后的程序計(jì)算流程如圖4所示。

圖4 DNB模型與子通道分析的耦合計(jì)算

特別地，為考慮通道間能量攪混，模型中熱流密度需由進(jìn)出口焓進(jìn)行修正，即：

(20)

式中：qj為軸向控制體j中熱流密度，W/m2；hout為軸向控制體j出口焓，J/kg；hin為軸向控制體j入口焓，J/kg；G為軸向控制體j質(zhì)量流量，kg/(m2·s)；A為軸向控制體j流通面積，m2；Det為軸向控制體j加熱直徑，m；l為軸向控制體j軸向長(zhǎng)度，m。

2.2 子通道模型與干涸型機(jī)理模型的耦合

液膜蒸干模型是基于單通道(如圓管、矩形通道)提出的，其與子通道模型的耦合需考慮子通道特殊的結(jié)構(gòu)形式，如圖5所示。典型的子通道包含角通道、邊通道和中心通道，一方面，機(jī)理模型需考慮通道間的質(zhì)量、能量攪混；另一方面，需考慮徑向功率的非均勻性，如圖3中的角通道和邊通道，由于存在非加熱壁面，因此壁面液膜軸向質(zhì)量方程中不存在蒸發(fā)項(xiàng)mv，對(duì)于中心通道，由于其周圍加熱棒功率的不一致性，導(dǎo)致液膜厚度不同，也需分開考慮。

圖5 棒束中液膜蒸干模型示意圖

因此，考慮干涸型機(jī)理模型在子通道中的應(yīng)用，式(10)液膜質(zhì)量流量的控制方程需改為以下形式：

(21)

式中：k為液膜標(biāo)識(shí)(k=1～N)；wl為子通道間液相橫向流質(zhì)量流量，kg/(m·s)；ζf為橫向流中液膜所占份額。

1個(gè)子通道中液膜可分為N份(對(duì)于角通道N=2，對(duì)于邊通道N=3，對(duì)于中心通道N=4)，它們共用同一區(qū)域的氣相和液相，對(duì)于液相橫向流，可由子通道程序計(jì)算得到，而液相又可分為壁面流動(dòng)液膜和中心夾帶的液滴，程序中兩者的比例ζf可認(rèn)為和通道中液膜與液滴的濃度比相同。

總地來說，由子通道程序計(jì)算得到各通道環(huán)狀流起始點(diǎn)，可認(rèn)為起始點(diǎn)處液膜均勻覆蓋在各壁面上，同時(shí)由子通道程序可確定通道間流體的橫向流動(dòng)，即可通過液膜質(zhì)量守恒方程及相關(guān)本構(gòu)方程確定各壁面軸向液膜質(zhì)量流量變化趨勢(shì)，當(dāng)某一壁面液膜質(zhì)量流量或空泡份額趨于0時(shí)即認(rèn)為發(fā)生干涸，耦合后的程序計(jì)算流程如圖6所示。

圖6 蒸干模型與子通道分析的耦合計(jì)算

3 驗(yàn)證和分析

子通道程序與臨界熱流密度機(jī)理模型耦合的驗(yàn)證基于EPRI的CHF實(shí)驗(yàn)數(shù)據(jù)庫(kù)[16-17]，實(shí)驗(yàn)工況包含壓水堆(PWR)和沸水堆(BWR)，可用于DNB型和干涸型臨界熱流密度機(jī)理模型的驗(yàn)證，通道類型主要為4×4棒束和5×5棒束，具體的通道幾何與熱工邊界參數(shù)列于表1。

表1 臨界熱流密度機(jī)理模型驗(yàn)證數(shù)據(jù)

本文中，首先由SACOS子通道程序建立棒束通道的幾何模型，給定邊界條件，包括出口壓力、入口流量、入口溫度和平均功率密度，通過子通道程序計(jì)算通道各處局部參數(shù)，進(jìn)而調(diào)用DNB或干涸型臨界熱流密度機(jī)理模型計(jì)算，逐漸增大平均功率密度，直至發(fā)生臨界。圖7～10示出了模型計(jì)算值與實(shí)驗(yàn)值的對(duì)比結(jié)果，其中圖7、8為DNB型臨界熱流密度的對(duì)比結(jié)果，圖9、10為干涸型臨界熱流密度的對(duì)比結(jié)果。

從圖7可看出，驗(yàn)證數(shù)據(jù)點(diǎn)為110個(gè)，94.4%的數(shù)據(jù)點(diǎn)誤差在20%之內(nèi)，機(jī)理模型對(duì)CHF的預(yù)測(cè)趨勢(shì)與實(shí)驗(yàn)結(jié)果一致。從圖8可看出，機(jī)理模型的預(yù)測(cè)結(jié)果滿足一般規(guī)律，且與實(shí)驗(yàn)結(jié)果相比有較好的一致性，臨界熱流密度(q)隨系統(tǒng)壓力的升高、入口過冷度的增大和入口流量的增大而增大。本文的DNB型臨界熱流密度機(jī)理模型是基于Lee和Mudawar[5]最早提出的微液層蒸干模型改進(jìn)而來，其中微液層的存在也被一些可視化實(shí)驗(yàn)所證實(shí)。但機(jī)理模型中仍存在經(jīng)驗(yàn)系數(shù)需人為確定，這也是模型計(jì)算誤差的主要原因。其主要表現(xiàn)在式(9)中的提升力系數(shù)CL，Lee和Mudawar提出了僅與雷諾數(shù)Re相關(guān)的CL表達(dá)式：

圖7 DNB型臨界熱流密度總體驗(yàn)證結(jié)果

圖8 DNB型臨界熱流密度的變化

CL=a1Rea2

(22)

而Beyerlein等[18]建議CL是平均空泡份額和Re的函數(shù)，即它與紊流波動(dòng)和當(dāng)?shù)仄莘植记闆r有關(guān)?；诖耍疚闹心Ｐ筒捎玫腃L表達(dá)式為：

CL=230Re-0.35-0.23exp(1.8α)

(23)

式中，CL的經(jīng)驗(yàn)處理使得模型計(jì)算值存在偏差，且具有一定的不確定性。可預(yù)料，當(dāng)不斷修正上式中的常數(shù)時(shí)，可獲得與實(shí)驗(yàn)值更為接近的預(yù)測(cè)結(jié)果。

從圖9可看出，驗(yàn)證數(shù)據(jù)點(diǎn)為696個(gè)，95.8%數(shù)據(jù)點(diǎn)誤差在20%之內(nèi)，機(jī)理模型對(duì)CHF的預(yù)測(cè)趨勢(shì)與實(shí)驗(yàn)結(jié)果一致，顯示了模型較好的預(yù)測(cè)能力。從圖10可看出，臨界熱流密度隨入口過冷度和流量的增大而增大，隨系統(tǒng)壓力的增大而減小，與實(shí)驗(yàn)趨勢(shì)相同?？傮w來看，干涸型機(jī)理模型可較為真實(shí)地反映棒束通道中壁面附著的液膜逐漸減薄至發(fā)生干涸的物理過程，兩者存在的偏差有兩個(gè)方面的原因：一方面，機(jī)理模型中液滴夾帶沉積模型存在計(jì)算偏差，目前關(guān)于液滴的夾帶沉積模型多基于有限的實(shí)驗(yàn)數(shù)據(jù)，且受限于特定的適用范圍，關(guān)系式的擴(kuò)展性較差，因而引入了一定的不確定性，另一方面，目前干涸型臨界熱流密度機(jī)理模型中并未考慮格架的影響，通常格架的存在會(huì)增強(qiáng)液滴的橫向運(yùn)動(dòng)，進(jìn)而改變液滴的沉積率，因此會(huì)影響臨界熱流密度。

圖9 干涸型臨界熱流密度總體驗(yàn)證結(jié)果

圖10 干涸型臨界熱流密度的變化

4 結(jié)論

本研究以子通道模型為基礎(chǔ)，通過耦合DNB型和干涸型臨界熱流密度機(jī)理模型，實(shí)現(xiàn)了對(duì)棒束通道內(nèi)臨界熱流密度的預(yù)測(cè)，其中借助于子通道程序計(jì)算通道中的焓分布和流量分布，為臨界熱流密度機(jī)理模型提供邊界條件，進(jìn)而判斷臨界的發(fā)生。通過與臨界熱流密度實(shí)驗(yàn)數(shù)據(jù)對(duì)比發(fā)現(xiàn)，耦合程序?qū)Π羰ǖ乐信R界熱流密度具有較好的預(yù)測(cè)精度，臨界熱流密度預(yù)測(cè)值隨相關(guān)熱工參數(shù)的變化規(guī)律與實(shí)驗(yàn)相比也有較好的一致性。

本研究中，基于棒束子通道的臨界熱流密度機(jī)理模型可有效彌補(bǔ)目前實(shí)驗(yàn)關(guān)系式適用范圍受限的不足，也可為后續(xù)臨界熱流密度實(shí)驗(yàn)工況設(shè)計(jì)提供參考，但需注意的是，目前的研究中臨界熱流密度模型的驗(yàn)證范圍仍有限，因此在后續(xù)的研究中有必要進(jìn)行更為充分的驗(yàn)證。同時(shí)，對(duì)于干涸現(xiàn)象，中心氣相、液膜和液滴的理論被普遍接受，除液滴的夾帶和沉積關(guān)系式外，模型中的不確定因素較少，而對(duì)于DNB現(xiàn)象，除微液層蒸干模型，還有近壁面汽泡壅塞、邊界層分離和界面提升等理論，模型中的不確定性參數(shù)也較多，因此后續(xù)的研究也需進(jìn)一步探討臨界熱流密度機(jī)理模型的可靠性。