一种用于随机介质中粒子输运的空间一致弦长采样方法

《Computer Physics Communications》：A spatially consistent chord length sampling method for particle transport in stochastic media

　　摘要：弦长采样（CLS）方法是一种高效的隐式几何建模技术，用于随机介质中的粒子输运模拟。在传统实现中，每个粒子独立采样弦长，导致同一空间位置处的物质分布因粒子而异。本研究开发了一种适用于随机介质的CLS输运方法，并在OpenMC中实现，该方法通过使用与位置相关的确定性种子，为弦长采样赋予空间一致性和可重复性。物理上一致的边界处理进一步消除了对传统填充分数校正参数的需求。该方法通过三阶渐进验证策略，利用中子临界性计算进行系统评估。对于四组白色边界球形基准案例，keff偏差均控制在100 pcm以内。含有28种参数组合的PWR单针基准案例以及具有多种厚度的MTR单板基准案例表明，CLS偏差主要受几何参数影响，偏差方向取决于球体与基体材料的相对中子特性。参数β定义为平均弦段周期与随机介质特征尺寸D的比值，被确定为有限几何结构分散燃料中CLS适用性的控制参数，它代表了粒子在随机介质中完成完整弦段循环次数的倒数。系统偏差呈Δkeff∝β的变化趋势。无论是棒状还是板状数据，在β偏差图中都表现出一致的衰减趋势，当β＜0.1时，偏差仍保持在100 pcm量级。值得注意的是，CLS方法的计算效率不会随分散球体数量的增加而下降，且在百万球体规模下的几何存储需求还不到显式方法所需存储量的十万分之一。该方法的几何加速框架与所输运的粒子类型无关，可扩展到更广泛的辐射输运场景。

引言：随机介质在核工程、辐射防护和医学物理领域极为常见，其中物质成分的空间分布具有固有的随机性，包括随机分布的粒子、孔隙和纤维等多种形态。辐射粒子（中子、光子等）在此类介质中的输运行为会受到物质不均匀分布的显著影响，因此准确模拟随机介质中的辐射输运对于核安全分析和辐射屏蔽设计至关重要。蒙特卡洛方法因其能适应复杂几何结构并具备连续能量处理能力，已成为随机介质辐射输运模拟的主要方法框架。然而，当随机介质的材料不均匀性表现为大量离散粒子的分布时，蒙特卡洛方法的几何建模效率就成了限制其实际应用的关键瓶颈。

在核反应堆工程中，由于其出色的裂变产物滞留能力和辐照稳定性，分散燃料颗粒已成为高温气冷反应堆（HTGR）燃料球、板状分散燃料元件以及各种实验性反应堆燃料设计的核心组成部分[1]。在这些系统中，大量燃料球随机分布在基体材料中，对这种随机介质系统进行精确的粒子输运模拟对于反应堆物理设计和安全分析而言极为重要。显式几何建模方法将每个球体视为独立的几何实体，直接记录和追踪球体与基体之间的物质界面。当球体数量达到百万级别时，界面追踪的计算成本就会变得极高[2]。

为了解决显式建模的效率瓶颈，研究人员在球体堆积、几何追踪加速以及蒙特卡洛代码的随机几何功能方面开展了大量研究。随机球体堆积是显式建模的基础步骤。Widom提出的随机顺序添加（RSA）算法是一种生成随机球体分布的经典方法，为后续研究奠定了理论基础。Brown和Martin率先在MCNP5中引入随机几何模块，实现了TRISO燃料颗粒的蒙特卡洛输运模拟。Liu等人将三种随机几何方法——随机格子法、弦长采样法（CLS）和显式建模法——应用于RMC代码中，而Rintala等人则将离散元方法（DEM）与Serpent相结合，实现了对砾石床反应堆的精确建模。Liu等人进一步将RMC的随机几何功能扩展到多类型球体和砾石床系统，并将其应用于HTR-10基准计算中。Wu等人在RMC中同时实现了显式建模法和CLS方法，为这两种策略的比较研究提供了统一的代码平台。随着砾石床反应堆工程需求的增长，Feng等人开发了ODR-VS方法，用于在高填充分数下快速生成球体，而Li等人则采用RMC与LAMMPS联合建模的方法，对HTR-PM首堆临界性进行了高精度模拟分析。Leppanen在Serpent中实现了基于Voronoi剖分的随机几何功能，将随机介质建模能力扩展到了碎片化燃料等新应用场景。Tan等人在RMC代码中实现了基于八叉树加速的并行RSA算法，以提高大规模堆积时的效率。在几何追踪加速方面，Liang等人开发了虚拟格子法，该方法利用哈希网格来管理随机球体的空间位置。Chen等人提出了一种多粒子格子建模方法，用于更精确地模拟板状分散燃料，而Liu等人则提出了一种基于自适应八叉树的射线追踪方法，该方法能够动态构建空间索引结构，从而加速大尺度随机介质中的射线追踪。不过，上述方法的核心局限在于需要存储所有球体的空间位置信息，内存消耗会随着球体数量的增加而呈线性增长。当球体数量达到百万级别，且每个球体都被视为独立的耗尽区时，内存需求就会成为实际工程应用中的严重瓶颈[2]、[10]。

弦长采样（CLS）方法提供了一种截然不同的技术途径：它通过在每个物质相中对粒子飞行距离进行实时采样，替代了显式几何追踪，从而避免了大规模球体界面的存储和计算。根据对粒子相的处理方式，现有的CLS方法可分为两类。第一类是完全隐式方法，这类方法对基体和粒子相都采用弦长采样，而不进行显式的球体几何建模，而是在输运过程中根据预定公式切换物质。Zimmerman和Adams提出了CLS的基本算法框架，为随机介质中中子输运的隐式模拟奠定了理论基础。Murata等人在MCNP中开发了MCNP-BALL程序，将CLS方法扩展到砾石床反应堆的连续能量蒙特卡洛计算中。Donovan和Danon将CLS方法扩展到二维随机介质，系统比较了全CLS方法和有限CLS方法的精度和效率，并分析了散射效应对CLS精度的影响。Liang和Ji系统分析了内存效应和边界效应对CLS精度的定量影响，阐明了在密集堆积系统中弦长独立性假设的局限性。Liang等人开发了一种连续能量CLS方法，并比较了三种填充分数校正方案，得出结论认为在HTGR基准问题中，Ji-Martin自洽校正效果最佳。Liang和Ji进一步提出了一种基于涂层均匀化的新的核间弦长分布模型，通过将涂层区域与基体材料均匀化，仅对燃料核之间的弦长进行采样，从而实现了TRISO燃料系统的高效建模。Li等人将CLS方法应用于一种氟盐冷却高温反应堆的燃耗分析，研究了在5%至30%的填充分数范围内CLS方法的精度表现。第二类是半隐式方法，这类方法仅对基体应用弦长采样，而对球体则保持显式几何建模。Liu等人在RMC中开发了一种改进的CLS方法，通过体积填充分数校准，消除了多球体混合配置下的系统偏差。Wu等人在RMC中同时实现了显式建模法和CLS方法，以支持多类型球体混合设计。

从研究问题的角度来看，CLS相关文献可分为两大研究方向：固定源问题和临界性问题。固定源问题侧重于传输和反射等积分量，以理想化的马尔可夫二元混合物作为研究对象，通过解析或半解析的基准解来评估CLS的基本精度。Adams等人建立了经典的棒状和板状基准配置，Brantley则比较了基本CLS算法和改进后的算法。Larmier等人首次为二维和三维泊松剖分提供了参考基准解，同年还研究了泊松、Voronoi和Box混合统计对三维随机介质中输运的影响[27]。Larmier等人将CLS方法扩展到d维马尔可夫二元混合物，系统研究了维度对CLS偏差的影响，并提出了泊松盒采样算法[29]，通过引入局部空间记忆来提高精度。Vu和Olson提出了条件点采样（CoPS）方法，利用完整的球体间记忆，在一维马尔可夫介质中实现了与参考解一致的输运精度。Tentori等人开发了一种具有记忆效应的CLS方法，用于空间非均匀的马尔可夫介质，还将记忆效应纳入LP方程中，建立了广义的理论框架[32]。临界性问题则关注特征值和裂变率分布，研究CLS在分散燃料工程中的应用精度，这类应用对弦长的空间一致性要求更为严格。Liang和Ji的一系列研究成果以及Liu等人改进的CLS方法都属于这一研究方向。Chen等人提出了针对边界效应和背散射效应的校正策略，Liang等人则基于等效传输概率构建了弦长校正因子，以抑制由体积排斥效应引起的偏差。Tan等人提出了具有局部记忆的半隐式弦长采样（SCLS）方法，通过引入局部显式几何结构，提高了强吸收体条件下的精度。Liu等人对CLS方法在分散燃料板中的应用进行了初步分析，评估了现有校正方案在板状几何结构中的适用性。

CLS方法的核心优势在于无需存储球体位置信息，内存消耗与球体数量无关，而且在输运过程中无需进行球体索引查询或边界搜索，因此其计算速度远远高于显式几何方法[20]。现有的CLS改进方案大致可分为两类：分布校正和空间记忆增强。前者通过调整有效填充分数等参数来减轻系统偏差，而后者则通过引入局部几何信息来恢复弦长之间的空间相关性。上述改进方案均为针对特定误差机制提出的有针对性的校正策略，目前尚未建立统一的CLS精度标度框架，用以预测该方法在不同几何结构下的适用范围。在本研究中，我们提出了一种基于确定性种子的弦长采样框架，该框架通过使用与位置相关的确定性种子，为弦长采样提供了可重复性和空间一致性。该方法属于半隐式CLS类别：球体保持显式几何建模，而基体的弦长则是通过确定性种子实时生成的。从方法论的角度来看，CLS本质上是一种几何加速方法：弦长的实时采样完全取决于物质的几何参数（球体半径、体积填充分数）和空间位置信息，与所输运粒子的物理类型无关；确定性种子机制的输入仅为空间坐标和飞行方向的量化哈希值，同样也不依赖于粒子类型。因此，本文提出的空间一致性框架既适用于中子输运，也适用于光子输运。本研究以中子输运临界性计算作为验证对象，而该方法本身具备直接扩展到其他辐射输运场景的能力。除了几何加速之外，CLS方法论还与蒙特卡罗输运领域的多个研究方向存在关联。第2节中关于指数分布假设的偏差分析，以及第7节所述的表格化PDF扩展方法，都与近期针对半解析输运方法中指数近似误差的修正工作具有相同的方法论基础[37]。第3节中的临界性计算及集合平均处理方法，也面临着与高级蒙特卡罗特征值算法相同的收敛性和统计波动问题[38]。全局方差降低方法旨在提升整个域内的计数效率，而CLS则着眼于几何级数的加速与偏差控制，为随机介质输运提供了另一种互补的技术路径。

OpenMC是由Romano、Forget及其合作者开发的开源连续能量蒙特卡罗粒子输运代码[39]，它支持基于DAGMC的CAD几何建模、实时耗尽计算以及多物理场耦合功能，已成为核工程领域重要的开源研究工具。目前，OpenMC尚未内置隐式随机介质输运模拟功能，分散燃料系统的建模仍依赖显式几何方法。为此，我们在OpenMC框架内开发了一种半隐式CLS方法：通过位置相关的确定性种子，利用空间坐标和飞行方向的量化哈希值作为随机数种子，使得弦长采样结果由空间位置和飞行方向唯一确定，从而为弦长采样提供了可重复的确定性框架。本文将阐述该方法的基本原理和OpenMC实现架构，通过三级逐步验证策略系统评估其计算精度，并基于跨几何配置的验证数据提出一个统一的CLS精度标度参数。

**基本原理**
考虑由两种材料构成的二元随机介质系统：颗粒相（材料1）和基体相（材料2）。颗粒相由半径为r的硬球组成，以下简称“球体”，它们以体积填充分数?随机分布在基体中。当颗粒沿随机方向穿过该介质时，其在颗粒相和基体相中的飞行路径长度——即弦长l1和l2——遵循特定的统计分布。

**验证案例**
为验证OpenMC中空间一致性CLS方法实现的正确性，设计了三级逐步验证策略。第一级采用球形白色边界基准验证，通过TRISO燃料颗粒和分散毒物这两种典型应用场景，检验在无限介质近似下的方法基本精度。第二级则是对PWR单棒有限几何结构的验证，采用UO2+Zr体系。

**通过β参数的统一分析**
上述验证结果表明，CLS偏差会随着几何参数的变化而呈现系统性衰减趋势，且这种趋势在棒状几何结构和板状几何结构中均一致。本节提出将参数β=(〈l1〉+〈l2〉)/D作为CLS精度的统一标度参数，其中D为随机介质区域的特征尺寸——对于球形几何结构，取为等效容器半径；对于棒状几何结构，则取为燃料区域半径。

**计算性能**
虽然空间一致性CLS方法在精度方面具有优势，但其隐式几何建模方法在计算效率和存储开销方面也可能带来好处。本节从三个角度定量评估了CLS方法相对于显式几何方法的性能特性：计算速率随球体数量的变化、不同材料间的计算效率比较，以及并行扩展能力。

**讨论**
跨材料比较表明，CLS偏差受几何参数?和r的影响，而不仅仅取决于吸收截面大小。对于同一材料而言，不同几何条件下的偏差变化远大于在相同几何条件下不同材料之间的偏差差异。这种几何因素的主导作用源于一个微观原因：用于描述基体弦长的指数分布忽略了连续弦长之间的空间相关性。

**结论**
本研究在开源蒙特卡罗代码OpenMC框架内开发并实现了一种半隐式CLS随机介质输运方法，通过位置相关的确定性种子，为弦长采样提供了空间一致性和可重复性。颗粒相仍采用显式几何建模，而基体弦长则通过确定性种子实时生成，从而保持了弦长分布的理论可控性。

**CRediT作者贡献声明**
刘鹏飞：撰写——初稿、可视化、软件、方法论、正式分析、概念设计。陈玉清：撰写——审稿与编辑、监督。李松：资源获取、资金筹集。陈彦伟：撰写——审稿与编辑。徐晓北：撰写——审稿与编辑。施学良：撰写——审稿与编辑。

**利益冲突声明**
作者声明不存在任何可能影响本文所述工作的已知财务利益或个人关系。

刘鹏飞|陈玉清|李松|陈彦伟|徐晓北|施学良
中国武汉430033，海军工程大学