钠冷快堆组件冷却剂沸腾子通道分析方法研究

方闻韬，佟立丽，曹学武

钠冷快堆组件冷却剂沸腾子通道分析方法研究

方闻韬，佟立丽*，曹学武

（上海交通大学机械与动力工程学院，上海 200240）

钠冷快堆发生超设计基准事故时，组件内冷却剂可能沸腾甚至干涸，准确预测其温度分布对钠冷快堆的安全评估具有重要意义。基于均相流模型构建守恒方程，采用Mikityuk对流传热模型以及Cheng-Todreas阻力模型等关系式，开发了适用于钠冷快堆两相流动模拟的子通道分析方法，与FFM-2A单相流动实验数据和KNS-37钠沸腾实验结果进行了对比验证，并与同类子通道分析程序的计算结果作比较，验证了方法的合理性。

钠冷快堆；钠沸腾；子通道分析

钠冷快堆发生超设计基准事故时，堆芯冷却剂可能发生沸腾，钠沸腾导致的反应性上升可能引发严重事故的风险[1]，相关研究指出为延缓甚至防止超设计基准事故发展为严重事故[2,3]，须采取合适的非能动安全措施。由于钠沸腾现象对事故发展进程影响大，因此，准确评估冷却剂沸腾甚至干涸下组件盒内温度分布对钠冷快堆的安全评估具有重要意义。

子通道分析方法具有精细且计算速度快的优点，常被用于模拟快堆组件盒内冷却剂流动传热过程。MATRA-LMR[4]使用分布式阻力模型，通过考虑流动阻力在竖直和水平方向上的分量，模拟了绕丝对冷却剂的影响，ANTEO+[5]内置多种不同的冷却剂物性模型，可用于液态钠、液态铅冷却剂的盒内流动过程模拟，SUBAC[6]耦合了较多近年来开发的各类热工水力模型，并通过模拟分析对相同类型的模型进行了比较评价。但是，这些子通道分析程序仅适用于单相流动传热模拟，未能扩展到快堆事故情况下两相冷却剂流动的模拟。SABENA[7]基于两流体模型，通过对冷却剂通道进行划分可对冷却剂沸腾过程进行一维或二维模拟；石晓波等[8]学者采用子通道分析方法对钠冷快堆完全堵流事故SCARABEE BE+1实验包壳熔化前的冷却剂沸腾进行了数值模拟。钠冷快堆低压条件下液态钠和钠蒸汽之间巨大的密度差使得钠的沸腾和冷凝往往伴随着较剧烈的压力波动，与轻水堆沸腾过程的现象有较大差异，钠沸腾子通道模拟研究仍需进一步完善。

现有子通道分析程序重点关注钠冷快堆稳态单相流动过程，针对瞬态沸腾流动的研究较少。而用于钠冷快堆沸腾模拟的系统级程序难以对通道内的流动传热进行精细的模拟。本文在传统单相子通道分析方法基础上耦合适用于液态金属冷却剂的传热与流动阻力关系式建立了快堆组件盒内两相子通道分析方法，通过模拟结果与实验数据的对比验证了分析方法的适用性，为钠冷快堆组件盒内沸腾现象研究提供理论方法。

1　计算模型

子通道分析方法将冷却剂流通面积划分为多个彼此相邻的区域，并沿流动方向多层划分，形成首尾连续的控制体。通过网格细化减小控制体内的各类热工参数差异，并取参数平均值用于构建守恒方程组。通过对守恒方程以及补充关系式进行迭代得到控制体内各参数的数值解，实现对流动与传热过程的模拟。

1.1　子通道守恒方程

连续性方程：

式中：——轴向流通面积；

——冷却剂密度；

——时间；

——冷却剂轴向质量流量；

——轴向流动长度；

——单位长度横向压差质量流量。

能量守恒方程为：

式中：——冷却剂比焓；

——燃料表面热流密度；

wall——燃料与冷却剂接触长度；

——冷却剂横向换热系数；

——相邻燃料棒的间隙宽度；

Δ——相邻通道冷却剂温度差；

Δ——相邻通道冷却剂比焓差。

轴向动量守恒方程为：

式中：——冷却剂轴向流动速度；

——冷却剂压力；

——重力加速度；

——轴向通道与重力方向的夹角；

——单相摩擦阻力系数；

φ——两相摩擦阻力因子；

h——通道水力直径；

——局部形阻系数；

Δ——相邻通道冷却剂轴向流动速度差。

径向动量守恒方程为：

式中：——相邻控制体质心距离；

G——横流形阻系数。

1.2　对流传热模型

计算燃料与冷却剂之间的传热时需要设定对流传热系数，通常采用努塞尔数表征。

采用Mikityuk模型[9]计算沸腾前液态钠流过冷却剂通道的努塞尔数：

式中：——冷却剂努塞尔数；

——相邻燃料棒栅距与燃料棒直径之比；

——佩克勒特数。

模型的适用范围为相邻燃料棒栅距与燃料棒直径之比介于1到1.95，佩克勒特数介于30到5 000的情况。

采用Dittus-Boelter模型[10]计算干涸下钠蒸汽流过冷却剂通道的努塞尔数：

式中：——雷诺数；

——普朗特数。

水在管道中沸腾流动时会产生复杂的流形变化，而法国GR19实验[11]指出钠沸腾开始后冷却剂流型在数秒内迅速由泡状流过渡到环状流，壁面不与钠蒸汽长时间接触。考虑到在这一过程中对流传热显著强于核态沸腾传热，参考了多孔介质程序THERMIT-6S[12]以及系统级程序TRACE[13]，将干涸前两相冷却剂与壁面传热视为液态钠的单相传热。当燃料热流密度达到临界热流密度时燃料表面将发生干涸，此时由于壁面与钠蒸汽接触，因此需要分别考虑燃料对液态钠和钠蒸汽的传热。1979年Autruffe学者开展的钠沸腾实验指出钠沸腾过程中，在空泡份额达到0.957之前不会出现蒸汽与壁面接触的现象[14]。基于该结论，国际上钠沸腾计算程序通常将空泡份额等于0.957视为发生干涸的临界条件。因此，本方法将发生干涸的判定依据设置为空泡份额达到0.957，此时燃料壁面与液态钠的相对接触面积通过对空泡份额差分得到。

式中：l——液态钠与壁面相对接触面积；

——空泡份额。

通过上述方法，将整个传热过程分为纯液态传热、沸腾传热、干涸后传热和纯蒸汽传热四个过程，通过空泡份额进行区分。

1.3　摩擦阻力模型

冷却剂与棒束、组件盒壁面之间的摩擦阻力使用达西公式进行计算，其中需要补充的参数包括单相摩擦阻力系数以及两相摩擦阻力因子φ。采用Cheng-Todreas模型[15]计算带绕丝组件的中央通道、边通道和角通道的摩擦阻力系数。Cheng学者与Todreas学者基于带绕丝棒束的流动压降实验数据，提出了带绕丝棒束摩擦阻力模型。各通道的带绕丝棒束流动阻力模型计算公式为：

w——通道湿周；

w——绕丝直径；

——绕丝螺距；

w——绕丝倾角。

由于该模型是在光滑棒束摩擦阻力模型的基础上扩展得到，因此模型中使用了光滑棒束摩擦阻力系数常数，该常数可以通过Cheng学者给出的半经验公式获得。d、s、r等参数没有明确物理含义，需要对通道形貌参数进行运算得到。为常数，在湍流情况下取0.18，层流情况下取1.0。过渡流动时的摩擦阻力系数计算分别使用湍流、层流公式进行计算并进行插值。该模型适用于组件盒内燃料棒数目介于7～271，相邻燃料棒栅距与燃料棒直径之比介于1～1.42，绕丝螺距与燃料棒直径之比介于8～42的情况。

在轻水堆两相流动研究中，两相流动摩擦阻力采用相同质量流量冷却剂的单相流动摩擦阻力与两相摩擦阻力因子相乘进行计算，而在液态金属两相流动压降研究时该方法也得到了广泛的使用，例如1974年的Kaiser钠沸腾实验[16]和1984年的Kottowski钠沸腾实验[17]，均基于该方法给出了两相摩擦阻力因子模型。本文采用了Kaiser模型[18]计算两相摩擦阻力因子2：

式中：LM——Lockhart和Martinelli提出的无量纲数，定义为：

式中：——混合物含气量；

l——液态钠密度；

g——钠蒸汽密度；

l——液态钠动力黏度；

g——钠蒸汽动力黏度。

1.4　计算流程

计算流程如图1所示。计算时从通道入口节点开始逐层求解热工水力参数，直到通道出口完成一次迭代。每层节点求解时首先计算燃料棒内导热以及对外传热，然后根据燃料棒传热求解能量守恒方程，根据冷却剂焓值判断流动状态并更新物性参数，随后根据前一节点相邻通道之间压力差求解横向动量守恒方程获得横向流量，再根据冷却剂物性以及横向流量求解连续性方程得到通道内流量，最后求解轴向动量守恒方程确定通道内压降，完成当前节点所有热工水力参数的求解。完成所有节点参数求解后对所有参数进行收敛性判断，决定是否继续迭代，若开始新迭代则继承上一次迭代保留的参数，否则执行输出。

图1 计算流程

2　计算结果验证

为了验证分析方法计算的准确性，选取美国FFM-2A实验[19]以及德国KNS-37实验[20]中的代表性工况进行模拟，并将模拟结果与实验结果进行比较分析。

2.1　FFM-2A实验验证

为了评估单相情况下分析方法的适用性，选取美国橡树岭国家实验室FFM-2A实验的两个工况进行模拟。FFM-2A实验使用由19根电加热棒组成的燃料组件模拟钠在带绕丝燃料棒组件盒内的流动过程。模拟时子通道与加热棒编号如图2所示，表1列出了燃料组件的具体几何参数。

图2 模拟FFM-2A实验所用子通道与加热棒编号

表1　FFM-2A实验组件几何参数

模拟选取的两个实验工况分别为：Test 2 Run 105冷却剂流量0.211 5 m3/s，总加热功率322.3 kW，入口冷却剂温度443.44 ℃以及Test 4 Run 101冷却剂流量0.002 8 m3/s，总加热功率5.17 kW，入口冷却剂温度308.43 ℃。图3～图4为两个实验工况的各通道出口冷却剂相对温升模拟值，同时将模拟结果与实验结果、其他子通道程序的模拟结果进行对比，可以看出模拟结果与实验结果吻合较好，且与其他程序模拟准确度相当。

2.2　KNS-37实验验证

为了评估两相情况下分析方法的适用性，选取德国卡尔斯鲁厄核研究中心开展的KNS实验的一个工况进行模拟。KNS实验使用由37根电加热棒组成的燃料组件模拟反应堆一回路失流情况下组件盒内钠从沸腾直到发生干涸的过程。考虑对称边界条件对1/12的组件盒进行模拟，子通道与加热棒编号如图5所示，表2列出了燃料组件的具体几何参数。

图3 Test 2 Run 105冷却剂相对温升

图4 Test 4 Run 101冷却剂相对温升

图5 模拟KNS-37实验所用子通道与加热棒编号

表2　KNS-37实验组件几何尺寸

模拟选取了编号L22的实验工况，其中加热棒平均加热功率215 W/cm，入口冷却剂温度380 ℃，运行压力1.0 bar，入口冷却剂流速初值3.4 m/s。用于与模拟结果对比的数据来自于位于一号通道且高度为870 mm的T555温度传感器、位于一号通道且高度为779 mm的T658温度传感器、位于一号通道且高度为779 mm的V858空泡份额传感器、高度为545 mm的P713压力传感器以及出口流量传感器。图6、图7为不同高度冷却剂温度变化模拟结果与实验测量结果对比情况，可看出程序对冷却剂从沸腾直到干涸情况的模拟结果较好，从12 s断电后温度模拟结果下降趋势低于实验结果，可能是由于模拟时未考虑组件盒对环境散热，尤其是沸腾后组件盒温度较高，对环境散热的热流密度随组件盒表面的温度上升而显著增强。图8为测点冷却剂空泡份额模拟结果与实验测量结果对比情况，可以看出实验时测点附近空泡份额振荡较为强烈，模拟结果与空泡份额平均值相近。图9、图10为入口压力变化以及出口流量变化模拟结果与实验测量结果对比情况，同时将模拟结果与子通道程序SABENA模拟结果进行了对比。从图9、图10可以看出模拟结果、SABENA模拟结果以及实验结果具有相同的趋势。压力模拟结果在数值上总体低于测量结果，产生偏差的原因可能是分析方法忽略了冷却剂在离开加热段后与通道的摩擦，在摩擦较大的两相情况下这种偏差更加明显。流量模拟结果存在偏差的原因可能是模拟时采用的空泡份额模型导致了通道内发生沸腾的区域在7～9 s时快速扩张，导致出口流量偏高，随后在9～10 s时收缩，导致出口流量提前下降。

图6 T555温度模拟结果与实验结果对比

图7 T658温度模拟结果与实验结果对比

图8 V858空泡份额模拟结果与实验结果对比

图9 P713压力模拟结果与实验结果对比

图10 出口流量模拟结果与实验结果对比

3　结论

本文开发了用于钠冷快堆组件两相热工水力分析的子通道分析方法。通过对FFM-2A实验的模拟研究，验证了分析方法在单相情况下的适用性，并将模拟结果与子通道程序MATRA-LMR、ANTEO+、SUBAC的模拟结果进行了对比，所开发的分析方法的模拟精度与同类程序相当。同时，通过对KNS-37实验的模拟研究，验证了分析方法在两相情况下的适用性，与子通道程序SABENA模拟结果的对比表明其模拟结果与同类子通道程序具有可比性。

A通道轴向面积，m2Pr普朗特数 A′通道实际流通面积，m2p燃料棒栅距，mm C冷却剂传热系数，W/m2*Kq燃料表面热流密度，W/m C’f光滑棒束阻力摩擦系数常数Re雷诺数 Cf液态钠与壁面相对接触面积s相邻燃料棒的间隙宽度，m d燃料棒直径，mmT温度，K DH水力直径，mt时间，s Dw绕丝直径，mmu冷却剂流速，m/s f摩擦阻力系数w单位长度横向压差质量流量，kg/（m·s） g重力加速度，（m/s2）w'单位长度湍流交混流量，kg/（m·s） H绕丝螺距，mmx混合物含气量 h冷却剂比焓，（J/kg）z轴向流动长度，m K局部形阻系数 l相邻控制体质心距离，m希腊符号 m冷却剂轴向质量流量，（kg/s）α空泡份额 Nu努塞尔数θ轴向通道与重力方向的夹角，（°） P冷却剂压力，Paθw绕丝倾角，（°） Pw通道湿周，mμ流体动力黏度，（Pa·s） Pwall燃料与冷却剂接触长度，mρ密度，（kg/m3） Pe佩克莱特数φ2两相摩擦阻力因子

［1］ Bachrata A，Bertrand F，Marie N，et al. A Comparative Study on Severe Accident Phenomena Related to Melt Progression in Sodium Fast Reactors and Pressurized Water Reactors［J］. Journal of Nuclear Engineering and Radiation Science，2021，7（3）：030801-1-030801-9

［2］ Schikorr M，Bubelis E，Carluec B，et al. Assessment of SFR reactor safety issues. Part I：Analysis of the unprotected ULOF，ULOHS and UTOP transients for the SFR（v2b-ST）reactor design and assessment of the efficiency of a passive safety system for prevention of severe accidents［J］. Nuclear Engineering and Design，2015，285：249-262.

［3］ Kruessmann R，Ponomarev A，Pfrang W，et al. Assessment of SFR reactor safety issues：Part II：Analysis results of ULOF transients imposed on a variety of different innovative core designs with SAS-SFR［J］. Nuclear Engineering and Design，2015，285：263-283.

［4］ Kim W S，Kim Y G，Kim Y J. A subchannel analysis code MATRA-LMR for wire wrapped LMR subassembly［J］. Annals of Nuclear Energy，2002，29（3）：303-321.

［5］ Lodi F，Grasso G，Mattioli D，et al. ANTEO+：A subchannel code for thermal-hydraulic analysis of liquid metal cooled systems［J］. Nuclear Engineering and Design，2016，301：128-152.

［6］Sun R L，Zhang D L，Liang Y，et al. Development of a subchannel analysis code for SFR wire-wrapped fuel assemblies［J］. Progress in Nuclear Energy，2018，104：327-341.

［7］ Ninokata H，Okano T. SABENA：Subassembly boiling evolution numerical analysis［J］. Nuclear Engineering and Design，1990，120（2-3）：349-367.

［8］ 石晓波，罗锐，王洲，等. 钠冷快堆单个燃料组件冷却剂沸腾的数值模拟［J］. 核科学与工程，2004，24（4）：312-317.

［9］ Mikityuk K. Heat transfer to liquid metal：review of data and correlations for tube bundles［J］. Nuclear Engineering and Design，2009，239（4）：680-687.

［10］Dittus F W，Boelter L M K. Heat transfer in automobile radiators of the tubular type［J］. International communications in heat and mass transfer，1985，12（1）：3-22.

［11］ Menant B. Preliminary results of sodium boiling through a 19 heating rod bundle［R］. CEA Centre d'Etudes Nucleaires de Grenoble，1975.

［12］ No H C，Kazimi M S. An investigation of the physical foundations of two-fluid representation of sodium boiling in the liquid-metal fast breeder reactor［J］. Nuclear Science and Engineering，1987，97（4）：327-343.

［13］ Chenu A. Single-and two-phase flow modeling for coupled neutronics/thermal-hydraulics transient analysis of Advanced Sodium-Cooled Fast Reactors［R］. École Polytechnique Fédérale de Lausanne，2011.

［14］ Autruffe M I，Wilson G J，Stewart B，et al. Proposed momentum exchange coefficient for two-phase modeling of sodium boiling［C］//International meeting on fast reactor safety technology. 1979.

［15］ Chen S K，Chen Y M，Todreas N E. The upgraded Cheng and Todreas correlation for pressure drop in hexagonal wire-wrapped rod bundles［J］. Nuclear Engineering and Design，2018，335：356-373.

［16］ Kaiser A，Peppler W，Vöröss L. Type of flow，pressure drop，and critical heat flux of a two-phase sodium flow［J］. Nuclear Engineering and Design，1974，30（3）：305-315.

［17］ Kottowski H M，Savatteri C. Fundamentals of liquid metal boiling thermohydraulics［J］. Nuclear engineering and design，1984，82（2/3）：281-304.

［18］ Kaiser A，Peppler W. Sodium boiling experiments in a seven-pin bundle：Flow patterns and two-phase pressure drop［J］. Nuclear Engineering and Design，1977，43（2）：285-293.

［19］ Fontana M H，MacPherson R E，Gnadt P A，et al. Temperature distribution in the duct wall and at the exit of a 19-rod simulated LMFBR fuel assembly（FFM Bundle 2A）［J］. Nuclear Technology，1974，24（2）：176-200.

［20］ Bottoni M，Dorr B，Homann C，et al. Experimental and numerical investigations of sodium boiling experiments in pin bundle geometry［J］. Nuclear technology，1990，89（1）：56-82.

Study on the Sub-channel Analysis Method of Coolant Boiling in the Sodium Cooled Fast Reactor Assembly

FANG Wentao，TONG Lili*，CAO Xuewu

（School of Mechanical Engineering，Shanghai Jiao Tong University，Shanghai 200240，China）

When the beyond-design-basis accident occurs for the sodium cooled fast reactor, the coolant in the assembly may boil or even dry up. Therefore, accurate prediction of coolant temperature distribution is important to the safety assessment of the sodium fast reactor. In this paper, the conservation equation is constructed based on the homogeneous flow model, and adopting the Mikityuk convective heat transfer model and the Cheng-Todreas resistance model, a sub-channel analysis method suitable for the two-phase flow simulation of sodium-cooled fast reactors is developed. The results are compared and verified with the data of FFM-2A steady state experiment and the KNS-37 loss-flow sodium boiling experiment and also compared with the calculation results of similar sub-channel analysis codes, which shows the rationality of the method.

Sodium cooled fast reactor; Sodium boiling; Sub-channel analysis

TL333

A

0258-0918（2023）03-0544-09

2022-08-15

国家自然科学基金资助项目（U1967202）

方闻韬（1997—），男，浙江衢州人，硕士研究生，现从事核反应堆热工水力方面研究

佟立丽，E-mail：lltong@sjtu.edu.cn