γ辐照装置多目标优化排源算法研究及其应用

杨 磊，刘 燕，罗志平，李文革，周毅吉，张美琴

(1.烟台大学，山东 烟台 264005；2.湖南省农业科学院，湖南 长沙 410300)

60Co γ辐照加工技术在农业、工业、医疗卫生、食品、环境保护等多个领域应用十分广泛[1],目前总体发展平稳，新建装置呈现大型化、自动化及泛用化等特点[2]。因60Co源半衰期较短，为维持和扩大加工能力，需定期补充和更换60Co源棒，而源价格昂贵，为此提高源的整体利用率十分重要。生产实践表明，优化60Co辐照装置排源设计是提高源利用率的核心技术因素，对优化装置剂量学性能、改善货物加工剂量、降低单位成本、提高装置加工量和经济效益等均有显著影响。

排源方法主要包括人工经验排源和计算机自动排源两大类，行业公司仍以使用前者为多，此类方法对工作人员专业素质要求高，人工和时间花费多[3]，研究快速、可靠、结果方案足够优化的自动排源方法成为必然需要。目前国内外有基于贪心算法、神经网络算法(深度学习算法)等的排源方法，贪心算法快速有效但全局搜索能力不够[4]。李磊等[5]结合经验方法提出了一种混合求解方法，提高了全局搜索水平且排源效果较好。神经网络算法潜力巨大，但其难点在于需搜集足够多的优秀排源方案对算法进行预训练[6]。文献[7-8]提出了一种全局随机搜索型方法，针对辐照装置剂量不均匀度(DUR)进行优化排源设计并证实其可行有效，文献[9-10]在此基础上对相关算法进行并行化改进并利用GPU加速，可满足任意装机规模辐照装置的排源工作。近年来通过与多家单位的交流合作发现，企业不仅对DUR，而且对射线能量利用率、源棒排布工作量等更多目标均有需求，因此本文对原方法的评价数学模型、求解算法及代码均进行改进，形成一种新型多目标优化排源方法。

1 排源评价数学模型

影响60Co γ辐照装置性能的主要因素及相互关系如图1所示[11-13]，优化选源、排源是保证和提高60Co γ辐照装置加工质量和加工效率的关键，因目前相关单位多直接购买成品源棒，源的活度分布选择余地小，为此排源成为工作核心。相关单位在评判排源设计方案时会综合考虑DUR、射线能量利用率及排布工作量等，而主流评价方法——参考平面法[14-15]难以直接评判货箱中剂量分布DUR和射线能量利用率。本文从真实辐照情景出发，通过计算货箱在辐照室中完成动态辐照后的累积剂量分布，可直接得到上述参量值。

1.1 辐射场分布和货箱累积剂量分布的计算方法

(1)

图1 影响60Co γ辐照装置性能的主要因素及相互关系Fig.1 Main factor affecting 60Co gamma irradiation facility and interrelation

常数，对60Co，Γ=2.503×10-18C·m2·kg-1·Bq-1·s-1；A为源棒活度，Bq；L为源棒长度，m；a为待计算空间点到线源垂直距离的绝对值，m；b为垂足到线源最近点的距离，m；C为校准因子(测量值/理论值)，以修订源棒包壳、剂量计、辐照装置附属机构(如导源管、护源板(罩))带来的射线散射和吸收影响，对已定装置而言C近似为常数且归一化处理后实际对当前研究无影响。

图2 源架、货箱及剂量计的空间布置Fig.2 Spatial arrangement of source rack, container and dosimeter

(2)

其中：oi为流水线层编号；si为流水线上货箱工位编号。

进一步考虑货箱换面照射，则货箱中参考点的最终累积剂量Di,j,k可表示为：

(3)

1.2 基于多目标的评价函数定义

(4)

2 排源模型的局部复合求解算法

从评价模型可知，排源问题为典型的非线性规划问题，相关研究表明模拟植物生长算法(PGSA)求解排源问题可行有效，该算法通用性强，原则上只需替换评价函数即可[16]。为减少现场布源工作量，本文改进了PGSA而形成局部复合求解算法。

2.1 局部复合求解的基本原理

PGSA潜在需至少1个在所有变量维数上尽量均匀分布的初始方案，文献[16]表明强鲁棒性虽使得初始方案对最终结果的优化影响不大，但计算量变化剧烈。为保证初始方案的合理性及降低计算量，在求解过程中将源架分割为若干区域，源棒按总活度进行比例分配，分别排源后组成1个次优化方案，再进行整体优化计算。

2.2 算法的求解流程

算法的求解流程及功能如下。

1) 读入计算所需的各种初始数据，包括源架、棒位、源棒、货箱及填充材料、相应空间位置及控制算法本身等的系列参数。

2) 源架可按子源板自然分为N0个子区域，或按源架对称轴划分；将源棒按总活度对应分为N0组，接近退役的源棒优先分配到一起。

3) 依次将源棒分组排布到各源板子区域，若全部排布完，则进行步骤11；若没有，则进行步骤4。

5) 以xB为起点，以λ(≥1的正整数)为步长，在各维变量上沿坐标轴的正负方向寻找满足约束条件[a,b]的生长点xTemp，若f(xTemp)

6) 求出X内目标函数值最小点xmin，若f(xmin)≤f(Xmin)且不重复，则把该点放入最小生长点集合Xmin。

7) 若生长点集合X容量太大，则以一定的修剪比例对X内的劣势生长点进行裁减，经验值取0.4～0.6为佳；本轮产生的xTemp应全部保留从而保持生长点集合X的生物多样性以防算法早熟而导致优化停滞。

8) 若当前生长周期内的f(xmin)相比f(Xmin)没有进步即视为生长停滞，且停滞次数超过设定值，或总的生长周期数大于设定值，那么本分组计算结束，保存Xmin、f(Xmin)等结果待合并，进入步骤3开始下一分组排源，若不满足则进入步骤9。

9) 计算X内各生长点的形态素浓度，即各自的状态空间概率[16]。

10) 产生随机数β并计算落入的状态空间，将所落区间对应的生长点记为新生长基点xB，转回步骤5开始新的生长周期。

11) 将所有分组的Xmin合并，完成所有源棒在整体源架上的排放。

12) 若不需进一步整体排源，则输出Xmin形成系列最终排源方案，进而精细计算得到这些方案在货箱中的DUR、能量利用率、平均累积剂量等评价用结果。若进一步需要整体计算，则将Xmin作为整体排源用初始生长点X0，依次调用步骤4～10完成整体优化计算，最终得到系列优化排源方案供客户选用。

2.3 基于GPU的并行加速计算

本文提出的新评价模型相比原模型，其理论复杂度更高而计算复杂度偏低，原模型所耗计算时间占总计算时间的90%以上[7]，新模型经初步测试计算占45%左右，需进一步提升其计算效率。文献[17]中基于GPU通用计算技术实现了PGSA的并行化，在费用增加不明显的前提下使得原方法的计算效率显著提高。基于同样思路，改进后的PGSA也进行了基于GPU的并行化计算改进。

限于篇幅，本文仅给出当前算法中最重要的改进部分，即目标函数f的GPU并行化实现原理，如图3所示，其中步骤1、2、4、5计算压力小，以GPU最大Block和Thread数执行即可，步骤3计算压力大，执行时设计为每个Block负责计算生长点集合1个生长点的目标函数值，Block中每个Thread计算1根源棒造成的剂量。

新评价模型对每个Thread的任务压力明显下降，求解代码不仅计算效率显著提升，且在GPU单精度计算模式下即可获得满意结果，可避免购买昂贵的专业双精度计算卡。

图3 PGSA中目标函数f基于GPU的并行化原理Fig.3 GPU parallelization principle of target function f by PGSA

3 实际应用及结果分析

3.1 相关装置参数

以不同大型辐照装置的实践排源设计结果为例，验证算法的可靠性和有效性，表1列出了辐照装置的有关参数。

表1 辐照装置参数Table 1 Parameter of irradiation facility

算法代码的测试平台主要参数为：CPU，Xeon E3 1230 V2(数量1)；GPU，GTX285(数量1，硬件规格1.3)；内存，DDR3 1333(8 GB)；软件平台，Windows 7 64bit、VS2008、CUDA4.2、C++。

算法的鲁棒性较好，参数限制少。计算时源板按每层板块数划分，源棒按总活度等比例划分并由程序随机抽取后均匀布放到棒位上，货箱中参考点的空间间隔为5 cm，算法修剪比例为50%，编写参数输入文件，利用不同版本的代码程序计算，其耗费时间情况列于表2。从表2可知，计算同等规模装置，本文提出的新方法对应代码的计算时间明显缩短，可显著缩短排源方案设计周期。

表2 不同版本排源代码的计算时间Table 2 Computation time for source pencil deployment of various codes

3.2 结果及分析

现场部署排源方案后，按国标采用重铬酸银剂量计布置测试用货箱，表3为测试货箱参数、装置剂量场的计算结果及实测结果。

从表3可知，计算采取均衡策略，在测试货箱充满模拟货物前提下，对于装置1和装置2，新模型的理论DUR相比原模型分别增大了7.0%和16.5%，实测DUR相比计算结果分别增加了4.195%和0.70%，均显著低于国标要求限值而实际影响不大，同时能量利用率显著提高5.72倍和3.02倍。实测DUR相比理论值升高，其主要原因是计算模型相比实际辐照装置做了合理简化，未考虑护源板、货箱壁，另外计算模型以空气为受照介质，而实测采用了密度更大的材料，射线的衰减和散射差别较大，IAEA证实这类问题可采用ρ-DUR关系曲线修订而不影响排源[18]，且排源工作核心点是寻找优化方案，DUR和能量利用率的值仅用于评判方案的优化程度。

优化的排源方案在装置1、2参考面(图2)上的等剂量线分布如图4所示。

表3 测试货箱参数、装置剂量场的计算结果及实测结果Table 3 Test container parameter and result of facility dose field calculation and measurement

注：( )内为模型改进后的数据；[ ]内为相对偏差

图4 装置1(a)、2(b)参考面上的剂量分布Fig.4 Dose distribution along reference plane in facility 1 (a) and facility 2 (b)

从图4可知，排源方案在装置1、2参考面上的剂量分布很均匀，基本以装置中心为圆心呈同心圆分布且梯度均匀，这对降低DUR有利，证实了排源方案设计的可靠性。

4 结论

本文提出了一种新型多目标优化排源方法，建立了综合考虑货物DUR、射线能量利用率、现场布源工作量等多目标参数的新数学模型，设计了新型的GPU并行PGSA算法用于模型求解，并对两座大型γ辐照装置进行了实践排源。计算和测量结果表明：新方法计算精确度得到保证，对昂贵专业计算卡的依赖度下降；计算效率提升显著，计算时间从原有方法的h量级下降到min量级；新方法的理论和实测DUR比原方法增大但仍显著低于国标控制标准，且能量利用率显著提升，同时分区域排源降低了现场布源工作量。检测结果验证了新方法的可行性和有效性。