改进混合遗传算法在MTVRPTW中的建模与优化

宋 强

(广东理工学院 信息工程系，广东 肇庆 526100)

0 引 言

众所周知的车辆路径问题(vehicle routing problem，VRP)是一个NP-Hard的组合优化问题，在一些城市配送系统中，货物被送到城市配送中心(city distribution center，CDC)，然后由带有容量限制的车辆交付给最终客户。这些车辆可能在工作日结束之前更早地返回仓库，然后用于另一个送货行程。这就引入了多行程车辆路径问题的研究思路。

在多级配送系统中，尤其是在城市环境下，物流卡车不能直接给客户送货，而是先把货物送到CDC，在CDC有现货的时间段内再进行最终的配送。通常情况下，考虑到城市路况的限制，只有采用电动三轮车或者小型货车才能执行送货的要求，当货物在客户点卸载后必须立即重新装载。同时，货物在工作日全天都可以被送到城市配送中心，这意味着在配送周期开始时这些货物在CDC的仓库是没有现货的。发货时间的概念是与每个货物相关联的，这里假定货物的发货时间在工作日开始前是已知的。

在国内，采用遗传算法求解带时间窗的VRP的研究并不少见。例如，葛显龙等[1]基于第三方带软时间窗约束的车辆路径模型，设计了一种改进遗传算法对模型进行求解；王永锋等[2]根据遗传算法容易产生早熟、收敛速度慢、局部寻优能力差等缺点，提出了将混沌搜索技术和遗传算法相耦合的混沌遗传算法来求解带时间窗的物流配送车辆路径问题；潘立军等[3]设计了基于时差的插入法，提出了一种求解带时间窗的取送货问题的遗传算法。

而针对带时间窗的多行程VRP，国内研究文献数量较少。例如许争争等[4]以车辆容量和绕行时间限制为约束参数，以面向一个城市的集中式顾客接送服务为背景，设计了一种基于协作的3阶段算法来求解多行程路径问题；宋强[5]为了同时解决多行程车辆路径问题和配送中心定位问题，采用模拟退火的算法来获得最优路线，并改善了当前温度控制的位置。

建立了车辆之间在时间上相互依赖的数学模型，引入混合遗传算法求解该问题，介绍了与货物有关的发货时间的概念。在相关的国外文献中较少有人采用改进的混合遗传算法对带时间窗的多行程车辆路径问题(multi-trip vehicle routing problem with time windows, MTVRPTW)进行研究，已有的研究大都采用精确算法来求解问题。

例如，N. AZI等[6]提出了一个精确算法来解决单一车辆的MTVRPTW，该解决方法开发了一个带资源限制的初级最短路径算法；另外N. AZI[7]还采用了一个嵌入到分枝定价算法中的列生成方法来求解MTVRPTW，由于算法的局限性，该算法只能解决客户数量50以内的实例；F. HERNANDEZ等[8]使用和N. AZI类似的方法，为每个问题给出一个集合覆盖公式，每一列表示一个行程而不是一个工作日；以上这3种精确方法都使用了F. DOMINIQUE等[9]提出的用于初级最短路径问题的算法。

除此以外，M. BATTARRA等[10]研究了MTVRPTW的扩展问题，不同的货物群集在一起，并且不能用相同的车辆运输，他们提出的方法是一个两级顺序启发式算法：独立地考虑每个货物并制定一组可行的路线，然后，再把路线分配给车辆；D. CATTARUZZA等[11]采用遗传算法来解决多行程车辆路径问题，并采用Adsplit算法来评估染色体，引入的局部搜索算子用于对车辆重新分配行程；M. DREXL[12]介绍了车辆路径调度在同歩性方面的问题，并引入了多行程的送货方法。

但是以上这些启发式算法均未能同时把时间窗、发货时间和多行程的因素全部纳入约束条件，而这3个约束条件正是笔者研究的重点。讨论了城市环境下带时间窗的多行程车辆路径问题(multi-trip vehicle routing problem with time windows，MTVRPTW)，结构如下文如述：第2节对这个问题进行了定义；第3节提出了一种基于局部搜索的混合遗传算法(hybrid genetic algorithm, HGA)来求解提出的MTVRPTW；第4节给出了仿真实验的研究结果；第5节进行了全文总结。

1 问题的定义和相关符号

MTVRPTW可以在一个完整的无向图G=(V,E)上进行定义，V={0,…,N}表示顶点的集合，E={(i,j)|i,j∈V,i

这里给定一个时间范围TH来定义工作日的持续时间，它可以被看作是与仓库相关的时间窗。因而假设[E0,L0]=[0,TH]，并且D0=0。MTVRPTW需要确定一组路线并把每条路线分配给每辆车，这样总的行驶时间就会最小并满足下列条件:

1)每条路线都从仓库开始和结束;

2)每个路线的开始行驶时间不早于与客户相关的每条路线本身的最大发货时间Ri;

3)每个客户只能被一个路线访问;

4)客户i的服务从相应的时间窗TW[Ei,Li]开始;

5)客户的需求总和在任何路线上都不超过车载容量Q;

6)每辆车行驶完最后的路线返回仓库的时间不迟于TH。

假定每个客户i可以被回程车辆服务，也就是说，Ri+T0i+Ti0≤TH和Di≤Q，否则不存在可行的解决方案。

符号σ表示一条路线。大写字母∑表示分配给一辆车的路线集合，也被称为旅程，旅程是由不同的路线(v1,…,vn)⊕(w1,…,wm)组成，符号⊕表示路径的连接。例如(v1,…,vn)⊕(w1,…,wm)表示两条路径的连接，如果v1和wm都表示仓库的话，就会产生一条闭合路线。σ1⊕σ2意味着同一辆车先执行路线σ1后执行σ2。在路线σi上，车辆访问客户的时间用τσi表示，车辆离开仓库后在路线σi上的行驶时间由Tσi表示。

符号“∈”表示“属于”。例如σ∈∑意味着路线σ属于旅程∑，v∈σ意味着客户v属于路线σ。但值得注意的是，在不考虑发货时间时，满足条件Tσi+1=Tσi+τσi。另一方面，在考虑到发货时间时，应满足条件Tσi+1=max{Tσi+τσi,maxvi∈σi(Rvi)}。

如果车辆在Li之前到达客户vi的位置，该路径规划就是可行的。

2 方 法

提出了一种改进的混合遗传算法 (hybrid genetic algorithm, HGA)用于求解MTVRPTW。该解决方案的主要思路可以用算法1表示。该算法是针对D. CATTARUZZA等[13]提出的基因算法解决方案的一个改进，主要改进之处在于选择双亲染色体产生子染色体后，基因算法对子染色体进行训练以达到可行的效果。而本算法是采用了局部搜索的方法对子染色体进行改进以组成新的种群。

算法1：改进的混合遗传算法

1:初始化种群

2:while不满足终端规则do

3:选择双亲染色体Ψ=(Ψ1,…,ΨN)和Ψ=(Ψ1,…,ΨN)

4:产生一个子染色体Ψ=(Ψ1,…,ΨN)

5:用LS改进Ψ=(Ψ1,…,ΨN)

6:IfΨ=(Ψ1,…,ΨN)是不可行 then

7:修复Ψ=(Ψ1,…,ΨN)

8:end if

9:把Ψ=(Ψ1,…,ΨN)插入到种群

10：if 种群大小Ψ=(Ψ1,…,ΨN) then

11:选择幸存者

12:end if

13:end while

首先，π随机染色体的生成是为了创建一个初始种群，这些个体采用局部搜索(local search，LS)进行了改进。采用传统的二进制锦标赛过程进行选择，通过对双亲进行交叉选择产生了下一代，修复了不可行染色体。当种群达到π+μ维度时，按照质量和对种群的多样化贡献相应地消除μ个体(正如T. VIDAL等[14]提出)。在后面的内容中提供了关于LS方法的细节，并采用辅助分割过程从一个染色体获得MTVRPTW的解决方案。

2.1 解决方案的表示和搜索空间

一个染色体是由N个客户节点组成的一个不带行程分隔符的有序排列Ψ=(Ψ1,…,ΨN)，Ψ可以视为一个旅行商问题(traveling salesman problem，TSP)的解决方案，通过分裂染色体(插入行程分隔符和分配行程给车辆)被转化为一个可行的MTVRPTW解决方案。从这一点来看，Ψ通常被称为一个巨网分割结构，即由车辆行驶路线和客户节点构成一个巨网，根据Ψ的分割方法来构建不同MTVRPTW的解决方案。

在搜索阶段，在适应度函数中过载和时间窗违规都是允许的，分别通过系数λ和θ进行惩罚。辅助分割过程用于从Ψ中获得一个MTVRPTW解决方案ξ。下面介绍两个符号:T∑(ξ)和TW∑(ξ)分别代表在方案ξ中旅程∑所需的行驶时间和可能产生的时间窗违规。Lσ(ξ)的是路线σ上的载货量，这里σ∈∑ 表示σ是旅程∑的一条路线。染色体Ψ的适应度函数F(Ψ)通过辅助分割可以找到求解最佳解决方案的成本ξ，定义如式(1)：

(1)

在以上表达式中，如果从Ψ中通过辅助分割获得一个可行的解决方案ξ，染色体Ψ就称为可行。

2.2 用于带时间窗VRP的局部搜索

局部搜索在带时间窗的VRP面前要比在经典VRP中变得更加复杂，不能直接在固定时间内检查其可行性。在本研究中，采用了T. VIDALET等[15]提出的方案(即Y. NAGATA等[16]提出的扩展计划)用于求解一大类带时间窗的问题。给定一个路线ρ，数值变量T(ρ)，TW(ρ)，E(ρ)，L(ρ)，D(ρ)分别表示最小持续时间、ρ的最小时间窗违规、从允许最小持续时间和时间窗违规的路线ρ上的第一个客户开始的最早和最晚服务时间，以及被服务客户的累积需求。给定两条路径ρ1和ρ2，持有的关系如式(2)～式(6)：

T(ρ1⊕ρ2)=T(ρ1)+T(ρ2)+Tv1,v2+ΔWT

(2)

TW(ρ1⊕ρ2)=TW(ρ1)+TW(ρ2)+ΔTW

(3)

E(ρ1⊕ρ2)=max{E(ρ2)-Δ,E(ρ1)}-ΔWT

(4)

L(ρ1⊕ρ2)=min{L(ρ2)-Δ,L(ρ1)}+ΔTW

(5)

D(ρ1⊕ρ2)=D(ρ1)+D(ρ2)

(6)

其中， Δ=T(ρ1)-TW(ρ1)+Tv1,v2

ΔWT=max{E(ρ2)-Δ-L(ρ1,0)}

ΔTW=max{E(ρ1)-Δ-L(ρ2,0)}

v1表示ρ1中的最后一个客户，v2表示ρ2中的第一个客户。

上面定义的这些数值变量可用于在连续时间内对典型LS进行评估，在预处理阶段可用于计算所有行驶路径及其返程路径。

2.3 发货时间的介绍和多行程实例的应用

该问题比带时间窗的VRP(vehicle routing problem with time windows, VRPTW)多出两个特点：车辆可以行驶多个行程；在整个工作时间范围内车辆可以完成货物的终端配送。由于车辆将完成的行程不止一个，或者由于车辆需要等待可以运输的货物，车辆在t>0时就可以开始一个行程。T. VIDAL等[15]证明了表达式(2)～式(6)的关系，并通过关联运算证明了式(7)、式(8)之间的关系：

T(ρ)(t)=T(ρ)+max{0,E(ρ)-t}

(7)

TW(ρ)(t)=TW(ρ)+max{0,t-L(ρ)}

(8)

表达式(8)用来计算由于离开仓库时间较晚而造成的时间窗违规。这里引入R(ρ)作为客户在路线ρ上的最大发货时间。我们定义R(vi)=Rvi并且具有下列关系：

R(ρ1⊕ρ2)=max{R(ρ1),R(ρ2)}

(9)

因而，发货时间不会给以后介绍的路径规划增加任何复杂性。

下一步需要考虑的仍然是多行程方面。针对每一条路线σi，行程都是从Tσi=0开始的，通过Tσi=max{R(σi),Tσi-1+τσi-1}，可以计算时间窗违规和完成时间的变化。给定一个旅程∑=(σ1⊕…⊕σn)，式(8)用于计算∑中所有的行程，而Tσi+1可以用式(10)～式(11)计算出来：

(10)

Tσi+1=max{t,R(σi+1)}

(11)

车辆在行程σi∈∑上行驶，可采用式(7)和式(8)，在连续的时间内对完成时间和时间窗违规的变化进行局部计算。而完整的计算需要O(∑max)，其中∑max表示在所有旅程中的最大行程数量。

2.4 修复过程

根据不可行性的特点，用惩罚因子乘以10，采用局部搜索可以修复不可行的染色体。如果由此产生的染色体仍不可行，λ和θ就会重新乘以10并重新应用LS。

2.5 用辅助分割算法求解MTVRPTW

2.5.1 构建辅助图

事实上，不保证一定存在这样一条可行的路径，并允许存在关于时间和负载约束的不可行解。下面用一条弧(i,j)来表示对应的路径成本。

(12)

在这里，路径成本不需考虑行程的具体路线，不考虑由于车辆延迟出发造成的时间窗违规。正如在文献[17]中用于求解VRP的原始过程，考虑到在H上形成最短路径的弧提供的解决方案的成本远高于最佳方案，下面通过一个标记过程来获得在染色体Ψ上构建的用图H表示的用于MTVRPTW的解决方案。

2.5.2 分配行程给车辆

在多行程VRP(multi-trip vehicle routing problem, MTVRP)的环境下，特别是在MTVRPTW环境中，一次可以分配多个行程给同一辆车。时间窗违规完全取决于离开仓库的时间和路线，并且不需要考虑图H上的和弧相关的固定成本。

标记过程可以定义如下：列举出所有在染色体Ψ中从节点0到节点N的可能路径，并构建在节点N上具有较低成本的标签的最佳解决方案。从节点0开始，标签沿着图进行逐步扩展。每个标签L都有M+3个字段域：前面的M个字段以降序存储车辆行驶次数，第(M+1)个字段记录总负载的不可行解，第(M+2)个字段存放前任节点，最后一个字段保存着采用表达式(1)评估的局部解决方案的成本，并作为标签L的成本c(L)。扩展一个标签时，需要构造M个新的标签，每个标签都用于标记每辆车的新行程的位置。当达到节点N时，选择和节点N相关的带有最低成本c(L)的标签L，可以形成相关的解决方案。

为了加快这个标记过程，按照下列规则淘汰支配标签: 让L1和L2成为和同一个节点i∈V′相关的两个标签。当且仅当满足下列条件时，称为L1强势支配L2：

(13)

在式(13)中，c(L)是和标签L相关的成本，δj(L1,L2)=max{0,Tj(L1)-Tj(L2)}，其中Tj(L)是和标签L相关的车辆j的行驶时间。这样，对于给定两个标签L1和L2，如果扩展L1而不以同样方式扩展L2的话，就会受到更多时间窗违规引起的处罚。如果不等式(13)成立，L2就不能采用比L1更好的方法进行扩展，那么L2就会被淘汰。

后面初步的计算实验显示这个过程的时间效率比较低。为此提出了一种启发式版本的支配规则，并引入一个参数γ≥1。当且仅当下列条件满足时，称为L1弱势支配L2：

(14)

c(L1)≤c(L2)

(15)

值得注意的是当γ=1时，弱势支配规则等价于强势的版本。γ>1时，不等式(14)很容易满足，可以删除大量的标签。添加条件(15) 是因为当γ>1时，即使在条件下c(L2)

γ的值是按照与每个节点相关的标签数量进行动态调整的，采用式(16)表示:

(16)

这里|Li|是与节点i相关联的标签数量，Lt是一个阈值参数，表示应该与每个节点关联的标签数量。应用关系表达式(16)后，如果γ小于1，就会被置回为1，这是由于γ<1的强势标签被保留。Lt越低，标记过程也就越快，获得解决方案的质量也就越差。

3 数值实验

3.1 自定义实例

Li计算如下:

1)如果si<240

_Li=240 概率为 0.7;

_Li=360 概率为0.2;

_Li=480概率为0.1.

2)如果240≤si<360

_Li=360概率为0.6;

_Li=480概率为0.4；

3)否则Li=480

发货时间Ri由下面决定:

Ri=0概率为0.5；否则Ri随机分布在 [0,Rσ]， 其中Rσ=minj∈σ{Tσ-sj,Lj-sj}和i∈σ。

最后Ei=Ri+S0+T0i，就可以创建带有若干客户和车辆数的多个样例。

3.2 确定Lt

为了确定Lt的值，随机生成一个带有60个染色体的实例，其中客户数N=20，车辆数M=8。然后，采用带有不同Lt值的强势支配规则和弱势支配规则分别通过辅助分割过程进行评估，结果如表1。

可以注意到，使用强势支配规则分割一个染色体平均需要超过3分钟，这证明了引入弱势支配规则的合理性。从表1可以看到获得解决方案的质量在有所下降，成本增加，但是Lt递减的过程下降更快。

表1 支配规则的比较Table 1 Comparison of dominated rules

注：对一个M=8和N=20的样例的60条染色体进行辅助分割(时间单位:秒)。

3.3 仿真结果

为了验证本算法的有效性，在MATLAB 2013a的环境下进行编程，并进行了相关测试和分析。其中实验电脑配置为Intel Core I7/8G/Win8.1操作系统。构建数学模型的相关参数指定如下。

客户数为20；配送中心坐标值为(20,30)，客户坐标值为(0,50)之间的随机值；配送车辆总数为8，车辆载货量为10；车辆行驶速度为60；假定每辆车可能产生的行程数为3， 客户需求采用均匀分布的随机数，范围为1～4；客户服务时间采用均匀分布的随机数，范围为5～10 min，3个时间间隔[6,9]、[12,15]和[18,21]分别用于生成3个不同行程的时间窗[Ei,Li]以及Ri的值；车辆的配置成本为200，由于迟到产生的违反时间窗的成本为300，单位里程的运输成本为30。

分别采用人工蜂群算法(artificial bee colony，ABC)[19]，杂草优化算法(invasive weed optimization，IWO)[20]，粒子群算法(particle swarm optimization，PSO)[21]，和笔者提出的改进混合遗传算法(hybrid genetic algorithm，HGA)进行了对比。其中ABC参数设置为：算法最大迭代次数800次，雇佣蜂数量100只，观察蜂数量100只，侦查蜂步数为50，加速度系数上限为1，惯性权重阻尼比为0.99。IWO参数设置为：算法最大迭代次数800次，初始种群数目为10，最大种群数目为30，最小种子数目为2，最大种子数目为10，方差缩减指数为2，标准差初始值为0.5，标准差最终值为0.001。PSO算法参数设置为：算法最大迭代次数为800次，种群数目为100，权重为1，惯性权重阻尼比为0.99，个人学习参数为1.5，全局学习参数为2。 HGA算法参数设置为：最大迭代次数800次，种群数目50，交叉概率0.8，变异概率0.2，减速系数为1，局部搜索概率为0.4，局部搜索歩长为10。

图1～图4为采用不同算法进行3个行程货物配送的优化路线图。

图1 HGA算法优化多行程路线Fig. 1 The multi-trip roadmap of HGA algorithm optimization

图2 ABC算法优化多行程路线Fig. 2 The multi-trip roadmap of ABC algorithm optimization

图3 IWO算法优化多行程路线Fig. 3 The multi-trip roadmap of IWO algorithm optimization

图4 PSO算法优化多行程路线Fig. 4 The multi-trip roadmap of PSO algorithm optimization

图5 算法收敛Fig. 5 The algorithm convergence graph

各个仿真程序运行获得的收敛曲线如图5。在图5的仿真结果中，在相同的数学模型条件下，采用HGA获得的最优成本为33 568；采用ABC获得的最优成本为37 450；采用PSO获得的最优成本为33 950；采用IWO获得的最优成本为35 065。通过比较证明，这四种算法中HGA获得的成本最低。

4 结 语

介绍了一种带时间窗和发货时间的多行程车辆路径问题。这个问题特别适用于城市物流环境，即物流卡车把货物不断地送到城市配送中心，再由配送中心仓库中的有容量限制的电动三轮车等小型车辆进行多次地配送。提出了一个基于辅助分割过程的改进混合遗传算法来处理这个问题，结合时间窗、发货时间和多行程等约束条件进行数学建模，引入了一组实例，通过和文献中提出的人工蜂群算法、杂草优化算法、粒子群算法进行对比，仿真实验结果证明，该改进的混合遗传算法具有较高的可行性和使用效率。

同时需要说明的是，遗传算法是寻优的一种计算方法，多行程的路径寻优问题的复杂性在于网络路径交通动态的本身，虽然考虑了路径时间问题，但路径的实际时间是随着交通流等多种因素的影响而变化的，比如交通阻塞、车辆故障、紧急状况等，如果把这些因素纳入约束条件，必将导致动态寻优的过程会更加复杂，这些可以作为以后的研究内容进一步研究。