时代的快速发展和人民生活的水平不断提高，使得电器电子产品的更新换代速率加快. 大量废旧的电器电子产品中含有许多可再利用的资源和危害材料，若大量闲置或随意丢弃，会对环境造成污染和带来资源的浪费. 为了有效地解决废旧电器电子产品的拆卸问题，Gungor等^[1]提出拆卸线平衡问题(disassembly line balancing problem, DLBP).

DLBP经过20年的发展，学者们从求解方法和拆卸布局方面不断深入研究. 求解方法主要包括启发式方法^[2]、数学规划法^[3]和智能算法等，而启发式方法过于依赖启发式规则，数学规划方法不适合求解大规模DLBP. 智能算法已广泛应用于求解DLBP，Ren等^[4]以拆卸利润最大为目标，提出重力搜索算法求解. 丁力平等^[5]结合Pareto解集思想与智能算法，更好地平衡多个优化目标，提出基于Pareto的蚁群算法求解多目标DLBP. 智能算法以其诸多优势已广泛应用于求解拆卸序列规划^[6-8]、车辆路径规划^[9]和车间调度^[10-11]等问题，均取得了良好的成果.

拆卸线布局的研究主要集中于直线型^{[1-5, 12-14]}. 与直线型拆卸线相比，并行拆卸线是将2条及以上的拆卸线平行放置，可同步拆卸多种不同的产品，有效地提高拆卸效率和生产线柔性. Hezer等^[15]提出并行拆线平衡问题(parallel disassembly line balancing problem, PDLBP)，采用基于最短路径的网络模型求解. Guo等^[16]设计可同时拆卸和装配的并行生产线，构建以拆卸成本和平滑度指标为优化目标的数学模型，应用变邻域和非支配排序遗传算法的混合算法求解. Zhu等^[17]构建并行不完全拆卸线平衡问题的数学模型，设计混合群邻域算法求解. Hezer等^[15-17]分别构建不同问题的PDLBP数学模型，但存在一定的不足，所建模型均属于概念模型，无法精确求解，即无法验证所有约束对于优化目标的正确性. 存在对2条拆卸线任务定义不明确的情况，例如提到并行线上具有相同编号任务时无法确定是哪条线上的任务. 本文考虑当前劳动力短缺的国情，拆卸线属于劳动密集型工作，需要大量的劳动. 替换部分工作站操作者由机器人拆卸，即在站间操作者不同的基础上与PDLBP结合，构建以工作站数目、机器人数、拆卸成本和空闲时间均衡指标为优化目标的PDLBP数学模型，应用CPLEX和LINGO对模型进行精确求解，验证该模型的正确性.

头脑风暴优化算法(brain strom optimization, BSO)^[18]是受头脑风暴会议中思维碰撞过程的启发而提出的群智能优化算法，现已成功应用于电动车路径优化^[19]和离散调度问题^[20]等离散型组合优化问题中，表现出优越的求解性能. DLBP属于离散型组合优化问题，为了解决本文的多目标PDLBP，提出改进头脑风暴优化算法(improved brain storm optimization, IBSO). 在算法中，将原始连续操作进行离散化处理，分别设计变异交叉操作，分别对应算法中的单个个体产生机制和2个个体产生机制，设计四点交叉操作策略. 引入Pareto解集思想和NSGA-Ⅱ拥挤距离^[21]的双精英策略，筛选非劣解.

并行拆卸线的拆卸示意图如图1所示，2条直线型拆卸线平行放置并同步运行，用于拆卸2种不同产品，拆卸线Ⅰ和Ⅱ之间通道用于安排工作站，并行拆卸线工作站安排存在2种不同情况. 第1种为共享工作站情况，如图1的工作站1和4所示，在该情况中拆卸机器人或工人需要在2条拆卸线中运行拆卸2种产品. 第2种为独立工作站情况，如图1的工作站2、3、5和6所示，该情况与直线型相同，只用于拆卸某一条线上任务. 为了在一定程度上减少劳动力，避免全自动化投入成本过高，并在一定程度上提高拆卸效率，考虑在站间采用不同的操作者.

采用零件的优先关系图，以直观地表示零件之间的拆卸约束关系. 如图2所示为分别包含5个拆卸任务和4个拆卸任务的2个产品任务优先关系图. 图中，圆圈和正方形内数字表示各产品的拆卸任务，各任务间由箭头相连，箭尾与紧前任务相连，箭头与紧后任务相连；P_l为拆卸线l上任务的优先关系矩阵，

${{\boldsymbol{P}}_l} = {\left[ {\begin{array}{*{20}{c}} {{e^l_{11}}}&{{e^l_{12}}}& \cdots &{{e^l_{1{n_l}}}} \\ {{e^l_{21}}}&{{e^l_{22}}}& \cdots &{{e^l_{2{n_l}}}} \\ \vdots & \vdots &{}& \vdots \\ {{e^l_{{n_l}1}}}&{{e^l_{{n_l}2}}}& \cdots &{{e^l_{{n_l}{n_l}}}} \end{array}} \right]_{{n_l} \times {n_l}}}. $

若拆卸线l中任务i是任务j的紧前任务，则表示 $e_{ij}^l$=1，否则 $e_{ij}^l $=0；l为拆卸线编号， $l \in \left\{ {1,2} \right\}$；n_l为拆卸线l上产品的任务总数.

将拆卸优先关系图转换为拆卸优先关系矩阵，用于问题求解，例如产品1中任务5的紧前任务有任务3和任务4，所以任务5须等所有紧前任务拆卸完成才可拆卸，即在P₁中P₁(3，5)和P₁(4，5)均等于1，如两任务即任务i和任务j间无优先关系约束，则P(i，j)=0. 为了更好地描述问题与模型，给出如下假设.

1)每个工作站只分配一个机器人或一名工人工作.

2)机器人和工人拆卸每种产品各项任务的时间已知.

3)忽略共享工作站时操作者的移动时间.

4)待拆卸产品完整，开展完全拆卸.

目标函数如下.

(1) $ {f_1} = \sum\limits_{w \in {{W}}} {{S_w}}, $

(2) ${f_2} = \sum\limits_{w \in W}{{S_w}{r_{2w}}} , $

(3) $ \begin{split} &{{f_3} = \sum\limits_{w \in {{W}}} {{S_w}} \left( {\sum\limits_{i \in {{{I}}_1} \cup {{{I}}_2}} {\sum\limits_{k \in {{K}}} {\left( {{C_{w1}} {x_{iwk}} {t_{ik}} + {C_{w2}} \left( {{C_{\rm{T}}} - {x_{iwk}} {t_{ik}}} \right)} \right)} } }+ \right.} \\ &{\left. { {C_{\rm{r}}} {r_{2w}} + \sum\limits_{i \in {{{I}}_1} \cup {{{I}}_2}} {\left( {{C_{{\rm{pt}}}} {x_{iw1}} {t_{i1}} + {C_{{\rm{rt}}}} {x_{iw2}} {t_{i2}}} \right)} } \right)}, \end{split} $

(4) $ {f_4} = \sum\limits_{w \in {{W}}} {{S_w}} {\left( {{C_{\rm{T}}} - \sum\limits_{i \in {{{I}}_1} \cup {{{I}}_2}} {\sum\limits_{k \in {{K}}} {{x_{iwk}} {t_{ik}}} } } \right)^2}, $

(5) $ F = \min \;\left\{ {{f_1},{f_2},{f_3},{f_4}} \right\} .$

式中：C_T为节拍时间；i为拆卸任务编号；I_l为拆卸线l上任务编号集合，拆卸线Ⅰ任务集为 ${{{I}}_1} = $ $ \left\{ {1,2,\cdots,{n_1}} \right\}$，拆卸线Ⅱ任务集为 ${{{I}}_2} = \{ {n_1} + 1, {n_1} + 2,\cdots, $ $ {n_1} + {n_2} \}$；W为工作站编号集合， ${{W}} = \{ 1,2,\cdots, {n_1} + {n_2} \}$；w为工作站编号，其中 $w \in {{W}}$；K为操作者编号集合， ${{K}} = \left\{ {1,2} \right\}$；操作者编号k=1为工人拆卸，k=2为机器拆卸，其中 $k \in {{K}}$；t_ik为任务i由操作者k拆卸所需的时间；C_w1为工作站工作时的单位时间成本；C_w2为工作站待机时的单位时间成本；C_r为投入一个机器人所需的成本；C_rt为机器人拆卸时的单位时间成本；C_pt为工人拆卸时的单位时间成本；S_w为0-1变量，若工作站w开启则S_w=1，否则为0；x_iwk为0-1变量，若任务i 分配给工作站w中的操作者k，则x_iwk=1，否则为0；r_kw为0-1变量，若操作者k分配给工作站w，则r_kw=1，否则为0.

式(1)表示优化目标f₁为开启工作站数目. 式(2)表示优化目标f₂为机器人使用数量. 式(3)表示优化目标f₃为拆卸过程所需成本，综合考虑成本，包括工作站工作成本、工作站待机成本、机器人数量成本、工人拆卸成本和机器人拆卸成本. 式(4)表示优化目标f₄为空闲时间均衡指标. 式(5)表示4个优化目标尽可能最小.

约束条件如下.

(6) $ \sum\limits_{w \in {{W}}} {\sum\limits_{k \in {{K}}} {{x_{iwk}}} } = 1,\; i \in {{{I}}_1} \cup {{{I}}_2}; $

(7) $ 1 \leqslant \sum\limits_{w \in {{W}}} {{S_w}} \leqslant {n_1} + {n_2}; $

(8) $ \sum\limits_{k \in {{{K}}}} {\sum\limits_{i \in {{{I}}_1} \cup {{{I}}_2}} {{x_{iwk}}} } \leqslant \left( {{n_1} + {n_2}} \right) {S_w},\; w \in {{W}}; $

(9) $ \sum\limits_{k \in {{{ K}}}} {\sum\limits_{i \in {{{I}}_1} \cup {{{I}}_2}} {{x_{iwk}}} } \geqslant {S_w},\; w \in {{W}}; $

(10) $ {S_{w - 1}} \geqslant {S_w},\; w \in \left\{ {2,3,\cdots,{n_1} + {n_2}} \right\}; $

(11) $ \sum\limits_{k \in {{K}}} {{r_{kw}}} = {S_w},\; w \in {{W}}; $

(12) $ \sum\limits_{i \in {{{I}}_1} \cup {{{I}}_2}} {\sum\limits_{k \in {{K}}} {{x_{iwk}}} } {t_{ik}} \leqslant {C_{\rm{T}}} {S_w},\; w \in {{W}}; $

(13) $\left. {\begin{array}{*{20}{l}} {\displaystyle\sum\limits_{k \in {{K}}} {\displaystyle\sum\limits_{w \in {{W}}} {\left( {{x_{iwk}}w} \right)} } \displaystyle\sum\limits_{k \in K} {\displaystyle\sum\limits_{w \in W} {\left( {{x_{jwk}}w} \right)} } ,} \\ {i,j \in \left\{ {\left( {i,j} \right)|e_{ij}^1 = 1||e^2_{(i - {n_1})(j - {n_1})} = 1} \right\};} \end{array} } \right\}$

(14) $ \sum\limits_{i \in {{{I}}_1} \cup {{{I}}_2}} {{x_{iwk}}} \leqslant {{M}} {r_{kw}},\; k \in {{K}},\; w \in {{W}},\;\; \\ {{M}} \geqslant {C_{\rm{T}}} ({n_1} + {n_2}). $

式中：i、j为拆卸任务；M为一个大的正数. 式(6)为任务分配约束，完全拆卸且每个拆卸任务仅能分配一个工作站. 式(7)为工作站开启数目约束，所开启工作站在1~(n₁+n₂). 式(8)、(9)为拆卸任务被分配时必有工位开启. 式(10)为工位约束，工作站须按顺序拆卸. 式(11)为工位分配操作者约束，每开启一个工作站须分配一种类型的操作者工作. 式(12)为节拍约束，每个开启工作站的工作时间不得超过节拍时间. 式(13)为优先关系约束，拆卸线l上任务拆卸须满足各自产品的优先关系. 式(14)表示若任务i分配至工作站w且分配给第k类型的操作者，则该类型的操作者须分配至工作站w.

BSO算法是以聚类中心为导向来产生新个体，对应于头脑风暴会议中以一个讨论得到的方案为参考生成新的方案过程. 结合PDLBP的特点，提出改进头脑风暴优化算法. 具体的改进措施如下. 离散化单个个体和2个个体产生机制的原始操作，应用L2范数方法在聚类前对目标矩阵进行归一化处理，设计基于交叉算子和变异算子的产生机制. 为了增加种群个体的多样性，设计四点交叉的操作策略.

1）分别从拆卸线Ⅰ和Ⅱ中的优先关系矩阵P_l中搜索无紧前工作的任务，将所有的无紧前工作的任务构成可分配任务集A.

2）从可分配任务集A中随机挑选一个任务分配到当前的拆卸序列中. 若该任务为拆卸线Ⅱ上的任务，则将该任务在P₂的位置索引+n₁后分配至当前的拆卸序列中.

3）若该任务为拆卸线Ⅰ中的任务，则解除拆卸线Ⅰ与该任务相关任务的紧前约束；若该任务为拆卸线Ⅱ中的任务，则解除拆卸线Ⅱ与该任务相关任务的紧前约束.

4）重复1）~3），直到拆卸线Ⅰ和拆卸线Ⅱ中所有的任务分配完成.

5）随机生成一组长度为n₁+n₂、只有0和1的随机数组，其中0表示工人，1表示机器人.

并行拆卸线的解码过程是按照拆卸序列的顺序，将各项拆卸任务合理地分配到工作站中. 由于本文是在节拍时间给定的情况下优化工作站数目，开启工作站数目未知. 为了将所有拆卸任务分配至工作站中，并为每个开启的工作站分配工人或机器人，在编码时随机生成一组长度为(n₁+n₂)的0-1随机数组作为工作站属性用于解码. 如图3所示为解码示意图. 图中，圆形为拆卸线Ⅰ任务，正方形为拆卸线Ⅱ任务. 具体的解码步骤如下.

1）输入带有双层编码的可行拆卸序列，开启工作站w=1. 令工作站剩余时间等于节拍时间，即R_T = C_T，拆卸序列指针d_n = 1.

2）提取拆卸序列第d_n位置的任务i，确定第2层0-1随机数组中第w位置对应的工作站属性. 若第w位置对应的0-1码为0，则该工作站为工人拆卸；若对应的0-1码为1，则该工作站为机器人拆卸，提取任务i对应拆卸操作者类型的拆卸时间t_ik.

3）判断t_ik是否大于R_T. 若不大于R_T，则转至4）；若大于R_T，则转至5）.

4）将任务i分配至第w个工作站，计算该工作站剩余时间R_T = R_T − t_ik.

5）开启新的工作站w=w+1，重新提取0-1随机数组中第w位置对应的工作站属性和任务i对应操作者类型的拆卸时间t_ik，R_T = C_T − t_ik.

6）判断d_n是否小于n₁+n₂. 若小于n₁+n₂，则d_n = d_n+1并返回2）；否则输出任务分配结果.

针对DLBP问题的特点，目标函数的各维特征尺度存在不同，采用K-mean聚类方法会造成各维特征对形成的距离的影响存在一定差异. 在聚类之前，需要对目标矩阵进行预处理. 形成聚类中心的步骤如下.

1）对适应度的目标矩阵，应用L2范数方法进行归一化预处理.

2）计算归一化得到的目标矩阵中每个个体间的欧氏距离，形成距离矩阵.

3）选取平均欧氏距离最大的个体作为第一个聚类中心，将距离聚类中心最近的(C_n/(m−1))个个体聚为一类.

4）计算每个个体到各聚类中心的欧氏距离，选取平均距离最大的个体作为下一个聚类中心，根据上述方法形成所有聚类中心并完成聚类操作.

基于单个个体的产生机制是一种无规则局部搜索操作，该过程相似于进化算法中的变异操作. 通过不同的选择概率选择聚类中心和非聚类中心，添加扰动，产生新个体. 离散化的扰动公式如下：

(15) $\left. {\begin{array}{*{20}{c}} {X_{{\text{new}}}^j = X_{{\text{section}}}^i + \xi N\left( {\mu ,\sigma } \right),} \\ {\begin{array}{*{20}{c}} {\xi \left( g \right) = \log {\text{sig}}\left( {\dfrac{{0.5 {M_{\rm{gen}}} - g}}{k}} \right) \cdot {\text{random}}(0,1),} \\ {\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\log \;{\text{sig}}\left( x \right) = {{\left[ {1 + \exp\; ( - x)} \right]^{-1}}} .} \end{array}} \end{array}} \right\} $

式中： $X_{\rm{new}}^j$为添加扰动后产生的新个体 ${X_{\rm{new}}}$的位置j； $X_{\rm{section}}^i$为随机选择需要添加扰动的聚类中心或非聚类中心的位置i； $N\left( {\mu ,\sigma } \right)$为均值为 $\;\mu $、方差为 $\sigma $的高斯随机函数； $\xi $为随机值，取值为0~1.0，是对于新个体的权重系数，即高斯随机函数的系数，用来控制扰动时的幅度. $\xi $由对数sigmoid传递函数 $ {\rm{log}}\;{\rm{sig}}(·) $和random(0，1)乘积组成. 其中，M_gen为设定的最大迭代次数；g为当前的迭代数；k为 $ {\rm{log}}\;{\rm{sig}}(·) $函数的斜率；random(0，1.0)为随机产生的(0，1.0)的随机数. 具体的基于单个个体的产生机制如图4所示. 图中，圆形为拆卸线Ⅰ任务，正方形为拆卸线Ⅱ任务. 随机选中的聚类中心或非聚类中心为[1-6-2-4-7-3-8-5-9]，0-1序列为[0-1-1-0-1-0-1-0-1]，随机选中添加扰动任务为任务3；根据优先关系矩阵确定任务3的紧前任务为任务2，紧后任务为任务5. 任务3紧前与紧后任务间的序列[4-7-3-8]即为可操作片段，随机生成与可操作片段长度相同的(0，1.0)递增数组{0.12，0.41，0.63，0.76}. 将任务3对应的随机数添加式(15)的扰动后，由0.63变为0.34，即重新递增排序后，任务3排在任务4与任务7之间，任务3下的0-1码随任务3移动后将0变为1，扰动完成.

基于2个个体的产生机制是头脑风暴过程中不同小组之间思维碰撞的过程，在多轮小组头脑风暴后，当前小组的思路通常会变窄，新个体的产生会变得越来越难，所以需要跳出当前的小组讨论，进行随机的小组间思想碰撞来产生新方案. 结合该过程与DLBP，引入进化算法中的两点交叉操作.

1）随机选取2个聚类中心或非聚类中心X₁和X₂，随机生成2个1~(n₁+n₂)的随机整数a和b，a < b.

2）对于拆卸序列，提取X₁序列位置中a~b的序列片段C，寻找X₂中与片段C相同的任务；根据X₂序列中任务的拆卸顺序，将片段C中的任务重新排序形成C_new. 在重新排序过程中，片段C之间的任务在X₂中满足产品各自的优先关系，所以重新排序后片段C_new任务序列均未破坏2个拆卸产品的优先关系约束.

3）在X₁中将片段C_new替换片段C.

4）对于0-1码，将X₂中与片段C中相同任务下的0-1码进行替换，将操作后的新的拆卸序列和0-1码重新组合，形成新个体X_new1.

为了增加种群个体的多样性，扩大个体可操作性片段和随机性，引入四点交叉操作策略. 具体如图5所示，对于第1层拆卸序列，随机生成4个1~n₁+n₂交叉点c₁<c₂<c₃<c₄. 将X_old1的拆卸序列分成5段，选择c₁~c₂和c₃~c₄ 2段序列作为信息片段，保持X_old1 2个信息片段外的序列顺序不变，找到X_old2中对应2个信息片段的任务顺序，将其序列重新排序后替换原始信息片段，重新排序时因各任务在X_old2即满足任务间优先关系. 结合图2的两产品优先关系图可知，重新排序后的新片段未违反零件间的优先关系. 第2层0-1码，将X_old2中与信息片段相同任务下的0-1码替换原始信用片段下的0-1码，将新的拆卸序列和0-1码组成新个体X_new.

DLBP是典型的多目标组合优化问题，这与单目标优化不同，由于优化目标有多个，不能同时求解出最优解，只能提供一组数量合适的非劣解. 引用Pareto解集思想处理所涉及的4个优化目标并筛选非劣解，即算法运行一次可以求得多个互不占优的非劣解. 假设决策空间的x¹和x²是满足约束条件式(6)~（14)的2个可行解，当含有M个目标函数的解x¹优于解x²时，则称解x¹支配于解x²，表示为x¹$ \prec $ x². 具体的支配关系如下：

(16) $ \left.{\begin{array}{*{20}{c}} {{f_i}\left( {{x^1}} \right) \leqslant {f_i}\left( {{x^2}} \right),\; i \in \left\{ {1,2,\cdots,M} \right\}}; \\ {{f_j}\left( {{x^1}} \right) < {f_j}\left( {{x^2}} \right),\; j \in \left\{ {1,2,\cdots,M} \right\}}. \end{array}} \right\}$

该算法经过多次迭代，将Pareto筛选出的非劣解更新并储存到外部档案Q中，当外部档案Q中的非劣解达到一定的数量后，算法的效率会下降. 为了避免该情况的产生，引入NSGA-Ⅱ拥挤距离机制为非劣解排序，有效地减少非劣解的个数，实现精英保留策略. 拥挤距离定义为

(17) $ {C_d} = \sum\limits_{m = 1}^M {\frac{{f_m^{d + 1} - f_m^{d - 1}}}{{f_m^{\max } - f_m^{\min }}}} .$

式中：d为非劣解的个数， $d \in \left\{ {2,3,\cdots,D - 1} \right\}$；m为目标函数的个数， $m \in \left\{ {1,2,\cdots,M} \right\}$； $f_m^{d - 1}$和 $f_m^{d + 1}$分别为第d−1个和第d+1个非劣解的第m个目标函数值； $f_m^{\min }$和 $f_m^{\max }$分别为所有非劣解中第m个目标函数的最小值与最大值.

Pareto解集评价：超体积指标^[22](hyper-volume)简称HV指标，表示算法运行求得的Pareto解集与参考点即各目标最差值所围成的目标空间中区域的超体积. 如图6所示，在最小化的双目标优化问题中，求解出包含5个非劣解的Pareto解集{x¹，x²，x³，x⁴，x⁵}，f₁和f₂的最差值构成参考点R的坐标点(a，b)，HV指标的值即为坐标点(a，b)与非劣解x¹、x²、x³、x⁴、x⁵坐标点在二维空间所围成的超体积，即为阴影部分面积. HV指标越大，说明该算法的综合性能越好. 引用该指标，评价对比实例应用中IBSO算法与GA、NSGA-Ⅱ算法求解实例的性能. 具体公式如下：

(18) $ {\rm{HV}} = \delta \left( { \cup _{i = 1}^{\left| P \right|}{v_i}} \right). $

式中： $\delta $表示 Lebesgue 测度，用来测量超体积； $\left| P \right|$为Pareto解集中非劣解的个数；v_i为参照点与Pareto解集中第i个解构成的超体积.

IBSO具体算法流程如下.

1）输入算法相应的初始参数：生成的种群规模为C_n，拆卸线Ⅰ和Ⅱ的优先关系矩阵为P_l，带有机器人拆卸时间、工人拆卸时间等的信息矩阵K_B，最大的迭代次数为M_gen，聚类个数为m，聚类参数为p₁，分组参数为p₂，单个体扰动参数为p₃，两个体扰动参数p₄和节拍时间C_T等.

2）生成初始种群，建立外部档案 $Q =\phi $.

3）计算C_n个个体的目标函数值，应用Pareto筛选并更新外部档案Q.

4）设置初始迭代次数g=1，开始迭代.

5）对初始种群进行聚类操作，将其聚类为m类.

6）生成一个随机数R_a，若R_a < p₁，则随机生成一个个体替换该聚类中心；否则执行7）.

7）随机生成C_n个随机数R_b. 对于R_b < p₂的个体转至8），进行单个体操作；其余个体转至9），进行2个个体操作.

8）依次选择个体，生成随机数R_c. 若R_c < p₃，则对该个体的聚类中心进行2.4节操作产生新个体，否则对选中个体操作产生新个体，并与原个体比较择优保留.

9）生成随机数R_d. 若R_d<p₄，则随机选取一个聚类中心与该个体聚类中心进行2.5节的操作，产生新个体；否则随机选一个非聚类中心与该个体操作产生新个体，并与原个体比较，择优保留.

10）将8）、9）产生的新个体进行随机匹配，进行四点交叉操作（见2.6节）.

11）应用Pareto和拥挤距离对8）~10）获得的结果进行筛选和排序，更新种群和外部档案Q，计算该迭代次数下外部档案中Pareto解集的HV指标并储存.

12）判断g是否大于M_gen. 若小于M_gen，则返回5）并令g=g+1；否则程序终止，输出最优拆卸方案与目标函数值.

为了验证所建的PDLBP模型的正确性，采用CPLEX和LINGO精确解求解器，对1.3节所建模型开发相关程序，求解P10和P8算例. 算例中，基础拆卸时间和各零件优先关系等信息均来自文献[1, 23]，具体的拆卸信息如表1所示. 表中，t_ip为工人拆卸时间，t_ir为机器人拆卸时间，c_ost为成本. 应用数学规划的方法对所提各优化目标求解时，各目标类型均不同，例如f₁为混合整数线性规划(mixed integer linear programming, MILP)问题，f₂和f₃为混合整数二次规划(mixed integer quadratic programming, MIQP)问题，上述目标均可用CPLEX求解，但是f₄为混合整数非线性规划(mixed integer nonlinear programming, MINLP)问题，无法应用CPLEX求解，应用LINGO对f₄进行精确求解. 为了验证IBSO算法的有效性，应用MATLAB运行IBSO算法求解该算例，与精确解对比. IBSO参数为：C_n=120，M_gen=200，p₁=0.2，p₂=0.5，p₃=0.5，p₄=0.6，精确解与IBSO算法的对比结果如表2所示. 表中，t_s为求解时间. 各目标最优值的任务分配图如图7所示. 图中，t_w为工作站工作时间；任务前为“−”的是拆卸线Ⅱ任务；工作站中为点划线，表示该工作站由机器人拆卸，拆卸线Ⅰ拆卸P10产品即n₁=10，拆卸线Ⅱ拆卸P8产品即n₂=8.

由表2对比可知，CPLEX与LINGO求得精确解与IBSO算法求得各目标的最优值均相同，由于各目标求解难度不同，求解时间具有较大差距. 求解目标函数f₂仅需要0.28 s，求解目标函数f₄需要183 342.7 s，求解时间大大增加. 应用IBSO算法求解该算例，算法运行一次就可以求解出包含全部精确解的目标函数，求解时间比较适中，为11.1 s，验证了所建PDLBP模型的正确性与IBSO算法的有效性.

为了验证IBSO算法求解性能的优越性，应用该算法求解包含25项拆卸任务的移动电话^[24]经典算例问题，与文献[23, 25~28]的求解结果进行对比. 利用IBSO算法求解P25经典算例问题时，根据问题的特征对解码操作进行修改，以适合求解该问题，但编码思想相同且算法主体结构和核心操作均未改变，因此应用P25经典算例问题验证IBSO算法的性能. P25问题的4个优化目标分别为：最小化工作站开启数f₁、空闲时间均衡指标f₂、危害指标f₃和需求指标f₄. IBSO求解P25问题的参数为：C_n=150，M_gen=300，p₁=0.2，p₂=0.4，p₃=0.5，p₄=0.6.

如表3所示为IBSO算法与文献[23, 25~28]的求解结果对比. 多目标角度如下：将所有算法结果进行一次Pareto筛选，筛选共得到12个Pareto前沿解，其中包含IBSO算法求解结果的有10个，可见IBSO所得的非劣解集更逼近真实Pareto前沿. 单目标角度如下：虽然所有对比算法的前2个目标值均求解为9，但对于目标f₃和f₄的最优值分别为73和811，均劣于IBSO算法所求的最优值70和802. 根据式(18)分别计算各算法所求结果的HV. 从表3可以看出，利用IBSO算法所求结果的HV最大为5.033×10⁸. 求解经典算例的HV指标参考点为[25，4 291，106，1 015]. 综上所述，以上对比有效验证了IBSO算法的求解性能.

为了验证IBSO算法求解本文问题的求解性能，应用27规模的电视机^[29]和25规模的电冰箱^[17]拆卸实例，构建2条并行拆卸线；运用所建模型和所提算法对2种废旧电器分配工人或机器人进行拆卸，应用2种经典算法GA^[30]和NSGA-Ⅱ^[21]和BSO求解实例与IBSO算法进行对比. 电视机和电冰箱的拆卸信息如表4所示，电视机和电冰箱的优先关系分别如图8、9所示.

如图10所示为利用IBSO、GA、NSGA-Ⅱ和BSO 4种算法迭代600代的HV指标迭代图. 可以看出，IBSO算法在119代时HV指标为1.061 9×10¹²，之后经过较长的一段平稳期；当达到305代时HV指标上升为1.062 4×10¹²，上升幅度仅为0.04%，之后直到算法终止达到收敛. NSGA-Ⅱ、GA和BSO分别于432代、349代和504代后算法达到收敛，IBSO算法的收敛速度均快于上述3种算法，且IBSO算法的HV指标均在3种对比算法上方. IBSO、GA、NSGA-Ⅱ和BSO算法运行10次的平均时间分别为117.8、2 926.25、3 133.21和93.82 s. 综上所述，以上对比验证了IBSO算法的求解性能. HV指标的参考点为[52，52，1000，622858]，图10中IBSO算法迭代600代时所求的非劣解集如表5所示. IBSO和BSO参数为：C_n =300，M_gen=600，p₁=0.2，p₂=0.4，p₃=0.4，p₄=0.6；GA和NSGA-Ⅱ参数为：交叉概率=0.9，变异概率=0.3.

如表5所示为利用IBSO算法求解实例的结果，共求得10个Pareto非劣解. f₁求解结果为6和7，f₂求解范围为0~6，f₃求解范围为281.88~366.66，f₄求解范围为7~5 040. f₄变化幅度较大的原因为f₄是用于衡量各工作站作业时间偏离节拍时间的程度，如式(4)所示，在计算时需要计算空闲时间的平方，因而目标值差异较大. 分析表5中4个优化目标的结果可知，随着某个指标的上升，存在其他指标的下降，各优化目标间存在一定冲突，因此无法同时求出4个优化目标的最优值. 若决策者更注重拆卸成本，则可选择方案8；若更注重空闲时间均衡指标，则可选择方案7. 若更注重各指标的均衡，则可选择其他方案. 表5中，拆卸任务前为“−”表示为拆卸线Ⅱ上任务，工作站中Wr表示该工作站分配机器人工作，拆卸线Ⅰ拆卸电视机即n₁=27，拆卸线Ⅱ拆卸电冰箱即n₂=25.

(1)设计站间采用不同操作者的可同步拆卸不同种废旧产品的并行拆卸线，2条并行放置的拆卸线可以分别拆卸不用类型废旧产品进行回收利用. 构建多目标PDLBP混合整数规划模型，应用CPLEX和LINGO对模型的小规模算例进行精确求解，验证了模型的正确性.

(2)提出IBSO算法，用于求解本文问题. 算法采用双层编码方式构造可行拆卸序列，设计基于单个个体产生机制和基于2个个体产生机制的变异交叉操作和四点交叉操作，引入Pareto解集思想和拥挤距离机制筛选多目标非劣解.

(3)应用IBSO算法求解25规模基准算例，与文献[23, 25~28]的求解结果进行对比，验证了IBSO算法的优越性. 将该算法引入电视机与电冰箱的大规模拆卸实例中，通过与GA、NSGA-Ⅱ和BSO算法求解结果和HV指标对比，验证了IBSO算法的求解性能.