高温气冷堆热工分析模型改进与堆芯温度场分析

孙世妍，张佑杰，郑艳华，夏 冰

(1.清华大学 核能与新能源技术研究院，先进核能技术协同创新中心，先进反应堆工程与安全教育部重点实验室，北京 100084; 2.中核战略规划研究总院，北京 100048)

球床式高温气冷堆具有良好的固有安全性和经济性，在电力生产和工艺热利用等方面具有广阔的应用前景[1]。超高温运行时，高温气冷堆的冷却剂出口温度能达到1 000 ℃以上，更高的热利用效率和发电效率得以实现，经济竞争力显著提升[2]。与此同时，在超高温环境下，反应堆的安全性更需时刻得到保证，除要开展更耐高温的燃料元件和反应堆结构材料研究外，如何更准确地模拟堆芯温度分布这一关键参数也是超高温运行研究中关注的首要问题之一。

在反应堆的设计和安全分析过程中，需建立热工分析模型，对正常运行和事故工况下的堆内温度做出预测，以选择合适的运行参数和评估反应堆的安全性。在球床式高温气冷堆安全分析程序THERMIX下建立的10 MW高温气冷实验堆(HTR-10)的热工分析模型，对HTR-10的设计和安全分析起到了重要作用[3]。试验验证和与各国模拟结果的横向比较表明，该模型能较准确地模拟反应堆在额定设计工况下满功率稳态运行以及多种安全特性试验和验证试验工况下的热工水力行为，为运行方案设计、安全评审和其他科学研究提供了有力的依据[4]。

由于堆芯球流的不确定性和堆芯散体结构的特点[5]，堆内冷却剂流动路径具有一定的不确定性。研究表明，由原有的HTR-10热工分析模型得出的堆芯外围结构材料温度监测点处的温度与试验测量值之间仍存在一定差异[4]。为更为精确地分析堆芯温度场，可从更精细的堆内气相对流行为模拟方面对该模型做出合理的改进，使其模拟堆芯温度分布的精确度进一步提升。

本文简要分析HTR-10的主要结构和传热情况，在更精细梳理堆内冷却剂流动路径的基础上，改进原有的热工分析模型，使堆芯冷却剂流动路径在满足总体热平衡的同时，更合理地模拟球床外围结构的局部传热过程，得到更准确的堆芯温度分布。

1 堆芯结构及传热概况

HTR-10是具有第4代反应堆特征的球床式高温气冷实验堆，采用内含多层包覆燃料颗粒的球形燃料元件，利用石墨作为慢化剂、氦气作为冷却剂。堆芯壳内布置了大量的石墨块和含硼碳砖，石墨反射层起到慢化和反射中子的作用，碳砖因其导热系数较小且含有热中子吸收材料，能有效减少堆芯的热损失，并减小热中子对堆芯壳和反应堆压力容器等金属结构的损伤。此外，石墨反射层的内部轮廓形成了球形元件的流动通道，石墨砖块内开有孔道，形成了冷却剂的流动通道和控制棒等结构的运动通道[6]。

堆内的绝大部分热量是在燃料球内产生的，其中大部分经堆内冷却剂对流带出堆芯，其余通过导热和辐射经堆芯外围结构材料传递至外界环境。因此，堆芯的温度分布是堆芯元件和多种结构材料发热、导热、辐射、对流的综合体现，需经反应堆物理、传热和冷却剂流动耦合计算得出。在堆芯热工分析模型中，将各结构准确地体现出来并合理地区分固相导热和气相对流区域，有利于更准确地模拟传热和流动过程，从而得到更精确的堆芯温度场。

2 堆芯冷却剂流动情况分析

球床堆堆内固相物质分布相对稳定，而冷却剂分布和流动情况则较为复杂。HTR-10的堆芯布置和冷却剂在反应堆压力容器(RPV)内主要的流动路径如图1所示[3,6]。冷却剂进入RPV后，经过RPV和堆芯壳之间的环形下降空间，汇集在RPV底部的球形腔室内。少量氦气流入堆底卸料管，冷却其中的燃料元件。此外的大部分氦气在堆芯壳底部进入碳砖和反射层中的20个冷氦上升孔道，向上流动至顶反射层的冷氦联箱中。其中的少量氦气流入控制棒孔道对控制棒进行有效冷却，其余的大部分冷却剂进入堆芯球床，将燃料产生的大部分热量带出堆芯。上述几部分冷却剂最终汇集在堆芯球床下方的热氦联箱，充分搅混后进入热气导管。此外，由于石墨块和碳砖块间存在大量窄缝，如图2所示，部分冷却剂会流入窄缝形成堆芯漏流。漏流的存在削弱了冷却剂对堆芯燃料的冷却作用[7]，并使得反射层附近的温度分布发生变化。

图1 HTR-10一回路系统Fig.1 Primary system of HTR-10

图2 HTR-10堆芯横截面Fig.2 Cross section of HTR-10 reactor core

所以最终汇入HTR-10热氦联箱的冷却剂包括4部分，流动路径分别为：1) 从顶反射层的冷氦联箱向下流经堆芯球床；2) 从RPV底部空腔向上流经卸料管；3) 流经石墨反射层内控制棒孔道；4) 流入石墨砖缝的堆芯冷却剂漏流。除漏流外，其他部分均具有确定的流动范围和方向。而由于砖缝数量众多，在反应堆长时间运行后，砖缝的大小可能有所变化，因此，漏流流量、位置和方向都难以准确估计，这对模拟冷却剂漏流造成了困难。

漏流分布与堆内构件缝隙分布情况关系密切。石墨块和碳砖块被叠放在反应堆堆芯壳内，构成石墨反射层和碳砖隔热层。在反射层高度范围内，周向上由内而外布置有20个石墨砖和20个含硼碳砖(图2)。在顶反射层中，砖块之间的20个纵向窄缝直接与冷氦联箱连通(图1)，冷却剂能较容易地通过这些窄缝从堆芯球床上方的冷氦联箱流入堆芯球床下方的热氦联箱。而堆芯下方，在冷却剂的流动路径上，由于堆芯壳等钢结构的存在，在一定程度上对砖缝起到封闭作用，砖缝未暴露于空腔或联箱中，氦气漏入窄缝中相对较为困难。

在原有的热工分析模型中，将冷却剂在反射层内的漏流等效地设置在压力容器底部冷氦空腔和堆芯球床下方热氦联箱之间[3-4]，如图3所示。该模型体现了漏流冷却剂不参与堆芯球床的换热，满足堆芯发热和冷却剂载出热量之间的总体平衡，得出的堆芯温度分布总体上较为准确，但在模拟漏流冷却剂在石墨反射层中的传热行为方面还可更为精确。另外，等效地设置在堆底的漏流流道会使当地的温度分布模拟产生一定偏差，可能进而影响堆芯球床的温度分布，需进行相应的改进。

图3 原有HTR-10热工水力模型中的冷却剂流动路径简图Fig.3 Brief flow process in original thermal hydraulic model for HTR-10

3 改进的堆芯热工分析模型

3.1 物理模型

反应堆温度分布问题是传热、流动、核裂变和中子输运等多个物理过程共同作用的结果，需建立不同层次、不同物理过程之间的耦合模型来完整描述反应堆的物理热工行为。第1个层次的耦合为反应堆热工水力过程内部的耦合，包括固相导热、气相流动及流固之间的耦合，以对流换热量作为耦合的桥梁；第2个层次的耦合是热工水力过程与中子输运过程之间的耦合，以核燃料或结构材料发热量作为耦合的桥梁。在物理模型中将这两部分热量内热源化，即将存在对流换热的区域的对流换热功率和具有内热源的区域的发热功率体积平均化，即可实现耦合[8]。

球床式高温气冷堆具有柱对称的堆体结构，考虑到固相的热惯性，堆内、外的复杂导热过程，可在柱坐标系下，采用二维瞬态导热方程(式(1))描述[9]。

(1)

(2)

(3)

燃料元件具有球对称结构，燃料球内的导热过程，可在球坐标系下，用一维瞬态导热方程(式(4))描述[9]，球心为坐标原点。

(4)

其中：T为燃料球温度，K；λfb为燃料球导热系数，W/(m·K)。

对于气相，忽略其容积惯性，冷却剂流场和温度场采用二维类稳态质量守恒方程(式(5))、动量守恒方程(式(6)、(7))和能量守恒方程(式(8))描述[9]。

(5)

(6)

(7)

(8)

其中：Gr、Gz分别为气相质量流量G的径向、轴向分量，kg/(s·m2)；p为气相压力，Pa；Ψ为球床摩擦阻力系数；d为燃料球直径，m；ε为球床孔隙率；ρf为气相密度，kg/m3；g为重力加速度，m/s2；λeff,r、λeff,z分别为球床的径向、轴向等效导热系数，W/(m·K)；cp为气相的比定压热容，J/(kg·K)。

3.2 计算模型

THERMIX是由德国于利希研究中心开发的球床式高温气冷堆安全分析程序[9]，清华大学核能与新能源技术研究院在原有程序的基础上做了完善和改进[10]，用于模拟反应堆在正常运行和事故条件下的热工水力行为。该程序中的模型和参数计算公式大部分在德国进行了试验论证，并实际用于AVR、THTR-300和HTR-Module等反应堆的设计，其分析结果具有较高的可信度和应用价值[4]。

利用THERMIX程序建立的模型包括反应堆的固相导热计算模型、压力容器内的气相流动模型、燃料球的一维导热模型以及一回路流网模型等，运用THERMIX程序根据反应堆的结构合理地构建堆芯几何模型、划分网格并给定边界条件，程序可结合堆芯物理分析软件VSOP给出的堆芯功率密度分布，对上述方程进行离散和迭代求解，得出相应工况下堆芯的固相温度场、燃料温度场，以及气相压力场、流场和温度场。模型中各成分区按照各自的材料属性求解相应的控制方程，当几何模型或材料属性设置改变时，控制方程形式不变，求解方程数量和物性参数按照修改情况发生相应变化。

本文采用THERMIX作为HTR-10的热工分析求解工具。HTR-10堆芯结构复杂，在模型建立过程中需对各结构进行适当简化，建模思路如下：1) 将堆芯球床简化为均匀多孔介质，该区域产生了堆芯的绝大部分热量；2) 将球床中的导热和辐射相耦合的传热问题简化为等效导热问题，其导热系数由燃料球和石墨球导热、球外壁面间的辐射、孔隙中的气体导热等多种效应叠加而成，由实验经验关系式综合给出；3) 从堆芯球床向外划分控制棒孔道和冷氦孔道两个主要流道，此外，按照冷却剂绕流路径设置若干竖直和水平流道；4) 将外界环境和堆腔表面冷却器作为传热边界；5) 将压力容器上下两端的半球形腔室内的结构和空腔按体积折算为圆柱体；6) 将除冷却剂流通路径外的其他气体空腔(反应堆压力容器内顶部氦气空腔及压力容器外所有空腔)均作为不流动区，仅做导热和辐射计算。

在以往的研究中，根据HTR-10本体及热气导管的结构特点、几何参数和系统布置，结合安全分析程序THERMIX的具体要求，在(r，z)柱坐标系下，建立了HTR-10的热工分析模型，坐标原点取在堆芯球床自由堆积面等效平面和堆芯中心线的交点处[4]。

1——堆芯球床；2——底反射层上部流道；3——底反射层下部流道；4——热氦联箱；5——堆芯球床上部空腔；6——不流动区；7——压力容器底部空腔；8——堆底环形冷却剂流道；9——堆芯壳底冷却剂绕流流道；10——冷氦流道；11——卸料管入口节流件；12——控制棒孔道出口节流件；13——控制棒孔道；14——冷氦联箱；15——底反射层内小联箱；16——压力容器内氦气入口；17——压力容器内环形流道；18——顶反射层流道；19——反射层内的冷却剂漏流流道；20——堆芯球床出口节流件图4 改进的HTR-10热工分析模型(气相流动部分)Fig.4 Gas convection model in improved thermal hydraulic model for HTR-10

本文对原有模型的冷却剂流动路径进行了改进，根据堆芯漏流冷却剂的分布特点，并结合堆芯热工计算模型中不同的冷却剂漏流模拟方法对堆芯温度分布的影响研究结果[11]，将漏流流道设置在石墨反射层内，流道在周向上位于堆芯球床和控制棒孔道之间，连接顶反射层内的冷氦联箱和底反射层内的热氦联箱。在改进模型中，固相导热模型模拟了从堆芯至反应堆混凝土舱室外的流体边界的固相导热和辐射，模型在径向被划分为35个网格，轴向被划分为61个网格，共2 135个栅元，分为42个不同的物质区，包括HTR-10的球床堆芯、石墨反射层、含硼碳砖、堆内金属构件、堆芯壳、反应堆压力容器、冷却剂流道、卸料管、余热排出系统的表面冷却器、保温层、混凝土层等结构；气相对流模型(图4)模拟了冷却剂在堆内的流动换热，反应堆压力容器内部的彩色区域为计算堆芯气相流动的区域，模型在径向被划分为18个网格，轴向被划分为39个网格，共702个栅元，分为20个不同的物质区，包括堆芯球床、侧反射层内的冷氦流道、控制棒孔道、石墨顶反射层和底反射层内的冷却剂流道、冷氦联箱、热氦联箱等，采用有气体流动的球床区、垂直导流管和气体沿所有方向流动的空腔区模型来模拟这些部件，不流动区则是气相流动计算模型的边界。

3.3 参数和边界条件

在热工分析中，氦气的热物性参数、燃料球的导热系数和表面换热系数、堆芯球床的摩擦压降等均根据德国安全导则KTA3102.1～3[12-14]中的公式计算。用于模拟堆芯球床的多孔介质等效导热系数可通过ZSB-R公式[15-16]求得，求解过程中需用到在放射性环境下石墨导热系数随温度变化的经验关系式。反射层石墨和碳砖的热物性参数由德国提供的经验关系式计算[9-10]。

模型中的冷却剂流量分配列于表1，流量分配方案的制定遵循HTR-10的热工设计准则[3]，以满足反应堆的运行安全要求，流量占比是指各部分冷却剂流量在冷却剂入口总流量中所占的份额。

表1 HTR-10热工分析模型中的冷却剂流量分配

本文运用该模型分别模拟10 MW额定设计工况和3 MW测温试验工况，其运行参数列于表2，并将其作为计算的参数和边界条件。设置对流模型中的冷却剂压力，将区域16设置为速度入口，按照表中所列数值分别设置其流量和温度，并设置导热模型中的外界环境温度。

表2 计算中的主要参数和边界条件Table 2 Main parameter and boundary condition in calculation

该热工水力模型中的燃料发热由反应堆物理分析程序VSOP算得的功率密度分布给出。功率密度分布由反应堆当前的燃料分布和中子通量分布确定，燃料分布情况与反应堆的换料方案关系密切，HTR-10采用燃料球5次通过堆芯的换料方案。上述两工况，反应堆分别处于初装堆芯和过渡堆芯状态，其功率密度分布如图5所示，总热功率分别为10.00 MW和3.20 MW。

图5 模拟工况的堆芯功率密度分布Fig.5 Power density profile of simulated condition

4 改进效果及堆芯温度场分析

4.1 模型改进效果

在HTR-10运行过程中，热工测量系统会实时监测各种运行工况下堆内构件及压力壳的温度分布情况[4]。为使模型更准确地反映堆芯温度分布，本文选取测量系统中与堆芯球床距离较近、位于侧反射层和侧碳砖内的12个热电偶测量数据作为模型改进效果的验证目标值，它们在堆内的分布情况如图6所示。

图6 堆芯侧向温度测点布置Fig.6 Layout of temperature measurement points around core

这12个测点分两层成对布置于堆芯球床圆柱区域外围的半高度处和底部，测点与轴线的距离两两相等。在模型坐标系中，T21～T26 6个测点的布置高度为z=78.3 cm，T31～T36 6个测点的布置高度为z=168.3 cm。测点处的堆内构件布置情况列于表3。

表3 温度测点位置Table 3 Location of temperature measurement point

运用原模型和改进后热工分析模型分别模拟10 MW额定设计工况和3 MW测温试验工况，得到反应堆侧向温度测点的温度对比结果，如图7所示。图中散点为堆芯侧向温度测点的温度实测值，曲线代表改进前、后的热工分析模型在测点高度的温度计算值。

由图7可看出，原模型对堆芯侧向构件的温度模拟存在一定偏差，偏差随着与堆芯距离的增大而减小。这是由于实际运行中堆芯漏流冷却剂对侧反射层石墨和碳砖等构件具有冷却作用，使得实际温度较模拟温度更低；与堆芯的距离越近，漏流效应越严重，偏差就越明显。体现在物理模型中即侧反射层中温度较高的固相向温度较低的气相传热，温度降低，最终达到热平衡。

图7 堆芯侧向构件温度对比Fig.7 Comparison of temperature of components around core

两个算例中测点温度的对比列于表4、5。改进前，10 MW额定设计工况下，温度模拟的最大偏差超过100 ℃，相对偏差约为30%，靠近堆芯球床的6个测点的平均相对偏差约为20%；3 MW测温试验工况下，温度模拟的最大偏差约为80 ℃，相对偏差约为20%，靠近堆芯球床的6个测点的平均相对偏差约为16%。可见，改进前的堆芯热工模型对堆芯球床外围构件的温度模拟偏差较大，这可能在一定程度上对堆芯球床温度的准确模拟产生影响。

表4 模型改进效果(10 MW额定设计工况)Table 4 Model improvement effectiveness (10 MW rated design condition)

改进后，测点的温度偏差整体上大幅缩小，对靠近堆芯球床的8个测点的模拟尤为准确。10 MW额定设计工况下，反射层内测点计算温度的最大相对偏差约为10%；3 MW测温试验工况下，最大相对偏差约为5%。可见，改进后的堆芯热工模型模拟堆芯冷却剂漏流对反射层的冷却效果更精确，模型对外围构件温度分布模拟准确度显著提高。

表5 模型改进效果(3 MW测温试验工况)Table 5 Model improvement effectiveness (3 MW temperature measurement experiment condition)

4.2 堆芯温度场分析

运用改进后的热工分析模型模拟HTR-10在初装堆芯状态下按照额定设计工况满功率稳态运行堆芯的温度分布情况，其中堆内固体温度场如图8所示。

图8 HTR-10额定设计工况满功率稳态运行堆内温度场Fig.8 Temperature distribution in HTR-10 core at full power and steady state under rated design condition

该工况下，燃料元件表面最高温度为911.08 ℃，燃料元件中心最高温度为941.80 ℃，该点位于堆芯轴线上，高度z=188.80 cm，底反射层最高温度为807.56 ℃，侧反射层最高温度为477.00 ℃。计算结果与CPR-5各国计算的HTR-10额定设计工况下满功率稳态运行时堆芯燃料元件及堆内构件的温度范围相符，说明改进后的热工分析模型对堆芯球床温度分布模拟的准确性良好。

HTR-10的燃料颗粒包覆材料能在不超过1 600 ℃的高温环境下保持其完整性，计算得到的燃料及反射层最高温度均未超过相应材料的温度限值，说明HTR-10能保证额定设计工况下的运行安全性，充分的温度裕量正是开展超高温运行研究所需的必要条件。

5 结论

本文分析了HTR-10冷却剂流动情况对堆芯温度场可能产生的影响，改进了原热工分析模型中的堆芯冷却剂漏流流道，并将计算结果与试验数据进行对比，而后简要分析了额定设计工况下满功率稳态运行时的堆芯温度场，得到如下主要结论。

1) 运用反应堆运行过程中堆内构件测点的实测值对模型进行了校验，改进后的热工计算模型对堆芯外围温度模拟的准确度有显著提升，在额定设计工况下的相对偏差由最大约30%缩小至最大约10%。

2) HTR-10(初装堆芯)额定设计工况满功率稳态运行，一回路氦气流量为4.32 kg/s，堆芯入口氦气温度为250 ℃时，出口氦气温度为700 ℃，燃料中心最高温度为941.80 ℃，反射层最高温度为807.56 ℃，均远低于材料的温度限值，安全裕度充足，因此HTR-10的冷却剂出口温度能进一步提高，具有实现超高温运行的潜力。

3) 改进后的HTR-10热工分析模型，对堆芯球床和堆芯外围构件的温度场模拟更加准确，可作为后续开展HTR-10最高温度及其不确定性研究以及超高温运行预设计的工具。