基于RELAP5的氦氙流动换热计算模块开发与验证

周 彪，孙 倩，孙 俊，孙玉良

(清华大学 核能与新能源技术研究院，先进反应堆工程与安全教育部重点实验室，北京 100084)

高效可靠的空间电源是空间任务顺利开展的重要保障。随着未来军民空间任务的不断拓展，以化学能、太阳能为主的常规能源将难以满足空间任务需求[1]。空间核反应堆电源具备工作寿命长、功率覆盖范围广、受空间环境影响小等特点，是未来开展多类别空间任务的理想电源形式[2]。自20世纪50年代，美国和俄罗斯就启动了空间核反应堆电源的相关研究工作，提出了多种技术路线。当系统功率达到百千瓦乃至兆瓦量级时，采用气冷堆结合布雷顿循环技术方案的优势更加明显，主要在于其能量转换效率高、技术成熟度高、系统比质量更小。为满足空间布雷顿循环系统紧凑性的设计要求，降低叶轮机械气动载荷，同时获得相对更佳的换热性能，通常选择氦氙混合气体(He-Xe)作为反应堆冷却剂和布雷顿循环工质[3-5]。

反应堆热工系统分析程序是开展反应堆设计与安全评价的重要工具。针对不同类型的堆型，国内外研究机构开发了多种热工系统分析程序，典型轻水堆系统分析程序有RELAP、CATHARE、TRACE等[6]；液态金属反应堆系统分析程序有SAS4A/SASSY-1、SIMMER、RELAP5-3D/ATHENA等[7]。尽管有THERMIX、TINTE等气冷堆热工系统分析程序[8-9]，但这些程序只能用于球床式反应堆的热工计算，应用范围有限。针对氦氙气冷反应堆，李杨柳等[10]开发了用于环形与圆形流道的反应堆单通道分析程序，但无法满足系统层面的热工分析。美国爱达荷国家实验室开发了RELAP5-3D/ATHENA程序[11]，可拓展用于超临界CO2、He-Xe等工质的高温气冷堆热工系统分析，但该程序对国内不开放。为更好地对氦氙气冷空间核反应堆系统开展热工分析与安全评价，开发相应热工系统分析程序十分必要。

考虑到RELAP5程序的物理模型和计算方法较为成熟，且广泛应用于核动力系统安全评审，本文选择对RELAP5/MOD4.0程序进行拓展，基于程序源码添加He-Xe物性计算模块与适用的传热关系式，使程序初步具备计算He-Xe流动换热的功能。

1 RELAP5/MOD4.0概述

RELAP5是由美国爱达荷国家实验室与美国核管会协同开发的一维瞬态热工分析程序，RELAP5/MOD4.0为当前最新的可用版本。程序采用自上而下的编程逻辑，如图1所示，顶层主要由INPUT、TRNCTL、STRIP 3个功能区组成，各功能区再分别由相应子程序执行特定功能。其中INPUT控制输入卡数据读取与计算初始化；TRNCTL控制水力学部件与热构件的步进计算与结果保存；STRIP提供数据提取、绘图等功能。本文基于RELAP5/MOD4.0开发He-Xe流动换热计算模块，编写He-Xe密度、动力黏度、热导率、比热容4个物性计算子程序，在RNEWP与HYDRO路径文件中添加物性子程序的调用接口，实现对He-Xe物性的计算。在HTADV路径下的dittus.for文件中新增适用的传热关系式，以满足He-Xe对流换热的计算需求。

图1 RELAP5/MOD4.0程序框架

2 程序功能拓展

2.1 He-Xe物性计算模块

RELAP5/MOD4.0程序中不凝结气体密度ρ与比定压热容cp采用理想气体状态方程进行求解。El-Genk等[5]研究表明，该方法求解稠密He-Xe物性时将引入较大误差。不凝结气体动力黏度μ与热导率λ分别采用式(1)、(2)进行计算，这些公式均是基于Sutherlands定律演变而来[12]。Sutherlands定律是基于理想气体分子动力理论与分子间势能理论提出的，在空气黏度等物性计算方面具有较好的适用性。

(1)

λ=λrTb

(2)

式中：S为Sutherlands温度，程序中S=114.0 K；μr、λr分别为参考温度Tr(Tr=273.13 K)下混合气体的动力黏度与热导率；b为与气体类型相关的常数；T为气体实际温度。用βr代表混合气体μr、λr、b，三者的数值均可由通式(3)计算得到。式(3)中下标1表示氙气，2表示氦气，x1为氙气摩尔分数。纯气体相关常数在RELAP5用户手册[13]中有详细说明，其数值列于表1。

表1 程序中计算物性的常数

βr=x1βr1+(1-x1)βr2

(3)

表2 混合稀有气体物性的半经验公式

基于Tournier等的物性公式，本文利用MATLAB开发了He-Xe物性计算程序PROHX，并在之前的研究[15]中验证了程序计算的准确性。为探究式(1)、(2)对He-Xe物性计算的适用性，将其计算值与PROHX程序计算结果进行对比，如图2所示，发现式(1)、(2)计算值与PROHX计算值差异明显，表明Sutherlands定律不适用于He-Xe物性计算。因此，需在程序中添加新的物性模块。本研究在源码中新增了rhohx.for(密度)、viscshx.for(黏度)、thcnhx.for(热导率)、cpphx.for(比热容)4个物性子程序文件以及hexeconstant.for库文件(定义物性计算时涉及的变量与相关函数)，构成完整的氦氙物性计算模块。上述物性均为温度、压力、混合比相关的函数。通过输入卡给定相应工况参数后，在程序初始化和步进计算过程中调用相应物性子程序，可实现对氦氙混合气体物性的计算。

图2 氦氙物性计算结果对比

2.2 He-Xe传热关系式

RELAP5程序中强迫对流换热计算关系式采用Dittus-Bolter(DB)公式。为探究DB公式对He-Xe流动换热计算的适用性，采用DB公式计算了不同Pr(空间堆推荐使用的He-Xe对应Pr为0.16

图3 DB公式计算值与实验值对比

图3中实验数据参考Taylor等[4]的圆管加热实验，实验装置示意图如图4所示，其中实验测试段为直径D=5.87 mm的圆管，由绝热段(l1=56D)和均匀热流加热段(l2=60D)组成，该实验探究了包括He-Xe在内的多种混合气体的流动换热过程。

图4 实验装置示意图

DB公式为常物性公式，而气体流经加热流道时，截面径向位置处的物性差异明显，此时气体变物性对流动换热的影响不可忽略。如表3所列，本研究首先选取4组变物性传热公式[16-18]进行适用性分析。式中：Nub为变物性Nu；Nu0为常物性Nu；Tw为壁面温度；Tb为He-Xe平均温度；n为变物性修正因子；ζ为中间变量。为覆盖目前空间堆设计方案中使用的工质范围，本文选择4个不同摩尔质量混合气体的实验组[4]，如表4所列，工况范围为：0.2130时，Petukhov公式与实验值符合较好，原因是其未考虑热入口段的影响(即x/D的影响)。综合考虑，本研究采用Taylor公式对程序原有传热关系式进行修正。

图5 现有Nu关系式计算值与实验值的对比

表3 现有的Nu关系式

表4 氦氙混合气体的实验参数

3 程序功能验证与结果分析

为验证拓展后的RELAP5程序对He-Xe流动换热计算的准确性，采用拓展后的程序对上述实验[4]的测试段进行建模。为探究程序拓展后的改进效果，引入程序默认的DB公式作对比。计算的边界条件与实验工况一致，将不同工况下的程序计算值与实验结果进行对比，结果如图6所示，发现程序中传热关系式的选择对壁面温度Tw的计算值有重要影响。采用DB公式计算得到的加热段Tw显著低于实验结果；采用Taylor公式得到的计算值与实验值符合较好；当Tw/Tb较大时，在接近热充分发展段处的Tw计算值略高于实验值，此时计算结果具有一定的保守性。

图6 RELAP5程序Tw计算值与实验值的对比

此外，计算了He-Xe平均温度Tb，如图7所示。由图7可看出，采用DB公式与Taylor公式的计算结果基本重合，且各工况下Tb计算值与实验值符合较好。这是因为RELAP5求解流体平均温度是直接通过能量守恒法来实现的，只要边界条件确定，各控制体的平均温度则确定，与传热关系式的选择无关。在流动阻力方面，计算了两组工况下第2个压力测点处的压力以及从起始加热点至第2个压力测点之间的压降，结果列于表5。由表5可看出，两个工况的压降计算相对误差均在10%以内。通过以上计算，验证了拓展后的RELAP5程序在计算He-Xe流动换热功能上的准确性。

图7 RELAP程序Tb计算值与实验值的对比

表5 拓展后的RELAP5程序计算He-Xe流动压降的误差

4 结论

为开发适用于氦氙气冷空间堆的热工系统分析程序，本文对RELAP5/MOD4.0程序开展了适用性分析，并对源码进行了二次开发。研究表明：程序默认的Sutherlands定律用于He-Xe物性计算时将引入较大误差；Dittus-Bolter公式对He-Xe对流换热时的Nu预测偏高，壁温计算结果不保守。在RELAP5源码实现了He-Xe物性计算模块和传热关系式的添加，对程序功能进行了拓展。通过与实验值的对比，发现拓展后的程序对He-Xe流动换热模拟结果与实验值吻合较好，验证了程序对He-Xe流动换热计算的准确性。在后续研究中，将开发布雷顿循环系统中相关部件的计算模块，基于具体空间堆设计方案，搭建系统计算模型，从而实现对系统层面的瞬态及稳态计算分析。