基于遗传算法的中子屏蔽材料组分优化研究

陈法国，李国栋，杨明明，韩 毅，梁润成

(中国辐射防护研究院，山西 太原 030006)

屏蔽设计本质上是一个带有约束条件的多目标寻优过程。一方面，屏蔽效果与射线类型和能谱分布密切相关；另一方面，在实际工程中，屏蔽设计还需要考虑重量、体积、成本、结构性能等多种设计限定条件。传统的屏蔽设计通常以设计者的个人经验为基础，优化设计效率相对较低，设计结果的全局最优性受人因等不确定因素影响较大。为此，遗传算法(GA)、人工神经网络算法(ANN)、特征统计算法(CSA)等多种优化计算方法[1-5]，被应用于屏蔽结构设计、复合屏蔽材料组分优化等问题求解中，以实现屏蔽设计方案的自动优化。

中子屏蔽设计需要考虑中子的慢化、吸收以及次级光子的屏蔽，需要同时采用快中子非弹性散射占优的重核材料、中能中子弹性散射截面大的轻核材料以及热中子俘获截面大的中子吸收材料，是一个典型的优化求解问题[6]。以中子复合屏蔽材料组分优化设计为例，以屏蔽剂量和材料密度最小化为目标，开展了遗传算法在多目标屏蔽优化设计中的应用研究。

1 中子复合屏蔽材料组分优化模型及遗传算法简介

1.1 中子复合屏蔽材料组分优化模型

利用优化算法自动设计中子复合屏蔽材料的组分比例，使其以最小的密度实现对特定能谱中子的最好屏蔽效果。中子源项设为单向面源，能谱分布为Watt裂变谱；中子平均能量为2.348 MeV，能谱分布如图1所示；复合材料由聚乙烯、铅、硼、锂、铁、铝6种常见材料组成。假设各组分均匀混合且混合前后总体积不变，用单位注量的中子穿过厚度为30 cm平板结构后的中子和次级光子剂量之和(定义为屏蔽剂量D)来表征屏蔽效果，屏蔽结构示意图如图2所示。

图1 入射中子的能谱分布Fig.1 Energy spectrum of incident neutron

图2 屏蔽结构示意图Fig.2 Schematic diagram of shielding structure

中子复合屏蔽材料组分优化的数学函数模型可表示为：

minF(X)=[D(X),ρ(X)]T

s.t. 0≤xi≤1，∑xi=1，1/ρ(X)=

∑xi/ρi，i=1,2, …6

(1)

式中，F(X)为两维设计目标，D(X)为组分为X屏蔽材料对应的屏蔽剂量，X=[x1,x2, …x6]分别为6种组分的质量分数，ρ为材料密度。

1.2 遗传算法简介

传统的屏蔽优化问题求解，采用加权法将多目标降维到单目标进行处理。加权法简化了问题求解，但子目标权重因子的设置具有一定的主观性，并且优化结果在理论上只是Pareto最优解集中的一个解，在权重因子变化后需要再次运行优化过程进行求解；而多目标优化算法可以在单次运行过程中，求解Pareto最优解集[7]。为此，采用遗传算法直接对式(1)中的双目标优化问题进行求解。

遗传算法是由Holland提出的模拟生物种群遗传和自然选择规律而演化的一种随机搜索方法。为提高优化效率、避免过早收敛，在基本遗传算法的基础上，发展了向量评估遗传算法(VEGA)、多目标遗传算法(MOGA)、非支配排序遗传算法(NSGA)、基于距离的Pareto遗传算法(DPGA)等[7]。其中，NSGA是学者Srinivas和Deb基于非支配排序思想提出的，通过分类排序可提高适应度较高的个体保留至下一代的概率，并采用适应度共享策略使个体分布更均匀，以保持种群的多样性，避免过早收敛；为优化NSGA的时间复杂度较高、未采用精英策略、优化过程和寻优效率对共享参数设置较为敏感的问题，Deb等人又提出了快速非支配排序遗传算法(NSGA-II)[8]。

2 屏蔽材料组分优化算法架构

利用NSGA-II建立式(1)数学模型所表示的中子复合屏蔽材料组分优化算法。优化算法主要包括屏蔽计算模块和遗传算法模块两部分，算法流程图如图3所示。遗传算法的主要构成要素有基本运行参数、染色体编码方式、个体适应度评价、遗传算子。

图3 中子屏蔽材料组分优化算法流程图Fig.3 Flow chart of optimization algorithm for neutron shielding material component

2.1 遗传算法参数设置

基本运行参数包括种群大小M(通常设为20～100)、终止进化代数T(通常设为50～200)、交叉因子Pc、变异因子Pm。交叉因子和变异因子决定了父代和子代个体的差异性，直接影响算法的收敛性和寻优性能；增加差异性可扩大搜索空间，但不利于收敛，反之亦然。

2.2 编码和解码

编码和解码是应用遗传算法的关键，用于建立优化问题的设计变量与遗传基因的映射关系。屏蔽材料组分优化是一个定和约束问题，即∑xi=1，0≤xi≤1；若直接用质量分数xi进行交叉或变异操作，子代种群中会出现不满足约束条件的个体。为此，采用球坐标转换方法，引入一个采用浮点数编码方式的5维向量θ=[θ1,θ2…θ5]T，0≤θi≤π/2，来表征设计变量X=[x1,x2…x6]T，其中：

xi=(sinθ1·sinθ2·…sinθ5)2，i=1

xi=(cosθi-1·sinθi·…sinθ5)2，i=2,3,4,5

(2)

xi=(cosθ5)2，i=6

在执行交叉或变异算子时，直接对θi进行操作，然后按照式(2)将其解码为质量分数xi，用于计算复合材料的密度和屏蔽剂量。通过球坐标转换，既能满足定和约束条件，又可保证交叉和变异操作的随机性。

2.3 屏蔽计算模块

遗传算法通过解码得到特定个体的材料组分后，调用屏蔽计算模块来计算对应的屏蔽剂量。在整个优化算法流程中，对每代种群中的新个体都需要调用一次屏蔽计算模块；若直接采用蒙特卡罗(MC)程序作为屏蔽计算模块，多次调用会增加时长。为了兼顾屏蔽计算精度与运行时间，屏蔽计算模块采用了基于BP(误差逆向传播算法)神经网络算法的屏蔽快速计算方法。

针对屏蔽优化问题的具体模型，建立可根据材料组分比例预测屏蔽剂量的BP神经网络，以便在优化过程中快速调用。首先，在整个搜索空间内对6种材料的组分比例进行均匀抽样，并根据抽样参数建立MC模型计算屏蔽后剂量D，形成初始训练样本。其次，利用训练样本对BP神经网络进行训练和学习，寻找组分比例与剂量之间的计算规律。最后，利用训练后的网络来预测任意组分比例所对应的屏蔽计算结果。

利用297个训练样本训练后，BP神经网络预测结果和MC计算结果的相对偏差在-5.2%～8.1%范围内；随机抽取33种新组分比例作为测试样本，神经网络预测结果与MC计算结果的相对偏差在-3.9%～5.1%范围内。在执行优化算法时，将个体编码进行解码后，可直接调用屏蔽快速计算模块，并将计算结果用于个体适应度评价。

2.4 个体适应度评价

适应度函数是用来衡量种群中的个体接近最优解的优良程度。NSGA-II算法中，利用非支配序rank和拥挤度两个属性来评价个体适应度。

对于两个变量个体X1和X2，若所有子目标函数满足fi(X1)≤fi(X2)，且存在子目标函数使fj(X1)

在对两个个体进行选择时，非支配序小的个体占优，非支配序相同则拥挤距离大的个体占优；非支配序小可保证被选择的个体较优，拥挤距离大可保证个体多样性，避免过早收敛于局部最优。

2.5 遗传算子

选择、交叉、变异是遗传算法中的基本算子。

NSGA-II中，选择算子采用锦标赛方法。从父代种群M(种群大小为N)中随机选择s1个个体，进行适应度评价，再从s1个个体中选择占优的s2个个体进入子代种群Ms(N>s1>s2)；重复上述过程，直到子代种群的个体数量与父代种群一致，然后按照设定的交叉因子和变异因子对子代种群中的个体进行交叉、变异操作。

交叉算子对两个父代个体中的编码向量进行部分重组，将其中的父代基因按照特定算法处理，生成两个新的子代个体；变异算子用一个基于特定概率分布的随机数替代父代个体中的某个基因位，从而产生新的子代个体。NSGA-II中采用模拟二进制的交叉算子(SBX)和多项式变异算子。为了防止种群中的较优个体因交叉或变异操作而丢失，NSGA-II引入了精英策略。在精英策略中，将父代种群与子代种群混合(种群大小为2N)，重新计算混合种群中所有个体的非支配序和拥挤度，并选择占优的N个个体作为新一代种群。

对每一代种群重复执行以上操作，当满足迭代条件后，将当前种群中的个体进行解码，就可输出屏蔽材料组分比例的优化设计结果。

3 优化结果及对比分析

3.1 优化结果

针对Watt谱分布的入射中子，在未屏蔽的条件下，单位注量对应的周围剂量当量为388.6 pSv·cm2。6种组分材料的密度取值列于表1，屏蔽计算时各组成元素的同位素设置为天然丰度。

利用Matlab软件编译了优化算法程序，并利用CPU主频为3.7 GHz的计算机运行。种群大小设为100，交叉因子设为20，变异因子设为20，以材料密度最小和屏蔽后剂量最小为目标，经过100代后中子复合屏蔽材料组分的优化方案分布如图4所示；图中显示了初始方案以及优化1代、10代、50代、100代后的方案(a包含了所有数据，b为局部放大图)。由图可知，初始方案在搜索区域内随机分布；随着迭代次数的增加，优化方案逐渐接近Pareto前沿区域，并且在该区域内较为均匀地分布。

表1 组分材料的密度(g/cm3)Tab.1 Densities of component materials

图4 中子屏蔽材料组分的优化结果Fig.4 Optimization results of neutron shielding material component

对于每一代种群，用该代种群中所有个体对应的密度和屏蔽后剂量的乘积之和，来表示种群个体与Pareto前沿的接近程度，即种群优良指数。如图5所示，随着进化代数的不断增加，种群优良指数逐渐降低；当进化代数超过40代时，优良指数最终趋于稳定。

图5 种群优良指数随进化代数的变化Fig.5 Figure of merit for population against evolutional generation number

3.2 结果对比分析

在最终优化方案的100个个体中，剂量最小方案、密度最小方案、典型中间最优方案的组分比列于表2。由表可知，密度最小方案中，复合材料的主要组分几乎全部为锂；剂量最小方案中，复合材料的主要组分是86.5%的聚乙烯、10.4%的天然硼以及少量的铝和铅；对于一系列中间方案，随着复合材料密度的增大，锂的含量逐渐减少，聚乙烯的含量逐渐增加，同时含有少量的硼。

Watt谱分布中子能量主要在2 MeV以下，屏蔽材料的有效原子量越小，慢化作用越好。6种材料中锂的密度最低、有效原子量最小且6Li与热中子的相互作用截面大，单位质量的屏蔽效果最好；对于密度最小方案的理想情况，优化结果几乎全部为锂，符合理论预期。剂量最小方案更符合实际情况，并且实际应用经验中含硼聚乙烯是对该能谱中子屏蔽效果最好的材料，与优化结果一致。上述定性结果说明了该屏蔽材料组分优化算法的可行性。

在得出优化方案后，根据优化的材料组分进行MC建模计算，并将屏蔽计算结果与优化结果进行对比。在最终的100个可行解中，遗传算法优化结果与MC计算结果的相对偏差在-14.9%～4.8%范围内。

4 讨论和结语

4.1 讨论

1)引入基于BP神经网络的快速屏蔽计算模块，避免了优化过程中对MC计算程序的大量调用，节省了计算时间。计算时间包括297个初始样本的MC计算时间、神经网络训练时间和遗传算法程序运行时间，分别为127.1 min、10.6 s和71.6 s。若直接调用MC程序，100个个体的种群进化100代，需要运行10 100次MC程序，则整个优化过程至少耗时72.4 h。但神经网络训练样本不足会导致预测结果与MC计算值的相对偏差增大，甚至误导遗传算法的寻优方向而陷入局部最优。最终的优化解集和训练样本相比，屏蔽后剂量的BP神经网络预测结果与MC计算结果的相对偏差从-5.2%～8.1%变到-14.9%～4.8%，整体呈现变大趋势。BP神经网络算法的初始训练样本在整个搜索空间内随机抽样得到，在整个变量区域内的分布是相对均匀的；随着种群的不断进化，优化解集逐渐向Pareto前沿区域靠近，初始训练样本在该区域的代表性不足，基于初始训练样本的神经网络预测结果的不确定度就会增大。为了降低神经网络预测的不确定度对优化结果的影响，可采用分阶段优化或部分调用MC计算程序的处理办法。

2)在适应度评价中引入不确定度因子，可提高种群的多样性，避免因屏蔽计算误差导致的过早收敛。对于屏蔽计算，不论是MC方法、点核积分方法，还是基于BP神经网络算法的快速预测，都存在不可忽视的计算误差。对个体适应度进行评价时，若直接以屏蔽后剂量的大小为判据，可能造成对个体优劣的误判，从而误导寻优过程；在算法调试时，多次独立执行优化程序，就出现了剂量最小方案中的组分比例随机陷入不同局部最优解的情况。为此，在适应度评价中引入了一个剂量不确定因子δD，若两个个体X1和X2对应的剂量子目标存在如下关系：

|D(X1)-D(X2)|/(D(X1)+D(X2))≤δD

则认为X1和X2剂量子目标的适应度相等；根据屏蔽计算误差设置不确定因子δD，在本工作中为计算误差的二分之一。通过引入不确定因子，在选择算子中对个体基于锦标赛进行选择时，可避免因屏蔽计算误差造成的个体优劣误判，从而保证种群的多样性。

3)NSGA-II算法在求解双目标和三目标问题时，具有很好的优化效率。但实际应用中，屏蔽设计通常还需要考虑成本、体积、材料性能等多种因素，优化子目标可能超过三个。NSGA-II算法在较多目标求解中，存在以下问题：①随着子目标数量增加，一个随机种群中非支配解的数量占比将成指数增长，经过几代进化后会很快占据整个种群，使搜索过程变慢；②随着设计目标维度增加，计算时间成本将大大增加；③对于三维以下的设计目标，可以借助于优化解集分布前沿的图形化表示，辅助对优化结果的评价和决策；但设计目标维度超过三个后，不容易用图形直观表示。对于三目标以上的优化问题，可以考虑参考自适应加权法，将多维设计目标降维后进行求解；也可以采用Deb等人最新提出的基于参考点的非支配遗传算法[9]。

4.2 结语

利用NSGA-II遗传算法，开展了屏蔽剂量和材料密度最小化的双目标中子复合材料组分优化研究，以Watt裂变中子谱以及聚乙烯、铅、硼、锂、铁、铝6种常见材料组成30 cm厚平板结构为例，进行了自动优化设计；优化结果定性说明了基于算法的屏蔽材料组分优化方法的可行性。本工作目前只关注遗传算法在多目标屏蔽设计中的应用，暂未关注屏蔽结构长期使用后的中子活化对屏蔽效果的影响及深穿透等问题，这些问题可在后续工作中继续关注。以设计者的经验为先验条件，结合该优化算法，可更好地进行辐射屏蔽设计，为更高效确定优化设计方案提供现实指导。