考虑碳排放的两阶段选址-路径问题及其算法

汤希峰 ，何 杰 ，张 浩

（1.河海大学土木与交通学院，江苏 南京 210098；2.东南大学交通学院，江苏 南京 211189）

城市是生产、生活的主要集中地，也是物流活动最为集中的地方.两层的物流网络作为一种可行的城市物流解决方案正被国内外许多大中城市所采用[1-4]，其核心问题之一是两阶段选址-路径问题 （2ELRP）.经典的2E-LRP 可以描述为：重型或中型货车负责将货物从一个位于城郊的物流园区运送至若干位于城市中心区周边的配送中心，货物在这些配送中心经过重新配装后再由轻型或微型货车送达终端客户，其目标是寻求最优的配送中心选址以及第一阶段和第二阶段货车的最优路径以实现总的物流成本最小化[5].

众所周知，货车相较于客车容易导致交通事故和交通拥堵，更是污染物和碳排放的主要贡献者.以我国为例，2019 年全国货车保有量虽然只占汽车总量的10.71%，货车CO、HC、NOx、PM 排放量却占到了汽车排放总量的29.7%、26.3%、83.5%、90.1%[6].随着生态环境保护问题日渐受到重视，如何优化2 层网络结构的物流运作以减少货车的碳排放引起了国内外学者越来越多的研究兴趣.Govindan 等[7]围绕生鲜配送，建立了以物流成本最小化和环境影响最小化为目标的2E-LRP 多目标优化模型，并设计了多目标粒子群优化算法（MOPSO）和改进的多目标可变邻域搜索算法（AMOVNS）2 种求解算法；胡大伟等[8]针对燃油车和电动车不同的行驶与排放特性，分别建立了燃油车和电动车的2 阶段开放式选址-路径问题模型，并提出了一种改进的模拟退火算法；刘成清等[9]以经典的2 阶段开放式选址-路径问题模型为基础，通过考虑车辆的燃油消耗和CO2排放以及车辆路径选择的灵活性，建立了基于路径灵活性的两阶段开放式低碳选址-路径问题模型，并运用Cplex进行了求解；Pitakaso 等[10]以燃油消耗最小化为目标，建立了考虑碳排放的绿色两阶段选址路径问题（G2ELRP）优化模型，并设计了相应的可变邻域策略的自适应搜索（VANSAS）算法.

一方面，鉴于车辆的碳排放受到车型、车速、载重、行驶距离以及驾驶行为等诸多因素的影响[11-12]，现有文献大多基于车辆的燃油消耗进行碳排放的计算，所需参数较多，实用性不强；另一方面，现有的用于求解各种2E-LRP 的算法多适用于求解小规模问题 （客户数量n≤100, 候选配送中心数量m<10），求解大规模问题 （n>100,m≥10） 时计算耗时较长，很难实际应用.为此，本文基于一种简单实用的以排放因子为主要参数的碳排放计算方法，建立了考虑碳排放的2E-LRP 优化模型，并针对所研究问题的计算复杂性，设计了一种可快速求解大规模问题的两阶段混合算法.

1 模型构建

1.1 问题描述

基于图论，考虑碳排放的2E-LRP 可以描述为G={V,E}，V=V0∪Vs∪Vc，V0为物流园区，Vs= （1,2,…,m）为配送中心候选集，Vc= （m+1,m+2,…,m+n）为客户集；E=E1∪E2，E1= {（i,j）:i,j∈V0∪Vs;i≠j}为物流园区与配送中心之间的路线集，E2= {（i,j）:i,j∈Vs∪Vc;i,j∉Vs×Vs;i≠j} 为配送中心与客户以及客户与客户之间的路线集；H为可供物流园区使用的重型或中型货车的集合，每台重型或中型货车具有相同的载重能力，记为QH，满载和空载时的碳排放系数分别为εvfH、εveH（单位： kg/km, 下同）；K为可供配送中心共用的轻型或微型货车的集合，每台轻型或微型货车也具有相同的载重能力，记为QK（QK

图1 2E-LRP 示意Fig.1 An illustration of the 2E-LRP

1.2 数学模型

若用qijh和qijk分别为货车h（h∈H）在路线 （i,j）∈E1、货车k（k∈K）在路线（i,j）∈E2上行驶时的载重量；决策变量xijh=1 表示货车h（h∈H） 经过路线 （i,j）∈E1，否则xijh=0；yijk=1 表示货车k（k∈K）经过路线 （i,j）∈E2，否则yijk=0；zi=1 表示配送中心i（i∈Vs） 被选用，否则zi=0；uij=1 表示客户j由配送中心i提供服务，否则uij=0；则以碳排放最小化为目标的2E-LRP 数学模型可表示为

式 （1） 为目标函数，表示从物流园区到配送中心以及配送中心到客户2 个阶段的货车碳排放量之和最小，与经典2E-LRP 模型的目标函数 （关于车辆通行成本或车辆行驶距离的函数） 不同，是一个关于车辆载重和车辆行驶距离的函数；式（2）、（3） 表示2 个阶段中某一货车要么不被使用，要么至多被使用一次；式（4）、（5） 保证2 个阶段中货车路线的连续性，同时保证货车完成任务后返回原先出发的物流园区或配送中心；式（6） 表示客户只能由一个配送中心提供服务；式（7） 表示客户若由某一配送中心提供服务，则必被从该配送中心出发的货车访问一次；式（8） 表示某一配送中心要么不被选用，若被选用，所服务的客户总需求量不能超过其作业能力；式（9）、（10） 表示2 个阶段中货车服务的客户需求均不能超过其载重能力；式（11）、（12） 为2 个阶段中路线的完整性约束；式（13）、（14）表示2 个阶段中货车载重量在相邻路段上的数量关系；式（15） ～（18）为决策变量的取值约束.

2 算法设计

由于设施选址问题和车辆路径问题都属于NPhard 问题，所以本文研究的考虑碳排放的2E-LRP是更为复杂的NP-hard 问题，求解难度很大，尤其是随着问题规模的增大，计算耗时呈指数增加.鉴于此，本文设计了一种可快速求解大规模问题的2 阶段混合算法（TSHA ），算法流程如图2 所示.第1 阶段，通过将2E-LRP 转化成不考虑车辆路径的2 阶段设施选址问题 （2E-FLP），直接求解配送中心选址方案和客户分配方案；第2 阶段，基于得到的配送中心选址和客户分配方案，物流园区到被选用的配送中心以及每个配送中心到所分配客户的车辆路径问题就变成了独立且较易求解的VRP 问题，进而加以解决.

图2 两阶段混合算法的算法流程Fig.2 Flowchart of the TSHA

2.1 第1 阶段算法

由于迭代算法在求解过程中非常低效，为提高所建模型的求解效率，TSHA 在第1 阶段先将考虑碳排放的2E-LRP 转化成不考虑车辆路径的2EFLP，模型如下：

式 （19） 为目标函数，表示客户到配送中心以及配送中心到物流园区的以需求量为权重的距离之和最小，约束条件对应原问题模型中的式 （6）、（8）、式（15）、（16）.不难看出，上述2E-FLP 模型为线性规划模型，可直接运用Cplex 求解器进行求解，从而得到配送中心选址和客户分配方案.

2.2 第2 阶段算法

基于得到的配送中心选址和客户分配方案，物流园区到被选用的配送中心以及每个被选用配送中心到所分配客户的车辆路径问题都变成了独立的VRP 问题.为求解以碳排放最小化为目标的VRP 问题，TSHA 在第2 阶段采用了一种改进的蚁群算法，如图3 所示.

图3 改进蚁群算法的伪代码Fig.3 Pseudocode of the improved ant colony algorithm

与经典的蚁群算法相同[13]，改进的蚁群算法也包含2 个关键步骤：1） 主要负责路线的构建，即基于当前节点逐步选择下一个将要访问的节点；2） 主要负责信息素浓度的更新，具体包括局部信息素更新和全局信息素更新.

从式 （1） 碳排放的计算式可以看出，货车的碳排放与其载重量和行驶距离密切相关.为此，改进的蚁群算法采用了一种新的路线构建规则，如式 （20） 、式（21） 所示.

式中：Ω为尚未被访问的节点集，Ω∈Vs或Vc；τij为当前节点i与未到访节点j之间路线上的信息素浓度；ηij为节点i与节点j之间距离的倒数；α、β为重要度系数；q0为给定的常数；q为 [0, 1] 上的随机数.

当随机产生的q≤q0时，选择使取值最大的j作为下一个将要访问的节点；当q>q0时，则根据的值所占的比例选择下一个将要访问的节点j，在TSHA 算法中采用轮盘赌的方式进行选择.

每只蚂蚁都按照上述的伪比例选择规则构建路线，当蚂蚁构建出一条路线（即达到货车的载重限制时，返回原先出发的物流园区或配送中心），然后重新出发构建下一条路线，如此反复，直至遍访所有客户.当蚂蚁构建出一个完整的路线方案后，对该路线方案中的每条路线进行局部信息素更新，第t+ 1 次信息素浓度为

式中：ρ（0 ≤ρ≤1） 为信息素的挥发系数；τ0为初始的信息素浓度，取值为1/（mEnns），Enns为基于最近邻搜索算法得到的货车碳排放量.

当蚁群中每只蚂蚁均构建出一个完整的路线方案后，从中选择一个碳排放量最小的路线方案作为当次最优的路线方案，为在当次最优路线方案和当前最优路线之间取得平衡，亦即鼓励随后的蚂蚁在当次最优路线和当前最优路线附近的解空间继续搜索，改进的蚁群算法采用Ting 等[14]推荐的策略对当次最优和当前最优路线方案中的每条路线进行全局信息素更新，如式（23）、（24）所示.

式中：EIW为当次最差路线方案产生的碳排放.

如此，通过若干次循环后即可得到所有以碳排放最小化为目标的VRP 的最优路线方案.考虑到由物流园区提供服务的配送中心数量较少而由每个配送中心提供服务的客户数量相对较多，对得到的配送中心到客户的VRP 最优路线方案，本文运用3 种局部搜索算法对其再进行改进：cross_relocation 算法将从某个配送中心出发路线上的客户逐一插入从其他配送中心出发的路线中，若得到的路线可行且碳排放量小于原路线，则用改进后的路线替代原路线；cross_swap 算法将从不同配送中心出发的两条路线上的客户两两交换并重新生成新的路线，若新路线可行且碳排放量较小，则替代原路线；最后，再对每一条路线运用3-opt 算法进行改进.至此，改进后的路线方案即可作为最终的最优路线方案.

3 算例测试

由于目前还没有可用于测试2E-LRP 碳排放的标准算例集，所以本文选取了Prodhon 标准算例集[15]中全部6 个最大规模的算例稍加修改，用以测试本文所提出的模型和算法.Prodhon 算例集中6 个最大规模的算例均包括200 个客户和10 个候选配送中心.算例的修改包括去掉原算例中的配送中心选址费用、增加不同类型货车的碳排放系数，其中重型或中型货车满载和空载时的碳排放系数分别为2.392、1.638 kg/km，轻型或微型货车满载和空载时的碳排放系数分别为1.096、0.772 kg/km，其他数据保持不变.

设计的TSHA 算法采用MATLAB 2018b 编程并嵌入Cplex V12.7 求解器，改进的蚁群算法控制参数包括：求解物流园区到被选用的配送中心VRP 问题时，循环次数设定为选用的配送中心个数的1/2，蚁群中蚂蚁的数量与选用配送中心个数相同；求解配送中心到被分配客户的VRP 问题时循环次数均设定为被分配客户数量的1/2，蚁群中蚂蚁的数量均等于被分配的客户的数量；重要度系数α= 2，β= 1；给定常数q0= 0.5；挥发系数ρ= 0.2.针对每个算例，独立运行程序20 次，并分别记下计算结果和时间.

3.1 对TSHA 算法思想的验证

为验证TSHA 的算法思想，本文先编写了基于同样算法思想的TSHA-Ⅱ算法用以求解Prodhon 算例集中6 个最大规模的原算例.TSHA-Ⅱ与TSHA略微不同的是：1） 用物流成本代替碳排放作为目标函数；2） 改进的蚁群算法中，路线构建不考虑客户的货运需求量，即式 （20） 、（21） 中去掉启发因子dj.TSHA-Ⅱ算法的求解结果与目前文献中两种最新算法的求解结果对比如表1 所示.表中：B为目前文献所给出的最优解；Cmin为基于某种算法得到的最优解；Tmin为某种算法得到最优解所需的时间；G为某个算法求得的最优解与B之间的差异，即G=(Cmin-B)/B.

表1 TSHA-II 算法与LNS-2E、Two-stage 算法求解2E-LRP 原算例的结果对比Tab.1 Comparison of the results obtained by LNS-2E, two-stage algorithm, and TSHA-II in solving original 2E-LRP instances

通过对比可以看出：TSHA-Ⅱ算法求解大规模的2E-LRP 算例时，求解质量稍逊于Breunig 等[16]提出的LNS-2E 算法，但明显优于Dai 等[17]提出的Twostage 算法；求解效率虽稍逊于后者，但大大优于前者.这表明，TSHA-Ⅱ算法能够在求解质量与求解效率之间取得较好的平衡，可以在较短的时间内得到较高质量的解，TSHA 的算法思想可行有效.

3.2 对考虑碳排放的2E-LRP 算例的求解

算法思想通过验证后，TSHA 再被用于求解考虑碳排放的2E-LRP 大规模算例，计算结果如表2所示.表中：Eavg为20 次运算得到的平均碳排放量；Emax、Emin分别为20 次运算得到的最大和最小碳排放量；Tavg为20 次运算的平均耗时.

表2 TSHA 算法求解考虑碳排放的2E-LRP 算例的结果Tab.2 Results obtained by TSHA in solving 2E-LRP instances considering carbon emissions

从表2 可以看出，TSHA 在求解考虑碳排放的2E-LRP 大规模算例时，求解结果的最大值、最小值和平均值以及计算时长的最小值和平均值之间相差很小，表现稳定.可见，TSHA 可以作为求解考虑碳排放的2E-LRP 大规模算例的一种有效算法.

4 结 论

1） 相较于经典的以物流成本最小化为目标的2E-LRP 模型，本文基于一种简单实用的以排放因子为主要参数的碳排放计算方法，建立了以碳排放最小化为目标的2E-LRP 优化模型.

2） 对Prodhon 标准算例集中全部6 个最大规模算例的测试，表明TSHA-Ⅱ算法比现有算法更具兼顾求解结果和求解效率的优势，也证明了TSHA算法思想的可行性和有效性.

3） TSHA 算法在求解考虑碳排放的2E-LRP 大规模算例时表现得高效而且稳定，能够满足实际应用场景的计算要求.

4） 本文提出的模型和算法有助于物流企业通过优化城市物流运作实现车辆的减排目标，服务国家“双碳”战略.