氟盐冷却球床高温堆燃料管理程序的开发与分析

2021-02-03贾国斌伍建辉陈金根顾国祥蔡翔舟

原子能科学技术 2021年2期

贾国斌，戴叶，伍建辉，陈金根，顾国祥，蔡翔舟,*

(1.中国科学院上海应用物理研究所，上海 201800；2.中国科学院先进核能创新研究院，上海 201800；3.中国科学院大学，北京 100049)

熔盐堆是第4代核能系统的6种堆型之一[1]，主要特征为使用高温熔盐作为冷却剂。按照运行时燃料物质形态，可分为固态燃料熔盐堆与液态燃料熔盐堆。固态燃料熔盐堆按照燃料组件的几何形状分为棱柱状反应堆[2]与球床堆[3-4]，由于球床堆可实现在线连续换料，因此更易实现深燃耗，从而提高燃料利用率与反应堆的经济性。钍基熔盐堆(TMSR)核能系统[5]中设计的堆型之一为氟盐冷却球床高温堆(PB-FHR)[6-7]。

目前研究PB-FHR连续换料的典型程序包括COBBLE[8]，其通过SRAC程序[9]完成燃料球双重非均匀性的求解，使用中子扩散程序CITATION与燃耗程序ORIGEN2[10]实现平衡态搜索的中子通量求解与燃耗核素求解。耦合程序PBRE[11]基于蒙特卡罗程序MCNP[12]与点燃耗分析程序ORIGEN2，也是平衡态搜索程序。以上程序未考虑从启堆到堆芯平衡态的过程，并不能为PB-FHR制定完整的燃料管理策略。

由于PB-FHR与高温气冷球床堆(HTGR)在堆芯结构与换料模式上相似，因此可参考其燃料管理程序的求解方法。目前用于求解HTGR的中子学程序包括VSOP[13-14]，其使用共振自屏因子的方法处理双重非均匀性，使用组件输运-堆芯扩散反复迭代的方法处理泄漏效应，在处理燃料管理时，燃料球每个燃耗步沿轴向移动到下一网格内，并通过追踪燃料球在堆芯内的通过次数，将达到最大通过次数的燃料球卸出，实现将最深燃耗燃料球卸出堆芯的目的。She等[15]使用两步法计算流程，使用ENDF/B7数据库代替原VSOP以ENDF/B4为主的微观截面数据库，使用细群数据库(其中包括温度、燃耗变量)，在线进行泄漏反馈计算。美国爱达荷国家实验室(INL)开发了两维R-z的HTGR换料程序PEBBED[16]，其将堆芯内燃耗的移动作为燃料球在堆芯内的移动。Liem[17]开发出HTGR两维R-z换料分析程序BATAN-MPASS，并针对HTGR进行了单次与多次通过的计算分析。以上程序主要针对HTGR进行开发。在HTGR中，冷却剂的吸收和散射截面很小，均匀化时可不用考虑，而PB-FHR的冷却剂熔盐必须考虑其吸收与散射作用。此外，HTGR与PB-FHR的燃料球运动方向不同，因此需专门开发适用于PB-FHR的燃料管理程序。

本文开发燃料管理分析程序PBMSR，其主要基于中子输运程序DRAGON5[18-19]进行均匀化，并使用扩散程序DONJON5[20]进行全堆计算。为验证PBMSR在换料过程中计算的正确性，使用换料模式相似的HTGR进行验证。

1 方法介绍

1.1 中子学分析程序需解决的问题

中子学分析程序求解PB-FHR的燃料管理需解决以下3个问题。

1) 双重非均匀性

燃料球具有双重非均匀性，具体表现在包含裂变核素的燃料颗粒triso随机分布在燃料球的石墨中，对中子的吸收导致triso颗粒内中子能谱与燃料球的石墨基底不同，此为第一重非均匀性；燃料区(包覆颗粒+基体)与非燃料区(包括石墨球壳区、冷却剂熔盐)吸收截面与裂变截面的不同引起中子能谱差异，带来第二重非均匀性。第一重非均匀性决定不能将triso内的裂变核素按照体积均匀分布在石墨基底中进行均匀化，第二重非均匀性决定不能将燃料区内的核素按照体积均匀分布在非燃料区中进行均匀化。

2) 泄漏效应与活性区底部熔盐的均匀化

在石墨慢化的球床堆中，典型的燃料球直径在3～6 cm之间，中子自由程与组件尺寸(燃料球尺寸)基本相当或更大，更容易受到周围能谱的影响，使用传统的两步法计算结果会造成较大误差[21]。此外，PB-FHR的冷却剂具有较大的中子吸收截面与散射截面，燃料球受到浮力作用堆芯活性区下部存在大量的冷却剂熔盐，其对堆芯的有效增殖因数及轴向功率分布有较大影响，均匀化时需予以考虑。

3) 连续换料

燃料球在堆芯内连续移动，首先需考虑由于轴向与径向的功率分布差异带来的燃耗分布影响，因此在燃料球的连续移动过程中径向与轴向的燃耗都不相同。其次，燃料球出堆后经过燃耗检测，未达到燃耗阈值的燃料球重新进入堆芯，超过燃耗阈值的燃料球卸出堆芯，同时加入等量的新燃料球。传统的轻水堆程序难以对以上过程进行模拟。

1.2 输运程序DRAGON5与扩散程序DONJON5

确定论分析程序DRAGON5与DONJON5是加拿大蒙特利尔综合理工学院(大学)开发的反应堆数值分析程序。其中DRAGON5是围绕中子输运方程的求解而设计的栅格程序，包含多种功能模块：精细群微观截面数据库处理的LIB模块；几何特征描述的GEO模块；空间离散模块，包括碰撞概率模块SYBILT、EXCELT、NXT，离散坐标模块SNT；共振处理模块，包括基于等价原理的共振处理模块SHI，子群共振模块USS，输运方程求解模块FLU，宏观群常数加工模块EDI，燃耗计算模块EVO等。以上模块通过软件包内的程序GAN[22]连接在一起，模块之间通过良好的数据结构交换数据。为对具有双重非均匀性的燃料球均匀化，使用GEO中的descDH模块进行几何描述。为找到合适的空间离散与共振处理模块，将本文计算结果与MCNP计算结果和基准题进行比较计算。

DONJON5是中子扩散方程求解程序，使用与DRAGON5相同的软件环境，也通过不同的独立功能模块及数据结构与GAN[22]连接在一起。本文主要用到三维几何模块RESINI、程序包内Trivac程序[23]的空间离散模块TRIVAT与中子扩散方程求解模块FLUD、群常数插值模块NCR、功率输出模块FLPOW。

1.3 泄漏效应处理

为考虑PB-FHR的泄漏效应，并解决堆芯活性区下部熔盐的均匀化问题，本文使用三步法的计算流程[21]。三步法可考虑由于石墨慢化球床堆中子自由程大于燃料球带来的泄漏效应问题。三步法计算流程如图1所示。第1步通过对燃料球与外层熔盐均匀化，得到堆芯活性区的精细群宏观群常数。第2步根据上一步得到的燃料球栅元的精细群宏观群常数，建立简单堆芯几何结构，以对其进行均匀化。在这一步的均匀化中，由于考虑了燃料球外部的信息，因此可解决由于石墨慢化的球床堆中子自由程大于燃料栅元组件带来的泄漏效应问题。第3步是在扩散计算中完成，首先建立三维精细的堆芯模型，根据第2步简单堆芯均匀化得到的少群常数进行全堆的扩散计算，得到有效增殖因数与全堆的功率分布。

图1 三步法均匀化计算流程[21]Fig.1 Three steps method homogenization calculation process[21]

1.4 固态熔盐堆燃料管理程序

为解决PB-FHR的换料问题，本文耦合DONJON5进行计算，理论模型参考PEBBED程序燃料球的移动模型[15]。(x，y)流道上的燃料球移动模型具体形式如下：

(1)

其中：B(x，y，z)为在堆芯(x，y，z)位置处燃料球的燃耗；B(x，y，0)为此流道内入口燃料球的燃耗；w(x，y)为在此流道内燃料球的移动速度；A为燃耗与裂变反应率的换算常数；Σf(x，y，z)为堆芯(x，y，z)位置处的裂变群常数；φ(x，y，z)为堆芯(x，y，z)位置处的中子通量。

式(1)可简化为如下形式：

(2)

其中：M为i网格内所有的重金属质量；ΔT为从z=0移动到z处所需时间；Pi为i网格内的功率。

燃料循环模式按燃料球通过堆芯的次数，分为单次通过与多次通过。其中单次通过指燃料球在堆芯内仅通过1次就卸出堆芯，多次通过指燃料球在堆芯内通过多次，当达到燃耗阈值后才卸出堆芯。单次通过模式下燃料球沿流道逐步向出口移动，程序主要解决燃料球在每个流道内的移动与跟踪问题。多次通过模式参考VSOP程序[12-13]，不同处在于PBMSR将通过堆芯不同次数的燃料球燃耗进行混合，再按照体积份额进行平均。之所以可将燃耗混合代替VSOP的核素混合，是因为燃料球在整个寿期内，中子能谱只发生很小变化，因此可忽略群常数的变化，此时核素浓度分布与燃耗深度为唯一对应关系。

为给出PB-FHR在多次通过模式下优化的卸料燃耗，需首先按照点燃耗模式计算全堆有效增殖因数keff达到1条件下的燃耗深度。PBMSR计算的全堆平均燃耗Bave为：

(3)

其中：Bi为i网格的燃耗；n为堆芯内网格的总数。

图2 PBMSR计算流程图Fig.2 Flow chart of PBMSR calculation

耦合程序PBMSR的计算流程如图2所示。DONJON5首先给出全堆的功率分布与有效增殖因数，PBMSR根据功率分布计算得到燃耗分布。当用户设置为连续移动模型时，进入换料模块；当用户设置为全堆燃耗混合时，进入全堆燃耗混合模块。进入换料模块后，再进行单次通过与多次通过的燃料循环模式判断。换料模块流程结束后，给出燃料球移动后的全堆燃耗分布；当选择全堆燃耗混合模块时，每个燃耗步堆芯内的燃耗深度相同。当达到用户定义的燃耗步，或堆芯有效增殖因数小于0.99时程序结束。

1.5 燃料管理模型

为实现堆芯达到平衡态卸出堆芯的燃料球达到最深燃耗，需给出合适的燃料管理模型参数。PB-FHR与HTGR燃料管理模型参数主要包括燃料球在堆芯内的移动速度分布及燃料球的卸料燃耗。本文计算时假设从启堆到平衡态加入的新燃料球核素浓度相同，暂时不考虑不同种类燃料球的燃料管理。

根据Tang等[24]的理论模型，燃料球在连续换料条件下能达到的卸料燃耗Boptimal约等于考虑堆芯中子泄漏率单个燃料球无限增殖因数kinf到达1时对应燃耗深度的两倍。但由于泄漏率与堆芯结构及燃料类型相关，使用单个燃料球难以给出准确的最深燃耗。为解决这一问题，PBMSR通过在每个燃耗步下将堆芯内燃料球燃耗混合，给出有效增殖因数等于1时对应的平均燃耗Bkeff，参考Tang的理论模型给出优化的卸料燃耗Boptimal为：

Boptimal=2Bkeff

(4)

卸料燃耗Boptimal确定后，根据燃耗的定义式可求出燃料球在堆芯内的总辐照时间Ttotal：

(5)

其中，P为每个燃料球受到的平均功率。

PB-FHR与HTGR中燃料球在堆芯内的流动速度w为：

(6)

其中：Δt为单个燃耗步长时间；Hnode为换料网格的高度。

为得到燃料球在堆芯内通过的次数Nshuffle，需求出燃料球在堆芯内移动的平均速度waverage：

(7)

其中，N为在堆芯中运动的燃料球总数。

根据总的辐照时间Ttotal及堆芯活性区轴向高度Htotal就可得到燃料球在堆芯内通过的循环次数Nshuffle：

(8)

2 验证

2.1 双重非均匀性验证

输运程序DRAGON5具备计算双重非均匀性的功能，但空间离散与共振处理有不同的功能模块，为选择更加精确的模块，使用MCNP进行比较。为验证DRAGON5计算的准确性，与HTGR的燃料球基准题[25]进行比较。

选用的燃料球参数参见文献[25]，当温度为1 000 K时的计算结果列于表1。由表1可看出，DRAGON5在求解双重非均匀性的燃料球时，空间离散模块使用一维碰撞概率模块SYBILT、共振自屏模块选取SHI模块时与MCNP计算结果最接近。下面计算均选用这两种模块对燃料球进行均匀化。

表1 DRAGON5不同计算模块与MCNP结果对比Table 1 Result comparison among different DRAGON5 modules and MCNP

为进一步验证DRAGON5处理具有双重非均匀性燃料球的正确性，使用基准题[25]中不同的分析程序进行对比，当温度为293.6 K时的计算结果列于表2。由表2可看出，使用DRAGON5与WIMS9和APLOLLO2计算得到的kinf相对误差很小，进一步证明DRAGON5处理双重非均匀性燃料球是可行的。

2.2 PB-FHR全堆计算验证

为验证程序PBMSR计算PB-FHR的正确性，与Li等[26]的计算结果进行对比。其计算使用MCNP与ORIGEN2的耦合程序，计算条件为燃料球不运动。由于使用点燃耗程序，堆芯内各区域的燃料球核素浓度随燃耗演化相同。堆芯与燃料球模型如图3、4所示，关键参数列于表3。

表2 DRAGON5计算结果与基准题结果对比Table 2 Result comparison between DRAGON5 and benchmark

图3 PB-FHR堆芯模型Fig.3 PB-FHR core model

为考虑均匀化的泄漏效应及得到堆芯活性区下方熔盐的群常数，使用三步法进行计算。

第2步针对简单的二维堆芯进行均匀化，堆芯模型如图5b所示。考虑到不同几何位置处石墨反射层能谱差异，将反射层分为4部分，分别为上、下反射层，侧反射层1与侧反射层2。反射层外侧使用强中子吸收体10B作为真空边界的近似，用来考虑中子泄漏效应。少群常数使用两群能群，能群边界为4.0 eV。

图4 燃料球模型Fig.4 Fuel pebble model

表3 堆芯与燃料球的关键参数Table 3 Key parameter of core and fuel pebble

图5 燃料球均匀化模型(a)和简单堆芯均匀化几何模型(b)Fig.5 Fuel pebble homogenization model (a) and simple core homogenization geometry model (b)

第3步使用扩散程序DONJON5对PB-FHR的堆芯进行三维精细建模，堆芯活性区径向模型如图6所示。

红色代表燃料球，黄色代表反射层，每个网格边长为25 cm图6 堆芯活性区模型Fig.6 Active core model

图7 有效增殖因数对比Fig.7 Comparison of keff

有效增殖因数随燃耗的变化与文献[26]中的数据进行对比，结果如图7所示。由图7可看出，keff的误差基本保持在300 pcm以内，两者符合较好。此外还可得到，当keff=1时，燃耗深度约为105 GW·d/tHM。

2.3 换料功能验证

1) 单次通过模式

为验证PBMSR单次通过模式的正确性，与BATAN-MPASS[16]的计算结果进行对比。堆芯关键参数列于表4。少群常数能群边界选取文献[16]中的5群能群结构。

单次通过模式下卸料燃耗和堆芯活性区的轴向功率分布如图8所示。由图8a可见，PBMSR计算结果的卸料燃耗为65.1 GW·d/tHM，BATAN-MPASS的计算结果为63.46 GW·d/tHM，两者符合较好。

堆芯达到平衡态后，堆芯活性区的轴向功率分布如图8b所示。由图8b可看出，功率峰靠近堆芯入口，且大部分集中在堆芯上半部分，而下半部分的功率基本为0。这主要是因为入口处的燃料球燃耗深度为0。尽管入口的燃耗最小，但功率峰不集中在堆芯的入口，这主要是因为靠近堆芯入口的中子泄漏较大。

随着物联网的快速发展，超低功耗设计已经成为了模拟集成电路领域研究者的重点研究方向[1]。在可预见的未来，超低功耗模拟集成电路的发展将给可穿戴智能设备、智能家居、无线传感器网络等低功耗应用带来革命性的变化[2]。基准电流源和基准电压源是模拟/混合信号集成电路的重要组成部分，它们为电路设计提供了与稳定的参考电流和参考电压。在超低功耗模拟/混合信号集成电路设计中，要求其拥有极低的功耗以及较高的电源电压和温度特性，以满足超低功耗的要求。

表4 单次通过模式下的堆芯关键参数Table 4 Core key parameter of once-pass mode

图8 单次通过模式的卸料燃耗(a)和轴向功率分布(b)Fig.8 Discharged burnup (a) and axial power distribution (b) of once-pass mode

2) 多次通过模式

为验证PBMSR计算燃料球在多次通过模式下的正确性，与PEBBED[15]的计算结果进行对比。堆芯关键参数列于表5。

多次通过模式下的卸料燃耗和堆芯活性区轴向功率分布如图9所示。由图9a可见，PBMSR计算出的最深卸料燃耗约为83.7 GW·d/tHM，与PEBBED的计算结果符合较好。

由图9b可看出，堆芯出口处(轴向高度等于940 cm)功率相差较大，这主要是因为PEBBED的计算模型未考虑下反射层的影响。此外还可看出，轴向功率峰靠近入口，主要原因是达到平衡态后，不断加入的新燃料球使得入口处的燃料球燃耗深度小于堆芯的下半部分。

表5 多次通过燃料循环模式下的堆芯关键参数Table 5 Core key parameter of multi-pass mode

图9 多次通过模式下的卸料燃耗(a)和轴向功率分布(b)Fig.9 Discharged burnup (a) and axial power distribution (b) of multi-pass mode

3 基于PBMSR的PB-FHR燃料管理初步分析

为研究燃料球通过堆芯不同次数下PB-FHR的中子学性能，使用2.2节中的PB-FHR堆芯参数进行计算。根据有效增殖因数到达1时的燃耗为105 GW·d/tHM，可得到燃料球的最深燃耗为210 GW·d/tHM，根据平均功率密度计算得到总的中子辐照时间为1 066.07 d。设置燃料球在堆芯内的移动速度为匀速流动，燃料球通过堆芯次数的参数列于表6。在纵向网格25 cm条件下，燃料球运动速度越慢，单个燃耗步长时间越长，从而导致通过次数越少的燃料循环模式，每个燃耗步加入的新燃料球个数越多。

表6 燃料循环模式的参数Table 6 Parameter of fuel cycle mode

由表6可看出，每天卸出的燃料球数随燃料球循环次数的增加而增加。根据HTR-PM燃耗测量装置[27]所述，目前燃耗测量装置对燃料球的检测能力只能选择燃料球通过次数为9次以下的燃料管理方案。

使用PBMSR计算卸出堆芯的燃料球燃耗和轴向功率分布，如图10所示。由图10a可看出，尽管燃料球通过堆芯的次数不同，但到达平衡态时燃料球的燃耗基本相同，约为209 GW·d/tHM。这是由于燃料球出堆后随机进入堆芯的流道内，相当于受到平均功率，在中子辐照时间不变、燃料球类型不变的情况下(重金属装载量相同)，卸料燃耗相同。

由图10b可看出，功率峰主要集中在堆芯底部。这主要是因为PB-FHR的燃料球从堆芯底部进入，靠浮力从下向上移动，因此堆芯底部的燃耗深度较浅。此外还可看出，燃料球通过堆芯的次数越多，轴向功率展平越好。这主要是因为燃料球通过堆芯的次数越少，每个燃耗步添加的新燃料球越多(表6)，导致入口处的燃耗深度越浅。