1. 引言
为了对轻水堆核电厂的各种热工水力瞬态过程进行分析，已经开发出了许多大型综合性的系统分析程序，如RELAP5、RETRAN、TRAC、TRACE、CATHARE、ATHLET等，这些程序经过多年的持续改进，版本层出不穷，模型日臻完善。其中，RELAP5和RETRAN是比较广泛应用的程序，但两个程序的壁面换热模型采用的是基于实验拟合的换热关系式，这种经验关系式必然对程序计算结果产生一定的不确定性。本文调研了两个程序最新的研究状况，比较了两个程序的壁面换热关系式和临界热流密度计算方法，重点比较了壁面换热关系式选择的依据、考虑的因素和差异性，从而有助于理解两个程序壁面换热模型的差异，为寻找更为准确的壁面换热关系式提供参考；可以为透过换热关系式的不确定性来量化分析计算结果的不确定性提供参考；同时也可为自主化反应堆热工水力系统程序的换热模型选择提供参考价值。
2. 程序研究现状
2.1. RETRAN-3D程序研究现状
RETRAN程序是美国电力研究院(EPRI)主持研制的用于压水堆或沸水堆核电站系统热工水力瞬态分析的系统程序[1] ，主要被用于：1) 核电厂许可证申请；2) 轻水堆的运行过程瞬态分析；3) 小破口失水事故短期瞬态；4) 未能紧急停堆的瞬态；5) 三维功率分布和反应性反馈效应瞬态；6) 沸水堆稳定性分析等。程序的维护以满足美国联邦法10CFR50附录B所要求的质量保证为前提。
RETRAN-3D程序是在RETRAN-02程序的基础上发展而来，它继承了RETRAN-02程序的模拟能力和计算特点，改进了RETRAN-02程序的一些局限性，拓展了程序的使用范围 [2] ，如：反向流特性、传热和压降模型的不连续性、直流蒸汽发生器的稳态初始化和沸水堆两相区域模型、过冷空泡份额计算模型、动态滑移模型、稳压器模型、气泡上升模型、低温低压下水物性参数延伸(压力小于0.6895 MPa，温度小于38℃)、五方程非平衡态调用选项、多维中子动力学调用选项等。该程序以一维均匀、平衡态(HEM)模型(三方程)为基础，同时，在此基础上增加了两种可选择性的、用来提高两相流动模拟精度的模型：一种是通过增加了滑移流(动态方程或代数方程形式)模型选项来模拟两相流动过程中速度的差异(四方程)；另一种是通过增加滑移流模型和蒸汽连续性方程来计算两相速度的不均匀特点和温度的不平衡特性(五方程)，从而在一定程度上(仅限于过冷沸腾和冷水注入两相混合物时汽液之间的温度不平衡差异)突破HEM模型的局限性；此外程序中还增加了不凝气体模型，可以模拟不凝气体对流动和换热的影响。
在RETRAN-3D程序中，堆芯物理(核功率计算)可以选择的模型有 [1] ：三维中子动力学模型、一维中子动力学模型以及点堆模型，热工水力采用一维流场模型，流场方程可以选择三方程、四方程或五方程模型。程序建立有核蒸汽供应系统通用部件模型：泵、阀门、汽轮机、换热器、分离器、安注箱、热结构部件和控制系统，还包括了一些特殊过程模型，比如：截面突变模型、临界流、临界热流密度、硼追踪、非冷凝气体传输等。RETRAN-3D程序的稳态求解和瞬态求解均采用隐式算法，而RETRAN-02则是采用了显式和半隐式混合的求解方法。
RETRAN-3D程序的换热关系式以经典的沸腾曲线作为划分换热机制的依据 [3] ，考虑了温度控制表面的池式沸腾换热机理模型，也考虑了热流密度控制表面的强制对流和自由对流沸腾换热模型，换热关系式中包含了冷凝换热模型。根据两相流动换热现象机理分为不同的换热区段：单相液、过冷沸腾、充分发展的核态沸腾、强制对流蒸发、过渡沸腾、膜态沸腾、单相蒸汽，每个区段可选择多种不同的换热关系式，这种丰富的换热关系式扩展了程序的应用范围，但却增加了用户对换热模式选择的难度。此外，除壁面与流体之间的换热关系式外，五方程模型中还提供了简单的气液界面换热关系式。
2.2. RELAP5程序研究现状
RELAP5程序是进行轻水堆冷却剂系统事故瞬态分析的最佳估算程序 [4] ，由爱德华国家工程实验室(INEL)为美国核管会(NRC)开发，用于规则制定、许可审查计算、事故减缓措施评估、操作员规程评价和实验计划分析，该程序是一个高度通用的热工水力瞬态计算程序，既可计算反应堆冷却系统的瞬态行为，还可用于核电和常规系统的各种不同类型的热工水力瞬态模拟，涉及到蒸汽、水、不冷凝气体和溶解物等混合行为。该程序也被用作核电厂分析器的基础以及进行热工水力实验设计的参考，并被世界很多国家的核能企业和高等院校的核专业领域所广泛验证和扩展应用。最新版本为：RELAP5-3D、SCDAP/ RELAP5/MOD4.0/RELAP-7。程序应用涵盖了整个轻水堆系统的瞬态分析，如：大、小破口的失水事故、未能紧急停堆的预计瞬态以及操作瞬态如：给水丧失、失去场外电源、全厂断电、汽轮机跳闸等。
RELAP5程序采用一维、瞬态、两流体六方程模型，堆芯物理采用点堆模型，程序包含许多通用的组件模型，用来模拟各种系统。组件模型包括：泵、阀门、管道、热构件、点堆动力学、电加热器、喷射泵、汽轮机、分离器、安注箱和控制系统组件。另外，还包括了一些特殊过程模型，比如：截面突变模型、横向交叉流动模型、临界流(壅塞流)、临界热负荷、硼追踪、非冷凝气体传输等。
RELAP5程序的换热模型源自沸腾和冷凝曲线，通过比较壁面温度和当地流体饱和温度、主流温度的差异来选择换热区段(单相对流区段、沸腾区段、过渡区段、液膜蒸干区段、蒸汽冷凝区段)。除壁面换热关系式外，程序中还提供气液界面换热关系式，根据6种流型(波状流、塞状流、弹状流、泡状流、环状流、滴状流)来选择相对应的换热关系式，此外，还包括辐射换热模型选项。
2.3. 两个程序计算模型与方法比较
RETRAN-3D程序热工水力模型基于两相均匀、平衡态三方程模型，同时包含有考虑两相不均匀特性的四方程模型，以及两相不均匀、不平衡特性的五方程模型。三方程模型即为均匀、平衡态模型，汽液两相按照混合物特性进行计算；四方程模型是在三方程模型的基础上增加汽液之间的滑移速度方程，用来反映汽液两相之间的速度不均匀差异；五方程则是在四方程的基础上增加了气相的质量方程，来考虑两相温度的不平衡差异。
RELAP5程序的六方程模型是针对汽液两相分别建立各自独立的质量、动量和能量守恒方程，同时考虑汽相与液相之间的质量、动量和能量交换，即充分考虑了相之间的速度不均匀特性和热力不平衡特性。与RETRAN-3D程序相比，RELAP5程序对两相流的热工水力特性描述更为精细，是目前最为广泛应用的两相流模型。
在计算方法上，RETRAN-3D程序采用全隐式求解方法，RELAP5程序可以选择近隐式或隐式计算方法，两种方法在单相水流动压降与换热计算中，热工水力参数的差别几乎可以忽略不计，但在汽液两相流动压降与换热计算过程中，RELAP5程序刻画两相过程比较详细，程序的应用范围较RETRAN-3D程序更为广泛，计算结果相对更为准确。
3. 壁面换热模型比较
3.1. RETRAN-3D壁面换热关系式选择依据
RETRAN-3D程序的壁面换热关系式依据流体受热过程，其中主要涉及到池式沸腾和流动沸腾。换热关系式的选择源自沸腾换热分区。
池式沸腾曲线如图1所示[3] ，横坐标是表面温度，纵坐标是表面热负荷，采用温度控制表面的沸腾曲线表现特征为：随着表面温度的增加，热流密度表现趋势为ABDEF。采用热流密度控制表面的沸腾曲线表现特征为：随着表面热流密度的增加，对于强制对流和自由对流沸腾过程，表面温度表现趋势为ABDD’F，D’表示临界热负荷qCHF值下的表面新平衡状态。
加热或冷却过程中所历经的沸腾过程不一定会遵循相同的路径曲线，图中实线表示加热过程，虚线表示冷却过程。对于自由对流过程可能有两条趋势：FEDBA和FEE’BA；对于非均匀加热强制对流换热，随着热流密度的增加，表面温度表现趋势是FD’DBA，这时候池式沸腾曲线已经不再适用。在稳态情况下，反应堆燃料棒与冷却剂之间属于热流密度控制表面(忽略空泡反馈影响)，并带有非均匀加热特征，实际试验表现趋势可能为FD’EE’。在过渡沸腾区域，还可能出现FECE’趋势，过渡过程中的最大热负荷(C点)可能低于初始过程的最大热负荷(D点)。
反应堆瞬态过程中，燃料与冷却剂之间会经历强制对流沸腾、池式沸腾或自由对流换热情形，存在热流密度控制表面沸腾情况，也可能存在温度控制表面沸腾情况，因此需要综合全部考虑。图2给出了棒束通道内低热负荷和高热负荷下垂直流动与换热机理分区图，值得注意的是：在低热流密度下，临界热流密度(CHF)最有可能出现在环状流区；在高热流密度下，CHF有可能在过冷沸腾区内发生。
在图2棒束入口单相液对流区，流道中的流体温度低于饱和温度，对应图1沸腾曲线的AB段，该区域的换热关系式相对比较成熟，需要注意的是在低流量情况下，重力或自由对流的影响比较显著，在选择换热关系式时要考虑这种影响。
在图2的过冷沸腾区，主流温度仍低于饱和温度，但贴近壁面附近的流体温度已经达到并超过主流的饱和温度，局部区域已经产生沸腾，有小气泡产生，但是很多的工程经验关系式中并没有考虑这种微观影响，只是确定出图1中的沸腾起始点(ONB)B点，之后沸腾沿着图1中的BCD趋势来进行计算。在核态沸腾区，主流温度等于饱和温度。在高热流密度下，CHF可能提前出现在过冷沸腾区域。图1中的D点是蒸干点(低热流密度)或烧毁点(高热流密度)，此处发生CHF，CHF是一个复杂的、多变量函数，与压力、质量流速、含气量、几何结构等有关，程序中采用多个经验关系式及插值法来包络各种压力范围。
图1中的过渡沸腾(不稳定的膜态沸腾)区DE或EC可能在瞬态过程的某一部分(流体经过非常高温度的壁面)发生，试验中会观察到该区段内有不稳定的蒸汽膜附着在加热表面，汽膜破裂后邻近的液体会快速润湿该加热表面，表面温度随时间变化很慢，但也存在表面温度升高的情形。
Figure 1. Curve: boiling curve representative at a fixed pressure
图1. 在特定压力下的沸腾曲线
Figure 2. Curve: flow and heat transfer regimes in rod array with vertical flow
图2. 垂直棒束间向上流动与换热分区
在CHF的下游区域，由于包壳表面已经被蒸汽所覆盖，热量交换主要是依靠蒸汽的对流换热，也存在与弥散的液滴进行热量交换，该区段对应图1沸腾曲线中的ED’F区，程序中采用的是多个经验关系式。最终，所有的液体都将会被蒸干，蒸汽就被加热成过热蒸汽，过热蒸汽区采用的是合适的单相汽换热关系式。
另外，蒸汽可能会被冷凝，液滴就会产生并积聚在固体表面，冷凝、再润湿和积聚可能发生在低于饱和温度的表面上，液滴会因为重力的影响而下落，公式中需要考虑。
RETRAN-3D程序基于水的沸腾过程换热机制不同而分区段选择合适的换热关系式，这些换热关系式源自宏观实验结果，使用程序时应与计算流体力学(CFD)程序建模方法区别对待。
3.2. RELAP5壁面换热关系式选择依据
RELAP5程序的壁面换热关系式的选择依据是沸腾和冷凝曲线，见图3所示[4] ，程序将壁面与流体之间的换热分为三个区段：冷凝区域、对流区域和沸腾区域，每个区域根据换热现象的不同而选择不同的换热关系式，共计考虑了12种换热现象，对应12种换热模式，同时将不同几何结构(28种类型)的影响加以考虑，程序考虑到的几何构造类型见表1所示[4] 。沸腾区域内所选择的公式与RETRAN-3D所选择的公式有些是一致的。
3.3. 壁面换热关系式比较
表2列出了RETRAN-3D程序在不同换热区段所采用的壁面换热关系式以及程序选择模式编号(不考虑不凝结气体)，表3列出了RETRAN-3D程序所用CHF关系式[3] 。
表4列出了RELAP5程序针对不同换热现象所采用的换热公式包以及程序选择模式编号，表5列出了不同几何类型、不同流动状态和沸腾类型下换热公式的具体选择[4] 。
从壁面换热关系式的比较可以看出，两种程序均是根据不同的换热区段类型(单相液、沸腾、单相汽、冷凝、超临界)而选择与之相对应的壁面换热经验关系式，换热关系式涵盖的范围均比较广，从单相水到饱和水，再从饱和水到过热蒸汽，直至超临界水，还包括冷凝换热和临界热流密度CHF计算。在某些换热区段两个程序所选用的换热关系式是一样的，如：单相段的Dittus-Boelter、沸腾段的Chen公式。
Figure 3. Curve: RELAP5 boiling and condensing curves
图3. RELAP5沸腾和冷凝曲线
Table 1. Hydraulic geometries of the RELAP5 code
表1. RELAP5程序所考虑的水力几何特征
Table 2. Heat transfer correlations of the RETRAN-3D code
表2. RETRAN-3D程序换热关系式包
Table 3. CHF correlations of the RETRAN-3D code
表3. RETRAN-3D程序所用CHF关系式
注：RELAP5程序的CHF计算采用表格形式并可以充分考虑8个乘子的影响，这与RETRAN-3D程序选用固定关系式方法不同。
Table 4. Heat transfer correlations of the RELAP5 code
表4. RELAP5程序换热关系式包
注：模式编码范围：0~71，可在现有模型基础上继续增加新模型。
Table 5. Available wall heat transfer model of RELAP5
表5. RELAP5程序用到的壁面换热模型
两个程序换热关系式的不同点如下：
1) 根据换热现象的不同，RETRAN-3D程序将换热区划分为9个类型：单相液、过冷核态沸腾、充分发展的核态沸腾、强制对流蒸发、过渡沸腾、稳定的膜态沸腾、单相汽、超临界水、冷凝。RELAP5程序将换热区段划分为13个类型：不凝气体–蒸汽–水、单相液或超临界、单相液或过冷壁面(空泡份额 < 0.1)、过冷核态沸腾、饱和核态沸腾、过冷过渡沸腾、饱和过渡沸腾、过冷膜态沸腾、饱和膜态沸腾、单相汽或超临界两相、膜状冷凝、蒸汽中的冷凝、核态沸腾，同时，将轻水堆核蒸汽供应系统流动几何特征划分为28种类型，比较这28种流动几何构造类型换热差异后，再来挑选合适的换热公式。
2) RETRAN-3D程序提供两种换热程序包，强制对流换热程序包和联合换热包。程序默认调用的换热程序包为强制对流换热包，当用户指定选用联合换热包时，程序才会按照联合换热包提供的换热公式来计算壁面与流体之间的换热，联合换热包中包含了强制对流和自然对流换热公式。由于某些特定几何结构的换热公式具有特定的适用性，因此程序提供了多个可供选择使用的公式，以满足不同几何结构下特殊的换热公式，如：在冷凝换热区，可以选择的公式就有5个类型：管型、两区非平衡态、水平管内侧2个、垂直表面。在联合选项的单相液换热区段，大流量采用Dittus-Boelter公式，小流量采用Collier公式。RELAP5程序是根据壁面与流体之间的几何构造特征和流动状态来选择壁面与流体之间的换热关系式。
3) 在CHF计算方法上，RETRAN-3D程序根据压力区间来选择不同的CHF经验公式，在公式之间的断点处采用插值方法使之连续，所选用的CHF计算公式为：B&W-2公式、Barnett公式和修正的Barnett公式；此外，程序还增加了附加的、备选CHF关系式，用于压水堆(PWR)的CHF关系式：B&W-2关系式、W-3关系式、MacBeth关系式、Bowring关系式，用于沸水堆(BWR)的CHF关系式：Barnett关系式、Janssen-Levy关系式、MacBeth关系式、Bowring关系式，这些公式依据棒束组件实验数据拟合而来，有各自的适用性条件，可见，RETRAN-3D可选择的CHF关系式是很丰富的，用户需要仔细甄别各个公式的适用性，比较公式差异后选择最佳适用的公式。RELAP5程序采用查找内嵌的CHF表来获得特定运行状态下的CHF值(除水平放置的棒束采用了Folkin公式来计算外)，通过对15,000个实验数据点的三维表格处理创建出4410个点，并进行插值延伸处理，从而大大扩展了CHF计算范围，CHF表的压力范围为0.1~20.0 MPa，质量流量范围为0.0~7500.0 kg/m2∙s，含汽率范围为−0.5~1.0，除此之外，程序还考虑8个乘子的影响，用户可以根据具体组件的实验情况来施加8个乘子来修正CHF值，这8个乘子是：水力因子、棒束因子、栅格因子、加热长度因子、轴向功率因子、水平因子、垂直流动因子、压力超限因子。
4) 在壁面换热关系式的适用性和应用范围方面，两种程序的换热包以解决轻水堆的冷却剂系统流动与换热问题为基本出发点，相比之下，RELAP5程序更多考虑了流道不同几何结构对流动与换热的影响，将水力控制体分为管型、环形、分支以及特殊部件模型(主要有：泵、阀、分离器、喷射泵、汽轮机、安注分支、安注箱)，RETRAN-3D程序简单将水力控制体统一按照控制体来处理，并考虑特殊部件模型(主要有：泵、阀、分离器、汽轮机、安注箱)。在结构材料的热传导方面，两种程序均是采用一维导热模型，RETRAN-3D程序可以进行矩形和圆柱形两种类型的热导体的一维瞬态导热计算，RELAP5程序可以进行矩形、圆柱形和球形三种类型的热导体一维瞬态导热计算。用户在使用程序时需要注意这些模型的选择，对不能模拟的特殊结构需要进行简化处理。
5) 除了壁面换热关系式外，RELAP5还更多考虑了气相与液相界面之间的换热计算，本文仅仅列出了壁面与流体之间的换热公式。表格中各个换热公式的具体适用性范围和误差将在后续研究中阐述。
6) 在水物性参数计算方面，RETRAN-3D采用了拟合公式法，即对1967 ASME水蒸气物性表进行公式拟合，出发点在于：简单、平滑、准确。其中，饱和水的焓值采用压力变量的函数进行计算，饱和蒸汽的焓值亦采用压力变量的函数进行计算，压力上限的最大值取临界压力；液体区或蒸汽区流体的温度采用压力与焓值两个变量的函数进行计算，压力上限为41.36 MPa，焓值上限约为4000 kJ/kg；液体区或蒸汽区流体的定压比容采用压力和焓值两个变量的函数进行计算，压力上限为41.36 MPa，焓值上限约为4000 kJ/kg；热导率采用温度与密度两个变量的函数进行计算；动力粘度采用温度与密度两个变量的函数进行计算；表面张力采用温度变量的函数进行计算。
7) RELAP5水蒸气表采用表格查询法[5] ，水物性包括：轻水(h2o)、重水(d2o)、1984轻水(h2on)和1995轻水(h2o95)，水物性范围包括：单相水、饱和水、饱和汽、单相汽、超临界水。其中，基本轻水物性(h2o)所调用的表格数据基于1967 ASME水蒸气表，即采用国际公式化委员会IFC提出的工业用水和水蒸气计算模型IFC-67公式；基本重水(d2o)物性所用表格数据基于加拿大原子能公司AECL和美国国家航空航天局WASP项目研究数据；1984轻水(h2on)物性基于1984美国国家政府标准局(NBS)和加拿大国家研究委员会(NRC)所提供的水蒸气表；1995轻水(h2o95)基于国际水和水蒸气性质学会(IAPWS)提出的工业用水和水蒸气热力学性质计算公式，即公认的著名公式IAPWS-95，IAPWS在1997年发布了比IAPWS-95计算效率更快的IAPWS-97公式，但比较而言，IAPWS-95公式在计算精度上优于IAPWS-97公式。因此，RELAP5程序可选择不同研究组织和机构的水表。
4. 结束语
本文对比了核电行业内广泛使用的两种最佳估算程序RETRAN-3D程序和RELAP5程序的壁面换热模型关系式，通过调研分析了各自挑选壁面换热关系式的依据和所考虑的因素，对比了两种程序的壁面换热关系式所模拟的范围，分析表明：
1) 两个程序的壁面换热模型涵盖了水被加热和冷凝的全过程，即：单相水、饱和水、过热蒸汽、超临界水、冷凝。各个换热现象区段所用换热模型的适用性和热工误差界限需要进一步详细研究。
2) 两个程序中水蒸汽表的压力范围均已涵盖超临界水物性，不同的是：RETRAN-3D程序只能选择一种标准的轻水物性参数，RELAP5程序可以选择四种不同标准水蒸气物性参数(h2o、d2o、h2on、h2o95)。
3) 两个程序计算堆芯CHF的选择不同，RETRAN-3D程序采用拟合的经验公式来计算CHF，内含多个可供选择调用的CHF经验公式，RELAP5程序主要采用表格法(LOOK UP TABLE)来获得CHF值，两个方法之间的差异和各自误差的上下限需要进一步深入比较研究。
4) 在使用程序建模之前，必须首先进行程序适用性评价，明确程序的使用限制条件和能力，在程序的有效使用范围内进行建模和分析。
5) 本文的研究工作对于理解最佳估算程序壁面换热模型的选择原则和公式特点以及开发新的壁面换热关系式具有实际意义，同时也有助于深入理解两个程序在壁面换热模型方面的差异。
基金项目
受国家863核安全研究主题项目《压水堆核电厂安全级冷却链改进研究》课题经费支持。