CN113553550B - 蒙特卡罗粒子输运计算方法及相关计算方法 - Google Patents

Abstract

本发明涉及一种蒙特卡罗粒子输运计算方法，用于计算点源发出的粒子经过屏蔽体后到达目标区域的概率值，该蒙特卡罗粒子输运计算方法包括如下步骤：计算点源发出的粒子的每一条标准轨迹在各个有关位置的粒子权重；基于每条标准轨迹进行变换生成对应的、经过所述目标区域的样本轨迹；基于每条标准轨迹与每条样本轨迹的对应关系确定每条样本轨迹在所述目标区域的权重w_j；根据事件概率分布函数来计算每条样本轨迹占有的概率宽度c_j；按照下列公式加权求和而得到粒子到达目标区域的概率值p：

其中，M是经过目标区域的所述样本轨迹的总数量。本发明还涉及一种采用上述蒙特卡罗粒子输运计算方法进行反应堆堆芯微扰计算的方法。

Description

蒙特卡罗粒子输运计算方法及相关计算方法

技术领域

本发明涉及一种蒙特卡罗粒子输运计算方法。本发明还涉及一种采用上述蒙特卡罗粒子输运计算方法进行深穿透屏蔽计算(包括屏蔽点通量计算和精细通量分布计算等)和反应堆堆芯微扰计算的方法。

背景技术

蒙特卡罗方法是在粒子输运计算领域被广泛使用的一种计算方法。它基于概率理论，从粒子与介质原子核之间的基本相互作用(碰撞、散射、吸收等)出发，通过随机抽样的过程模拟粒子在介质内的飞行轨迹和核反应事件。利用大量粒子跟踪所形成的统计规律，得到粒子输运问题的近似解。当所模拟的粒子数量足够多时，该统计近似解可以逼近粒子输运问题的精确解。

蒙特卡罗方法用来求解粒子输运问题和确定论方法相比有很多优势：物理概念简单，计算流程清晰明了；不用划分空间网格；粒子飞行角度无需离散；粒子能量无需分群。这些特点减少了几何形状、飞行方向、粒子能量分布上的近似，使它适用的问题范围比一般的确定论方法更广。

但是蒙特卡罗方法也有其自身的不足。当粒子输运问题的计算规模不断变大时，粒子行为的事件树规模会爆发式增长，要想得到较为精确的模拟结果，所需要的粒子数量也会爆发式增长，从而导致其计算效率变得极低。例如复杂几何条件下的粒子深穿透屏蔽问题计算、大规模问题的粒子精细通量分布计算等，都需要消耗大量的计算资源。另外，由于蒙特卡罗方法中每个粒子行为的随机性，使得它难以直接用来分析一些微扰过程(例如反应堆物理领域的温度系数计算、密度系数计算、控制棒微分价值计算等)，而需要采用另外的理论方法来间接计算。

为了解决蒙特卡罗方法针对深穿透屏蔽等大规模问题的计算效率问题，前人研究了各种方法来解决，其中最主要的方法就是偏倚抽样。通过设置粒子的位置重要性、能量权窗等方法，引导粒子尽量多地向目标区域飞行。但是这些已有的方法对偏倚的设置要么是基于物理经验知识(依靠经验，改进效果有限)，要么是需要先采用确定论方法求解一遍输运方程得到价值函数(确定论计算量大)，并不能很好地解决这些问题，因此深穿透屏蔽计算、精细通量分布计算、微扰计算等问题仍然是当今蒙特卡罗方法的研究热点。

图1给出了一个典型的深穿透屏蔽计算问题示意图。图中S是一个放射源，I是一个点通量探测器，放射源和探测器之间是由各种屏蔽材料组合形成的屏蔽体，屏蔽体的屏蔽效果是使粒子通量衰减约13个量级。现在已知放射源源强和屏蔽体内的材料分布，需要计算探测器位置的粒子通量。

图2示出了传统蒙特卡罗方法的计算流程。首先根据屏蔽体内的材料分布，结合基础截面数据库得到粒子在每个材料区域的事件概率分布(例如首次碰撞自由程、碰撞后的散射和吸收概率、散射后的飞行角度概率分布、能量概率分布等)，然后根据事件概率分布函数进行随机采样，使粒子轨迹按照所采样的事件向后发展，直到粒子到达目标区域或者消亡为止。所有粒子跟踪完毕后，可根据目标区域的统计计数得到粒子到达该区域的概率值。因此，传统蒙特卡罗方法求解深穿透屏蔽问题的核心困难在于总样本空间非常巨大，而实际有效统计样本所占比例极少。

因此，这类问题如果采用蒙特卡罗方法来计算，将会面临很大的挑战。由于S位置处的源粒子穿透屏蔽层到达右侧表面的概率只有10^-13左右，而这些粒子恰好从探测器位置I附近穿透的概率又要比10^-13再降低几个量级，因此在S位置处需要模拟极大数量的源粒子发射才能在I位置得到有统计意义的通量。即使采用偏倚抽样的方法引导源粒子尽量向I位置飞行，由于飞行轨迹抽样的随机性，而且粒子飞行方向的偏倚难以依靠经验或者确定论计算来设置，也无法保证大部分轨迹能够有效地到达I点附近，从而导致求解该问题过程较为复杂而且所需的计算量仍然很大。

发明内容

针对上述问题，本发明提出了一种基于新的计算流程的蒙特卡罗粒子输运计算方法，该方法可以保证所跟踪的每条轨迹均能够从源粒子处到达目标区域，去除无效轨迹的跟踪，可大幅减少蒙特卡罗方法求解此类问题的计算量，适用于深穿透屏蔽计算、精细通量分布计算以及微扰问题的计算。

为解决上述问题，本发明提出了一种蒙特卡罗粒子输运计算方法，用于计算点源发出的粒子经过屏蔽体后到达目标区域的概率值，该蒙特卡罗粒子输运计算方法包括如下步骤：计算点源发出的粒子的每一条标准轨迹在各个有关位置的粒子权重，其中，所述标准轨迹为所述点源发出的粒子在无限均匀介质内的N条飞行轨迹，每条飞行轨迹跟踪到粒子权重低于预设统计下限为止，其中，所述无限均匀介质的宏观散射截面和宏观吸收截面均近似于所述屏蔽体的整体平均截面；基于每条标准轨迹生成对应的、经过所述目标区域的样本轨迹；基于每条标准轨迹与每条样本轨迹的对应关系确定每条所述样本轨迹在所述目标区域的权重w_j；根据事件概率分布函数来计算每条所述样本轨迹占有的概率宽度c_j；在考虑每条所述样本轨迹在所述目标区域的权重w_j的情况下对全部所述样本轨迹占有的概率宽度c_j加权求和而得到粒子到达目标区域的概率值p，公式如下：

其中，M是经过所述目标区域的所述样本轨迹的总数量。

图3示出了本发明提出的蒙特卡罗粒子输运计算方法的计算流程。与传统蒙特卡罗方法的不同之处在于，在本发明的蒙特卡罗粒子输运计算方法中，粒子的样本轨迹不是按照事件概率分布函数来采样得到，而是先根据最终统计目标的需求生成样本轨迹，然后再根据事件概率分布函数来计算每条样本所占有的概率宽度，对所有样本轨迹的概率宽度加权求和，即可得到粒子到达目标区域的概率值。在本发明的蒙特卡罗粒子输运计算方法中，所生成的每一个样本轨迹均对目标区域的计数有贡献，所有样本轨迹均为有效样本轨迹。

根据本发明的蒙特卡罗粒子输运计算方法的一种优选实施形式，样本轨迹采用对标准轨迹进行变换操作的方式得到，其中，所述变换操作使每条标准轨迹成为经过所述目标区域的所述样本轨迹。优选的是，所述变换操作包括对所述标准轨迹的旋转、伸缩。

根据本发明的蒙特卡罗粒子输运计算方法的一种优选实施形式，对所述标准轨迹的旋转是使该标准轨迹保持自身轨迹形状不变绕点源进行旋转，直至穿过目标区域，其中，在下列两种情况下舍弃所述标准轨迹，不对所述标准轨迹进行旋转：(i)如果经过旋转得到的穿过目标区域的轨迹的所有段并未都在所述屏蔽体内；或者(ii)如果所述标准轨迹的总长度短于所述点源至所述目标区域的距离。有的标准轨迹的总长度短于所述点源至所述目标区域的距离的原因是，该标准轨迹在无限均匀介质中较早衰减完毕，总位移还没有达到时其权重已经非常低，低于预设的下限。

根据本发明的蒙特卡罗粒子输运计算方法的一种优选实施形式，在所述标准轨迹有多个旋转角度可穿过目标区域的情况下，则将经过这些旋转角度得到的所有未被舍弃的样本轨迹都用于上述的根据事件概率分布函数来计算每条样本轨迹占有的概率宽度的步骤，其中，

每条样本轨迹的概率宽度为

其中，s是目标区域的面积，R是目标位置与粒子源的距离。

根据本发明的蒙特卡罗粒子输运计算方法的一种优选实施形式，在上述的根据事件概率分布函数来计算每条样本轨迹占有的概率宽度的步骤中，包括：计算所述屏蔽体中每条样本轨迹占有的概率宽度与所述无限均匀介质材料中的每条标准轨迹的概率宽度的比值c，以及计算上述比值c与每条标准轨迹的概率宽度1/N的乘积作为每条样本轨迹占有的概率宽度，其中，N为所述点源发出的粒子的所述标准轨迹的数量。优选的是，通过下列公式计算所述屏蔽体中每条样本轨迹占有的概率宽度与所述无限均匀介质材料中的每条标准轨迹的概率宽度的比值c：

其中，k为该条样本轨迹的总段数，其中，第i段轨迹为第i-1次和第i次碰撞的碰撞点之间的线段；∑_s1是第i段末端碰撞点处所述屏蔽体材料的散射截面；∑_s0是所述无限均匀介质材料在该点的散射截面；D_Ω1是屏蔽体材料中该点散射后飞向第i+1段对应方向Ω的概率密度；D_Ω0是所述无限均匀介质材料中散射向Ω方向的概率密度；∑_t1是第i段轨迹所经过的所述屏蔽体材料的总截面，L_i为第i段轨迹的累计飞行路程；∑_t0是所述无限均匀介质材料内的总截面。

在第i段轨迹经过了多种不同材料的情况下，则∑_t1·L_i为所经过的每一段屏蔽体材料的总截面∑_t与在该段屏蔽体材料中的累计飞行路程L的乘积∑_t·L之总和。

附图说明

下面结合附图阐释本发明的实施例。在附图中：

图1是深穿透屏蔽计算问题示意图

图2是传统蒙特卡罗方法的计算流程

图3是根据本发明的蒙特卡罗粒子输运计算方法的计算流程

图4是用根据本发明的蒙特卡罗粒子输运计算方法求解点通量问题时的旋转操作的示意图。

具体实施方式

下文结合附图详细解释本发明的蒙特卡罗粒子输运计算方法的主要步骤。

步骤1：生成有效轨迹样本

理论上，可以有很多种方法生成从源粒子到目标区域的样本轨迹，只要是根据确定论方法可以方便直接计算得出变换后的概率宽度的，都是可行的变换，例如进行镜面对称操作等。但是为了便于在后续步骤中计算每条样本所占有的采样概率宽度，在下文对于根据本发明的蒙特卡罗粒子输运计算方法的介绍中采用从标准轨迹进行旋转、拉伸变换的方式得到有效轨迹。

首先，采用传统蒙特卡罗方法生成点源在无限均匀介质内的N条飞行轨迹，每条轨迹跟踪到粒子权重低于预设的统计下限为止。这样每条轨迹在该均匀介质模型下的概率宽度均为1/N。这些轨迹我们称之为标准轨迹。在此应当理解，上述无限均匀介质的宏观散射截面和宏观吸收截面与屏蔽体的整体平均截面大致接近即可。对粒子的吸收采用改进了蒙特卡罗计算方法中的常规隐俘获的方式进行处理，即采用散射截面Σs来代替宏观碰撞截面Σt，平均碰撞自由程为1/Σs，粒子的权重为e^-Σa·L，其中Σa为介质的宏观吸收截面，L为源粒子的累计飞行路程。

然后保持上述N条标准轨迹的形状和位置不变，将上述无限均匀介质材料替换成真实屏蔽问题中的材料分布。此时由于材料截面发生变化，每条轨迹占有的概率宽度也会随之发生变化，需要重新计算每条轨迹的概率宽度，具体计算方法在步骤2中介绍。

为了使每条轨迹都能够到达目标区域(例如图一中的探测器位置I点)，还需要对这N条轨迹进行旋转、拉伸等变换操作，使得每条经过变换操作的轨迹都能够经过目标区域。由于进行了变换操作，每条轨迹的概率宽度再次发生变化，需要重新计算。

步骤2：计算每条样本所占有的概率宽度

假设用来产生标准轨迹的均匀介质截面为∑₀，真实屏蔽材料内的截面为∑₁(随位置有变化)。由于已知N条标准轨迹样本中每条样本的概率宽度均为1/N，当屏蔽材料截面发生变化后，只需要计算真实屏蔽材料中每条样本轨迹占有的概率宽度与所述无限均匀介质材料中的每条标准轨迹的概率宽度的比值，即可得到真实屏蔽材料下每条样本轨迹的概率宽度。

经过理论推导，得到真实屏蔽材料中每条样本轨迹概率宽度与无限均匀介质材料中的每条标准轨迹的概率宽度比值

其中，k为该条轨迹的总段数，其中，每两个相邻碰撞点之间的线段为一段，第i段轨迹为第i-1次和第i次碰撞的碰撞点之间的线段；∑_s1是第i段末端碰撞点处真实材料的散射截面，∑_s1是原均匀材料在该点的散射截面；D_Ω1是真实材料中该点散射后飞向第i+1段对应方向Ω的概率密度，D_Ω0是原均匀材料中散射向Ω方向的概率密度；∑_t1是第i段轨迹所经过的真实材料的总截面(如果第i段轨迹经过了多种不同材料，则∑_t1.L_i应变为所经过的每一段材料的∑_t.L之和)，∑_t0是原均匀材料内的总截面。

这样，真实屏蔽材料下每条样本轨迹的概率宽度即可通过上述计算公式得到，结果为c/N。上述理论方法已经通过数值测试实例进行了验证。使用某均匀材料截面产生标准轨迹，然后直接将轨迹用于由三种不同材料组合形成的屏蔽体，每条轨迹的概率宽度采用上述计算公式处理后，得到的屏蔽计算结果与直接计算组合屏蔽体的结果完全一致。

如果需要计算屏蔽体某个目标区域点的点通量(即，上述的目标区域现在变为一个点)，则还可以通过将标准轨迹进行旋转、拉伸等几何变换操作，使每条样本轨迹均到达目标区域点。进行这些操作时，轨迹的概率宽度也会随之发生变化，仍然可以基于上述理论思想进行计算。以旋转操作为例，如果希望每条标准轨迹都能到达目标点区域，则可将无限均匀介质内的每一条标准轨迹保持其自身轨迹形状完全不变(即视每条轨迹为一个刚体)，以点源为中心进行旋转，直到该轨迹穿过目标区域(如果该轨迹有多个不同旋转角度可穿过目标区域，则进行多次旋转；如果该轨迹总长度比目标位置和源之间的距离短，则舍弃该轨迹)。这样我们就获得了N1条左右从源粒子处到达探测器位置的样本轨迹，这些轨迹的概率宽度需要进行相应的调整，每条轨迹的概率宽度变为

其中s是目标区域的面积，R是目标位置与粒子源的距离。图4是用根据本发明的蒙特卡罗粒子输运计算方法求解点通量问题时的旋转操作的示意图，S为放射源，I为探测器，阴影区域为屏蔽体。

得到这些经过目标区域的样本轨迹后，将均匀材料替换为真实屏蔽材料，即可采用前面公式同样的方法计算得到真实屏蔽材料内这些样本轨迹的概率宽度。

步骤3：计算源粒子到达目标区域的概率

有了步骤1和步骤2的标准轨迹和概率宽度计算公式后，可以很容易地对源粒子到达目标区域的概率进行计算。首先采用传统方法计算得到点源在均匀介质内N条标准轨迹样本中每一条在各个位置的粒子权重w(由于隐俘获操作导致的权重)，然后将均匀介质替换为真实屏蔽材料，采用步骤2的方法计算每条从源粒子点到达目标区域的轨迹所对应的概率宽度c_j，基于每条标准轨迹与每条样本轨迹的对应关系确定每条所述样本轨迹在所述目标区域的权重w_j，然后将每条经过目标区域的轨迹在目标点的权重w_j与该轨迹的概率宽度c_j相乘，求和即得到源粒子到达目标区域的概率

其中M是经过目标区域的轨迹总数量。

下面以3个不同的具体问题为例来说明本方法的具体实施方式。

屏蔽点通量计算问题

以图1中的点通量计算问题为例。假设左侧源S为中子源，中间的以阴影表示的屏蔽体为混凝土墙，需要求右侧I点的中子通量。I点到S点的距离为R。

首先采用传统蒙特卡罗方法建立点中子源S在无限均匀混凝土材料中的计算模型，中子吸收截面采用步骤1中的隐俘获方式处理。从S处产生了100000个初始中子并进行跟踪，最终有5000条轨迹到达了以S为中心R为半径的球面之外(其它轨迹在到达球面之前权重已降低至预设下限)，那么这5000条轨迹可以作为我们后续计算所使用的标准轨迹，每条轨迹的概率宽度为10^-5。每条轨迹在球面位置处的权重为w_j。

然后用真实的屏蔽体几何形状替代原来的无限均匀混凝土，按照I点与S点的相对位置将这5000条轨迹进行旋转，使每条轨迹从S点出发，最终都穿过I点。同时判断每条轨迹是否所有段都在真实屏蔽体空间内部，穿出屏蔽体的轨迹被舍弃，剩余的全部线段都在屏蔽体空间内的有3000条。

针对这3000条轨迹求和得粒子穿透概率

再进行I点位置的角度处理，即可得到I点的中子通量。

如果屏蔽体不是同一种均匀材料，而是多种材料的组合，则前期计算也是采用上述过程，只是最后一步需要根据材料截面差别根据步骤2中的公式来计算每条轨迹的概率宽度c_j。对于屏蔽体内存在空腔的情况，则可以直接在轨迹旋转至空腔位置时进行相应的拉伸即可。

精细通量分布计算问题

同样以图1中的模型为例。假设左侧源S为中子源，中间的屏蔽体为混凝土墙，需要求屏蔽体右侧表面的中子通量精细分布。S到屏蔽体右侧表面的垂直距离为D。

前面的过程与上面计算点通量的过程相同，也是先构造点源S在无限均匀介质内的100000条轨迹，最终有8000条轨迹到达了以S为中心D为半径的球面之外。

将这8000条轨迹的每一条以S点为中心在4π立体角内进行旋转，当该轨迹扫过右侧表面某个位置点y时，即可通过计算w_j.c_j得到这条轨迹对y点的通量值贡献。这样将8000条轨迹操作完毕，将每条轨迹对y点的通量值贡献相加，即可得到y点的通量，同理可以得到右侧表面上任意一点的通量值。

微扰计算问题

假设需要计算某个反应堆堆芯的等温温度系数。现在采用传统蒙特卡罗方法进行临界计算，得到温度T1下的堆芯有效增殖系数(k_eff)k1，统计标准差为w1；继续计算温度T2下的堆芯有效增殖系数k2，统计标准差为w2。如果采用(k2-k1)/(T2-T1)来直接计算温度系数，则由于统计涨落的存在，而k2-k1的值很小，因此k2与k1的细微差别容易淹没在蒙特卡罗统计方差里，导致结果产生很大的计算误差。

采用本发明中的蒙特卡罗粒子输运计算方法则可以保持两次计算方差的一致性，从而能够采用上面的直接法来计算温度系数。

具体步骤是：首先采用传统蒙特卡罗方法建立温度T1下的堆芯模型，产生粒子轨迹进行计算模拟；然后在跟踪每一条粒子轨迹时，同时采用类似本发明步骤2中的公式来计算T2温度下(截面发生变化)该轨迹的概率宽度变化系数c_j：

其中，k为该条样本轨迹的总段数，其中，第i段轨迹为第i-1次和第i次碰撞的碰撞点之间的线段；∑_s2是第i段末端碰撞点处T2温度下的散射截面；∑_s1是T1温度下在该点的散射截面；D_Ω2是T2温度下该点散射后飞向第i+1段对应方向Ω的概率密度；D_Ω1是T1温度下散射向Ω方向的概率密度；∑_t2是第i段轨迹所经过的T2温度下的总截面，L_i为第i段轨迹的累计飞行路程；∑_t1是T1温度下的总截面。

在计算原模型堆芯有效增殖系数的同时，采用概率宽度系数c_j修正后的轨迹权重来计算T2温度下的堆芯有效增殖系数k2，最后将得到的堆芯有效增殖系数之差k2-k1除以温度差T2-T1，即可得到准确的温度系数k_eff。

可以用根据本发明的蒙特卡罗粒子输运计算方法来直接计算温度系数、密度系数、控制棒微分价值等各种微扰问题。对于控制棒微分价值的计算过程，仅需将上述公式中的∑_s1、∑_s2、∑_t1、∑_t2等替换为分别表示控制棒动作前和控制棒动作后的截面。

以上记载了本发明的优选实施例，但是本发明的精神和范围不限于这里所公开的具体内容。本领域技术人员能够根据本发明的教导而做出更多的实施方式和应用，包括对上述实施例的各个技术特征的任意组合。这些实施方式和应用都在本发明的精神和范围内。本发明的精神和范围不由具体实施例来限定，而由权利要求来限定。

Claims (11)

1.一种蒙特卡罗粒子输运计算方法，用于计算点源发出的粒子经过屏蔽体后到达目标区域的概率值，该蒙特卡罗粒子输运计算方法包括如下步骤：

计算点源发出的粒子的每一条标准轨迹在各个有关位置的粒子权重，其中，所述标准轨迹为所述点源发出的粒子在无限均匀介质内的N条飞行轨迹，每条飞行轨迹跟踪到粒子权重低于预设统计下限为止，其中，所述无限均匀介质的宏观散射截面和宏观吸收截面均近似于所述屏蔽体的整体平均截面；

基于每条标准轨迹进行变换生成对应的、经过所述目标区域的样本轨迹；

基于每条标准轨迹与每条样本轨迹的对应关系确定每条所述样本轨迹在所述目标区域的权重w_j；

根据事件概率分布函数来计算每条所述样本轨迹占有的概率宽度c_j；

在考虑每条所述样本轨迹在所述目标区域的权重w_j的情况下对全部所述样本轨迹占有的概率宽度c_j加权求和而得到粒子到达目标区域的概率值p，公式如下：

其中，M是经过所述目标区域的所述样本轨迹的总数量。

2.根据权利要求1所述的蒙特卡罗粒子输运计算方法，其特征在于，

所述样本轨迹采用对标准轨迹进行变换操作的方式得到，其中，所述变换操作使每条标准轨迹成为经过所述目标区域的所述样本轨迹。

3.根据权利要求2所述的蒙特卡罗粒子输运计算方法，其特征在于，

所述变换操作包括对所述标准轨迹的旋转、拉伸，

其中，对所述标准轨迹的拉伸是指在屏蔽体内存在空腔的情况下，使所述标准轨迹当前变换得到的样本轨迹的长度直接伸长一个该空腔在该样本轨迹方向上的宽度。

4.根据权利要求3所述的蒙特卡罗粒子输运计算方法，其特征在于，

对所述标准轨迹的旋转是使该标准轨迹保持自身轨迹形状不变绕点源进行旋转，直至穿过目标区域，

其中，在下列两种情况下舍弃所述标准轨迹，不对所述标准轨迹进行旋转：

(i)如果经过旋转得到的穿过目标区域的轨迹的所有段并未都在所述屏蔽体内；或者

(ii)如果所述标准轨迹的总长度短于所述点源至所述目标区域的距离。

5.根据权利要求4所述的蒙特卡罗粒子输运计算方法，其特征在于，

在所述标准轨迹有多个旋转角度可穿过目标区域的情况下，则将经过这些旋转角度得到的所有未被舍弃的样本轨迹都用于上述的根据事件概率分布函数来计算每条样本轨迹占有的概率宽度的步骤，其中，

每条样本轨迹的概率宽度为

其中，s是目标区域的面积，R是目标位置与粒子源的距离。

6.根据权利要求2至5中任一项所述的蒙特卡罗粒子输运计算方法，其特征在于，

在根据事件概率分布函数来计算每条样本轨迹占有的概率宽度的步骤中，包括：

计算所述屏蔽体中每条样本轨迹占有的概率宽度与所述无限均匀介质材料中的每条标准轨迹的概率宽度的比值c，以及

计算上述比值c与每条标准轨迹的概率宽度1/N的乘积作为每条样本轨迹占有的概率宽度，其中，N为所述点源发出的粒子的所述标准轨迹的数量。

7.根据权利要求6所述的蒙特卡罗粒子输运计算方法，其特征在于，

通过下列公式计算所述屏蔽体中每条样本轨迹占有的概率宽度与所述无限均匀介质材料中的每条标准轨迹的概率宽度的比值c：

其中，k为该条样本轨迹的总段数，其中，第i段轨迹为第i-1次和第i次碰撞的碰撞点之间的线段；∑_s1是第i段末端碰撞点处所述屏蔽体材料的散射截面；∑_s0是所述无限均匀介质材料在该点的散射截面；D_Ω1是屏蔽体材料中该点散射后飞向第i+1段对应方向Ω的概率密度；D_Ω0是所述无限均匀介质材料中散射向Ω方向的概率密度；∑_t1是第i段轨迹所经过的所述屏蔽体材料的总截面，L_i为第i段轨迹的累计飞行路程；∑_t0是所述无限均匀介质材料内的总截面，在第i段轨迹经过了多种不同材料的情况下，则∑_t1·L_i为所经过的每一段屏蔽体材料的总截面∑_t与在该段屏蔽体材料中的累计飞行路程L的乘积∑t·L之总和。

8.根据权利要求1至7中任一项所述的蒙特卡罗粒子输运计算方法，其特征在于，对粒子的吸收采用这样的隐俘获方式进行处理，即，采用散射截面Σs来代替宏观碰撞截面Σt，平均碰撞自由程为1/Σs，粒子的权重为e^-Σa·L，其中，Σa为介质的宏观吸收截面，L为源粒子的累计飞行路程。

9.一种屏蔽点通量计算方法，包括：

根据权利要求1至8中任一项所述的蒙特卡罗粒子输运计算方法计算粒子到达目标区域点的概率值，

对所得到的概率值进行相应于目标区域点的角度处理，即可得到该点的中子通量。

10.一种精细通量分布计算方法，用于计算屏蔽体背离中子源一侧的表面或屏蔽体内的一个区域内的精细通量分布，包括：

根据权利要求1至8中任一项所述的蒙特卡罗粒子输运计算方法计算每条样本轨迹以源点为中心在4π立体角内旋转到所述屏蔽体背离中子源一侧的表面内或屏蔽体内的所述区域内的一点时的概率值p，并进而经过角度处理后得到该条样本轨迹对于该点的通量值贡献，累加所有样本轨迹对该点的通量贡献值得到该点的通量值，

重复上述步骤计算所述屏蔽体背离中子源一侧的表面内或所述屏蔽体内的所述区域内的其它各点的通量值而得出所述屏蔽体背离中子源一侧的表面内的精细通量分布。

11.一种反应堆堆芯微扰计算的方法，包括：

采用根据权利要求1至8中任一项所述的蒙特卡罗粒子输运计算方法建立温度T1下的堆芯模型，产生粒子轨迹进行临界计算模拟，得到温度T1下的堆芯有效增殖系数k1；

对于每一条粒子轨迹，采用下列公式来计算T2温度下相对于T1温度下该轨迹的概率宽度变化系数c_j：

其中，k为该条样本轨迹的总段数，其中，第i段轨迹为第i-1次和第i次碰撞的碰撞点之间的线段；∑_s2是第i段末端碰撞点处T2温度下的散射截面；∑_s1是T1温度下在该点的散射截面；D_Ω2是T2温度下该点散射后飞向第i+1段对应方向Ω的概率密度；D_Ω1是T1温度下散射向Ω方向的概率密度；∑_t2是第i段轨迹所经过的T2温度下的总截面，L_i为第i段轨迹的累计飞行路程；∑_t1是T1温度下的总截面；

对于每一条粒子轨迹，采用系数c_j修正后的轨迹权重来计算T2温度下的堆芯有效增殖系数，得到温度T2下的堆芯有效增殖系数k2；

将得到的堆芯有效增殖系数之差k2-k1除以温度差T2-T1，得到温度系数。

Images