壓水堆燃料棒束通道內過冷沸騰分析

張 蕊，干富軍，左巧林，田文喜，蘇光輝，秋穗正

（1.西安交通大學動力工程多相流國家重點實驗室，陜西西安 710049；2.上海核工程研究設計院，上海 200233）

壓水堆燃料棒束通道內過冷沸騰分析

張 蕊1，干富軍2，左巧林2，田文喜1，蘇光輝1，秋穗正1

（1.西安交通大學動力工程多相流國家重點實驗室，陜西西安 710049；2.上海核工程研究設計院，上海 200233）

使用Fluent14.5兩流體模型中的RPI（Rensselaer Polytechnic Institute）壁面沸騰模型，對堆芯燃料棒束通道內過冷沸騰現象進行數值模擬，得到了通道內的流場、溫度場以及空泡份額的分布，分析了定位格架和攪混翼的存在對熱工水力特性的影響。數值結果表明，格架的存在會造成很大的壓降，而攪混翼會對流場、溫度場和空泡份額分布產生顯著影響；RPI壁面沸騰模型的模擬結果與Bartolemei試驗數據符合很好。

過冷沸騰；燃料棒束；定位格架

過冷沸騰指壁面溫度高于飽和溫度，而流道內的液體仍處于過冷狀態的流動。核電廠的工程應用設備（如反應堆堆芯、蒸汽發生器）中也有可能出現過冷沸騰。過冷沸騰可極大地提高換熱表面的換熱系數，但由于換熱表面熱流密度較高，易發生傳熱惡化，當過冷沸騰發生偏離核態沸騰（DNB）現象時，由于加熱表面被氣膜覆蓋，加熱表面溫度瞬間飛升并發生燒毀。因此，進行堆芯流道內過冷沸騰分析很有必要。目前常規的單通道和子通道分析，不能詳細地描述燃料組件內的局部流場細節特征。隨著計算機技術的飛速發展，大型計算流體動力學商用軟件Fluent在諸多工程領域獲得了廣泛應用。因此，采用CFD方法對反應堆燃料組件過冷沸騰模擬顯得十分可行而且必要，同時也具有重要的工程價值和學術意義。

目前燃料棒束通道內CFD模擬主要集中在單相，Navarro等［1］用商用軟件CFX11.0模擬了5×5棒束帶定位格架結構內的單相流動。Liu等［2］采用商用軟件Fluent12.0分析了不同湍流模型（RNG k-ε模型、realizable k-ε模型、低Re k-ε模型、標準k-ε模型、SST k-w模型和雷諾應力模型）對5×5棒束通道計算結果的影響。Li等［3］對含攪混翼燃料組件特征通道內的單相流動進行了模擬研究。

而對燃料棒束通道內兩相流動換熱現象的研究還十分少見。因此，本文應用Fluent14.5計算軟件，采用Kurul等［4］所提出的RPI壁面沸騰模型，對燃料棒束單通道帶攪混翼和定位格架過冷沸騰現象進行模擬分析。

1 數學物理模型

過冷沸騰的數值分析基于歐拉兩流體模型進行，對氣相和液相分別求解質量、動量、能量守恒方程，同時求解相間質量、動量、能量傳遞模型。

1.1 歐拉兩流體模型

質量守恒方程：

動量守恒方程：

能量守恒方程：

式中：αi、ρi、vi、Si、pi、=τi、hi和qi分別為i相的體積份額、密度、速度、源項、壓力、應力張量、比焓和熱流量；˙mij和Qij分別為i相至j相的質量和能量傳遞。

1.2 相間動量傳遞模型

氣相、液相間作用力包括相間拽力、升力、壁面潤滑力、湍流耗散力，拽力FD為單位體積內氣泡施加在液相上的力：

其中：CD為升力系數，由Ishii關系式［5］計算；μf為液相黏度；Aif為相界面密度：

Re為相間相對雷諾數：

式中：dg為氣泡直徑；vg和vf分別為氣相、液相速度；αg為氣相體積份額；ρf為液相密度。

相間升力為液相速度梯度施加在氣相上的力，由下式計算：

其中，CL為升力系數，由Moraga關系式［6］計算。

壁面潤滑力作用在氣泡上，推動氣泡向主流方向運動，由下式計算：

其中：Cwl為壁面潤滑力系數，由Antal模型計算［7］；nw為壁面法向矢量。

湍流耗散力用來描述兩相間的湍流相互作用，由下式計算：

其中：CTD＝1；σfg＝0.9。

1.3 相間能量傳遞

氣相脫離壁面后，進入過冷流體時向液體釋熱，單位體積內的傳熱量為：

其中：kf為液相導熱系數；Nu采用Ranz-Marshall公式［8］計算。

1.4 相間質量傳遞-RPI壁面沸騰模型

采用RPI壁面沸騰模型計算壁面上的沸騰現象及氣泡的冷凝等傳質過程。RPI壁面沸騰模型包括基本熱流密度分配模型和輔助模型。

1）壁面熱流密度分配模型

RPI模型將壁面傳遞給液體的熱流密度qw分為3個部分，即單相液體對流帶走的熱流密度qc、蒸發熱流密度qe和淬滅熱流密度qq。

式中：hc為單相對流換熱系數；Af為壁面上液相潤濕面積；Tw為壁面溫度；Tf為液相溫度；Vd為基于氣泡脫離直徑的氣泡體積；Nw為氣泡核化密度；hfg為汽化潛熱；f為氣泡脫離頻率。

2）輔助模型

RPI模型的輔助模型包括液相潤濕面積Af、氣泡脫離頻率模型f、氣泡核化密度Nw、氣泡脫離直徑Dw等模型，計算關系式分別為：

式中，Ja為過冷度的Jacob數。

2 模型驗證

Bartolemei等［9－10］對不同直徑的豎直圓管內向上流動的過冷沸騰進行試驗研究，試驗獲得不同壓力、熱流密度、進口質量流密度、進口過冷度下的空泡份額沿軸向分布以及壁面溫度沿軸向分布。選取其中1個工況（d＝15.4 mm，G＝900 kg／（m2·s），q＝570 k W／m2，p＝4.5 MPa）用于模型的驗證。圖1a為根據試驗段參數建立的幾何模型及計算工況，圖1b為該試驗工況下Fluent模擬的沿軸向平均空泡份額、壁面平均溫度和液相平均溫度結果與試驗測量數據的對比。從圖中可看出，計算結果和試驗數據符合良好，因此RPI模型可用于過冷沸騰的模擬。

圖1 模型驗證Fig.1 Model validation

3 幾何模型及網格劃分

本文選取的計算區域為帶有格架和攪混翼、由4根1／4棒束圍成的冷卻劑通道。計算區域的幾何模型及網格劃分如圖2所示。冷卻劑從幾何模型的底部流入，在燃料棒表面熱流的加熱作用下，從單相對流換熱到發生過冷沸騰。幾何模型整體結構簡單，在定位格架和攪混翼局部位置幾何復雜，因此針對該幾何模型采用基于STAR-CCM＋的trimmer網格，整個幾何域內為六面體網格，在格架和攪混翼幾何較為復雜的區域進行網格劈分和加密。冷卻劑通道的4個邊界設置為兩對周期性邊界，允許流體在該邊界上流入流出，來考慮臨近流道對該流道的影響。

圖2 簡化幾何模型和網格劃分Fig.2 Simplified computational geometry and mesh

計算幾何的進口邊界設為速度進口邊界，出口邊界為壓力出口邊界，定位格架和攪混翼設置為無滑移壁面邊界條件，燃料棒表面設為給定熱流密度邊界。具體邊界條件設置列于表1。

表1 邊界條件設置Table 1 Calculation condition

為分析網格數量對計算結果的影響，分別用5套不同數量的網格對計算工況進行分析計算，網格參數列于表2，采用軸向平均空泡份額作為網格無關性的考量依據。圖3為沿軸向平均空泡份額的分布?？煽闯?，第4套網格和第5套網格的結果基本一致，可滿足網格無關性要求，因此以下的計算中均采用第4套網格。

表2 網格無關性分析Table 2 Grid independence analysis

圖3 沿軸向平均空泡份額的分布Fig.3 Void fraction profile along flow direction

4 結果分析

計算工況中系統的壓力選取一回路壓力15.8 MPa，燃料棒束表面熱流密度根據堆芯功率計算得到，冷卻劑的物性查表得到。

圖4為流道內的氣液兩相流線圖。從圖中可明顯看出，在攪混翼上游的流道，流線排列整齊近乎平行分布，經過攪混翼的導向之后，明顯觀察到通道內出現了渦旋。而氣相在攪混翼的下游會向流道中心集中，這是由于攪混翼的存在會在下游形成漩渦，由漩渦產生的離心力會使液相發生旋轉，同時使氣相積聚在漩渦中心，如圖4b所示。

圖4 氣液兩相流線圖Fig.4 Streamline presentation of gas-liquid two-phase flow

圖5為不同高度橫截面上平均壓力沿軸向變化曲線。從圖中可看出，格架位置壓降劇增，這是因為格架的存在會使該位置處的流通面積減小，從而導致很大的壓降，而攪混翼的存在對壓降影響不大。

以格架條帶上表面攪混翼根部所在平面為基準面，圖6示出攪混翼下游5、10、20、40 mm高度截面處的液相橫向速度矢量圖。從圖中可看出，流體經過攪混翼的導流作用會在下游形成渦旋，燃料棒束通道的中心會形成渦旋，在攪混翼兩側也會形成與中心渦旋方向相反的兩個小漩渦。這是由于攪混翼的存在，部分流體在攪混翼的導向下由周期性邊界流出；另一部分流體在攪混翼的導向下流動接觸燃料棒流向另一側邊界（圖6a中箭頭所示）。1對周期性邊界上的參數相同，從而形成如圖6a所示的渦旋。隨著軸向高度的增大，渦旋的強度逐漸減小、渦旋影響面積逐漸增大。在攪混翼下游40 mm處，中心渦旋和兩個方向相反的小渦旋仍存在，但強度明顯減小。

圖5 壓力沿軸向變化曲線Fig.5 Pressure change along flow direction

圖6 液相橫向速度矢量圖Fig.6 Liquid lateral velocity vector field distribution

圖7為不同高度截面處液相溫度的分布。5 mm高度截面上液相溫度隨著到燃料棒表面的距離增大而減小，距離燃料棒表面越近溫度越高。受攪混翼的影響，該截面上的流體按圖示方向流進或流出邊界。如x正y負象限內的攪混翼會將流體導向x＋和y－方向，致使冷流體流向邊界，將熱流體推擠出邊界。因此，會在4個周期性邊界上形成如圖所示的分布。在高度方向，隨著高度的增加，流體不斷被加熱，截面上的平均溫度不斷上升。由圖中明顯看出，經過攪混翼后流體在流道中會發生旋轉，使得溫度場的分布也逐漸逆時針方向旋轉。

圖8為燃料棒束表面的壁面溫度和熱流密度分布。圖中曲線1表示壁面溫度隨高度的變化，從圖中可看出，壁面溫度整體呈上升趨勢，但在格架和攪混翼位置會出現明顯的降低，這是由于格架和攪混翼的存在會對流場起到攪混作用，尤其是攪混翼會將冷流體導向燃料棒，使燃料棒表面溫度出現明顯降低。曲線2～5分別表示過冷沸騰RPI模型中壁面總熱流、蒸發熱流、單相對流熱流以及淬火熱流。過冷流體進入流道，隨著過冷度的降低在某一高度處發生過冷沸騰，單相對流換熱逐漸減小，蒸發熱流逐漸增加。從圖中可看出，蒸發熱流曲線和單相對流曲線在格架和攪混翼處會分別出現降低和升高，這是由于流場被攪混后會使冷流體流向或直接沖刷壁面，該位置由過冷沸騰轉變為單相強迫對流換熱，壁面總熱流不變的情況下，單相對流換熱量上升，從而導致蒸發熱流下降。

圖7 液相溫度的分布Fig.7 Distribution of liquid temperature

圖8 壁面溫度和熱流密度分布Fig.8 Distributions of wall temperature and heat flux

圖9為軸向5、10、20、40 mm高度截面處空泡份額分布?？张莘蓊~隨著與燃料棒表面距離增大而減小。由于攪混翼的存在，會將流道中心部分冷流體導向燃料棒表面，引起該處空泡份額降低（圖9a中1位置）；攪混翼會將冷流體導向邊界，而周期性邊界的存在，致使該邊界上有冷流體的流入（圖9a中2位置）?？张莘蓊~的分布隨高度變化發生不同程度的旋轉，攪混翼會撥動過冷流體流向燃料棒表面。20、40 mm截面上空泡份額的最大值出現在y方向的兩個邊界中心處，這正是兩個小渦旋出現的位置，這是由于流道中該邊界處的流通面積最小，該處的流體溫度和空泡份額較高，而氣泡在渦旋的作用下向渦旋中心聚集所致。

5 結論

本文采用經過試驗驗證的RPI過冷沸騰模型，對含定位格架和攪混翼的燃料棒束單通道內過冷沸騰現象進行了數值模擬及分析，得到以下結論：

1）通過與過冷沸騰基準試驗Bartolemei試驗數據進行對比，RPI模型可準確地預測過冷沸騰現象；

2）準確地模擬了含定位格架和攪混翼的單通道內的過冷沸騰現象，得到了流道內的流場、溫度場以及空泡份額的分布；

3）定位格架的存在會造成很大壓降，而攪混翼的存在會使流道形成渦旋，導致溫度和空泡份額的分布沿流道發生旋轉。

圖9 空泡份額分布Fig.9 Distribution of vapor volume fraction

［1］ NAVARRO M A，SANTOS A A C.Evaluation of a numeric procedure for flow simulation of a 5×5 PWR rod bundle with a mixing vane spacer［J］.Progress in Nuclear Energy，2011，53（8）：1 190-1 196.

［2］ LIU C C，FERNG Y M，SHIH C K.CFD evaluation of turbulence models for flow simulation of the fuel rod bundle with a spacer assembly［J］.Applied Thermal Engineering，2012，40：389-396.

［3］ LI X C，GAO Y.Methods of simulating largescale rod bundle and application to a 17×17 fuel assembly with mixing vane spacer grid［J］.Nuclear Engineering and Design，2014，267：10-22.

［4］ KURUL N，PODOWSKI M Z.On the modeling of multidimensional effects in boiling channels［C］∥ANS Proceedings of 27th National Heat Transfer Conference.Minnesota，USA：［s.n.］，1991：28-31.

［5］ ISHII M.Two-fluid model for two-phase flow［J］.Multiphase Science and Technology，1990，5（1-4）：1-63.

［6］ MORAGA F，BONETTO F，LAHEY R.Lateral forces on spheres in turbulent uniform shear flow［J］.International Journal of Multiphase Flow，1999，25（6-7）：1 321-1 372.

［7］ ANTAL S，LAHEY R，Jr，FLAHERTY J.Analysis of phase distribution in fully developed laminar bubbly two-phase flow［J］.International Journal of Multiphase Flow，1991，17（5）：635-652.

［8］ RANZ W，MARSHALL W.Evaporation from drops［J］.Chem Eng Prog，1952，48（3）：141-146.

［9］ BARTOLEMEI G G，BRANTOV V G，MOLOCHNIKOV Y S，et al.An experimental investigation of true volumetric vapor content with subcooled boiling in tubes［J］.Thermal Engineering，1982，29（3）：132-135.

［10］BARTOLEMEI G G，CHANTURIYA V M.Experimental study of true void fraction when boiling subcooled water in vertical tubes［J］.Thermal Engineering，1967，14（2）：123-128.

Analysis of Subcooled Boiling in PWR Rod Bundle Channel

ZHANG Rui1，GAN Fu-jun2，ZUO Qiao-lin2，TIAN Wen-xi1，SU Guang-hui1，QIU Sui-zheng1
（1.State Key Laboratory of Multiphase Flow in Power Engineering，Xi’an Jiaotong University，Xi’an 710049，China；
2.Shanghai Nuclear Engineering Research ＆Design Institute，Shanghai 200233，China）

The RPI（Rensselaer Polytechnic Institute）wall boiling model，embedded in the context of two-fluid model of Fluent14.5，was applied to simulate subcooled boiling in core fuel bundle channel with mixing vanes and space grid.The subcooled boiling key parameters including flow field，temperature field and void fraction distribution were obtained.The analysis of the effect of mixing vanes and space grid on flow channel thermal-hydraulic characteristics was conducted.It is shown that the existence of space grid can cause a large pressure drop，and mixing vanes have a great influence on the key parameters of flow field and temperature field.The results of RPI model show a good agreement with Bartolemei experiment data.

subcooled boiling；rod bundle；space grid

TL334

：A

1000-6931（2015）09-1579-07

10.7538／yzk.2015.49.09.1579

2014-05-04；

2014-12-08

國家杰出青年科學基金資助項目（11125522）

張 蕊（1990—），女，陜西咸陽人，博士研究生，核能科學與工程專業