预处理Krylov子空间方法在三维全芯实际反应堆数值模拟中的应用研究

《Annals of Nuclear Energy》：Application study of preconditioned Krylov subspace method on three-dimensional full-core practical reactor numerical simulation

高效预处理Krylov子空间方法在中子扩散方程中的应用研究，通过Bi-CGSTAB和GMRES结合ILU等预处理器，在M310和VVER-1000反应堆的三维稳态与瞬态模拟中验证了高精度（keff偏差<10pcm）和高效性（CPU加速比1.10-3.93）。

　　

Jiahe Bai|Hanbo Gao|Chenghui Wan|Hongchun Wu

西安交通大学核科学与技术学院，中国

摘要

在高梯度中子通量空间分布下进行高效且稳健的数值模拟仍然是反应堆物理学中的一个关键挑战。传统的扩散求解器在遇到强扰动（如控制棒移动）导致通量分布严重扭曲时，往往会出现收敛速度慢或数值不稳定等问题。为了解决这一问题，本研究开发并评估了一系列基于预处理Krylov子空间方法的中子扩散求解器，包括Bi-CGSTAB和GMRES，并结合了ILU、MILU、BILU和MBILU预处理器。这些求解器被应用于我们自主研发的核心分析代码SPARK中，用于对中国M310和VVER-1000反应堆的稳态和瞬态三维全核心模拟。从准确性和效率两个方面对它们进行了定量评估。关键核心参数（如有效增殖因子$k_{eff}$、反应性和相对功率分布）与原始求解器的结果进行了对比。数值结果表明，这些求解器的精度很高：$k_{eff}$和反应性的偏差分别在10pcm和±3 × 10^-5范围内，相对功率误差低于0.05%。在所开发的求解器中，使用ILU预处理的求解器表现最佳，Bi-CGSTAB的加速比介于1.10到2.10之间，GMRES的加速比介于2.20到3.93之间。此外，两种使用ILU预处理的求解器在CPU时间上没有显著差异，表明它们的效率相当。这些结果表明，特别是那些使用ILU预处理的子空间求解器，为大型复杂反应堆模拟中的中子扩散计算提供了一种高效的方法。

引言

三维（3-D）全核心多组中子扩散方程的数值解是反应堆物理分析中的关键任务，因为它直接影响核心分析代码的可靠性和适用性。通常，扩散方程使用有限差分法（FDM）、有限元法（FEM）或粗网格有限差分法（CMFD）等方法离散化，这些方法通常会产生大型、稀疏且非对称的线性系统。在反应堆核心计算中，计算负担尤为沉重，不仅因为线性系统的规模庞大，还因为极端的中子通量扭曲（例如由快速控制棒移动或强局部扰动引起）会导致传统求解器的收敛速度慢或数值不稳定。在核设计、安全分析和重新装载操作中，经常需要执行大规模且时间敏感的模拟，这种情况很常见。

先前的研究表明，在中子通量在空间中均匀变化的正常运行条件下，基于扩散的节点方法表现良好。然而，在VVER型反应堆中，核心/反射器界面附近以及插入控制棒的区域可能会出现非常陡峭的通量梯度，这使得传统的六边形几何形状的NEM公式难以保持准确性和数值稳定性（Ivanov等人，2002年）。Graham等人进一步证明，在传输-扩散耦合计算中，如部分插入的控制棒这样的强局部吸收体会导致明显的通量下降和尖锐的轴向和径向梯度，从而产生显著的功率分布误差（Graham等人，2016年）。因此，开发在高度扭曲的中子通量分布下仍然有效的高效且稳健的数值方法对于推进核心物理分析非常重要。

通常使用两种主要方法来求解线性系统：直接方法和迭代方法。然而，由于固有的局限性，直接方法（如高斯消元法）和经典迭代方法（如逐次松弛法）在实际应用中效率低且收敛速度慢。为了克服这些限制，Krylov子空间方法作为一种稳健且高效的方法应运而生，适用于大规模科学计算。一般来说，Krylov子空间方法根据其底层的正交化过程分为两大类：Arnoldi正交化和Bi-Lanczos正交化（Saad，2003年）。基于Arnoldi过程的方法包括共轭梯度（CG）方法和广义最小残差（GMRES）方法（Saad和Schultz，1986年）。而基于Bi-Lanczos正交化过程的方法包括共轭梯度平方（CGS）方法、准最小残差（QMR）方法和双共轭梯度稳定（Bi-CGSTAB）方法（Van der Vorst，1992年）等。实际上，GMRES和Bi-CGSTAB在数值求解器领域占据主导地位，特别是在大规模问题上，因为它们在非对称线性系统上的稳健性和效率较高。然而，在安全分析过程中，有时控制棒以极端方式插入反应堆核心，使得条件变得非常恶劣，纯Krylov子空间方法难以收敛。

为了降低条件数并加速收敛，预处理Krylov子空间方法受到了关注，它们可以显著提高计算效率并节省时间成本。许多研究调查了某些预处理Krylov子空间方法在反应堆物理分析中的应用（Yang等人，1993年；Chen等人，1997年；Zou等人，2017年；Xu和Hao，2020年；Cherezov等人，2024年），这些研究表明，基于预处理Krylov子空间方法的线性求解器在计算效率方面比直接方法和未预处理的子空间方法具有显著优势。应当注意的是，现有研究主要集中在评估核反应堆在相对简单和标准运行条件下的性能。然而，在实际应用中应模拟更复杂的情景，例如从反应堆活动区域中插入除一根棒以外的所有控制棒。在这种情况下，中子通量会呈现高梯度空间分布。

因此，本研究探讨了基于不同预处理子空间方法的中子扩散求解器在复杂和高梯度条件下的性能。具体来说，开发了基于Bi-CGSTAB和GMRES的中子扩散求解器，并结合了不完全分解（ILU）及其变体（包括块状不完全分解BILU、修改后的不完全分解MILU、修改后的块状不完全分解MBILU），并将其集成到我们自主研发的核心分析代码SPARK中。数值实验涵盖了稳态和瞬态情景，其中控制棒的位置和移动被设计为形成复杂的中子通量空间分布。此外，上述研究是在矩形组装和六边形组装的核心上进行的，分别对应M310和VVER-1000反应堆。此外，实现了无矩阵架构以最小化内存使用，从而加快功率迭代速度。

本文的其余部分组织如下：第2节基于CMFD推导了中子扩散方程，并讨论了基于GMRES和Bi-CGSTAB的求解器及其使用ILU及其变体作为预处理器的方案。第3节提供了矩形组装和六边形组装核心在稳态和瞬态条件下的数值验证和分析。第4节总结了本文的内容。

节选内容

理论和方法

在本研究中，使用粗网格有限差分（CMFD）方法和节点展开方法（NEM）的非线性耦合方案，在SPARK代码中解决了稳态和瞬态条件下的三维（3-D）多组中子扩散方程（Wan等人，2025年）。CMFD作为主要求解器，将中子扩散方程离散化为非对称线性系统，然后需要用Krylov子空间方法求解。因此，本节首先提供了

数值实验

本节提供了全面的数值实验，以确定子空间方法和预处理器的最佳组合。为了确保通用性，进行了广泛的应用测试，包括矩形组装和六边形组装反应堆、稳态和瞬态三维全核心模拟，以及不同类型的棒插入状态。研究基于稳定性和效率评估了求解器的性能，重点关注关键核心参数，如$k_{eff}$、反应性等

结论

本研究调查了基于预处理Krylov子空间方法的多种中子扩散求解器的性能，重点研究了结合了各种ILU类型预处理器的Bi-CGSTAB和GMRES算法。这些求解器被集成到SPARK核心分析平台中，并应用于M310和VVER-1000反应堆在稳态和瞬态条件下的全核心模拟。结果表明，所有求解器都保持了较高的精度，$k_{eff}$和反应性的偏差都在可接受范围内

利益冲突声明

作者声明他们没有已知的可能会影响本文报告工作的财务利益或个人关系。

致谢

本工作得到了国家自然科学基金（编号：12375175）的财政支持。