范厚明，孫秀娜，張躍光，任曉雪，田攀俊

大連海事大學(xué) 交通運(yùn)輸工程學(xué)院，遼寧 大連 116026

時(shí)變路網(wǎng)下帶混合時(shí)間窗的車輛路徑問題（vehicle routing problem with mixed time windows under timedependent network，VRPMTWTDN）是綜合考慮多中心聯(lián)合配送、混合時(shí)間窗及車輛行駛速度隨時(shí)間連續(xù)變化的VRP（vehicle routing problem）拓展問題。現(xiàn)實(shí)中，物流企業(yè)的多個(gè)配送中心共同為客戶提供配送服務(wù)，并根據(jù)客戶對(duì)貨物送達(dá)時(shí)間的不同要求，將客戶收貨時(shí)間約束分為硬時(shí)間窗約束和軟時(shí)間窗約束。如京東推出的“京準(zhǔn)達(dá)”服務(wù)、天貓推出的基于門店的“定時(shí)送”服務(wù)及哈根達(dá)斯的哈根速遞等，在客戶下單時(shí)能預(yù)約精準(zhǔn)收貨時(shí)間段，企業(yè)須在該時(shí)間窗內(nèi)將貨物送到達(dá)；未預(yù)約精準(zhǔn)收貨時(shí)間段的客戶允許配送時(shí)間存在一定偏差。同時(shí)，路網(wǎng)交通流量的時(shí)空差異使得配送車輛的行駛速度不斷變化，對(duì)貨物的送達(dá)時(shí)間產(chǎn)生一定影響。因此，時(shí)變路網(wǎng)下帶混合時(shí)間窗的多中心車輛路徑問題貼近實(shí)際配送生產(chǎn)活動(dòng)，具有重要的理論與應(yīng)用價(jià)值。

VRPMTWTDN 是對(duì)多中心車輛路徑問題（multidepot vehicle routing problem，MDVRP）、時(shí)間依賴型車輛路徑問題（time-dependent vehicle routing problem，TDVRP）、帶混合時(shí)間窗的車輛路徑問題（vehicle routing problem with mixed time windows，VRPMTW）的集成問題，已有很多學(xué)者對(duì)MDVRP、TDVRP 及VRPMTW展開了研究。

針對(duì)MDVRP，范厚明等[1]針對(duì)多中心聯(lián)合配送模式下集貨需求隨機(jī)的同時(shí)配集貨VRP建立兩階段優(yōu)化模型，并改進(jìn)變鄰域文化基因算法進(jìn)行求解。Zhang等[2]針對(duì)MDVRP，以碳排放量最小為目標(biāo)建立兩階段優(yōu)化模型，并設(shè)計(jì)改進(jìn)蟻群算法進(jìn)行求解。Zhen等[3]針對(duì)帶有時(shí)間窗的多行程MDVRP，以配送時(shí)長(zhǎng)最小為目標(biāo)建立優(yōu)化模型，并設(shè)計(jì)混合粒子群優(yōu)化算法和混合遺傳算法進(jìn)行求解。Li 等[4]針對(duì)MDVRP，以收益最大和配送成本、時(shí)間及碳排放量最小為目標(biāo)建立優(yōu)化模型，并設(shè)計(jì)改進(jìn)蟻群優(yōu)化算法進(jìn)行求解。

針對(duì)TDVRP，范厚明等[5]考慮客戶需求模糊及車輛行駛時(shí)間依賴路網(wǎng)速度變化特征，提出時(shí)變路網(wǎng)下多車型動(dòng)態(tài)VRP，建立兩階段優(yōu)化模型，并設(shè)計(jì)自適應(yīng)大規(guī)模鄰域搜索算法進(jìn)行求解。Ehmke等[6]考慮車輛行駛速度隨時(shí)間變化及車輛行駛速度等對(duì)油耗的影響，以碳排放量最小為目標(biāo)建立優(yōu)化模型，并設(shè)計(jì)禁忌搜索算法進(jìn)行求解。Alinaghian 等[7]考慮車輛行駛速度的時(shí)間依賴性及車輛行駛速度、載重量、道路坡度對(duì)車輛油耗的影響，以油耗量最小為目標(biāo)建立TDVRP優(yōu)化模型，并設(shè)計(jì)改進(jìn)高斯螢火蟲算法進(jìn)行求解。張凱慶等[8]考慮軟時(shí)間窗約束和車輛速度隨時(shí)間變化情況，建立以平均客戶滿意度最大、配送距離及成本最小為目標(biāo)的優(yōu)化模型，并設(shè)計(jì)兩階段算法進(jìn)行求解。Maha 等[9]考慮車輛行駛速度在路網(wǎng)中不同路段的差異性，建立帶時(shí)間窗的VRP優(yōu)化模型，并設(shè)計(jì)禁忌搜索算法進(jìn)行求解。付朝暉等[10]針對(duì)生鮮電商末端配送問題，考慮生鮮農(nóng)產(chǎn)品送達(dá)客戶時(shí)的鮮活度及城市路網(wǎng)交通狀況的時(shí)變性等因素，以總配送成本最小為目標(biāo)建立TDVRP 優(yōu)化模型，并設(shè)計(jì)改進(jìn)蟻群算法進(jìn)行求解。

針對(duì)VRPMTW，Wang 等[11]針對(duì)冷鏈車輛配送問題，考慮混合時(shí)間窗及碳交易政策中的碳價(jià)格和交易配額設(shè)置對(duì)模型求解的影響，以成本最小、碳排放量最小和客戶總滿意度最大為目標(biāo)建立優(yōu)化模型，并設(shè)計(jì)自適應(yīng)遺傳算法進(jìn)行求解。Wang 等[12]考慮混合時(shí)間窗約束，以配送成本和派遣車輛數(shù)最小為目標(biāo)建立多目標(biāo)優(yōu)化模型，并將智能水滴算法與遺傳算法相結(jié)合進(jìn)行求解。王勇等[13]針對(duì)VRPMTW，將客戶重要度與時(shí)間窗相結(jié)合，建立雙層數(shù)學(xué)模型，并將遺傳算法與禁忌搜索算法相結(jié)合進(jìn)行求解。鄒宗峰等[14]以運(yùn)輸風(fēng)險(xiǎn)、時(shí)間窗懲罰成本、交通狀況、道路能力、運(yùn)輸時(shí)間最優(yōu)為目標(biāo)建立VRPMTW優(yōu)化模型，并設(shè)計(jì)多目標(biāo)遺傳算法進(jìn)行求解。

通過梳理以上相關(guān)文獻(xiàn)可見，目前還沒有針對(duì)VRPMTWTDN 的研究成果，且MDVRP、TDVRP 及VRPMTW的研究還存在以下不足：（1）現(xiàn)有針對(duì)TDVRP的研究多以階梯分段函數(shù)表示全天道路車輛速度變化情況，未考慮車輛行駛速度的連續(xù)變化。（2）現(xiàn)有考慮車輛行駛速度連續(xù)變化的研究多假設(shè)車輛單位距離的油耗是固定的，鮮有考慮車速、載重量等對(duì)油耗的影響。（3）現(xiàn)有關(guān)于帶時(shí)間窗的車輛路徑問題的文獻(xiàn)，多數(shù)單獨(dú)考慮軟時(shí)間窗或硬時(shí)間窗，較少研究VRPMTW。針對(duì)以上不足，本文綜合考慮多中心聯(lián)合配送，不同客戶對(duì)服務(wù)時(shí)間的需求差異，車輛行駛速度隨時(shí)間連續(xù)變化以及車速、載重量對(duì)油耗的影響，對(duì)VRPMTWTDN展開研究。

1 問題分析

1.1 問題描述

由配送區(qū)域路網(wǎng)構(gòu)成完備的有向圖G=(V,E)中，節(jié)點(diǎn)集合為V=V0?V1?V2，其中V0={1,2,…,m} 為配送中心集合，V1為未預(yù)約精準(zhǔn)收貨時(shí)間段的客戶集合，V2為預(yù)約精準(zhǔn)收貨時(shí)間段的客戶集合；E={(i,j)|i,j∈V,i≠j}為邊集合，各配送中心間無路徑連接，lij為i、j兩節(jié)點(diǎn)之間的距離，tij為車輛在i、j兩節(jié)點(diǎn)之間的行駛時(shí)間；K為所有配送中心的車輛集合，車型相同，最大載重量為Q,Ki(i∈V0)為各配送中心i(i∈V0)的車輛集合，k為配送車輛集合K中任意一輛車，Sk為車輛k服務(wù)的客戶集合，F(xiàn)ijk為車輛k在i、j兩節(jié)點(diǎn)之間的燃油消耗量，燃油價(jià)格為c1，單位車輛的派遣成本為c2；車輛行駛速度v={v1,v2,…,vl}是連續(xù)變化的；車輛k從配送中心i出發(fā)的時(shí)刻為Tik，到達(dá)節(jié)點(diǎn)i的時(shí)刻為，對(duì)客戶i的開始服務(wù)時(shí)刻為，車輛完成所有配送任務(wù)后在各配送中心車輛數(shù)守恒原則下就近返回任一配送中心；客戶i的需求量為di，單位貨物裝卸時(shí)間為δ,配送中心有統(tǒng)一的工作時(shí)間窗[Ts,Tf],[ETi,LTi]為預(yù)約精準(zhǔn)收貨時(shí)間段客戶i(i∈V2)的時(shí)間窗，車輛須在該時(shí)間窗內(nèi)為客戶i提供服務(wù)[ETi,LTi]為未預(yù)約精準(zhǔn)收貨時(shí)間段客戶i(i∈V1)的期望服務(wù)時(shí)間窗，[EETi,ELTi](EETi≤ETi≤LTi≤ELTi)為其最大容忍時(shí)間窗，車輛須在[EETi,ELTi]內(nèi)為客戶i提供服務(wù)，在[EETi,ETi] 到達(dá)的單位懲罰成本為c3,在[LTi,ELTi] 到達(dá)的單位懲罰成本為；決策變量xijk表示車輛k是否從點(diǎn)i到達(dá)點(diǎn)j，是為1，否為0。

圖1為研究問題的示意圖，配送區(qū)域有2 個(gè)配送中心、25個(gè)客戶，每個(gè)配送中心擁有若干車輛。根據(jù)現(xiàn)有客戶位置及時(shí)間窗進(jìn)行優(yōu)化得到圖1所示的5條配送路徑，其中車輛1、車輛3 和車輛5 均返回原配送中心，車輛2和車輛4返回就近的配送中心。

圖1 研究問題的示意圖Fig.1 Schematic diagram of research problem

1.2 速度時(shí)間依賴函數(shù)的確定

現(xiàn)有對(duì)TDVRP 的研究，多以階梯分段函數(shù)表示路段的車輛行駛速度，在某一時(shí)間節(jié)點(diǎn)速度突變，不能很好地刻畫高峰時(shí)段配送區(qū)域路網(wǎng)對(duì)車輛行駛速度的限制，不符合現(xiàn)實(shí)中車輛行駛速度連續(xù)變化的特點(diǎn)。鑒于此，本文參考文獻(xiàn)[15]，以多項(xiàng)式函數(shù)表示車輛行駛速度隨時(shí)間的連續(xù)變化，如圖2所示。

圖2 車輛速度時(shí)間依賴函數(shù)Fig.2 Vehicle speed-time dependence function

1.3 車輛油耗計(jì)算

車輛油耗量取決于車輛行駛時(shí)間、行駛速度、車型和道路狀況等因素，本文借鑒文獻(xiàn)[16]研究結(jié)果，考慮車速、載重量對(duì)油耗的影響，建立油耗計(jì)算模型。

行駛時(shí)間dt內(nèi)的燃油消耗量dfij為：

式中，U(v(t))是空載車輛的燃油消耗率，單位是L/km。

U(v(t))與空載車速之間的關(guān)系為：

式中，α、β、θ和ω由道路和車輛狀況決定。

空載車輛k在路段(i,j)上的燃油消耗量fijk為：

式中，Tl ik為車輛k離開i點(diǎn)時(shí)刻，Ta jk為車輛k到達(dá)j點(diǎn)時(shí)刻。

假設(shè)油耗增量與車輛載重增量呈線性關(guān)系[21]，載重每增加L，油耗增加s（百分比），則載重量為Qik的車輛k在路段(i,j)上的油耗Fijk為：

2 模型建立

式（1）為目標(biāo)函數(shù)式，第一部分為油耗成本，第二部分為車輛派遣成本，第三部分為時(shí)間窗懲罰成本；式（2）表示車輛的載重限制；式（3）表示各客戶點(diǎn)車輛進(jìn)出平衡，且客戶只能被服務(wù)一次；式（4）表示各配送中心出發(fā)和返回的車輛數(shù)一致，且不超過其最大車輛數(shù)；式（5）表示各車輛須在配送中心工作截止時(shí)間前返回配送中心；式（6）和式（7）為車輛對(duì)客戶i開始服務(wù)時(shí)刻的約束；式（8）為車輛從節(jié)點(diǎn)i出發(fā)到達(dá)節(jié)點(diǎn)j時(shí)刻的計(jì)算方法，其中M為無窮大值；式（9）為消除子回路約束；式（10）為決策變量屬性。

3 算法設(shè)計(jì)

本文所提問題考慮多中心聯(lián)合配送、混合時(shí)間窗及車輛行駛速度連續(xù)變化等因素對(duì)配送方案制定的影響，屬于NP-hard問題，求解過程復(fù)雜，用精確算法求解較為困難，因此考慮采用啟發(fā)式方法進(jìn)行求解。遺傳算法（genetic algorithm，GA）具有快速隨機(jī)的全局搜索能力，可同時(shí)對(duì)搜索空間中的多個(gè)可行解進(jìn)行評(píng)估，但搜索具有盲目性，在求解到一定范圍時(shí)常做大量無為冗余迭代，使求解效率低；大規(guī)模鄰域搜索算法既能保持大規(guī)模鄰域的尋優(yōu)能力，又能克服低效的弱點(diǎn)，從而提高算法的求解能力。本文根據(jù)VRPMTWTDN 模型特點(diǎn)，將傳統(tǒng)遺傳算法與大規(guī)模鄰域搜索算法相結(jié)合，設(shè)計(jì)自適應(yīng)遺傳-大鄰域搜索算法（adaptive genetic algorithm with large neighborhood search，AGA-LNS）進(jìn)行求解。

自適應(yīng)遺傳-大鄰域搜索算法相比于傳統(tǒng)遺傳算法，通過大規(guī)模鄰域搜索提高尋優(yōu)能力，改進(jìn)的交叉變異算子對(duì)可行解不斷進(jìn)行優(yōu)化，且以最優(yōu)個(gè)體持續(xù)未改變代數(shù)作為迭代停止的判斷依據(jù)，避免做大量無為冗余迭代，提高求解效率；相比于傳統(tǒng)鄰域搜索，既能保持大規(guī)模鄰域的尋優(yōu)能力，又能克服低效的弱點(diǎn)，局部最優(yōu)解的質(zhì)量更好，但迭代耗時(shí)較長(zhǎng)。自適應(yīng)遺傳-大鄰域搜索算法早期尋優(yōu)主要依靠遺傳算法的全局搜索能力和大鄰域算法的局部搜索能力，通過改進(jìn)的交叉操作和變異操作可以快速接近全局最優(yōu)解；后期尋優(yōu)主要依靠大鄰域搜索算法的局部搜索能力，通過移除-插入操作對(duì)當(dāng)前最優(yōu)解進(jìn)行摧毀和重建，有利于算法跳出局部最優(yōu)解，彌補(bǔ)了遺傳算法后期收斂速度慢，容易造成冗余迭代的缺點(diǎn)。遺傳算法與大鄰域搜索算法相結(jié)合可以有效平衡算法的局部尋優(yōu)能力和全局尋優(yōu)能力，具有較強(qiáng)的求解性能。

自適應(yīng)遺傳-大鄰域搜索算法流程如圖3所示。

圖3 AGA-LNS流程圖Fig.3 Flow chart of AGA-LNS

3.1 初始解的構(gòu)造

初始解的質(zhì)量很大程度上影響著算法的求解質(zhì)量和求解效率，若隨機(jī)生成初始種群，極易產(chǎn)生大量無效個(gè)體，增加收斂代數(shù)。為提高可行解的生成概率，采用整數(shù)編碼方式，設(shè)計(jì)“先客戶后中心”的兩階段方法生成初始解。首先進(jìn)行客戶點(diǎn)分組并確定服務(wù)順序，然后確定發(fā)車時(shí)間和返回的配送中心及到達(dá)各節(jié)點(diǎn)的具體時(shí)刻。此方式可保證生成的初始種群均為可行解，且具有一定隨機(jī)性，增加了初始種群的多樣性和有效性，避免過早陷入局部最優(yōu)。具體步驟如下：

階段1 客戶點(diǎn)分組并確定服務(wù)順序。

步驟1 隨機(jī)派遣一輛車，用虛擬配送中心0表示其出發(fā)的配送中心，根據(jù)客戶允許的最早開始服務(wù)時(shí)刻由早到晚對(duì)所有未服務(wù)的客戶進(jìn)行排序，并從前n個(gè)客戶中隨機(jī)選取一個(gè)作為新派車輛服務(wù)的第一個(gè)客戶并將其插入到當(dāng)前路徑。

步驟2 求出從當(dāng)前客戶出發(fā)到所有未服務(wù)客戶的到達(dá)時(shí)刻并判斷其是否在該客戶時(shí)間窗內(nèi)，從到達(dá)時(shí)刻在客戶時(shí)間窗內(nèi)的所有客戶集合中選擇開始服務(wù)時(shí)刻最小的客戶插入到路徑，若所有客戶的到達(dá)時(shí)刻都早于其時(shí)間窗開始時(shí)刻，則選擇到達(dá)時(shí)刻距其時(shí)間窗開始時(shí)刻最近的客戶插入到當(dāng)前路徑。重復(fù)以上操作，直至車載量、未服務(wù)客戶的時(shí)間窗或返回配送中心時(shí)刻中任一約束不滿足，則新派一輛車。

步驟3 重復(fù)步驟1 和步驟2，直至所有客戶的配送任務(wù)安排完畢。

階段2 確定車輛出發(fā)和返回的配送中心及到達(dá)各節(jié)點(diǎn)的具體時(shí)刻。

步驟4 替換車輛出發(fā)和返回的虛擬配送中心0：計(jì)算每條路徑首尾客戶距離各配送中心的距離，依次選擇與各路徑首個(gè)客戶距離最小的配送中心替換車輛出發(fā)的虛擬配送中心0，若所選配送中心無可派遣車輛，則選擇距離次近的配送中心，直至替換完畢；然后依次選擇與各路徑末尾客戶距離最小的配送中心替換車輛返回的虛擬配送中心0，若所選配送中心不滿足車輛守恒，則選擇距離次近的配送中心，直至替換完畢。以各車輛服務(wù)的首個(gè)客戶的時(shí)間窗開始時(shí)刻作為其到達(dá)該客戶點(diǎn)的時(shí)刻，計(jì)算車輛從配送中心的出發(fā)時(shí)刻、到各客戶點(diǎn)的服務(wù)時(shí)刻及返回配送中心的時(shí)刻。

圖4 為染色體編解碼過程示意圖，以10 個(gè)客戶（編號(hào)為1～10）和2個(gè)配送中心（編號(hào)為11～12）為例，按照上述步驟生成初始解。首先在染色體起始位置插入車輛起始虛擬配送中心0，再將所有客戶按最早開始服務(wù)時(shí)刻由早到晚排序，從前5 個(gè)客戶中隨機(jī)選取客戶3 作為車輛服務(wù)的首個(gè)客戶，再計(jì)算從客戶3到所有未服務(wù)客戶的到達(dá)時(shí)刻，選出最佳后續(xù)客戶1 插入到路徑中，按照相同方法選出后續(xù)客戶2、5，直至不滿足任一約束時(shí)，車輛返回虛擬配送中心0。以此類推，劃分后的路徑如圖4（a）所示，再根據(jù)步驟4 所述操作先后確定車輛出發(fā)、返回的配送中心替換虛擬配送中心0，得到初始解如圖4（c）所示。

圖4 染色體解碼示意圖Fig.4 Chromosome decoding diagram

3.2 進(jìn)化操作

進(jìn)化操作首先進(jìn)行基于時(shí)差插入法的自適應(yīng)交叉變異，然后通過成對(duì)的移除、插入算子將路徑中的客戶按照一定規(guī)則移除，再將被移除的客戶按照一定規(guī)則插入到路徑。由于每一次交叉、變異操作或插入-移除操作都可能減少染色體中的車輛數(shù)，每次進(jìn)化操作后按照3.1節(jié)中步驟4描述的方法重新確定車輛出發(fā)和返回的配送中心。

3.2.1 時(shí)差插入法

客戶m插入到客戶i和客戶j之間的時(shí)間約束條件為：車輛在客戶m的開始服務(wù)時(shí)刻不晚于客戶m容忍的最晚接受時(shí)刻Tm；插入客戶m后，車輛到達(dá)客戶j的時(shí)刻不晚于客戶j的最晚開始時(shí)刻。插入的時(shí)間約束如下：

式中，車輛在客戶m的最早完成時(shí)刻EFm和最晚開始時(shí)刻LSm計(jì)算方法如式（12）、（13）所示，其中客戶i為客戶m的前序客戶，客戶j為客戶m的后序客戶。

時(shí)差插入法的具體步驟如下：

步驟1 計(jì)算所有客戶的最早完成時(shí)刻和最晚開始時(shí)刻。

步驟2 按照時(shí)間約束確定客戶m可插入的位置，并計(jì)算假設(shè)插入客戶m后各位置的目標(biāo)增量。

步驟3 選擇目標(biāo)增量最小的位置插入客戶m，若所有位置均不滿足插入條件，則新派遣一輛車。

步驟4 更新所有客戶的最早完成時(shí)刻和最晚開始時(shí)刻。

步驟5 重復(fù)上述操作，直至所有待插入客戶全部插入路徑中。

3.2.2 交叉算子和變異算子

交叉變異概率值的設(shè)定對(duì)解的尋優(yōu)過程有很大影響，傳統(tǒng)遺傳算法的交叉概率和變異概率根據(jù)經(jīng)驗(yàn)和主觀判斷設(shè)置為定值，若設(shè)置較小，則會(huì)使個(gè)體進(jìn)化速度較慢，增加算法的求解時(shí)間，若設(shè)置較大，則容易在算法求解后期破壞優(yōu)秀個(gè)體，使算法陷入局部最優(yōu)解。為提升算法的性能，參考文獻(xiàn)[17]設(shè)置自適應(yīng)交叉變異概率，相比于常規(guī)交叉變異概率，自適應(yīng)交叉、變異可以根據(jù)每次迭代所得新個(gè)體的適應(yīng)度值適應(yīng)性地調(diào)整交叉、變異概率，在進(jìn)化前期以較大的交叉變異概率加速個(gè)體進(jìn)化，增加解的多樣性，避免過早陷入局部最優(yōu)，在進(jìn)化后期交叉變異概率逐漸減小，避免破壞優(yōu)秀個(gè)體。自適應(yīng)交叉、變異概率計(jì)算方法如式（14）、（15）所示。

式中，Pc_max、Pm_max分別是交叉概率區(qū)間和變異概率區(qū)間的最大值，Pc_min、Pm_min分別是交叉概率區(qū)間和變異概率區(qū)間的最小值，f為個(gè)體適應(yīng)度值，其大小為目標(biāo)函數(shù)值的倒數(shù)，f′為進(jìn)行交叉父本中最優(yōu)個(gè)體的適應(yīng)度值，fmax為所有個(gè)體適應(yīng)度值中的最大值，favg為所有個(gè)體的平均適應(yīng)度值。

考慮車輛的服務(wù)時(shí)刻不能違反硬時(shí)間窗客戶的時(shí)間要求，若按傳統(tǒng)交叉變異進(jìn)行操作，則會(huì)破壞原有可行解，產(chǎn)生大量不可行解，對(duì)此，在交叉算子和變異算子中引入時(shí)差插入法，以客戶的時(shí)間窗約束作為種群優(yōu)化的前提條件，保證每次迭代后均生成可行解，提高可行解的生成率，進(jìn)而提高求解質(zhì)量。其中交叉操作對(duì)兩個(gè)父代個(gè)體進(jìn)行部分交換后保留較優(yōu)解，提高種群多樣性，變異操作可減少個(gè)體中子路徑的數(shù)量，即減少配送方案中派遣的車輛數(shù)，有利于提高求解質(zhì)量，降低配送成本。

（1）交叉操作：在兩個(gè)父代染色體中分別隨機(jī)選取一條子路徑并進(jìn)行交換，再對(duì)交換后的染色體進(jìn)行修復(fù)。首先刪除因交換導(dǎo)致重復(fù)的客戶，然后在滿足車輛載重等約束的基礎(chǔ)上，通過時(shí)差插入法插入因交換路徑導(dǎo)致染色體中缺失的客戶。以圖5為例，將父代染色體a中子路徑11—1—7—4—10—11與b中子路徑11—9—8—10—12 進(jìn)行交換得到a1 和b1，刪除a1 中重復(fù)客戶8、9及b1中重復(fù)客戶1、7、4得到a2和b2，再通過時(shí)差插入法將a2 中缺失的客戶1、7、4 插入a2 中、將b2 中缺失的客戶8、9插入b2中得到a3、b3。

（2）變異操作：在父代染色體中隨機(jī)刪除客戶數(shù)最少的子路徑，然后在滿足所有約束的前提下采用時(shí)差插入法重新插入刪除的客戶。以圖6 為例，刪除染色體a中包含客戶數(shù)最少的一條子路徑13—6—8—12得到a1，再采用時(shí)差插入法將客戶節(jié)點(diǎn)6、8重新插入得到a2。

3.2.3 移除算子和插入算子

每次迭代時(shí)，在交叉變異后引入基于時(shí)差插入法的移除算子和插入算子進(jìn)行種群優(yōu)化。先通過移除算子刪除當(dāng)前解中的μ 個(gè)客戶，再通過插入算子將刪除的客戶重新插入構(gòu)造新解。在每次迭代中多次使用移除-插入組合算子對(duì)可行解進(jìn)行摧毀和重建，可以增加解的多樣性，擴(kuò)大搜索范圍，提高算法的全局尋優(yōu)能力。移除算子和插入算子如下：

（1）移除算子

①隨機(jī)移除算子：隨機(jī)移除當(dāng)前解中的μ 個(gè)客戶。

②最遠(yuǎn)距離移除算子：計(jì)算當(dāng)前解中所有客戶的距離成本，移除距離成本D(m) 最大的客戶m ，其中m 的前序客戶，j 為m 的后序客戶。具體操作如圖7所示，圖中D1、D2 代表配送中心，其余為客戶點(diǎn)，箭頭上的數(shù)值為節(jié)點(diǎn)間的距離。

③最遠(yuǎn)偏離客戶中心移除算子：移除每輛車服務(wù)的客戶集合中與其他所有客戶的中心點(diǎn)距離最遠(yuǎn)的客戶，即移除客戶位置偏差D(j)最大的客戶?？蛻鬸 的位置偏差為為客戶j的位置坐標(biāo)，(X,Y)為一輛車服務(wù)的除客戶j 外其他所有客戶的中心點(diǎn)坐標(biāo)。

（2）插入算子

①貪婪插入算子：在滿足所有約束的基礎(chǔ)上采用時(shí)差插入法找出客戶m 的可插入位置，計(jì)算插入距離成本C(m)，并將客戶m 插入C(m)最小的位置，其中C(m)=dim+dmj-dij,i 為插入位置的前序客戶，j 為插入位置的后序客戶。

②最好時(shí)間插入算子：在滿足所有約束的基礎(chǔ)上采用時(shí)差插入法找出客戶m 的可插入位置，并計(jì)算其時(shí)間偏差量BT(m)，并將客戶m 插入BT(m)最小的位置，其中

3.2.4 迭代終止條件

當(dāng)最優(yōu)個(gè)體持續(xù)未改變的代數(shù)達(dá)到最優(yōu)個(gè)體持續(xù)未改變的最大代數(shù)MAXGEN 時(shí)停止迭代，輸出結(jié)果。

4 算法驗(yàn)證與算例分析

4.1 算子的選擇

移除和插入算子是成對(duì)出現(xiàn)的，為確定不同算子求解的有效性，采用Cordeau 標(biāo)準(zhǔn)算例庫(kù)中MDVRPTW（multi-depot vehicle routing problem with time windows）標(biāo)準(zhǔn)算例對(duì)交叉、變異算子和6組移除-插入算子進(jìn)行測(cè)試，各算子求解各算例所得結(jié)果的平均值如表1 所示。其中Ai(i=1,2,3)表示三種移除算子，Bi(i=1,2)表示兩種插入算子，BKS為標(biāo)準(zhǔn)算例庫(kù)中已知最優(yōu)解，Avg為各算子求解所有算例所得結(jié)果的平均值，GAP為Avg與BKS的誤差。由表1 可知，交叉算子的求解性能最好，平均誤差為4.1%，其次為變異算子，平均誤差為4.6%，移除-插入算子組合中，隨機(jī)移除算子和貪婪插入算子組合求解性能最好，平均誤差為12.3%，其次是最遠(yuǎn)距離移除算子和貪婪插入算子的組合，平均誤差為32.1%。為保證求解質(zhì)量及求解速度，選取以上兩種移除-插入算子組合應(yīng)用于算法中進(jìn)行可行解的摧毀和重建。

表1 算子測(cè)試結(jié)果Table 1 Test results of operators

4.2 算法驗(yàn)證

4.2.1 算法測(cè)試

目前沒有VRPMTWTDN 的標(biāo)準(zhǔn)算例庫(kù)，選擇MDVRPTW測(cè)試AGA-LNS的性能。采用MATLAB R2018b進(jìn)行算法編程，操作系統(tǒng)為Windows10，電腦內(nèi)存為8 GB，CPU 為Intel i7-7700M，主頻3.60 GHz。設(shè)置種群規(guī)模NIND為100，MAXGEN為30～50。

選擇Cordeau 標(biāo)準(zhǔn)算例庫(kù)中MDVRPTW 的20 個(gè)算例測(cè)試AGA-LNS，表2給出了禁忌搜索算法（tabu search，TS）[18]、變鄰域搜索算法（variable neighborhood search，VNS）[19]、基于復(fù)合鄰域的離散螢火蟲算法（discrete firefly algorithm with compound neighborhoods，DFACN）[20]、混合遺傳算法（hybrid genetic algorithm，HGA）[21]及AGA-LNS 的求解結(jié)果。其中N為客戶規(guī)模，M為配送中心數(shù)，K為配送中心最大車輛數(shù)，Q為車輛最大載量，BKS為標(biāo)準(zhǔn)算例庫(kù)中已知最優(yōu)解，C為算法求解的最優(yōu)值，GAP為算法求得結(jié)果與BKS的相對(duì)誤差，t為算法運(yùn)行10 次的平均耗時(shí)。由表2 可知，在求解質(zhì)量方面，AGA-LNS 求得的最優(yōu)值與BKS的最小誤差為-3.36%，最大誤差為4.39%，當(dāng)客戶規(guī)模增大時(shí)，相對(duì)誤差有所增大。AGA-LNS的求解結(jié)果優(yōu)于TS和VNS，略差于DFACN和HGA，求解耗時(shí)小于HGA。通過對(duì)比分析運(yùn)算結(jié)果，可驗(yàn)證AGA-LNS的有效性。

表2 算法測(cè)試結(jié)果Table 2 Algorithm test results

4.2.2 交叉變異概率分析

為驗(yàn)證自適應(yīng)交叉變異概率相比于常規(guī)交叉變異概率的優(yōu)越性，設(shè)置三種交叉概率0.40、0.65、0.90，三種變異概率0.02、0.11、0.20，采用Cordeau標(biāo)準(zhǔn)算例庫(kù)中前10 個(gè)MDVRPTW 標(biāo)準(zhǔn)算例對(duì)九種概率組合進(jìn)行測(cè)試，并與本文自適應(yīng)交叉變異概率下算例的求解結(jié)果進(jìn)行對(duì)比，如表3所示。其中C為運(yùn)用自適應(yīng)交叉變異概率求得的最優(yōu)解，Avg為所有算例最優(yōu)解的平均值，GAP為Avg與C平均值的相對(duì)誤差。由表3 可知，九種交叉變異概率組合中，交叉、變異概率分別為0.65、0.20時(shí)求解結(jié)果最優(yōu)，但與自適應(yīng)交叉變異概率下的求解結(jié)果相比仍相差1.20%，可見運(yùn)用自適應(yīng)交叉變異概率進(jìn)行求解相比于常規(guī)交叉變異概率更有優(yōu)越性。

表3 不同交叉變異概率下算例求解結(jié)果Table 3 Example results under different crossover and mutation probability

4.3 算例分析

VRPMTWTDN問題復(fù)雜，沒有通用的算例集，設(shè)計(jì)兩個(gè)算例并用AGA-LNS求解其配送方案。

4.3.1 算例1

假設(shè)在配送路網(wǎng)中有4 個(gè)配送中心和48 個(gè)客戶，配送中心坐標(biāo)分別為（4.163，13.559）、（21.387，17.105）、（-36.118，49.097）、（-31.201，0.235），隨機(jī)生成48 個(gè)客戶坐標(biāo)并將其作為客戶節(jié)點(diǎn)，其中橫坐標(biāo)與縱坐標(biāo)分別在[-75，50]和[-30，100]范圍內(nèi)隨機(jī)生成，客戶需求量在[0.015，0.375]內(nèi)隨機(jī)生成，客戶時(shí)間窗在[05：00，17：00]內(nèi)隨機(jī)生成。每個(gè)配送中心的工作時(shí)間窗均為[05：00，17：00]，最大車輛數(shù)均為4，車輛載重量為3 t，單位車輛的派遣成本為500元，油耗成本c1=5.5 元/升，早到單位懲罰成本c3=30 元/h，晚到單位懲罰成本c4=60 元/h，每噸貨物裝卸時(shí)間為0.5 h，交叉概率Pc、變異概率Pm的取值區(qū)間分別為[0.4，0.9]、[0.02，0.2]。根據(jù)文獻(xiàn)[20]相關(guān)研究，參數(shù)α、β、θ和ω設(shè)為30、100、1 和-5。配送中心根據(jù)各車輛第一個(gè)待服務(wù)客戶的時(shí)間窗合理安排出發(fā)時(shí)間，確保其在首個(gè)待服務(wù)客戶的時(shí)間窗內(nèi)到達(dá)。用AGA-LNS求解得到如表4所示的配送方案，其中YC代表油耗成本，TC代表時(shí)間窗懲罰成本，ZC代表配送方案總成本。

表4 算例1配送方案Table 4 Distribution scheme of Example 1

4.3.2 算例2

假設(shè)在配送路網(wǎng)中有4個(gè)配送中心和96個(gè)客戶，配送中心坐標(biāo)分別為（6.229，10.590）、（32.663，44.730）、（48.807，49.792）、（33.179，-4.968），隨機(jī)生成96 個(gè)客戶坐標(biāo)并將其作為客戶節(jié)點(diǎn)，其中橫、縱坐標(biāo)分別在[-40，87]和[-42，80]范圍內(nèi)隨機(jī)生成，客戶需求量在[0.015，0.375]內(nèi)隨機(jī)生成，客戶時(shí)間窗在[05：00，17：00]內(nèi)隨機(jī)生成，其他參數(shù)設(shè)置與算例1相同，用AGA-LNS求解得到如表5所示的配送方案。

表5 算例2配送方案Table 5 Distribution scheme of Example 2

4.3.3 實(shí)驗(yàn)結(jié)果分析

選用不同客戶規(guī)模的算例（4.3.1、4.3.2 小節(jié)中的算例及MDVRPTW 標(biāo)準(zhǔn)算例pr03、pr04）對(duì)AGA-LNS 進(jìn)行算法分析及敏感性分析，未給出的其他參數(shù)與算例1設(shè)置相同。

（1）算法分析

為驗(yàn)證AGA-LNS 的有效性，采用自適應(yīng)遺傳算法（AGA）、大鄰域搜索算法（LNS）求解不同客戶規(guī)模算例的配送方案，并與AGA-LNS 的求解結(jié)果進(jìn)行對(duì)比分析。其中AGA 采用本文設(shè)計(jì)的交叉變異操作，LNS 選用隨機(jī)移除算子和貪婪插入算子組合、最遠(yuǎn)距離移除算子和貪婪插入算子組合。表6 給出了3 種算法運(yùn)行10次的結(jié)果，其中C為10次運(yùn)算結(jié)果中的最優(yōu)解，t為算法平均耗時(shí)，GAP1 為AGA-LNS 與AGA、LNS 所求結(jié)果的誤差，GAP2 為AGA-LNS 與AGA、LNS 求解耗時(shí)的誤差。由表6 可知，AGA-LNS 與AGA、LNS 相比，在求解時(shí)長(zhǎng)和求解質(zhì)量上都有顯著優(yōu)勢(shì)，因此AGA-LNS能較好地解決VRPMTWTDN；對(duì)比AGA 和AGA-LNS的求解結(jié)果可知，客戶規(guī)模相同時(shí)，AGA-LNS的求解結(jié)果相對(duì)于AGA 的求解結(jié)果均有所改善，平均改善了2.73%，可見在自適應(yīng)遺傳算法的交叉變異操作后加入大鄰域搜索可改善自適應(yīng)遺傳算法陷入局部最優(yōu)解的問題，提高求解質(zhì)量。

表6 不同算法求解結(jié)果對(duì)比Table 6 Comparison of solution results of different algorithms

（2）敏感性分析

①混合時(shí)間窗比例的影響

為驗(yàn)證軟硬時(shí)間窗客戶數(shù)量的比例變化對(duì)配送方案制定的影響，本文針對(duì)不同客戶規(guī)模，分別設(shè)置了6種軟硬時(shí)間窗客戶數(shù)量比例（軟時(shí)間窗客戶數(shù)∶硬時(shí)間窗客戶數(shù)=S∶H）。表7給出了算法運(yùn)行10次的結(jié)果，其中N為客戶規(guī)模，C為10次運(yùn)算結(jié)果中的最優(yōu)解，GAP為不同軟硬時(shí)間窗客戶數(shù)量比例下10次運(yùn)算結(jié)果的最優(yōu)解與S∶H=1∶1 時(shí)求得最優(yōu)解的相對(duì)誤差。由表7 可知，相同客戶規(guī)模下，最優(yōu)解隨著硬時(shí)間窗客戶數(shù)量比例的增大而增大，即軟時(shí)間窗和硬時(shí)間窗客戶總數(shù)一定，硬時(shí)間窗客戶越多，其最優(yōu)配送方案總成本越高。

表7 軟硬時(shí)間窗客戶不同比例下最優(yōu)配送方案總成本Table 7 Total cost of optimal distribution scheme under different ratios of customers with hard and soft time windows

②車輛行駛速度的影響

為驗(yàn)證車輛行駛速度對(duì)配送方案制定的影響，本文對(duì)不同客戶規(guī)模算例在車輛勻速行駛情況下進(jìn)行求解。表8給出了算法運(yùn)行10次運(yùn)算結(jié)果，其中N為客戶規(guī)模，Case1 代表速度連續(xù)變化，Case2 代表速度恒為40 km/h，Case3代表速度恒為50 km/h，Case4代表速度恒為60 km/h，C為10次運(yùn)算結(jié)果中的最優(yōu)解，GAP為車輛勻速行駛時(shí)10次運(yùn)算結(jié)果中的最優(yōu)解與車速時(shí)間依賴時(shí)求得最優(yōu)解的相對(duì)誤差。由表8可知，相同客戶規(guī)模下，不同行駛速度下求得的最優(yōu)解與速度時(shí)間依賴時(shí)求得的最優(yōu)解具有一定誤差，如客戶規(guī)模為48時(shí)，誤差分別為1.05%、-0.79%、-1.68%。其中車輛行駛速度為50 km/h時(shí)的求解結(jié)果與本文設(shè)置的速度時(shí)間依賴函數(shù)下求解結(jié)果的誤差最小，車輛行駛速度增大或減小，與本文設(shè)置的速度時(shí)間依賴函數(shù)下求解結(jié)果的誤差均增大。綜合分析，在制定配送方案時(shí)考慮車輛行駛速度時(shí)間依賴性具有重要意義。

表8 不同車速下最優(yōu)配送方案總成本Table 8 Total cost of optimal distribution scheme under different vehicle speeds

表9給出了速度連續(xù)變化時(shí)與恒定速度時(shí)各自運(yùn)行10 次的平均求解耗時(shí)。其中Case1 代表速度連續(xù)變化，Case2代表速度恒為50 km/h，GAP為速度連續(xù)變化時(shí)算例求解的平均耗時(shí)與恒定速度時(shí)算例求解的平均耗時(shí)的相對(duì)誤差。由表9可知，速度連續(xù)變化時(shí)算例求解平均耗時(shí)比速度恒定時(shí)算例求解平均耗時(shí)平均多63.82%，因此考慮車輛速度連續(xù)變化增加了算法實(shí)現(xiàn)的時(shí)間成本約63.82%，但配送中心通常在配送前一天根據(jù)已有訂單制定配送方案，因此速度時(shí)間依賴下的算法時(shí)間成本是可行的。

表9 不同車速下算例求解耗時(shí)Table 9 Example solving time under different vehicle speeds

5 結(jié)束語

本文考慮配送路網(wǎng)中車輛速度的連續(xù)變化，并按照客戶對(duì)配送時(shí)間的不同要求對(duì)其服務(wù)時(shí)間窗進(jìn)行分類，對(duì)VRPMTWTDN 進(jìn)行研究，建立了以配送成本最小化為目標(biāo)的優(yōu)化模型，并設(shè)計(jì)了自適應(yīng)遺傳-大鄰域搜索算法求解配送方案，并實(shí)驗(yàn)驗(yàn)證了其有效性，結(jié)論如下：

（1）模型中不僅考慮了配送中心間配送資源的共享及混合時(shí)間窗，還考慮了車輛速度連續(xù)變化及車速、載重量對(duì)油耗的影響，雖增加了問題求解難度，但更貼近配送生產(chǎn)活動(dòng)的實(shí)際。

（2）AGA-LNS 引入時(shí)差插入法改進(jìn)交叉、變異算子，并嵌入移除、插入算子。針對(duì)4 個(gè)不同客戶規(guī)模的算例，AGA-LNS 與AGA 和LNS 相比，求解質(zhì)量平均改進(jìn)了2.73%和21.27%，求解時(shí)長(zhǎng)平均縮短了7.63%和2.13%，求解質(zhì)量高，收斂速度快，因此AGA-LNS 能較好地求解帶混合時(shí)間窗約束的車輛路徑問題。

（3）不同混合時(shí)間窗比例敏感性分析表明，相同客戶規(guī)模下，硬時(shí)間窗客戶數(shù)量比例越小，配送成本越低；不同車輛行駛速度敏感性分析表明，不同行駛速度對(duì)配送計(jì)劃的制定及配送成本有較大影響，因此在優(yōu)化配送方案時(shí)所設(shè)車輛行駛速度應(yīng)盡量貼近現(xiàn)實(shí)。

未來研究將考慮實(shí)時(shí)交通信息、電動(dòng)車配送等因素，使問題更貼近實(shí)際，并改進(jìn)、開發(fā)更加有效的求解算法。