《1前言》

1前言

硼中子俘获治疗(boron neutron capture therapy, BNCT)是一种针对脑部神经胶质瘤和黑色素瘤等恶性肿瘤较为理想的治疗方法。该方法的优点在于,肿瘤细胞中引入含有10B 的化合物后[1],利用低能中子束照射头部,形成低能中子对富集10B的肿瘤细胞的靶向照射,通过10B( n, α )7 Li 反应,放出 α 粒子、7 Li以及足以杀死肿瘤细胞的 2. 79 MeV能量。从而达到最大限度地杀死肿瘤细胞、保护正常组织的目的。

医院中子照射器 I 型堆(IHNI -1) B堆芯设计有 BNCT热中子束流孔道[2,3],治疗时利用从热中子束流孔道引出的中子束照射于头部病变部位,几何示意图见图 1。热中子束流孔道的模拟计算属于粒子深穿透问题,需要特殊的抽样技巧,如果从堆芯开始跟踪粒子、模拟计算人体头颅等效模型内的剂量分布将耗费大量机时,且计算结果方差较大;而从 BNCT 热中子束流孔道出口处跟踪粒子、模拟计算人体头颅等效模型内的剂量分布只需较少机时,其收敛速度较快、计算结果方差较小,因此,利用 MCNP / 4B程序在 BNCT 束流孔道出口处建立等效平面源的方案,为人体头颅等效模型剂量分布的快速计算提供了一个较为可靠的中子、平面源。

《图1》

图1 IHNI -1 堆 B 堆芯热中子束流孔道几何结构示意图

Fig. 1 Transverse sectional chart of the thermal neutron duct of IHNI-1 reactor with B-core

《2 IHNI -1 堆热中子束流孔道模拟计算》

2 IHNI -1 堆热中子束流孔道模拟计算

热中子束流孔道的数值计算采用蒙特卡罗耦合抽样技巧和分裂与轮盘赌技巧(或权窗游戏),计算模式为中子和耦合输运,中子和粒子抽样重要性(IMP : N ,P)在堆芯内置 1,沿束流孔道方向依次增大。堆芯核功率为 30 kW,堆芯单位时间内产生,堆芯单位时间内产生的中子数为每秒 2.279 48 × 1015

堆芯临界计算是中子束流孔道参数计算的前提,全堆的有效倍增系数 不同,孔道出口的束流参数也将不同,因此,需要进行临界搜索计算,得到 =1的堆芯装载模式,在此基础上进行的束流孔道模拟计算结果才可靠。

《2.1 堆芯临界搜索计算》

2.1 堆芯临界搜索计算

B堆芯由中心控制棒栅元、燃料栅元、水栅元组成,共计 11 圈,340 个栅位。冷态临界搜索计算时,根据临界计算结果,调整最外圈燃料栅元和水栅元的数量。

IHNI -1 堆B堆芯冷态临界搜索计算采用了 3 种堆芯布置,计算结果见表 1。由表 1 可知,堆芯装载 319.16 根燃料元件、235U 装载量为1 176.626 9 g时,堆芯达到冷态临界。

《表1》

表1 3 种堆芯布置下有效倍增系数和冷态临界235U装料量蒙特卡罗模拟计算

Table 1 The Monte Carlo calculation results of  and load of 235U for 3 core schemes

《2.2 热中子束流孔道参数计算》

2.2 热中子束流孔道参数计算

采用蒙特卡罗耦合抽样法进行了热中子束流孔道模拟计算,图 2 给出了热中子束流孔道轴线上不同位置处中子通量密度的变化曲线。由图 2 可知,堆芯核功率为 30 kW 时,热中子通量密度在铍反射层内达到最大值 1.599 99 × 1012 cm-2 · s-1,该点到堆芯中心点的距离为 12.5 cm 。同时,热中子束流孔道出口中心处的热中子通量密度达到 2. 377 95 × 109 cm-2 · s-1

《图2》

图2 热中子束流孔道中心轴线上中子通量密度衰减曲线

Fig.2 Attenuation curve of neutron flux density along axial line of thermal neutron duct 

《3中子、等效平面源计算模型》

3中子、等效平面源计算模型

理论上,置于反应堆 BNCT 中子源外的人体头颅等效模型内的通量密度分布可进行一次蒙特卡罗中子—耦合输运计算得到,但该方法计算速度慢、计算精度较差,为了加快计算速度,提高计算精度,将上述一次计算过程划分为两步:a. 采用蒙特卡罗耦合抽样方法计算反应堆中子束流孔道出口处垂直于束流孔道轴线的整个平面上的中子、平面源;b. 采用 MCNP /4B 程序分别计算中子、平面源在人体头颅等效模型内所产生的中子、通量密度。

上述中子、平面源 是空间、能量和方向的函数,其中,为距离,单位为 cm;为能量,单位为 MeV;μ 为散射角余弦,数值范围 -1~1。由于 MCNP /4B 程序在平面源问题的输运计算时,只需考虑 μ>0的源项,为了简化平面源的计算模型,设正向平面源

式(1)中,)为粒子正向流密度的空间分布;为与空间有关的粒子正向流密度的能谱分布;为与空间有关的粒子正向流密度的角分布;为归一化系数。

为了得到平面源的空间分布、能量分布和角分布,采用 MCNP /4B 程序中相应的分段计数卡(FSn)、分段除数卡(SDn)、计数能量卡(En)和计数余弦卡(Cn)对平面源的空间、能量和角度进行网格划分。

针对束流孔径为 10 cm和 14 cm的孔道设计方案,将平面源的空间变量 分别划分为 5 个网格,网格外径=1 , … ,5 )  分别为 5 cm、 10 cm、 20 cm、30 cm、 40 cm 和 7 cm、 10 cm、 20 cm、 30 cm、 40 cm。

等效中子平面源划分为 37 个能群,能群划分见表 2;等效 平面源划分为 18 个能群,能群划分见表 3;在平面源的角度上划分为 17 个网格,网格划分见表 4 。

《表2》

表2 等效中子平面源能群划分

Table 2 Energy groups of equivalent neutron surface source

《表3》

表3 等效平面源能群划分

Table 3 Energy groups of equivalent  surface source

《表4》

表4 等效平面源角分布网格划分

Table 4 Grids of angular distribution for equivalent surface source

采用蒙特卡罗耦合抽样方法[4]计算了各空间网格内的多群中子、正向流密度的能谱分布,并计算了穿过平面源表面上各空间网格内的正向流密度的角分布

进行方向归并,得到各空间网格上的正向流密度

式(2)即为平面源的空间分布。根据式(2),可得平面源的总源强

式(3)中,为各空间网格的面积,cm2为人体头颅剂量分布计算时的源强归一化常数。

由于 MCNP /4B 程序的输入卡片不能模拟上述一维空间几何各网格内的平面源,因此还需将此一维几何平面源转换为孔道设计的二维 xy 等效几何平面源,其平面源转换的等效几何见图 3,图 3 中一维同心圆平面和二维 x平面对应网格的面积相等。

《图3》

图3 平面源等效转换

Fig.3 Equivalent surface source mutation

《4中子、等效平面源计算结果》

4中子、等效平面源计算结果

采用等效平面源计算模型对B堆芯热中子束流孔道出口处的中子、等效平面源进行了模拟,计算了中子、正向流密度角分布、能谱分布等束流参数。

中子束流孔径分别采用10 cm 和14 cm。当中子束流孔径为10 cm时,等效中子源的归一常数为2. 282 246 47 ×1011 n/s,等效 源的归一常数为3. 334 001 535 ×1010 ;当中子束流孔径为14 cm时,等效中子源的归一常数为3. 461 319 5 ×1011 n/s,等效 源的归一常数为3. 273 953 14 ×1010 

图 4 和图 5 分别给出了中子束流孔径为 10 cm 时,热中子束流孔道出口处中子流密度的角分布和中子正向流密度的能谱分布。

《图4》

图4 热中子束流孔道出口处群中子流密度的角分布

Fig.4 Angular distribution of neutron flow density at exit of the thermal neutron duct

《图5》

图5 热中子束流孔道出口处正向中子流密度的能谱分布

Fig.5 Energy spectra distribution of neutron flow density at exit of the thermal neutron duct

《5结语》

5结语

采用蒙特卡罗程序 MCNP 对医院中子照射器 I 型堆(IHNI -1)B堆芯进行了临界搜索计算,在此基础上进行了热中子束流孔道模拟,计算了束流参数,采用等效平面源计算模型,建立了两种出口孔径的等效中子、平面源。BNCT 热中子束流孔道出口处平面源的建立,为人体头颅等效模型吸收剂量分布的快速计算提供了一个较为可靠的中子、平面源[5]