熔盐堆热管式非能动余热排出系统建模及程序开发

刘明皓，张玉龙，王成龙，张大林，秋穗正，李 毅，尹莎莎，刘 航，杨红发

(1.中国核动力研究设计院 核反应堆系统设计技术重点实验室，四川 成都 610041； 2.西安交通大学 能源与动力工程学院，陕西 西安 710049)

熔盐堆作为第4代核反应堆中唯一一种液态燃料堆型，已得到全世界范围内的广泛关注[1-3]。福岛核事故后，人们对核反应堆系统提出了更高的非能动安全要求。为提升熔盐堆非能动安全性，同时使系统满足小型化、模块化的设计需求，研究者结合高温热管技术提出了熔盐堆热管式非能动余热排出系统(HP-PRHRS)概念设计，相关理论和实验研究已先行展开，初步验证了高温热管用于熔盐堆系统的可行性[4-8]。HP-PRHRS主要由卸料罐、排热烟囱和高温热管组成。当反应堆发生一回路大破口、失流等事故时，压力容器中燃料盐温度迅速上升使冷冻阀熔断，冷冻阀开启，燃料盐依靠重力作用快速下泄到卸料罐中，插在卸料罐中的高温热管迅速启动导出其热量。HP-PRHRS运行由3个自然循环(对流)耦合完成：1) 卸料罐内燃料盐与高温热管间的自然对流传热；2) 高温热管内部工质的气液两相自然循环；3) 排热烟囱内空气的自然循环。HP-PRHRS的运行不需额外的外界驱动力，仅依靠系统自身的自然循环(对流)即可实现燃料盐有效冷却，使系统具备良好的非能动安全特性。同时系统结构得到大幅简化，有效避免了中间环节故障，对未来熔盐堆模块化、小型化设计非常有利。

本文基于HP-PRHRS结构和运行特点，建立一套较为完整的数学物理模型，耦合熔盐堆堆芯物理热工模型、高温热管模型和非能动余热排出系统模型等，开发熔盐堆HP-PRHRS分析程序PRAC，并采用MSRE基准题和瞬态实验数据对程序进行验证。

1 主要数学物理模型

1.1 堆芯物理热工模型

堆芯物理热工模型采用点堆模型与燃料盐流动换热模型耦合，将堆芯功率分解为时间和空间的函数，采用反应堆中子动力学模型计算反应堆瞬态功率，考虑了燃料多普勒效应、熔盐密度等反应性反馈和燃料盐流动带来的缓发中子先驱核迁移机制影响。方程组如下：

(1)

(2)

(3)

燃料盐在石墨通道内流动，热量由熔盐带出堆芯，同时燃料盐同石墨间存在传热，进而得出下列传热方程：

(4)

(5)

Pout(t)=w(t)cpf(Tout-Tin)

(6)

其中：Mf、Mg分别为燃料盐和石墨的质量，kg；cpf、cpg分别为燃料盐和石墨的比定压热容，J·K-1·kg-1；Pf、Pg分别为核释热分别进入燃料盐和石墨的功率，W；Pout为燃料盐流动传出堆芯的功率，W；h为石墨通道与燃料盐之间的对流换热系数，W·K-1·m-2；A为石墨通道与燃料盐的对流面积，m2；Tg、Tf分别为石墨和燃料盐的平均温度，K；Tout、Tin分别为堆芯出口和进口温度，K。

燃料盐温度和石墨温度对反应性的温度反馈可按下式计算：

ρr,T(t)=αf(Tf-Tf0)+αg(Tg-Tg0)

(7)

其中：ρr，T(t)为温度反应性反馈；αf为燃料盐温度反馈系数，K-1；αg为石墨温度反馈系数，K-1；Tg0、Tf0分别为初始时刻石墨和燃料盐的平均温度，K。

1.2 高温热管模型

图1 高温热管原理图Fig.1 Schematic of high temperature heat pipe

热管是一较为复杂的系统，原理图如图1所示。热管依靠由吸液芯毛细力驱动的气液两相自然循环完成热量传递，其工作过程包含了多种不同的传热形式和连续的相变过程，因此国际上尚未有准确获得热管内部瞬态特性分析解的理论模型。而实际应用关注的是高温热管的整体宏观性能，主要是热管传热过程中的功率和温度等关键参数。基于上述理念，本文采用热管热阻网络模型，进行高温热管的求解分析，大幅提升数值稳定性和计算效率。基于热管传热流程，可将热管的传热过程以热阻网络的形式表示。经合理简化后，最终的高温热管热阻网络模型如图2所示。图2中：热阻1和2分别为蒸发段管壁径向热阻和吸液芯径向热阻；热阻3和4分别为蒸发段和冷凝段的液环热阻；热阻5和6分别为冷凝段管壁径向热阻和吸液芯径向热阻；热阻7和8分别为绝热段管壁轴向热阻和吸液芯轴向热阻；Q为输入功率，W；T∞为环境温度，K。

图2 热管热阻网络模型Fig.2 Network model of high temperature heat pipe

对于网络中的单个热阻控制体，有如下控制方程：

(8)

其中：δi为热阻导热厚度，m；αi为热扩散系数，m2·s-1；Ti，1和Ti，2为热阻两侧边界温度，K；Ti为热阻中心温度，K。

进而联立各热阻单元，可获得如下热管热阻网络模型控制方程，通过求解可获得热管各区温度，进而计算出热流密度、温度变化率、热管效率等参数。

(9)

εn(n+1)-2)Tn+ε(n+1)nTn+1)

n=2,3,4,5

(10)

(11)

(12)

(13)

(14)

(15)

(16)

其中：ki为热导率，W·m-1·K-1；h∞，c为冷凝段与环境换热系数，W·K-1·m-2；Sc为冷却换热面积，m2；T∞，c为环境冷却温度，K。

热管的正常工作依赖其内部工质的两相循环，因此需补充工质循环判别模型[9]，Ψ>1表明热管建立内部工质循环。

(17)

其中：Ψ为循环判别无量纲数；μv为蒸汽动力黏度，Pa·s；hfg为工质汽化潜热，J·kg-1；ρv为蒸汽工质密度，kg·m-3；Rv为蒸汽腔当量半径，m；Le、La、Lc分别为热管蒸发段、绝热段、冷凝段长度，m；δ为吸液芯厚度，m；D为热管直径，m；THP为热管运行温度，K；keff为吸液芯有效热导率，W·m-1·K，考虑吸液芯材料热导率kw和吸液芯内工质热导率kl，keff按照下式计算：

(18)

其中，φ为吸液芯孔隙率。

然而在实际运行中，热管的传热能力受其传热极限的限制。当对热管的输入功率超过热管传热极限时，热管内部的循环便会遭到破坏，导致热管失效。传热极限与热管尺寸、工质、吸液芯结构、工作温度等因素有关[10-11]。因此在热管的计算分析中，需结合热管的传热极限(主要是高温热管，尤其是碱金属热管，一般为毛细极限、声速极限和携带极限等)，对热管传热功率进行判别校对，以确定热管处于正常工作状态。

1.3 非能动余热排出系统模型

非能动余热排出系统模型分为卸料罐内燃料盐流动换热和卸料罐外排热烟囱内空气流动换热两部分，二者通过高温热管模型实现参数交互。

1) 卸料罐内燃料盐流动换热

卸料罐内燃料盐控制方程为：

(19)

其中：m为燃料盐质量，kg；Tsalt为燃料盐温度，K；τ为系统运行时间，s；cp为燃料盐比定压热容，J·K-1·kg-1；Qd为衰变热功率，W；Qc为系统冷却功率，W。

对于燃料盐的衰变热功率，若其衰变曲线未知，可采用Shure修正预测模型进行燃料盐衰变热功率计算[12]。

(20)

其中：Q0为停堆前反应堆的运行功率，W；Qd为燃料盐排出时间τ后的衰变热功率，W；τ0为停堆前反应堆已连续运行时间，s；a和b为时间相关系数。

卸料罐内燃料盐自然循环。罐内热管蒸发段与燃料盐间的自然对流传热可采用下列模型计算。

对于单根独立的热管[6]：

(21)

其中：Hhp为热管高度，m；Hmax为最大液位高度，m。

对于热管竖列管束，热管间距P与管径D比值为4.33≤P/D≤8.67时[7]：

Nu=3.710 9Ra0.116 06

(22)

为使高温热管处于良好的工作状态，热管倾斜插入，罐外热管冷凝段高于罐内蒸发段。对于倾斜管，浮升力(重力)可分解为平行和垂直于表面的两方向作用力，因此圆管的长度和管径对于自然对流都有影响[13]。因此对于倾斜圆管的自然对流传热，其特征长度需采用倾斜角进行修正：

(23)

其中：Lc为倾斜管特征长度，m；L为圆管长度，m；d为圆管直径，m；θ为圆管倾斜角。当倾角θ=0°时即为水平管，此时Lc=d；当θ=90°时即为竖直管，此时Lc=L。

2) 卸料罐外排热烟囱内空气流动换热

空气由排热烟囱底部流入，经过卸料罐区域对热管冷凝段进行冷却，然后从顶部出口流出。由于空气的流速远低于声速，且空气压缩性影响很小，因此采用一维不可压缩模型进行空气的流动和换热计算。空气流动换热整体控制方程如下。

连续性方程：

wi=wair

(24)

动量守恒方程：

(25)

(26)

其中：wair、wi分别为空气流量和第i个空气流动换热控制体内流量，kg·s-1；Δpi为排热烟囱内第i个空气流动换热控制体内的压降，Pa；Δpstack为排热烟囱内总压降，Pa；Δp0为空气从排热烟囱出口到入口高度所对应的环境压降，Pa；fi为传热区摩擦阻力系数；Kj为局部阻力系数；ρi为第i个控制体内空气密度；Ai为第i个控制体换热面积；g为重力加速度；Δzi为高度差，m；De为当量直径，m。

能量守恒方程：

(27)

其中：Vi为排热烟囱内第i个空气流动换热控制体的容积，m3；Si为控制体i的换热面积，m2；qi为控制体i换热面上的热流密度，W·m-2；hi-1、hi分别为控制体i的进、出口焓，J·kg-1。

空气与热管冷凝段间通过对流传递热量。当空气流动较弱、流速较低时，可用Corcione模型进行空气外掠竖列管束的自然对流换热计算[14]。

Nu=Ra0.235(0.292ln((P/d)0.4·

N-0.2)+0.447)

(28)

其中：5×102≤Ra≤5×105，P/d≤10-lgRa；N为竖列管束中的管数量。

此时，空气流动可视为外掠光滑圆管，阻力计算可采用如下关系式：

(29)

其中：Δp为管外空气流动的压力损失，Pa；Gmax为最大质量流速，kg·m-2·s-1；ξ为考虑管束排列方式的修正系数；gc为重力加速度，m·s-2。

2 PRAC程序开发

2.1 程序结构

HP-PRHRS分析程序PRAC采用标准FORTRAN 90程序设计语言编写。基于模块化设计思路，根据功能进行模块划分，使每个模块间具有较高的独立性。PRAC程序的模块结构如图3所示。其中堆芯物理热工模块由中子动力学模块和流体动力学模块耦合；余热排出系统模块由高温热管模块、卸料罐模块和排热烟囱模块3部分耦合；辅助模块、物性模块和数值计算模块为公用模块，可根据需要受到其他模块调用。

2.2 程序流程

PRAC程序流程图如图4所示。程序逻辑与HP-PRHRS实际运行流程一致，主要包括堆芯计算和非能动余热排出系统计算两大部分。首先依据输入参数开展熔盐堆堆芯的热工物理耦合计算，获得温度、功率、流量和反应性反馈等堆芯参数。在计算过程中程序进行判断，如果冷冻阀开启，燃料盐进入余热排出系统，余热排出系统随即开始运作，然后开始非能动余热排出系统计算，获得相关系统参数。相对地，如果无需冷冻阀开启，在堆芯模块计算至设定时间时程序便会终止，余热排出系统计算模块不会启动。

2.3 数值算法

PRAC程序在运行求解中会同时涉及点堆中子动力学模型和系统热工水力模型，结合其余辅助模型，构成一套封闭的方程组。由于不同模型间时间量级差距较大，构成的方程组具有很大的刚性。针对这一问题，程序采用吉尔(Gear)算法进行方程组求解，确保数值计算的精确稳定。对于热管热阻网络模型建立的一阶常微分方程组，则采用高精度单步算法龙格库塔法(R-K法)，选取区间上若干点的斜率进行加权平均，通过基于泰勒级数展开的待定系数法实现微分方程组的步进求解。

图3 PRAC程序模块结构Fig.3 Module structure of PRAC code

3 程序及模型验证

3.1 堆芯物理热工模块验证

本文采用MSRE的启泵基准题和停泵基准题数据来验证其堆芯物理热工模块的可靠性和准确性[15]。MSRE堆芯主要参数列于表1。MSRE在启泵、停泵变流量的过程中，通过调节控制棒插入的深度来保持堆芯临界，并记录控制棒的位置，通过计算控制棒移动所引入的反应性就可得出熔盐流量变化所导致的堆芯反应性的变化。

MSRE启泵和停泵基准题的计算结果如图5所示。由图5a可看出：启泵时反应性降低，通过控制棒提升来补偿损失的反应性，10 s后流量达到稳定，控制棒停止运动，控制棒补偿的反应性与熔盐流动造成的反应性损失平衡。由图5b可看出：熔盐流动造成的反应性损失与实验结果吻合较好；停泵时熔盐流量减少，导致堆芯反应性增加，控制棒插入堆芯来减少多余反应性，20 s后回路中的流量降为0，且最终控制棒减少的反应性与流量减少所增加的反应性相平衡。两个基准题的反应性计算值均与实验数据符合较好，证明了程序与模型的准确性和可靠性。

3.2 余热排出系统模块验证

余热排出系统模块采用前期研究获得的瞬态实验数据进行对比验证[8]。瞬态实验在熔盐堆HP-PRHRS实验系统上进行，如图6所示。瞬态实验系统和高温热管主要参数列于表2。卸料罐总高为0.92 m，内径为0.6 m，侧面安装6根高温钾热管。卸料罐放置于排热烟囱内，烟囱高度在0～4 m内可调。

图4 PRAC程序流程图Fig.4 Flow chart of PRAC code

表1 MSRE基准题堆芯参数Table 1 Core parameter in MSRE benchmark

图5 MSRE启泵(a)和停泵(b)基准题计算结果Fig.5 Benchmark calculation result of MSRE pump strat-up (a) and pump stop (b)

本文针对系统瞬态实验中的不同工况进行了计算，将程序计算结果与实验数据进行对比，如图7所示。对比实验值与PRAC计算值，在不同实验工况下PRAC计算值与实验值均吻合良好。PRAC燃料盐温度计算值与实验值相对偏差小于1.2%，热管壁温计算相对偏差小于2.8%。PRAC计算结果能准确反映出系统冷却过程中氟盐温度的变化情况和热管壁温的变化情况，证明了程序与模型的准确性和可靠性。

图6 熔盐堆HP-PRHRS实验系统Fig.6 HP-PRHRS experimental system for MSR

表2 瞬态实验系统和高温热管的主要参数Table 2 Main parameter of transient experimental system and high temperature heat pipe

图7 瞬态实验燃料盐温度和热管壁温对比Fig.7 Comparison salt temperature and heat pipe temperature for transient experiment

4 结论

本文根据熔盐堆堆芯和HP-PRHRS的结构及运行特点建立了一套合理完整的系统模型，主要由堆芯物理热工模型、高温热管模型和非能动余热排出系统模型3部分耦合。采用模块化编程设计，开发了熔盐堆HP-PRHRS分析程序PRAC。采用MSRE基准题和瞬态实验数据进行对比验证，PRAC计算结果与基准值吻合良好，且与燃料盐温度实验值的相对偏差小于1.2%，与热管壁温实验值的相对偏差小于2.8%，验证了模型与程序的合理性与准确性。本文模型和程序能为后续开展熔盐堆HP-PRHRS的深入设计提供模型和软件基础。