压水堆堆芯Pin-by-pin均匀化计算中非均匀泄漏修正模型研究

张 斌，李云召，吴宏春

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

20世纪70年代以后，得益于粗网节块法[1-2]的迅速发展，以组件均匀化理论[3-4]和粗网节块方法为理论框架的两步法计算方案逐渐成为压水堆工程计算中普遍采用的燃料管理中子学数值计算方法。两步法中，堆芯计算所用燃料组件均匀化截面由全反射边界条件下的单组件计算产生。栅格计算下的全反射边界条件无法考虑燃料组件在实际反应堆中径向与轴向的中子泄漏，求得的中子能谱实质上是无限介质能谱。

无限介质能谱的修正方法通常是采用泄漏修正模型[5-6]通过对堆芯几何曲率的计算近似考虑中子泄漏对全反射边界条件下单组件计算的影响，从而修正用于求解等效均匀化常数的中子通量密度分布。在泄漏修正模型中，通过对中子泄漏率的计算可得到下游堆芯低阶输运方程计算所用的扩散系数，即中子泄漏系数。

组件均匀化栅格计算中通常采用均匀泄漏修正方法对无限介质能谱进行修正(基模修正)，这是由于组件均匀化计算把整个组件均匀化成一种“等效”均匀介质，此特点完美契合均匀泄漏修正方法中组件各处中子泄漏率相等的设定。径向及轴向的中子泄漏导致燃料组件中子“消失”的部分可由中子泄漏系数和几何曲率表示。然而，在组件各处中子泄漏有较大差异且均匀化成多区“等效”均匀材料的情况下，均匀泄漏修正方法是不适用的。这是因为在此类情况下，中子射流效应和中子能谱干涉效应会对等效均匀化少群常数带来较大的影响，而均匀泄漏修正模型无法描述多个均匀化区域之间的中子学差异。

在压水堆燃料组件内各处中子泄漏率相等的假设对于Pin-by-pin 均匀化计算[7-8]显然是不合理的，尤其是当燃料组件中存在可燃毒物棒、控制棒、水洞等特殊栅元时。非均匀泄漏修正模型[7]假设在燃料组件内部中子通量密度和中子泄漏率在空间上是变化的，是一种更为精确的泄漏修正模型。由于栅元均匀化后的燃料组件与均匀化前的燃料组件之间存在一定的差异，因此非均匀泄漏修正模型在Pin-by-pin均匀化计算中应用时必须与栅元均匀化技术搭配使用，以保证均匀化前后目标问题的一致性。

本文在压水堆堆芯Pin-by-pin均匀化计算中分别采用均匀泄漏修正模型及非均匀泄漏修正模型对组件计算的无限能谱进行修正，研究栅元均匀化方法和非均匀泄漏修正模型在Pin-by-pin均匀化计算中的联合实现方式，并分析压水堆堆芯Pin-by-pin均匀化计算中不同扩散系数来源对堆芯计算精度的影响。

1 方法

中子泄漏率主要受中子各向异性散射和中子射流效应的影响。在压水堆中，慢化剂内大量轻核的存在导致中子各向异性散射非常严重，因此压水堆泄漏修正模型中会精确考虑一阶散射截面的影响；中子射流效应是由于堆芯内强非均匀性所引发的中子流动现象，在传统组件均匀化计算中，整个组件被均匀化成一种等效均匀介质，即假设组件内不存在中子流动现象，这种泄漏修正方式即为均匀泄漏修正方法。在Pin-by-pin均匀化计算中，由于燃料栅元所处位置不同，水洞、导向管甚至可燃毒物棒等特殊栅元的存在，导致组件内部中子射流效应的影响无法再被忽略。区别于均匀泄漏修正方法的非均匀泄漏修正方法能有效处理空间各处不同的中子泄漏率并提供空间相关的中子泄漏系数。

在堆芯任意位置，中子通量密度可由其空间分布和基态分布构成[9]：

φ(r，E，Ω)=Ψ(r)φ(r，E，Ω)

(1)

式中：r、E、Ω为空间、能量、角度；φ(r，E，Ω)为中子角通量密度分布；Ψ(r) 为中子通量密度空间分布；φ(r，E，Ω)为中子角通量密度基态分布。

由反应堆物理基础理论，在堆芯处于稳态时，堆芯内中子通量密度的空间分布满足波动方程：

(2)

式中，B2为堆芯几何曲率。

式(1)中，中子通量密度的基态分布是否与空间变量相关决定了泄漏修正模型的种类。如果与空间变量无关，称为均匀泄漏修正方法；反之，则称为非均匀泄漏修正方法。

1.1 均匀泄漏修正

均匀泄漏修正方法假设中子通量密度对空间变量与能量变量可进行变量分离，且中子通量密度基态分布呈现与空间无关的均匀分布，即式(1)可写成如下形式：

φ(r，E，Ω)=φ(E，Ω)Ψ(r)

(3)

在有限均匀介质内，中子输运方程如式(4)所示。考虑到压水堆中各向异性散射的重要性，对散射截面进行勒让德展开并舍弃2阶及以上的高阶散射截面，由式(3)、(4)可推导得到均匀泄漏修正模型的方程，如式(5)所示。

电子文档是记录图书馆生产的良好凭证，通过对电子文档的分析和整理，容易发现图书馆业务和读者使用的规律，从而改进服务方式和手段，为图书馆管理人员提供决策支撑材料。电子文档存贮在计算机中，通过检索软件进行检索，经过网络进行传播，方便用户对文档资源的搜索，节省查阅资料的时间。

(4)

(5)

(6)

式中：Σt(E)为宏观总截面；Σs(E←E′，Ω←Ω′)为宏观散射截面；χ(E)为中子裂变谱；φ(r，E)为中子标通量密度分布；keff为反应堆特征值；νΣf(E)为宏观中子产生截面。

均匀泄漏修正方法在栅格程序中通常有3种模型：B1模型、P1模型和扩散模型。

1) B1模型：假设散射截面是线性各向异性散射，角通量基态分布没有进行任何近似，方程组为式(7)。

(7)

2) P1模型：假设散射截面是线性各向异性散射，且角通量基态分布通过球谐函数展开后高于1阶的展开项均设为0。方程组形式除系数γ(B，Σt(E))=1.0外，与B1模型相同。

3) 扩散模型：假设散射是各向同性散射，总截面采用输运修正后的输运截面，且角通量基态分布通过球谐函数展开后高于1阶的展开项均为0。方程组形式除1阶散射截面Σs1(E←E′)=0外，与P1模型相同。

均匀泄漏修正方法应用于Pin-by-pin均匀化计算时，采用B1模型进行能谱修正，由于B1模型中除散射截面是线性各向异性散射假设外，未引入其他假设条件，能谱修正精度高于P1模型和扩散模型。计算中采用统一的中子能谱进行各栅元中子通量密度的能谱修正，堆芯所用扩散系数基于菲克定律由各栅元输运截面计算得到。

1.2 非均匀泄漏修正

非均匀泄漏修正方法假设中子通量密度分布由其空间分布与基态分布构成，且中子通量密度基态分布在空间上呈非均匀布置，如式(1)所示。在有限非均匀介质内，其中子输运方程如式(7)所示。

对散射截面进行勒让德展开并舍弃2阶及以上的高阶散射截面，由式(1)、(7)出发，通过源项各项同性近似、角度积分可推导得到非均匀泄漏修正模型的方程组，如式(8)所示。

(8)

(9)

(10)

非均匀泄漏修正方法中对各栅元的中子泄漏系数定义如下：

(11)

在Pin-by-pin均匀化计算中采用非均匀泄漏修正方法时，通过临界曲率搜索计算求解式(10)，保证燃料组件处于临界状态。非均匀泄漏修正方法能计算得到单组件内与空间相关的临界能谱，并提供各栅元中子泄漏系数作为下游堆芯低阶方程求解时的扩散系数。

根据Pin-by-pin均匀化理论模型，仅靠均匀化少群截面与扩散系数无法保证燃料组件均匀化前后反应率、中子泄漏率以及特征值守恒，因此非均匀泄漏修正方法无法直接应用于栅元均匀化后的组件问题。在非均匀泄漏修正计算前，栅格程序需计算包含不连续因子[10-11]或超级均匀化因子[12-13]在内的等效均匀化常数以保证栅元均匀化组件问题与非均匀组件问题保持守恒。本文对非均匀泄漏修正采用离散纵标方法(SN方法)进行求解，因此需提前产生针对SN方法的不连续因子。此处等效均匀化少群常数必须采用不连续因子而非超级均匀化因子，这是因为超级均匀化因子在均匀化理论中通过改变截面而保证反应率守恒，间接影响非均匀泄漏修正计算得到的中子通量密度。中子通量密度的人为改变决定了包含超级均匀化因子在内的等效均匀化常数无法用于临界曲率搜索中的渐进能谱计算。利用SN方法对含有不连续因子的栅元均匀化组件问题进行非均匀泄漏修正方程组求解，产生各栅元的渐进能谱和泄漏系数。求解非均匀泄漏修正方程组后，采用渐进能谱对中子通量密度进行能谱修正，并重新产生Pin-by-pin均匀化少群常数。由于各栅元中用于能谱修正的渐进能谱各不相同，因此在非均匀泄漏修正后，等效均匀化常数的计算必须使用超级均匀化方法。这是因为栅元不连续因子的计算与栅元交界面上的净中子流密度相关，交界面上净中子流密度受相邻栅元不同渐进能谱的影响而无法被唯一确定的渐进能谱修正，具体计算流程如图1所示。本文基于广义等效均匀化方法及超级均匀化方法进行栅元不连续因子[14]及超级均匀化因子[15]计算，并采用基于中子平衡的快速临界曲率搜索方法进行非均匀泄漏修正计算。

图1 非均匀泄漏修正计算流程Fig.1 Flow of heterogeneous leakage model

2 数值计算与分析

本文基于二维临界堆芯问题进行数值计算，验证与分析不同泄漏修正方法在压水堆堆芯Pin-by-pin均匀化计算中的应用效果。二维堆芯布置如图 2所示，堆芯中各类型燃料组件的几何结构采用商用压水堆中常用的AFA3G燃料组件，其富集度及可燃毒物棒数量示于图2。

图2 1/4压水堆堆芯布置Fig.2 Layout of 1/4 PWR core

采用4种产生等效均匀化少群常数的计算方案，进而在全堆芯Pin-by-pin均匀化计算中加以分析和比较，在全反射边界条件下，对单组件问题进行栅格计算，产生Pin-by-pin等效均匀化少群常数：1) 不带泄漏修正模型(W/O)的栅格计算，其扩散系数根据菲克定律计算得到(DC)；2) 带均匀泄漏修正模型(Hom)的栅格计算，其扩散系数根据菲克定律计算得到(DC)；3) 带非均匀泄漏修正模型(Het)的栅格计算，其扩散系数根据菲克定律计算得到(DC)；4) 带非均匀泄漏修正模型(Het)的栅格计算，其扩散系数由各栅元的中子泄漏系数提供(LC)。

二维临界堆芯参考解由MOC一步法计算获得，特征值keff参考解为1.0，棒功率分布如图3所示。全堆芯Pin-by-pin均匀化计算的特征值误差和棒功率相对误差列于表1。

图3 1/4压水堆堆芯棒功率分布Fig.3 Pin-power distribution of 1/4 PWR core

表1 特征值误差和棒功率相对误差Table 1 Errors of keff and pin-power

由表1可知：1) 全堆芯Pin-by-pin均匀化计算采用相同扩散系数时，均匀泄漏修正方法与非均匀泄漏修正方法的修正效果相近。2) 均匀泄漏修正方法和非均匀泄漏修正方法采用由菲克定律计算的扩散系数时，统计结果的计算精度与不带泄漏修正模型的相当。为进一步分析比较不同计算方案下的计算效果，本文统计了全堆芯的棒功率误差分布，如图4所示。可发现，采用菲克定律计算的扩散系数时，泄漏修正模型能改善堆芯Pin-by-pin均匀化计算的整体功率倾斜程度。3) 扩散系数由非均匀泄漏修正方法中各栅元的泄漏系数决定时，全堆芯Pin-by-pin均匀化计算的精度显著提高，特征值误差保持在50 pcm以下的同时，棒功率相对误差的均方根偏差与平均误差都降低到1%以下。

由图4可发现，带各栅元中子泄漏系数的泄漏修正模型能降低最大棒功率相对误差及其均方根误差的同时，全面降低堆芯各处Pin-by-pin均匀化计算的棒功率相对误差。最大棒功率相对误差出现在燃料组件与反射层接触的角点上，这是因为堆芯外围燃料组件的真实边界条件最不符合栅格计算采用的全反射边界条件。

图4 1/4压水堆堆芯棒功率相对误差分布Fig.4 Pin-power relative error distribution of 1/4 PWR core

3 结论

本文在压水堆堆芯Pin-by-pin均匀化计算中采用非均匀泄漏修正模型对组件计算的中子能谱进行修正，研究了非均匀泄漏修正模型和栅元均匀化方法在Pin-by-pin均匀化计算中的联合实现方式，并采用非均匀泄漏修正模型产生的中子泄漏系数作为各栅元扩散系数进行堆芯Pin-by-pin均匀化计算。数值结果表明，非均匀泄漏修正模型能有效提高压水堆堆芯Pin-by-pin均匀化计算的精度。