宋士雄 魏 泉 蔡翔舟 郭 威

基于CFD方法的球床式高溫氣冷堆穩(wěn)態(tài)熱工水力分析

宋士雄1,2魏 泉1,2蔡翔舟1,3郭 威1,3

1（中國科學(xué)院上海應(yīng)用物理研究所 嘉定園區(qū) 上海 201800）
2（中國科學(xué)院大學(xué) 北京 100049）
3（中國科學(xué)院核輻射與核能技術(shù)重點(diǎn)實(shí)驗(yàn)室 上海 201800）

基于CFD通用計(jì)算程序Fluent，開發(fā)了多孔介質(zhì)流固兩相局域非熱平衡模型，對PBMR-400滿功率名義工況下穩(wěn)態(tài)熱工水力行為進(jìn)行了分析。高溫氣冷堆堆芯的熱工水力計(jì)算模型通常為帶有強(qiáng)內(nèi)熱源宏觀的多孔介質(zhì)模型，F(xiàn)luent程序本身的多孔介質(zhì)模型為局域熱平衡模型，無法恰當(dāng)描述堆芯熱工水力行為。通過開啟Fluent的自定義標(biāo)量場(UDS)功能并在動(dòng)量方程上附加阻力源項(xiàng)，實(shí)現(xiàn)了局域非熱平衡模型，能夠精確地計(jì)算球床堆堆芯熱工水力行為。CFD程序使用了非常精細(xì)的網(wǎng)格，可以計(jì)算更加精細(xì)的流場和溫場分布。計(jì)算表明，F(xiàn)luent的結(jié)果與國際上著名的球床堆熱工水力程序THERMIX，TINTE符合較好，可以作為球床堆熱工水力分析工具。

CFD (Computational Fluid Dynamics)，高溫氣冷堆，多孔介質(zhì)

準(zhǔn)確可靠的反應(yīng)堆熱工水力分析對于球床式高溫反應(yīng)堆的設(shè)計(jì)是非常重要的。目前國際上球床式高溫氣冷堆熱工水力行為分析程序THERMIX[1?3]和TINTE[4–5]，均采用二維軸對稱局域非熱平衡多孔介質(zhì)模型，流場方程為附加多孔介質(zhì)阻力源項(xiàng)的準(zhǔn)靜態(tài)方程，冷卻劑溫場和固相球床的溫場分別獨(dú)立求解，各個(gè)模塊之間進(jìn)行迭代并收斂。隨著計(jì)算機(jī)性能的大幅度提高以及算法的改進(jìn)，在反應(yīng)堆設(shè)計(jì)及研究中CFD (Computational Fluid Dynamics)方法開始被采用[6–7]，并且實(shí)現(xiàn)了CFD程序與中子學(xué)程序、系統(tǒng)熱工程序的耦合計(jì)算[8?10]。Becker等[11]基于CFD開發(fā)了球床反應(yīng)堆熱工程序，通過與SANA裝置的實(shí)驗(yàn)結(jié)果比較，驗(yàn)證了程序的適用性。

CFD程序在計(jì)算過程中可以使用更加精確的流體力學(xué)方程與很精細(xì)的網(wǎng)格，因此可以對局部流場、溫場進(jìn)行更加精確的分析。本文基于通用CFD程序Fluent[12–13]開發(fā)了局域非熱平衡多孔介質(zhì)模型。流場求解采用Fluent的SIMPLE算法解算器，氣體冷卻劑的溫場方程求解使用Fluent自帶的能量場方程解算器，固相溫場方程通過附加標(biāo)量場方程UDS (User Defined Scalar)來求解，流固兩相的熱耦合通過交換源項(xiàng)來實(shí)現(xiàn)。使用Fluent的UDF (User Defined Function)功能，實(shí)現(xiàn)了球床邊界與反射層邊界的熱耦合，冷卻劑邊界與反射層邊界的熱耦合。

PBMR-400 (Pebble Bed Module Reactor)是南非國營電力設(shè)計(jì)的熱功率為400 MW、電功率為165MWe的模塊式高溫氣冷球床式反應(yīng)堆[14–15]。反應(yīng)堆堆芯直接連接氦氣輪機(jī)，設(shè)置了具有閉路直接循環(huán)特征的能量轉(zhuǎn)換單元，堆芯燃料球采用環(huán)形布置，具有很強(qiáng)的非能動(dòng)余熱排出能力。反應(yīng)堆堆芯活性區(qū)共有大約452000個(gè)燃料球，內(nèi)外半徑分別為1.0 m、1.85 m，活性區(qū)高度為11 m。堆芯運(yùn)行壓力為9 MPa，入口冷卻劑溫度約為500 °C，出口氦氣溫度約為900 °C，堆芯的平均功率密度為4.8MW·m?3。本文基于Fluent的自定義標(biāo)量場及自定義函數(shù)功能開發(fā)了局域非熱平衡多孔介質(zhì)模型，對PBMR-400滿功率名義工況的熱工水力行為進(jìn)行了分析，驗(yàn)證了程序的適用性。

1 計(jì)算模型和方法

1.1理論模型

本文采用了簡化的二維軸對稱模型，以堆芯中心反射層圓心為對稱軸線，模型包括：(1)堆芯多孔介質(zhì)區(qū)；(2)中心反射層區(qū)；(3)外部反射層。邊界條件：Fluent對邊界條件有幾類方式，通過測試不同的邊界條件對計(jì)算結(jié)果的影響，并考慮程序的收斂性。本文中采取速度入口，壓力出口邊界，在計(jì)算過程中通過調(diào)節(jié)出口表壓，使入口壓力達(dá)到反應(yīng)堆的操作壓力，程序自身的操作壓力設(shè)置為0 MPa。在壓力9 MPa，溫度773 K工況下，入口質(zhì)量流193kg·s?1折算入口速度為4.5875 m·s?1。出口條件為壓力出口，表壓為8.73 MPa。外部反射層熱邊界條件為固定溫度邊界773 K，反射層與球床熱邊界設(shè)置為非耦合邊界，使用UDF自定義邊界條件。簡化后的模型經(jīng)網(wǎng)格無關(guān)性驗(yàn)證后使用48400個(gè)計(jì)算網(wǎng)格如圖1所示。

圖1 計(jì)算網(wǎng)格Fig.1 Calculation mesh.

二維軸對稱模型的控制方程如下，其中物理量均為體積平均的表觀物理量，

連續(xù)性方程：

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

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

動(dòng)量源項(xiàng)：

式中，Ψ是球床的阻力因子，ε為球床的孔隙率，μ為氦氣的動(dòng)力學(xué)粘度。本文使用德國Nuclear Safety Standards Commission KTA[16]推薦的阻力因子：

式中，雷諾數(shù)Re=ρvDp/μ。

冷卻劑能量場方程：

固相球床的能量場方程：

式中，Q.為單位體積內(nèi)的功率密度，W·m?3。為流固兩相之間的熱傳導(dǎo)項(xiàng)。

主要采用KTA[17–18]推薦的物性和熱工參數(shù)：

密度：ρ=48.14P/T(1+0.4446P/T1.2)?1

定壓比熱容：Cp=519 J/kg·K

定容比熱容：Cv=3177 J/kg·K

動(dòng)力學(xué)粘度：μ=3.674×10?7×T0.7Pa·s

氦氣冷卻劑的熱傳導(dǎo)系數(shù)：kf=2.682×10?3×(1+ 1.123×P)×T0.711?0.0002×PW·m?1·K?1

球床的熱傳導(dǎo)系數(shù)的計(jì)算采用Zehner-Schlunder[19]模型。

1.2計(jì)算方法簡介

Fluent帶有一個(gè)基于局域熱平衡假設(shè)的多孔介質(zhì)計(jì)算模型，流固兩相的熱傳導(dǎo)方程合并為一個(gè)總的能量場方程，熱傳導(dǎo)系數(shù)通過流固兩相體積權(quán)重平均得到。對于不含內(nèi)熱源的，或流固兩相之間溫度梯度不大的多孔介質(zhì)模型，局域熱平衡模型的理論計(jì)算能給出和實(shí)驗(yàn)符合較好的結(jié)果。高溫氣冷堆堆芯具有強(qiáng)內(nèi)熱源，流固兩相的溫度梯度較大，局域熱平衡假設(shè)的多孔介質(zhì)模型在這種工況下就不再適用。本文不使用Fluent自帶的多孔介質(zhì)模型，在基于表觀速度的動(dòng)量方程上附加新的阻力源項(xiàng)式(4)表征多孔介質(zhì)區(qū)域的流動(dòng)行為。氦氣冷卻劑的能量場方程(6)使用Fluent的流體能量場方程，開啟標(biāo)準(zhǔn)的湍流模式。固相球床的溫場方程(7)的求解通過開啟Fluent中的UDS來實(shí)現(xiàn)，流固兩相在迭代過程中交換源項(xiàng)，實(shí)現(xiàn)溫場的耦合計(jì)算。堆芯多孔介質(zhì)區(qū)域與中心反射層和外部反射層的熱傳導(dǎo)包括堆芯流體和反射層的熱傳導(dǎo)、堆芯固相球床和反射層的熱傳導(dǎo)，熱邊界的耦合通過在邊界上設(shè)置自定義函數(shù)根據(jù)壁面溫差調(diào)節(jié)熱通量來實(shí)現(xiàn)。

程序計(jì)算流程圖見圖2。

圖2 計(jì)算流程圖Fig.2 Computation program flow diagram.

2 PBMR-400滿功率名義工況熱工水力分析

名義工況下，PBMR-400堆芯冷卻劑入口質(zhì)量流率為193 kg·s?1，入口溫度為500 °C(773 K)，出口溫度約為900 °C(1173 K)，運(yùn)行壓力9 MPa。基準(zhǔn)參考堆芯的功率密度(MW·m?3)分布如表1所示，堆芯徑向1.85 m共5個(gè)分區(qū)，軸向11 m共22個(gè)分區(qū)[20](R/Z為堆芯徑向/軸向坐標(biāo)位置)。

表1 PBMR-400基準(zhǔn)參考功率分布Table 1 PBMR-400 Benchmark power density.

圖3是堆芯冷卻劑溫場及反射層溫場分布，堆芯冷卻劑溫場最熱的區(qū)域位于堆芯最下方靠近中心反射層的出口位置，最高溫度達(dá)到1210 K。出口處的平均溫度約為1172 K。中心反射層距頂部3 m左右，靠近堆芯球床一側(cè)的溫度高于相同高度上的堆芯流場溫度，從該處開始堆芯功率密度的提高使得球床固相溫度升高很快，堆芯固相球床向中心反射層的熱傳導(dǎo)相對于流相熱傳導(dǎo)比例加大，這一區(qū)域的冷卻劑被中心反射層和球床加熱(圖3、圖4數(shù)據(jù)相對于中心軸線做了對稱處理)。

堆芯固相溫場(UDS0標(biāo)量場)分布距頂部從3 m左右開始升高到1100 K以上(圖4)，堆芯固相最熱的區(qū)域位于堆芯中下部靠近中心反射層的區(qū)域，固相溫場的最高溫度達(dá)到1220 K。

圖3 堆芯冷卻劑及反射層溫場分布Fig.3 Coolant and reflector temperature distribution.

圖5 (a)顯示，堆芯不同軸向位置上冷卻劑溫度的徑向分布較為均勻。球床固相在中心反射層附近區(qū)域出現(xiàn)了較大的溫降(圖5b)，溫度降低的原因在于此處的空隙率較大，流場的流速相對較高(圖6)，熱對流傳導(dǎo)效率較高。外反射層附近區(qū)域溫度下降較快的原因是此區(qū)域孔隙率較高，流速較快，雖然此處的功率密度略有上升，但是由于流體與固相球床對流傳熱性能的提高，固相的溫場仍下降較顯著。

圖4 球床固相溫度分布Fig.4 Solid phase temperature.

堆芯內(nèi)總壓降為2.76×105Pa(圖7)，略小于TINTE及THERMIX的計(jì)算結(jié)果，部分原因在于本文中采取的簡化堆芯模型，堆芯底部的出流區(qū)域引起的壓降沒有考慮在內(nèi)。

圖5 不同軸線位置處冷卻劑溫度分布(a)和固相溫度分布(b)Fig.5 Fluid temperature (a) and solid phase temperature (b) distribution.■ z=10.95m, ● z=7.85m, ▲ z=5.45m, ▼ z=2.45m, ? z=0.75m

圖6 冷卻劑速度分布Fig.6 Fluid velocity distribution.

圖7 堆芯壓降分布Fig.7 Reactor core pressure drop.

3 程序的比較和驗(yàn)證

表2中的結(jié)果表明，本文的計(jì)算結(jié)果與球床式高溫氣冷堆熱工流體程序THERMIX與TINTE的結(jié)果符合較好。堆芯氦氣冷卻劑的平均溫度，固相球床的平均溫度，以及氦氣出口溫度的相對偏差如表2所示，與國際通用的高溫氣冷球床堆熱工程序主要計(jì)算結(jié)果的偏差均在可以接受的范圍之內(nèi)。

表2 Fluent結(jié)果和THERMIX、TINTE[21?23]的比較Table 2 Calculation result comparison between CFD and THERMIX, TITE.

4 結(jié)語

本文基于CFD通用計(jì)算程序Fluent，研究了高溫氣冷堆PBMR-400穩(wěn)態(tài)熱工水力行為。通過Fluent的自定義函數(shù)功能，開發(fā)了局域非熱平衡的多孔介質(zhì)模型，計(jì)算了反應(yīng)堆內(nèi)冷卻劑的溫場分布和固相球床的溫場分布，分析了堆芯流場及不同軸線位置處的溫場分布。通過和國際上應(yīng)用廣泛的球床式高溫氣冷堆熱工水力程序計(jì)算結(jié)果比較，表明基于CFD的高溫氣冷堆熱工計(jì)算結(jié)果與THERMIX與TINTE符合的較好，CFD程序可用于球床式反應(yīng)堆堆芯熱工水力行為的分析。本文分析了PBMR-400穩(wěn)態(tài)熱工水力行為，在此工作的基礎(chǔ)上可以繼續(xù)開發(fā)瞬態(tài)熱工水力計(jì)算程序并和中子學(xué)程序耦合，研究反應(yīng)堆堆芯反應(yīng)性引入以及失壓失冷等事故工況下堆芯的瞬態(tài)熱工水力行為。

1 Petersen K. Zur Sicherheitskonzeption des Hochtemperaturreaktors mit natiirlicher Warmeableitung aus dem Kern im Storfall[R]. FZ Jülich, Germany, 1972

2 Struth S. Thermix-Direkt: Ein Rechenprogramm zur instation?ren zweidimensionalen simulation thermohydraulischer transienten[R]. FZ Jülich, Germany, 1985

3 Cleveland J C, Greene S R. Application of THERMIX-KONVEK code to accident analyses of modular pebble bed high temperature reactors (HTRs)[R]. NUREG/CR-469 ORNL/TM-9905, 1986

4 Gerwin H, Scherer W. The two-dimensional reactor dynamics programme TINTE. Part1: basic principles & methods of solution[R]. Jül-2167: TINTE Part: 1, 1987

5 Gerwin H, Scherer W, Teuchert E. The TINTE modular code system for computational simulation of transient processes in the primary circuit of a pebble-bed high-temperature gas-cooled reactor[J]. Nuclear Science Engineering, 1989,103: 302?312

6 Smith B L. Assessment of CFD codes used in nuclear reactor safety simulations[J]. Nuclear Engineering and Technology, 2010,42(4): 339?364

7 李林森, 王侃, 宋小明. CFD在核能系統(tǒng)分析中應(yīng)用的最新進(jìn)展[J]. 核動(dòng)力工程, 2009,5(增刊): 28?33

LI Linsen, WANG Kan, SONG Xiaoming. International research progress of CFD application in analysis of nuclear power system[J]. Nuclear Power Engineering, 2009,5(Supp): 28?33

8 桂學(xué)文, 蔡琦, 陳玉清. 基于CFD的反應(yīng)堆局部三維流動(dòng)模型與時(shí)空中子動(dòng)力學(xué)模型耦合研究[J]. 核科學(xué)與工程, 2010,3: 216?222

GUI Xuewen, CAI Qi, CHEN Yuqing. Study on coupling of local three-dimension flow model based on CFD method and space-time neutron kinetics model[J]. Nuclear Science and Engineering, 2010,3: 216?222

9 劉余, 張虹, 賈寶山. 核反應(yīng)堆熱工水力多尺度耦合模擬初步研究[J]. 核動(dòng)力工程, 2010, (s1): 11?15

LIU Yu, ZHANG Hong, JIA Baoshan. Preliminary research on nuclear reactor thermal-hydraulic multi-scale coupled simulation[J]. Nuclear Power Engineering, 2010, (s1): 11?15

10 劉余, 張虹, 賈寶山. RELAP5與CFX程序耦合研究[J].原子能科學(xué)技術(shù), 2010, 44(3): 304?308

LIU Yu, ZHANG Hong, JIA Baoshan. Research on coupling between RELAP5 and CFX codes[J]. Atomic Energy Science and Technology, 2010,44(3): 304?308

11 Becker S, Laurien E. Three-dimensional numerical simulation of flow and heat transport in high-temperature nuclear reactors[J]. Nuclear Engineering and Design, 2003,222: 189?201

12 Fluent 6.3 User’s Guide. Fluent Inc.http://aerojet.engr. ucdavis.edu/fluenthelp/html/ug/main_pre.htm, 2006

13 Fluent 6.3 UDF Manual. Fluent Inc.http://aerojet.engr. ucdavis.edu/fluenthelp/html/udf/node2.htm, 2006

14 Mulder E, Teuchert E. Plutonium disposition in the PBMR-400 high temperature gas-cooled reactor. PHYSOR, Chicago, 2004

15 Mulder E. Core basic design report for the PBMR nuclear power plant project[R]. PBMR Report, Doc. No.02979-34

16 KTA 3102.3. Reactor core design of high-temperature gas-cooled reactors Part 3: loss of pressure through friction in pebble bed cores[S].1981

17 KTA 3102.1. Reactor core design for high-temperature gas-cooled reactor part 1: calculation of the material properties of helium[S]. 1978

18 KTA 3102.2. Reactor core design of high-temperature gas-cooled reactors part 2: heat transfer in spherical fuel elements[S]. 1983

19 Zehner P, Schlunder E U. Thermal conductivity of granular materials at moderate temperatures[J]. Chemie Ingenieur Technik, 1970,42: 933?941

20 Gougar H D. The application of the PEBBED code suite to the PBMR-400 coupled code benchmark. Idaho National Laboratory Idaho Falls, Idaho 83415, 2006

21 Boer B, Lathouwers D, Kloosterman J L, et al. Validation of the DALTON-THERMIX code system with transient analyses of the HTR-10 and application to the PBMR. 2010,170(2): 306?321

22 Dudley T, Bouwer W, Piet de Villiers, et al. The thermal-hydraulic model for the pebble bed modular reactor (PBMR) plant operator training simulator system[J]. Nuclear Engineering and Design, 2008,238: 3102?3113

23 Boer B, Kloosterman J L, Lathouwers D, et al. In-core fuel management optimization of pebble-bed reactors[J]. Annals of Nuclear Energy, 2009,36: 1049?1058

CLCTL349

High temperature gas-cooled pebble bed reactor steady state thermal-hydraulics analyses based on CFD method

SONG Shixiong1,2WEI Quan1,2CAI Xiangzhou1,3GUO Wei1,3
1(Shanghai Institute of Applied Physics, Chinese Academy of Sciences, Jiading Campus, Shanghai 201800, China)
2(University of Chinese Academy of Sciences, Beijing 100049, China)
3(Key Laborotary of Nuclear Radiation and Nuclear Energy Technology, Chinese Academy of Sciences, Shanghai 201800, China)

Background:Based on general purpose CFD code Fluent, the PBMR-400 full load nominal condition thermal-hydraulics performance was studied by applying local thermal non-equilibrium porous media model.Purpose:In thermal hydraulics study of the gas cooled pebble bed reactor, the core of the reactor can be treated as macroscopic porous media with strong inner heat source, and the original Fluent code can not handle it properly.Methods:By introducing a UDS in the calculation domain of the reactor core and subjoining a new resistance term, we develop a non-equilibrium porous media model which can give an accurate description of the core of the pebble bed. The mesh of CFD code is finer than that of the traditional pebble bed reactor thermal hydraulics analysis code such as THERMIX and TINTE, thus more information about coolant velocity fields, temperature field and solid phase temperature field can be acquired.Results:The nominal condition calculation results of the CFD code are compared to those of the well-established thermal-hydraulic code THERMIX and TINTE, and show a good consistency.Conclusion:The extended local thermal non-equilibrium model can be used to analyse thermal-hydraulics of high temperature pebble bed type reactor.

CFD (Computational Fluid Dynamics), Gas-cooled pebble bed reactor, Porous media

TL349

10.11889/j.0253-3219.2013.hjs.36.120601

中國科學(xué)院戰(zhàn)略性先導(dǎo)科技專項(xiàng)資助項(xiàng)目(XDA02010200)資助

宋士雄，男，1983年出生，2011年于蘭州大學(xué)獲碩士學(xué)位，現(xiàn)為中國科學(xué)院上海應(yīng)用物理研究所博士研究生，反應(yīng)堆熱工水力

2013-09-05，

2013-11-05