核数据敏感性与不确定性分析及其在目标精度评估中的应用

刘 勇，曹良志，吴宏春，郑友琦，万承辉

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

反应堆是一个极端复杂的系统，其数值模拟过程往往被分解成若干个计算环节。而每个计算环节都存在不确定性的引入与传递。核数据作为反应堆物理计算最重要的输入参数之一，其不确定性将被引入并传递到后续计算当中。量化核数据引入的不确定性是反应堆不确定性分析的重要内容。另一方面，不确定性分析不仅能给出计算结果的不确定性，还能为数值计算的精度改善提供足够的数据支持，识别不确定性来源，给出能有效降低计算结果不确定性的着手点。

出于经济性和安全性的考虑，为减少反应堆设计或操作的安全裕量，高精度的数值模拟计算结果是重要前提，因此设计目标精度的概念被提出。具有较高精度的输入参数是保证计算结果不确定性满足目标精度的重要条件。目前，核数据的精度水平对反应堆计算结果引入的不确定性是显著的[1-4]，因此，有必要对核数据的不确定性水平提出定量要求，为与核素、反应和能量相关的核数据的改进需求提供参考标准，以达到设计目标精度，一般称这个过程为目标精度评估[5]。与不确定性计算不同的是，不确定性计算是根据输入参数的不确定性获取输出参数的不确定性，而目标精度评估是根据目标参数的不确定性限制，反求出输入参数的不确定性水平要求。因此相对来讲，二者是一个互逆过程。目前，针对4代堆的不确定性分析及目标精度评估正受到国际上的广泛关注[5]。

对于以上问题，本文以快堆为对象进行研究。首先，为能有效获取堆芯参数对微观核数据的敏感性与不确定性信息，基于堆芯两步法的计算方案，提出基于敏感性传递的堆芯参数敏感性与不确定性计算方法。在西安交通大学自主研发的快堆程序SARAX[6]的基础上，进行敏感性与不确定性分析程序COLEUS的开发。以快堆基准题BN-600[7]为敏感性与不确定性分析对象，量化堆芯有效增殖因数keff由核数据引入的不确定性及不确定性的主要贡献者。其次，针对目标精度评估问题，建立相应的优化问题数学模型。为克服传统优化算法收敛速度慢、易陷入局部最优解以及变量初值依赖性较强等问题，在COLEUS中实现差分进化算法计算功能，进行优化问题的求解，并基于BN-600问题的敏感性与不确定性分析结果，进行目标精度评估的研究。

1 理论模型

1.1 敏感性与不确定性计算方法

设反应堆物理计算某积分响应R为核数据的函数，根据一阶不确定性传递公式可得R的相对方差[8]为：

rcov(R)≈SCST

(1)

式中：rcov为响应的相对方差；S为相对敏感性系数矩阵；C为核数据相对协方差矩阵。

核数据的相对协方差矩阵一般可从评价核数据库获取，因此不确定性量化的关键问题在于敏感性系数的计算。本文以keff为研究对象，根据微扰理论[9-10]，其敏感性系数计算式为：

(2)

式中：α为某种核数据；Φ、Φ*分别为前向和共轭中子通量；L为中子泄漏、吸收和散射算子；F为裂变算子；λ=1/keff。

针对快堆堆芯两步法的计算方案，对于堆芯keff的敏感性系数计算，目前大多直接计算堆芯对核素少群截面的敏感性系数，然后利用核素微观截面的协方差数据进行不确定性量化。这种计算过程通常会带来两个问题：1) 在敏感性系数计算中，针对的是组件计算得到的少群截面，未考虑截面扰动对共振计算、截面归并权重谱的影响；2) 在不确定性量化中，敏感性系数与协方差数据存在不一致，前者针对经过共振计算问题的相关截面，后者针对与问题无关的截面。

对于这两个问题，本文提出如下计算方法：1) 在组件计算中，通过直接数值扰动方法[11-12]逐能量段扰动各截面再进行组件计算，考虑截面扰动对共振计算和权重谱的影响，获取组件少群常数对各核素微观截面的敏感性系数；2) 在堆芯计算中，采用微扰理论，获取堆芯keff对组件少群常数的敏感性系数；3) 根据式(3)描述的敏感性系数传递法则，结合1、2步获取的敏感性系数，得到堆芯参数对各核素截面的敏感性系数。最后根据一阶不确定性传递公式进行不确定性量化。

(3)

式中：J为少群常数总数，包括类型和能群；Σ为组件少群常数。

基于以上理论模型，COLEUS进行堆芯参数敏感性与不确定性的计算流程如图1所示。

1.2 目标精度评估方法

目标精度评估是根据目标参数的不确定性限制，给出输入参数的不确定性水平要求，它需要平衡如下问题：一方面，为使核数据对目标参数的不确定性贡献尽量小，需要其自身标准偏差(di)尽量小；另一方面，由于核数据在测量、计算等方面存在难度，实际操作期望di尽量大。因此提出如下优化问题：

(4)

使得，

(5)

式中：Q为目标函数；i为核数据标识；I为核数据总数量；λ为代价因子，表征核数据评价的难易程度；d为核数据标准偏差；n为目标参数标识；N为目标参数总数量；R0为目标参数标准偏差要求。

图1 COLEUS堆芯参数敏感性与不确定性计算流程Fig.1 Flow chart of sensitivity and uncertainty calculation for core parameter in COLEUS

核数据分为核素、反应类型和能群，如果考虑所有核截面数据，上述优化问题的规模将是巨大的。考虑到有些核数据的不确定性贡献较小，或有些核数据不在关注范围内，可将其排除在优化问题外。为考虑协方差的影响，同时排除掉不参与求解的核数据的不确定性贡献，式(5)将约束条件分为如下两项，第1项表示参与求解的核数据的不确定性贡献项G(包括协方差项)，第2项表示不参与求解的核数据的不确定性贡献项P。

1.3 差分进化算法

针对1.2节建立的优化问题，可发现目标精度评估优化问题是1个带约束条件的、非线性的优化问题，具备较大规模。数学上这是一较难求解的问题。此时对优化算法的选择显得特别重要。目前国际上对于目标精度问题的求解研究较少，大多采用序列二次规划算法进行求解[13-14]。该方法属于确定性方法，该算法将实际的优化问题转化为一系列的二次规划问题，通过求解二次规划问题产生搜索方向进行迭代求解。该算法的研究一直是非线性优化问题的研究热点。传统的序列二次规划算法需大量求解二次规划问题，计算量和存储量均较大，当变量维度增大，问题计算规模变得相当可观。此外该方法还存在对二次规划子问题可能不相容以及Maratos效应等问题[15]。因此，目前该方法处于长期的发展和改善中，以使得该算法具备全局收敛性和超线性的收敛速度，同时克服算法的固有问题。

差分进化算法属于随机性方法，也属于智能优化算法，经典的测试函数实验表明差分进化算法是目前诸多智能优化算法中优化性能最好的算法[16-17]。该算法的主要优点是执行简单、收敛速度较快、全局优化性能突出等。本文拟采用该方法对目标精度评估问题进行求解。由于本文的研究重点不在算法本身，因此本节将对差分进化算法作简要描述，参考已有的研究实现该算法，并进行程序开发。

在差分进化算法中，一般称维度与目标函数决策变量数相同的单个向量为个体，由NP个个体构成的集合称为种群，NP称为种群规模。差分进化算法在算法执行过程中，保持种群规模不变，通过进化的方式改善种群中个体的质量，在可行域中搜索最优解。

典型的算法操作包括变异、修复、交叉和选择。变异操作釆用变异策略来生成变异个体，一般通过对上一代(父代)的差分操作，生成下一代(子代)的个体；修复操作的目的是使新生成的个体处于可行域内，保证生成的变异个体的可行性；交叉操作通过父代个体和子代个体的交叉策略来产生新的尝试个体；选择操作根据尝试个体和父代个体的评价函数来决定进入下一代种群的优秀个体。算法通过以上操作进行迭代计算，保证在每次迭代中淘汰劣质个体，保留优良个体，使得种群整体上向全局最优解区域靠近，优秀个体逼近全局最优解。

本文通过采用文献[18]提出的Oracle罚函数方法，改进典型的差分进化算法，使其能处理约束条件。最后形成的程序流程如图2所示，其中g为种群进化代数。

图2 COLEUS差分进化算法计算流程Fig.2 Flow chart of differential evolution algorithm in COLEUS

2 数值结果

BN-600快堆基准问题是由IAEA发布的一套全MOX燃料钠冷快堆问题，该问题的堆芯径向和轴向几何布置如图3所示。该问题包含低富集度钚燃料(LEZ)、中富集度钚燃料(MEZ)和高富集度钚燃料(HEZ)，位于LEZ中心的内部增殖区(IBZ)厚5.1 cm，3组控制棒(SCR)和1组补偿棒(SHR)插于LEZ区，HEZ区外有两个钢屏蔽区(SSA1(第1层)和SSA2(第2、3层))，最外层是径向反射层区(REF)。组件中心距为9.902 cm。外部边界为真空边界。

本文采用均匀计算模型，对控制棒全插问题进行计算。首先，SARAX对各组件进行零维超细群共振计算，获取1 968群有效自屏截面，然后进行零维组件计算获取1 968群能谱，对1 968群微观截面进行归并得到各组件26群宏观截面。其次，在堆芯层面，采用三角形节块SN方法进行输运计算，轴向分为26层，径向上每个六角形分为6个三角形网格。

a——径向几何结构及组件布置；b——轴向几何尺寸及材料布置图3 BN-600堆芯结构Fig.3 Core structure of BN-600

参考解采用蒙特卡罗输运程序MCNP建立均匀模型计算。在几何上严格描述各区域，但不考虑各区域内部详细几何结构而进行打混处理。计算条件为执行500代计算，每代投入100 000个粒子，前100代不参与最终统计，得到keff为1.018 02±0.000 07。SARAX的keff计算结果为1.012 08，与MCNP计算结果的相对偏差为-0.58%，具有较好的计算精度。

2.1 敏感性与不确定性分析

本文对BN-600中燃料、结构材料及冷却剂材料中最为关注的11种核素进行敏感性与不确定性分析，包括238U、235U、239Pu、240Pu、241Pu、56Fe、58Ni、52Cr、23Na、16O、10B。其中各核素的协方差数据来自于ENDF/B-Ⅶ.1评价核数据库，采用NJOY[19]程序制作。

图4 主要核素截面敏感性系数Fig.4 Sensitivity coefficient of main nuclide cross section

图4示出了该问题中keff对扰动最敏感的核素的截面(图中ν为平均裂变中子数，σf为裂变截面，σγ为(n,γ)反应截面，σinel为非弹性散射截面，σelas为弹性散射截面)，可发现，keff对作为主要裂变核素的239Pu的裂变截面和平均裂变中子数具有显著的敏感性。此外，238U的俘获截面、裂变截面和非弹性散射截面和反应堆中的结构材料的重要成分56Fe以及冷却剂23Na也具有较大的敏感性系数。

图5 主要核素截面不确定性贡献Fig.5 Uncertainty contribution of main nuclide cross section

敏感性系数表征反应截面扰动对反应堆数值计算结果的相对影响，为确定各反应截面自身带有的不确定性对计算结果引入的实际的不确定性，图5示出了核数据对keff引入的不确定性以及主要核素截面的不确定性贡献(图中σα、σn2n分别为(n,α)、(n,2n)反应截面)。可发现，由于核数据导致的keff总的不确定性约为0.935%，其中238U和239Pu占主导地位。238U是堆芯keff不确定性的主要贡献者，特别是238U的非弹性散射截面和俘获截面，占有相当大的比例。位于第2位的不确定性贡献者为239Pu，其俘获截面和裂变截面是不确定性的主要来源。56Fe和23Na也对keff不确定性同样具有可观的贡献。

2.2 目标精度评估

首先，为测试COLEUS在优化问题求解上的数值可靠性，本文采用4组30维的标准测试函数[20]进行测试(表1)。程序运行环境为32位Windows7操作系统，英特尔I5处理器，内存为4G。参考解来自文献[16]，为公平比较，采用与文献[16]相同的计算条件，将种群规模设为100，最大函数评价次数设为50万，统计30次独立运行的结果，以排除随机效应带来的影响，比较结果列于表2。可发现，COLEUS在这些测试问题上，较文献具有更靠近理论最优解的结果，证明了COLEUS在优化问题求解上的可靠性。

表1 优化问题标准测试函数Table 1 Standard test function for optimization problem

从2.1节分析可知，对于BN-600问题，由核数据对keff引入的总的不确定性约0.935%，而根据OECD对未来快堆keff不确定性的要求来看，期望其不确定性控制在0.3%以下。基于此，本文将keff的目标精度设置为0.3%。本文暂不考虑工程实际经验，将代价因子λ设为1.0。为控制变量数，排除掉对问题不确定性贡献很小的截面，提高计算效率，通过设置不确定性贡献占总不确定性的最小比例来对需参与优化求解的核数据进行筛选。本文将该比例设置为0.5%，参与问题求解的核素截面包括238Uσelas、238Uσinel、238Uσγ、238Uν、239Puσinel、239Puσf、239Puσγ、239Puν、56Feσelas、56Feσinel、56Feσγ、23Naσinel，此时以上核素总的不确定性贡献为0.929%，剩余不确定性(式(5)的Pn项)为0.101%，变量维度为250。

为保证COLEUS结果的可靠性，本问题种群数设为750，进化代数为10 000，独立执行30次计算进行统计，目标函数最小值为5.447 7×106±6.652 8×102，计算时间为790 s。30次独立运算均能使优化后的目标精度达到期望的0.3%，30次独立运算结果波动仅0.01%，证明了COLEUS的差分进化算法计算功能具有较好的鲁棒性和精度。单核平均每代计算时间为790 s，计算速度可接受。

鉴于数据量较大，对于截面优化结果，本文只列出各核素截面需进行精度改善的能量段及其初始和优化后的相对标准偏差(表3)。可看到，优化结果对不同的核素和截面提出了不同的不确定性要求，对于不少截面数据，其相对标准偏差需在原来的基础上降低到一半以下，甚至降低1个数量级。要求较高的截面包括238U的非弹性散射截面、239Pu的非弹性散射截面、辐射俘获截面和56Fe、23Na的相关截面。以238U的非弹性散射截面(图6)为例，可看到其不确定性从原来的10%～30%左右整体上需降低到10%以下。

图6 238U非弹性散射截面优化结果Fig.6 Optimization result of inelastic cross section of 238U

3 小结

通过提出针对堆芯两步法计算方案的敏感性与不确定性计算方法，有效获取了堆芯参数对核素截面数据的敏感性与不确定性信息。对BN-600快堆基准题进行了敏感性与不确定性分析，发现堆芯keff由核数据引入的不确定性达0.94%，其中238U的非弹性散射截面的贡献是主要的。

建立了目标精度评估问题的数学模型，并采用差分进化算法进行优化问题的求解，得到在BN-600堆芯keff由核数据引入的不确定性(目标精度)为0.3%的要求下，各核素截面的不确定性需要。数值结果表明，要满足该要求，目前大多核数据的不确定性需降低一半以上，对于不确定性贡献最大的238U非弹性散射截面，需从原来的10%～30%左右，整体上降低到10%以下。值得指出的是，本文尚未引入工程经验，代价因子暂设为1.0，后续工作可考虑各截面在工程实际中的评价难易度，获取更贴近工程应用的结果。