求解带时间窗车辆路径优化问题的改进细菌觅食算法

李 珺，郝丽艳，何奕涛，段钰蓉

（兰州交通大学电子与信息工程学院，兰州 730070）

0 概述

车辆路径问题（Vehicle Routing Problem，VRP）一直是物流配送活动中的基本问题，受到国内外学者的广泛关注。随着电子商务的发展以及客户需求的提高，企业在设计物流系统时需要在规定时间窗内完成服务，从而提高服务水平。带时间窗的车辆路径问题（Vehicle Routing Problem with Time Windows，VRPTW）是经典VRP 问题的一个重要扩展，每个顾客有预先设定的时间窗，车辆不能早于顾客允许的最早服务时间提供服务，也不能晚于顾客允许的最晚服务时间提供服务。

求解VRPTW 问题的算法大致可以分为精确算法和启发式算法两类。精确算法在理论上可以找到问题的最优解，但由于在实际应用中消耗的空间和时间成本较大，因此计算机内存要求较高，仅适用于求解较小规模的路径优化问题。启发式算法不管是求解小规模的问题还是大规模的问题，都能够在一定范围和较短的时间内给出满意解或最优解。因此，目前相关领域的学者专注于设计不同的启发式算法寻找该问题的近似最优解，特别是对现代启发式算法的研究。文献［1］利用分布式多agent系统实现分布式求解VRPTW问题。文献［2］介绍基于量子蚁群算法的VRPTW 问题。文献［3］介绍一种求解VRPTW 问题的离散粒子群算法。文献［4］介绍改进的蚁群算法求解带时间窗的车辆路径问题。文献［5］介绍一种新的求解带时间窗车辆路径问题的遗传算法。文献［6］介绍考虑载荷约束的求解时间窗车辆路径问题的遗传算法。

目前，还有一些学者采用深度强化学习方法求解车辆路径问题。文献［7］提出一个基于深度强化学习框架的改进启发式规划方法来求解车辆路径问题。通过设计基于自我关注的深层架构策略网络，指导下一步解决方案的选择。学习到的策略能很好地概括不同的初始解和问题规模，并在实际数据集上给出了较好的解决方案。文献［8］提出一种基于紧急程度的插入启发式算法来构造初始解，在局部搜索改进阶段使用强化学习来指导搜索。实验结果表明，基于单个解的算法中的抽样方案并没有显著提高解的质量，但是可以大幅降低搜索过程中发现的不可行解的比例，在大规模搜索空间内具有很强的搜索能力。文献［9］提出一种具有动态编译码结构的动态关注模型，该模型能够在不同的构建步骤中动态探索节点特征并有效挖掘隐藏的结构信息，但是训练阶段比较耗时，仅可针对中小规模问题实例进行训练。文献［10］使用基于蒙特卡罗预测及深度策略梯度学习的智能车辆路径规划方法，设计一种基于环境感知预测、行为决策和控制序列生成的深度强化学习框架，并将多个单一的强化学习算法应用到智能车路径规划系统的不同子任务中。实验结果验证了车辆具有的预估风险能力，并以尽量大的安全距离调整车辆状态，行驶时有更加连续的转角控制序列，实现快速、安全、平稳的路径规划。文献［11］提出一种基于强化学习的超启发式算法，设计算法的高层启发式策略，包括选择策略和解的接收准则，算法在实验求解过程中整体性能较好，但在客户规模较大的车辆路径问题上仍具有改进空间。目前，关于使用深度强化学习方法求解VRPTW 问题的研究相对较少，可进行更深入的探讨。

细菌觅食算法（Bacterial Foraging Algorithm，BFA）作为现代启发式算法之一，是由PASSINO［12］于2002 年提出的一种仿生随机搜索算法，具有局部搜索能力强、并行搜索、实现简单等优点，目前已在优化计算、图像处理、模式识别、系统仿真、人脸识别等领域得到广泛应用［13-15］。本文提出一种改进的细菌觅食算法来求解VRPTW 问题，使用K-means 算法根据地理位置对待配送点进行聚类，依据聚类结果采用贪婪插入法生成初始解。将细菌觅食的趋化操作与大邻域搜索相结合，设计4 种removal 算子进行寻优，扩大搜素范围，在不增加时间复杂度的情况下，提高收敛精度。

1 考虑时间窗的车辆路径优化问题

VRPTW 问题描述：有一定数量的相同型号的车辆，从中心仓库出发，对若干个客户节点进行服务，要求满足约束条件，合理安排配送路径，使总配送成本最小。基本假设如下：

式（1）表示最小化车辆总行程；式（2）表示每个客户仅被一辆车服务一次；式（3）表示每辆车的行驶距离都不得超过最大行驶距离；式（4）表示每个客户只被同一辆车服务；式（5）表示车辆从配送中心出发最终必须再回到配送中心；式（6）表示每辆车配送的客户总需求量不得超过车辆最大载重；式（7）表示客户i的开始服务时间必须满足给定的时间窗；式（8）表示消除子回环。

2 改进的细菌觅食算法

改进的细菌觅食算法（Improved Bacterial Foraging Algorithm，IBFA）是在细菌觅食算法的基础上结合大邻域搜索（Large Neighborhood Search，LNS）的优化算法［16-18］。该算法主要通过趋化操作、复制操作和迁徙操作的迭代计算来搜索问题的最优解。趋化操作是算法的核心操作，因此本文在趋化操作中结合了大邻域搜索算法中的removal 算子，运用4 种removal 算子进行寻优，扩大细菌的搜索范围。大邻域搜索算法由SHAW［19］于1998年提出，克服局部搜索的缺陷，提高BFA优化能力，能够更快地找到最优解，对于求解VRPTW问题十分有效。

2.1 编码

本文采用非负整数编码方式表示解空间，配送中心编号为0，客户节点编号为1，2，…，N，令S=［S1，S2，…，SK］（K为车辆数目）为VRPTW 问题的一组可行解，假设车辆路径规划为S=［S1，S2］=［｛0，3，2，8，7，6，0｝，｛0，1，4，5，9，0｝］，第1 条路径从配送中心出发，服务客户3、2、8、7、6 后回到配送中心，第2 条路径从配送中心出发，服务客户节点1、4、5、9 后回到配送中心，如图1 所示。

图1 配送示意图Fig.1 Schematic diagram of distribution

2.2 初始解构造

本文采用基于K-means 算法的贪婪插入法构造VRPTW 问题的初始解。在运用贪婪插入法构造初始解时，由于前期客户的插入次序对解的构造有较大影响，因此使用K-means 算法对客户节点进行预处理，其中V为配送中心所拥有的最大车辆数目，k为随机选取该范围中的一个整数，利用客户节点的位置信息对客户节点聚类，根据聚类结果产生客户的插入序列，使用贪婪插入法生成路径。贪婪插入法首先构建一条从配送中心出发再返回配送中心的路径，即［0，0］，然后运用K-means 算法排列好的客户节点依次插入路径中的最佳位置，当路径中无法再插入新客户时，重新构建一条［0，0］的路径，重复上述操作，直到所有客户都插入为止。

为验证采用K-means 构造初始解对算法的影响，本文选取Solomon 测试集中的R211 算例进行验证。测试集共100 个客户节点，车辆载重为1 000，使用K-means 算法聚类后的客户排列顺序和原始客户集中的客户排列顺序分别生成初始解，共进行10 次实验，采用K-means 聚类后得到的最终解平均值为782.43，使用初始客户集生成最终解的平均值为791.08。图2 给出了使用K-means 算法得到的最优解构造的初始路径。图3 给出了使用原始序列得到的最优解构造的初始路径。实验结果表明，使用K-means算法对客户集进行预处理得到的解更优。

图2 基于K-means 算法最优解的初始路径Fig.2 Initial path based on the optimal solution obtained by the K-means algorithm

图3 基于原始序列最优解的初始路径Fig.3 Initial path based on the optimal solution obtained by the original sequence

2.3 趋化操作改进

在细菌觅食优化算法中，对问题的解空间进行寻优的核心操作是趋化操作，包括旋转和游动。细菌向任意方向移动单位步长定义为旋转。旋转后计算适应度值，若得到的路径适应度值更小则继续沿着该方向游动，直到达到最大步长或者适应度值不发生改变，此过程称为游动。由于细菌的步长和方向会影响到寻优效果，因此本文结合VRPTW 问题的特性，将大邻域搜索中的removal 算子结合到趋化操作中。为达到更好的搜索效果，共设计4 种removal算子，分别为Random removal算子、Least profit removal 算 子、Route removal 算子和Relatedness removal 算子。这4 种算子作为细菌游动的方向，游动步长设置为m。用轮盘赌方法选择一种removal算子，用来模拟细菌在旋转过程中任意选择的方向，确定方向后进行移除插入操作。如果该算子取得的适应度值变好，则继续采用，直到达到细菌的最大游动次数，此过程模拟了细菌在该方向上的游动。若旋转后适应度值没有改善，则选择另一个removal 算子，即另一个方向，继续进行游动。每个细菌都从旋转操作开始，在达到趋化算子次数前不断反复交替执行旋转和游动操作，更新路径。将removal 算子应用到趋化操作中，扩大搜索范围，有利于局部搜索。在移除了m个客户节点后，运用之前的贪婪插入法重新将移除的客户进行插入。

下面以一个Solomom 数据集中的实例说明利用removal 算子寻优。假设VRPTW 的配送中心和客户节点信息如表1 所示，编号0 表示配送中心，编号1～9表示客户节点，每辆车的载重为100。采用贪心策略插入法生成的细菌个体为S=［｛0，8，3，9，5，0｝，｛0，1，7，2，4，0｝，｛0，6，0｝］，即由三辆车配送，适应度值为328.883 0。removal 算子的移除个数m=3，最大游动次数为4。

表1 配送中心及客户节点信息Table 1 Information of distribution center and customer nodes

假设选定Random removal算子第一次游动移除的客户节点为7、2、8，在之前生成的解S中找到7、2、8 并将其移除，移除后的路径Snew=［｛0，3，9，5，0｝，｛0，1，4，0｝，｛0，6，0｝］。首先，插入客户节点7，在满足约束条件后，在每条路径中插入的具有最小配送成本的最佳位置如图4 所示。客户节点7 插入3 条路径最佳位置的适应度值分别为354.033 4、301.987 5、273.882 6，通过比较得出最小适应度值为273.882 6，因此将客户节点插入路径3，此时Snew=［｛0，3，9，5，0｝，｛0，1，4，0｝，｛0，7，6，0｝］。然后，插入客户节点2，在满足约束条件后，在每条路径中插入的具有最小配送成本的最佳位置如图5所示。客户节点2 插入3 条路径最佳位置的适应度值分别355.298 4、293.456 7、374.964 9，通过比较得出最小适应度值为293.456 7，因此将客户节点2插入路径2，此时Snew=［｛0，3，9，5，0｝，｛0，1，2，4，0｝，｛0，7，6，0｝］。最后，插入客户节点8，在满足约束条件后，在每条路径中插入的具有最小配送成本的最佳位置如图6 所示。客户节点8 插入3 条路径最佳位置的适应度值分别296.257 3、319.856 7、326.907 1，通过比较得出最小适应度值为296.257 3，因此将客户节点8 插入路径1，此时Snew=［｛0，8，3，9，5，0｝，｛0，1，2，4，0｝，｛0，7，6，0｝］。Snew的适应度值为296.257 3，S的适应度值为328.883 0，经过比较，Snew的适应度值较小，则用Snew的位置更新S，继续采用removal算子进行局部寻优，直到达到最大游动次数。如果该removal算子没有改善适应度值，则维持S的位置不变，结束此次游动，细菌进行翻转，用另一个removal 算子进行第二次游动，重复以上操作，直到适应度值不再发生改变。

图4 客户节点7 插入每条路径的最佳位置Fig.4 The best position for inserting client node 7 into each path

图5 客户节点2 插入每条路径的最佳位置Fig.5 The best position for inserting client node 2 into each path

图6 客户节点8 插入每条路径的最佳位置Fig.6 The best position for inserting client node 8 into each path

图7 和图8 分别给出了初始路径和经过趋化操作的改进路径。

图7 初始路径Fig.7 Initial path

图8 经过趋化操作的改进路径Fig.8 Improved path after chemotaxis

4 种removal 算子的具体描述如下：

1）Random removal 算子。该算子是随机选择q个客户节点进行移除，再将其进行插入。该算子的随机性增加了搜索的多样性。

2）Least profit removal 算子。该算子主要移除成本较高的客户，给定一个客户i和一个解s，定义客户节点i的成本为Ci(s)=f(s)-f-i(s)，其中f-i(s)表示移除客户节点i之后的适应度值，对C(s)值进行降序排列，选择成本较高的q个客户节点进行移除。

3）Route removal 算子。该算子主要是减少寻优过程中车辆的使用数目，选取配送客户节点较少的车辆中的客户节点进行移除。在采用Route removal算子进行局部搜索过程中，可以发现某些车辆所服务的客户节点数量非常少，这样对车辆资源造成较大的浪费，也不利于总配送路径的减少。因此，在Route removal 算子中应该优先移除这些车辆配送数目较少的客户节点，再将移除的客户节点用贪婪插入策略重新插入路径中。

4）Relatedness removal 算子。该算子移除相关性较高的客户节点，根据两个客户节点i和j的距离Cij、送货需求的差异|Di-Dj|以及两个客户的最早开始时间|Ei-Ej|的差异来衡量两个客户节点i和j的相关性。每个局部相关性度量分别使用α、β、γ进行加权，并对问题实例给出的所有客户S集合的各个极值进行归一化。因此，两个客户节点i和j的相关性度量Rij表示如下：

其中：max（Cij）是任意两个客户之间的最大距离；max（Di）是交付的最大需求值；min（Di）是交付的最小需求值；max（Ei）是客户开始服务时间的最大值，min（Ei）是客户开始服务时间的最小值。对任意两个客户节点间的相关性进行排序，移除相关性较大的q个客户节点，重新将它们插入到路径中，从而得到新的解决方案。

2.4 算法流程

改进的细菌觅食算法步骤具体如下：

1）初始化参数。初始化细菌种群数目N、趋化操作次数Nc、最大游动次数Ns、复制次数Nre、迁徙次数Ned、迁徙概率Ped、removal 算子中的移除个数q，根据基于K-means 的贪婪插入法构造初始解。

2）对种群中的所有细菌个体进行趋化操作。

3）对细菌的适应度值进行评价并按照升序排列。将前面N/2 个细菌个体的位置复制给排列在后面的N/2 个细菌个体。

4）为每个细菌个体随机生成一个概率rand。如果rand 小于Ped，则随机生成Nnum到2 的全排列数据，将生成的数据按照贪心策略插入法生成细菌个体；如果rand大于等于Ped，则维持原细菌个体位置不变。

5）重复步骤2～步骤4，直至达到最大迁徙次数Ned。

3 实验与结果分析

3.1 实验环境

实验中使用的计算机配置为2 GHz CPU、8 GB RAM、64 位操作系统。编程语言使用Matlab 语言，版本为Matlab r2017b。

3.2 算例测试与比较

为评价改进算法的性能，使用Solomon 数据集进行实验，将所得最短配送距离与对比算法进行比较。Solomon 数据集的VRPTW 共分为C1、C2、R1、R2、RC1、RC2 等6 类。C1 和C2 类的客户呈 现分块聚集分布；R1 和R2 类的客户是随机分布的；RC1 和RC2 类的客户既有随机分布，又有分块聚集分布。R1、C1、RC1 数据集的车辆容量小、时间窗口窄，每条路线只允许少数客户；R2、C2 和RC2 数据集的车辆容量大，且在较长的计划周期内允许许多客户使用同一辆车进行服务。本文算法在每个算例上独立运行10 次。将实验结果与文献［20］中的具有交叉操作的人工蜂群（Artificial Bee Colony with Crossover operation，ABC-C）算法、具有扫描策略的人工蜂群（Artificial Bee Colony with Scanning strategy，ABC-S）算法和改进人工蜂群（Improved Artificial Bee Colony，IABC）算法所给出的算例进行比较，将全部算例分别与文献［21］中的混合遗传算法（Hybrid Genetic Algorithm，HGA）和文献［22］中的自适应记忆算法（Adaptive Memetic Algorithm，AMA）进行比较。具体比较情况如表2～表4 所示，其中“—”部分表示对比算法未计算相应数据。改进的细菌觅食算法中的参数设置为：细菌种群数目N=30，趋化操作次数Nc=50，最大游动次数Ns=3，复制次数Nre=5，迁徙次数Ned=2，邻域搜索长度W=10。

表2 C 类算例测试结果比较Table 2 Comparison of test results in C class example

表3 R 类算例测试结果比较Table 3 Comparison of test results in R class example

表4 RC 类算例测试结果比较Table 4 Comparison of test results in RC class example

在与ABC-C、ABC-S 和IABC 算法的比较中：由表2 可以看出，对于C1 类算例，IBFA 算法除了C101求解性能较弱，其他算例均优于对比算法；由表3 可以看出，对于R1 类算例，IBFA 算法除了R101、R102、R103、R109、R110 和R112 的求解性能较弱以外，其他算例都优于对比算法；由表4 可以看出，对于RC类算例，IBFA 算法除了RC201 算例求解性能较弱以外，其他算例都优于对比算法。以上结论说明IBFA算法是有效的。

在与HGA 算法的比较中：由表2 可以看出，IBFA 算法对于C 类算例的求解性能均优于对比算法；由表3 可以看出，IBFA 算法对于R 类算例的求解效果也优于对比算法；由表4 可以看出，在RC 类算例 中，IBFA 算法除 了RC103、RC104 和RC107 求 解性能弱于对比算法以外，其他算例均优于对比算法。整体而言，IBFA 算法的改进效果较好。

在与AMA 算法的比较中：由表2 可以看出，在C 类算例中，IBFA 算法与对比算法取得了相同的效果；由表3 可以看出，在R 类算例中，IBFA 算法仅有R106 的求解性能劣于对比算法，其他算例均优于对比算法；由表4 可以看出，在RC 类算例中，IBFA 算法仅有RC103、RC104、RC107 和RC108 的求解性能劣于对比算法，其他算例均优于对比算法。整体而言，IBFA 算法的寻优效果较好。

为进一步验证IBFA 算法的有效性，同样使用Solomon 数据集进行测试，将其与已知最优的启发式结果进行比较，已知最优启发式结果中的数据从http：//w.cba.neu.edu/～msolomon/heuristi.htm 中选取。在表5～表7 中，TD 为已知启发式最优解，IBFA 为本文算法得到的最优解，GAP 表示本文算法最优解与已知启发式最优解之间的百分比偏差，用于评估算法的质量。

表5 C 类算例中IBFA 与已知最优解的比较Table 5 Comparison of IBFA and known optimal solutions in C class example

表6 R 类算例中IBFA 与已知最优解的比较Table 6 Comparison of IBFA and known optimal solutions in R class example

表7 RC 类算例中IBFA 与已知最优解的比较Table 7 Comparison of IBFA and known optimal solutions in RC class example

由表5 可以看出，对于C1 和C2 类算例，IBFA算法与已知启发式最优解之间没有百分比偏差。由表6 可以看出，对于R 类算例，IBFA 算法除了R106 的百分比偏差是正值以外，其余均是负值，范围为-15.03% 到-0.18%，证明求解质量较好。由表7 可以看出，对于RC 类算例，IBFA 算法中RC103、RC104、RC107 和RC108 的百分比偏差为正值，其余均是负值，范围为-19.59%到-0.43%，证明求解效果较好。

图9 给出了IBFA 算法与其他算法在不同算例上的最短配送距离比较结果。由图9（a）和图9（b）可以看出，在C 类算例上，IBFA 算法基本取得了最佳效果。由图9（c）和图9（d）可以看出，改进细菌觅食算法整体上优于其他算法。由图9（e）和图9（f）可以看出，RC1 类算例RC103、RC104、RC107 和RC108 效果比较差，其他算例效果整体较好。以上结论证明了IBFA 算法改进效果较好。图10 给出了IBFA 算法的部分算例的路径规划结果。

图9 6 类算例上的算法最短配送距离比较Fig.9 Comparison of the shortest distribution distances of algorithms on six class examples

图10 部分算例的IBFA 路径规划结果Fig.10 IBFA path planning results of some examples

3.3 算法运行时间比较

为验证算法性能，将本文IBFA 算法与ABC-C、ABC-S 和IABC 算法进行运行时间对比，依据文献［14］中的数据经过计算得到：IBFA 算法在C 类、R 类以及RC 类算例的平均运行时间分别为61.33 s、119.42 s 和114.63 s；ABC-C 算法在C 类、R 类以及RC 类算例的平均运行时间分别为77.67 s、218.33 s和185.75 s；ABC-S 算法在C 类、R 类 以及RC 类算例的平均运行时间分别为74.11 s、213.5 s 和178.13 s；IABC 算法在C 类、R 类以及RC 类算例的平均运行时间分别为68.44 s、199.92 s 和168.75 s。从平均运行时间上看，IBFA 算法在各类算例上的平均运行时间均优于ABC-C、ABC-S 和IABC 算法，证明IBFA算法运行效率较高。

4 结束语

本文提出一种改进的细菌觅食算法，用于求解带时间窗的车辆路径优化问题。将大邻域搜索方法与趋化操作相结合，设计4 种removal 算子，扩大搜索范围，提高求解质量。实验结果证明，改进算法具有较强的寻优能力，在不增加时间复杂度的基础上，能够得到更优的配送路径。下一步将从环境保护角度出发，设计减少碳排放、提升客户满意度、降低配送成本的多目标函数，研究低碳环保的带时间窗车辆路径规划问题。