开放科学（资源服务）标志码（OSID）：

0 概述

近年来，车辆与无人机协同配送成为应急救援的重要手段，在新冠肺炎疫情下的医疗物资配送^[1]、西昌泸山森林火灾救援^[2]等事件中均有应用。其中，合理规划车辆与无人机路径对于提高应急救援的时效性至关重要，这也使得相应的路径规划问题成为一个备受关注的研究热点。

车辆与无人机协同路径规划的研究尚处于起步阶段，国内外学者对该问题的研究主要集中在带无人机的旅行商问题（Travelling Salesman Problem with Drone，TSPD）、带无人机的车辆路径问题（Vehicle Routing Problem with Drones，VRPD）及其衍生问题的建模与求解算法上。在对TSPD的拓展研究中，文献[3-4]分别探讨了无人机与车辆在途结合问题以及单车辆多无人机的TSPD变体。文献[5]建立高寒山地环境下的车机协同物资配送模型，并提出相应的启发式求解算法。文献[6]提出一种车辆-无人机串联交付系统优化模型，其利用K-means算法寻找无人机发射点并利用遗传算法完成车辆路径构造。文献[7]提出带无人机站点的车辆路径规划问题（TSP-DS），文献[8]考虑多车辆与异构无人机并提出多车辆带无人机旅行商问题（mFSTSP）。在考虑车辆容量限制的情况下，TSPD会转化为更加一般化的VRPD。文献[9-10]在最坏情况下分析VRPD中无人机速度与数量这2个关键因素。文献[11-12]考虑VRPD中的单车多无人机变体以及车机在途结合变体。文献[13]介绍带无人机和时间窗的车辆路径问题。文献[14-15]分别考虑无人机补给的车辆路径问题与无人机取货的车辆路径问题。此外，文献[16]在其研究中创新性地提出用强化学习方法来解决车机组合路径规划问题，这为利用人工智能技术^[17-18]解决复杂车机协同问题提供了新思路，在该问题中，车辆为无人机提供起降与补给服务，但不负责配送任务。

综上，车机协同路径规划及其衍生问题目前已取得一定的进展，且现有研究主要关注无人机每次配送单个用户的情形。然而，随着无人机技术的不断发展，支持单次多用户配送的高负载无人机开始出现，这使得更加灵活的车机协同方式成为可能，亟待开展相应路径规划问题的建模与求解工作。本文对VRPD进行拓展，提出一种车辆与无人机协同路径规划问题，即考虑多投递的带无人机车辆路径规划问题（Vehicle Routing Problem with Drones Considering Multi-Delivery，VRPD-MD）。在VRPD-MD中，无人机可以在任意架次中为多个受灾区域配送物资，而车辆在配送物资的同时支持无人机完成多次起飞和降落。本文建立以最小化最大完工时间为目标函数的混合整数规划模型，并设计一种基于遗传方法的自适应算法（Adaptive Algorithm Based on Genetic Method，AAGM）对该模型进行求解。

1 数学模型 1.1 问题描述

在VRPD-MD中，车辆与无人机协同完成多个受灾点配送任务，无人机在任意架次中服务多个受灾点，车辆在配送的同时支持无人机多次起降，目的是寻找耗时最短的配送路径。本文所提出的车机协同模式更加贴近实际配送需求，可以充分发挥无人机的速度优势，同时车机并行交付方式也能有效减少总体交付时间。为了建立数学模型，本文作出如下合理假设：

1）每辆车搭载一架无人机，且无人机必须从车辆起飞并返回起飞车辆。考虑车辆与无人机的容量约束，无人机电量限制其续航时长。

2）车辆与无人机行驶在欧氏空间，且只能在受灾点交会。

3）车辆与无人机必须在交会位置进行时间同步，即较早到达交会节点的一方需要等待另一方。

4）有足够的电池可供无人机替换，并且忽略电池更换和电池充电的时间。

1.2 问题建模

为便于模型描述，表 1给出建模过程中所应用的符号及其含义。

基于上述符号和变量，本文建立如下混合整数规划模型：

$ {\rm{Min\;}}z=\sum\limits_{k\in {K}_{T}{}^{}}{T}_{0\left(r\right)}^{k} $

(1)

s.t.

$ \sum\limits_{k\in {K}_{T}}{v}_{i}^{k}+\sum\limits_{k\text{'}\in {K}_{D}}{d}_{i}^{k\text{'}}=1, \forall i\in {V}_{C} $

(2)

$ \sum\limits_{j\in {V}_{C}\bigcup 0\left(r\right)}{x}_{0\left(s\right), j}^{k}=\sum\limits_{j\in {V}_{C}\bigcup 0\left(s\right)}{x}_{j, 0\left(r\right)}^{k}=1, \forall k\in {K}_{T} $

(3)

$ \sum\limits_{i\in {V}_{C}\bigcup 0\left(s\right)}{x}_{i, j}^{k}=\sum\limits_{h\in {V}_{C}\bigcup 0\left(r\right)}{x}_{j, h}^{k}={v}_{j}^{k}, \forall k\in {K}_{T}, \forall j\in {V}_{C} $

(4)

$ \sum\limits_{k\text{'}\in {K}_{D}}\sum\limits_{j\in {V}_{C}}{y}_{i, j}^{k\text{'}}=1, \forall i\in {V}_{L} $

(5)

$ \sum\limits_{k\text{'}\in {K}_{D}}\sum\limits_{j\in {V}_{C}}{y}_{j, i}^{k\text{'}}=1, \forall i\in {V}_{R} $

(6)

$ {d}_{j}^{k\text{'}}\left(\sum\limits_{i\in {V}_{C}{}_{}}{y}_{i, j}^{k\text{'}}-1\right)=0, \forall k\text{'}\in {K}_{D}, \forall j\in {V}_{C} $

(7)

$ {d}_{j}^{k\text{'}}\left(\sum\limits_{h\in {V}_{C}}{y}_{j, h}^{k\text{'}}-1\right)=0, \forall k\text{'}\in {K}_{D}, \forall j\in {V}_{C} $

(8)

$ \sum\limits_{i\in {V}_{C}}{D}_{i}{v}_{i}^{k}+\sum\limits_{i\in {V}_{C}}{D}_{i}{d}_{i}^{k\text{'}}\le Q, \forall k\in {K}_{T}, \forall k\text{'}\in {K}_{D} $

(9)

$ \sum\limits_{i\in {V}_{C}}{D}_{i}{d}_{i}^{k\text{'}}{y}_{ipk\text{'}} < q, \forall k\text{'}\in {K}_{D}, p\in {R}_{k\text{'}} $

(10)

$ \sum\limits_{i\in {V}_{C}}\sum\limits_{j\in {V}_{C}}{y}_{i, j}^{k\text{'}}{\tau }_{i, j}^{k\text{'}}{y}_{ipk\text{'}}{y}_{jpk\text{'}}\le B, \forall k\text{'}\in {K}_{D}, p\in {R}_{k\text{'}} $

(11)

$ \sum\limits_{h\in {V}_{C}}{y}_{i, h}^{k\text{'}}({T}_{i}^{k}-{t}_{i}^{k\text{'}})=0, \forall k\in {K}_{T}, \forall k\text{'}\in {K}_{D}, \forall i\in {V}_{L} $

(12)

$ \sum\limits_{l\in {V}_{C}}{y}_{l, j}^{k\text{'}}({T}_{j}^{k}-{t}_{j}^{k\text{'}})=0, \forall k\in {K}_{T}, \forall k\text{'}\in {K}_{D}, \forall j\in {V}_{R} $

(13)

$ {T}_{j}^{k}\ge {T}_{i}^{k}+{\Gamma }_{i, j}^{k}-M(1-{x}_{i, j}^{k}), \forall i, j\in V, \forall k\in {K}_{T} $

(14)

$ {t}_{j}^{k\text{'}}\ge {t}_{i}^{k\text{'}}+{\tau }_{i, j}^{k\text{'}}-M(1-{y}_{i, j}^{k\text{'}}), \forall i, j\in {V}_{D}, k\text{'}\in {K}_{D} $

(15)

$ {v}_{i}^{k}, {d}_{i}^{k\text{'}}, {x}_{i, j}^{k}, {y}_{i, j}^{k\text{'}}, {y}_{ipk\text{'}}\in \left\{\mathrm{0, 1}\right\} $

(16)

$ {T}_{i}^{k}, {t}_{j}^{k\text{'}}\ge 0\text{，}\forall i\in V\text{，}\forall j\in {V}_{T}\bigcup {V}_{D}\text{，}\forall k\in {K}_{T}\text{，}\forall k\text{'}\in {K}_{D} $

(17)

其中：目标函数式（1）表示最小化所有车辆返回仓库的运行时间；约束条件式（2）确保每一个受灾点必须被车辆或无人机访问；约束条件式（3）为车辆在仓库节点的流量约束；约束条件式（4）保障车辆在受灾点的流量平衡；约束条件式（5）、式（6）分别保障无人机在起飞节点与降落节点的流量平衡；约束条件式（7）、式（8）保障无人机在提供服务的受灾点流量平衡；约束条件式（9）、式（10）分别为车辆容量约束与无人机架次容量约束；约束条件式（11）为无人机架次续航约束；约束条件式（12）、式（13）将车辆与无人机在起飞节点和降落节点的时间调整为一致；约束条件式（14）、式（15）为车辆与无人机的运行时间不等式，其中，无人机的起飞和降落消耗一个单位时间，该时间分别被计算到车辆与无人机的运行时间内；约束条件式（16）、式（17）表示所有变量的取值约束。

2 基于遗传思想的自适应算法设计

遗传算法^[19-20]是求解无人机路径规划问题的常用方法之一，本文考虑VRPD-MD的特点，设计一种基于遗传思想的自适应算法AAGM。

2.1 染色体编码

在AAGM算法中，一个解被表示为由多条车辆与无人机组合路线组成的染色体。本文采用自然数编码方式，染色体表示卡车和无人机的访问序列，每个基因位置存储一个受灾节点。如图 1所示，车机组合1表示一个车辆与无人机组合路径，整条染色体包括多条组合路径，当需要执行算子操作时对其进行解码操作。

解码操作分为两步：第一步获取车机组合路径；第二步确定节点属性。以图 1为例，首先得到2条车机组合路径，组合路线中第一个重复出现的节点（如车机组合2中的节点1）是无人机第一架次的起飞节点，第一架次的降落节点是第二个重复出现的点（如车机组合2中的节点3）。同样，第三和第四个重复出现的节点（节点9和节点13）是无人机第二架次的起飞点与降落点。通过该对比选择规则，可以识别出每个车机组合路线上的所有交会节点与非交会节点。

2.2 初始种群生成

本文采用先聚类后构造路径的思想生成初始种群中的个体，具体步骤如下：

步骤1  利用K-means聚类得到满足车辆载重约束的多个受灾点簇，结合旅行商问题，利用随机规则生成各簇相应的车辆路径。

步骤2  选择一条车辆访问路径，随机选择一个节点作为无人机第一架次的起飞节点。

步骤3  判断当前起飞节点是否为当前车辆路径最后的访问点，如果是，则进入步骤8；否则，进入下一步。

步骤4  选择当前车辆路径上2个连续受灾节点加入无人机路径，将后续节点设置为无人机降落节点，若所选节点包括最后的访问节点，则将该最后访问节点设置为降落节点。

步骤5  判断式（10）与式（11）是否均满足，若是，则进入下一步；否则，撤销当前所构建的无人机路径，选择当前起飞点的相邻后续节点作为新的起飞点。

步骤6  生成一条有效的无人机路径，并将当前无人机路径降落点的相邻后续节点设置为下一架次的起飞节点。

步骤7  判断起飞点是否为当前车辆路径的最终访问节点，若是，则进入下一步；否则，返回步骤4。

步骤8  返回当前车辆路径相应的车机协同配送路径。

步骤9  判断所有车辆路径是否完成重构，若是，则结束；否则，返回步骤2。

2.3 遗传算子

本文针对VRPD-MD的特点，设计2类内部搜索算子以生成邻域解：访问节点交叉算子用于对车辆访问节点与无人机访问节点进行交换；交会节点变异算子用于对无人机起飞节点与降落节点进行调整。

2.3.1 访问节点交叉算子

访问节点交叉算子基于同一车机组合路线进行内部变换，通过从车辆路径与无人机路径中随机选择不同数量的访问节点进行交换以生成邻域解。图 2所示为Truck（1）-Drone（1）算子的一种操作过程，从无人机路径中选择节点2、车辆路径中选择节点5进行交换，得到一个新的邻域解。基于此规则，根据交换节点数量的不同，设计Truck（2）-Drone（2）、Truck（0）-Drone（1）、Truck（1）-Drone（0），加上Truck（1）-Drone（1）共4个访问节点交叉算子，其中数字对应选择交换的节点数目。

2.3.2 交会节点变异算子

交会节点变异算子用来对无人机起飞节点与降落节点进行调整以获得邻域解，包括起飞节点变异（Off-Change）与降落节点变异（Down-Change）2类算子。图 3所示为Off-Change的一种操作过程。

在图 3中，无人机起飞点由节点1调整至节点2，类似地，Down-Change可以实现降落点调整操作。需要注意的是，算子每次仅执行一个移动操作，即调整起飞点（或降落点）为其相邻节点。

2.4 算子自适应选择机制

为了提升算法的收敛速度与求解质量，本文设计一种算子自适应选择机制。通过为不同算子设置不同的权重，在调用算子时根据其权重选择当前最优算子，并根据优化结果对算子权重进行动态更新。

算子权重动态更新为自适应选择机制的关键，根据所得邻域解的质量更新算子权重：若某算子所得邻域解优于原解，则将算子的权重增加0.7（算子初始权重为1）；若所得邻域解为劣解，以Metropolis规则接受该解，并将算子的权重增加0.5；如果所得邻域解被舍弃，则算子权重保持不变。

2.5 算法整体步骤

AAGM算法的整体步骤如下：

步骤1  设置控制参数，如种群规模$ {P}_{\mathrm{n}\mathrm{u}\mathrm{m}} $、最大迭代次数$ {N}_{\mathrm{m}\mathrm{a}\mathrm{x}} $、当前温度$ {T}_{s} $、降温系数△、终止温度$ {T}_{e} $等。

步骤2  生成算法初始化种群，令种群最优解为全局最优解。

步骤3  判断$ {T}_{s} $是否小于$ {T}_{e} $，若是，进入下一步；否则，输出全局最优解。

步骤4  判断是否达到$ {N}_{\mathrm{m}\mathrm{a}\mathrm{x}} $，若是，进入步骤10；否则，进入下一步。

步骤5  执行第$ n $（$ n=\mathrm{1, 2}, \cdots $）次进化，初始化算子权重为1。利用轮盘赌策略结合个体适应度（目标函数值的倒数），并根据优化个体占种群的比例挑选出进化种群。

步骤6  选择进化种群中的第$ j $（$ j=\mathrm{1, 2}, \cdots $）个个体，利用轮盘赌策略结合算子权重选择执行算子，利用算子生成当前个体的邻域解。

步骤7  若邻域解优于全局最优解，则令邻域解为全局最优解；若邻域解优于当前个体，则将进化种群中的当前个体替换为邻域解；若邻域解为劣解，则利用Metropolis规则接受该解，根据算子自适应选择机制动态更新算子权重。

步骤8  判断进化种群中所有个体是否都完成进化，若是，进入下一步；否则，$ j=j+1 $，返回步骤6。

步骤9  更新进化迭代次数：$ n=n+1 $。更新初始种群：从当前初始种群中随机选择一定数量的个体，将当前进化种群规模扩充至$ {P}_{\mathrm{n}\mathrm{u}\mathrm{m}} $，并设置为新的初始种群，返回步骤4。

步骤10  更新当前温度：$ {T}_{s}={T}_{s}·{\mathit{\Delta}} $，返回步骤3。

3 实验结果与分析

本文利用不同的实例集进行实验，这些实例从考虑容量限制的车辆路径规划标准实例集转化而来^[21]。将实例中的客户点视为受灾点，保留其原始的车辆和客户信息，在此基础上增加无人机信息，其中，无人机的速度是卡车的1.5倍，其他参数设为一个合理值。算法运行环境为2.9 GHz、Intel Core i10、8 GB RAM、Windows 10、64位系统的计算机，通过Java编程实现。

3.1 配送模式对比实验

本组实验一方面验证车机协同配送模式相比纯车辆配送模式的优势，另一方面通过比较几种不同的车机协同模式证明本文所提多架次多投递无人机配送模式的有效性。通过修改初始解生成策略，本文增加无人机多架次单投递、无人机单架次多投递2种常见车机协同模式。从CVRP算例中选择不同问题规模的15个实例，在不同模式下基于每个实例运行10次AAGM算法，并记录最优解、平均解、平均求解时间、平均解与CVRP最优解的GAP值，具体结果如表 2所示，其中，n表示节点数，k表示所需车辆数。从表 2可以看出，3种车机协同模式下平均目标值在大多数实例中均优于CVRP最优解，其中，多架次多投递模式下的车机协同取得最佳结果。对比多架次单投递与多架次多投递模式，可以看出，具有多点投递能力的无人机配送模式配送效率提升明显。多架次单投递与单架次多投递模式的对比结果显示，在本文实验案例规模下，增加无人机单架次投递节点数优于增加无人机架次数。当扩展到更大规模的问题时，多架次多投递模式相较单架次多投递可能具有更大的优势。

3.2 算法性能对比实验

为进一步验证AAGM算法的性能优势，将其与已有算法进行对比。由于无人机多架次多投递模式下的车机协同问题尚未有相关求解算法，因此无法直接进行算法比较。文献[22]提出了求解无人机单架次多投递模式下车机协同路径规划问题的大规模邻域搜索（LNS）算法，因此，本文采用与LNS算法相同的测试用例，在无人机单架次多投递模式下开展算法对比实验。基于每个实例运行10次AAGM算法，并记录最优解、平均求解时间、最优解与LNS最优解的GAP值，具体结果如表 3所示。从表 3可以看出，AAGM在大多数算例（8/10）中的最优解优于对比算法LNS，平均差距为-3.51%。在求解速度上，AAGM同样具有明显优势，对于大多数算例，其能够在10 s内完成计算，而LNS需要的计算时间则超过1 min。

3.3 自适应策略性能验证

本组实验的目的：一是通过灵敏度分析进行算法参数调优；二是研究自适应策略对于AAGM算法性能的影响。对比最优参数配置的非自适应算法（None-Adaptive Algorithm Based on Genetic Method，NAAGM）与非最优参数配置的AAGM算法在求解质量与求解速度上的差异。在参数调优环节，首先根据经验确定参数的测试区间，通过选取不同规模算例更改参数组合以进行多组测试，得到最佳参数配置如表 4所示。分别在15个测试用例上运行AAGM与NAAGM各10次，记录平均运行时间与最优解。图 4（a）展示AAGM与NAAGM在求解时间上的对比，可以看出，在所有测试算例上，增加自适应策略后的算法求解时效性明显提升，求解时间平均节省30%。图 4（b）展示自适应策略对于算法求解质量的影响，从中可以看出，在大多数算例（11/15）中，AAGM算法的最优解优于NAAGM算法，求解质量平均提升1.83%。

4 结束语

本文研究一种考虑多点投递的带无人机车辆路径规划问题，针对该问题建立相应的混合整数规划模型。为对模型进行求解，设计一种高效的AAGM求解算法，通过集成算子自适应选择机制与基于Metropolis规则的劣解接受机制，以提高算法的全局寻优性能。利用改进的CVRP算例集进行3组验证实验，结果表明，多架次多投递车机协同模式在带无人机的车辆路径规划任务中可以发挥优势，且AAGM算法在求解VRPD-MD时具有可行性和有效性。多点投递的车机协同路径规划研究目前尚处于起步阶段，本文下一步将从2个方面展开研究：一是探索更高效的求解算法，如设计基于车辆与无人机路径同步生成规则的启发式算法，利用该算法生成初始解；二是研究多种变体问题，如无人机跨车辆停靠的车机协同配送问题。