徐林浩，錢 斌，胡 蓉，于乃康

（1. 昆明理工大學(xué) 信息工程與自動化學(xué)院，云南 昆明 650500；2. 昆明理工大學(xué) 機電工程學(xué)院，云南 昆明 650500）

車輛路徑問題(Vehicle Routing Problem,VRP)由Dantzig和Ramser于1959年首次提出[1]，VRP作為經(jīng)典的組合優(yōu)化問題，得到了廣大學(xué)者的研究。帶容量的車輛路徑問題(Capacitated Vehicle Routing Problem,CVRP)是VRP的一類，廣泛存在交通運輸、物流配送等方面，同時，當今社會推行綠色低碳發(fā)展，綠色車輛路徑問題受到了重視。在上述背景下，研究綠色帶容量的車輛路徑問題(Green Capacitated Vehicle Routing Problem, GCVRP)具有重要的現(xiàn)實意義。由于VRP為NP-hard問題，而VRP可歸約為GCVRP，故GCVRP也屬于NP-hard問題，對其展開研究具有較大的理論價值。

VRP的求解方法可以分為兩類，一類為智能優(yōu)化算法，如遺傳算法、蟻群算法和禁忌搜索算法等[2-4]，該類算法通過建立問題的排序模型進行求解，一般不能保證得到全局最優(yōu)解以及無法判斷所求解的優(yōu)劣程度；另一類為精確算法，如分支定界、列生成和拉格朗日松弛算法等[5-10]，該類算法通過建立問題嚴格的數(shù)學(xué)模型，能在合理時間內(nèi)求得小規(guī)模問題最優(yōu)解，針對大規(guī)模問題，該類算法能在合理時間內(nèi)求得問題的精確下界，由于該類算法所求解為問題最優(yōu)解或逼近最優(yōu)解的下界，能較大程度保證解的求解質(zhì)量，同時求解時間也基本合理。

近年來，相關(guān)學(xué)者對VRP及其拓展問題的精確算法或基于精確算法的混合算法進行了研究。Bouzid等[11]針對CVRP，以最小化路徑成本為優(yōu)化目標，建立了整數(shù)規(guī)劃模型，并提出一種拉格朗日分裂(Lagrangian Split, LS)與變鄰域搜索(Variable Neighbourhood Search, VNS)結(jié)合的混合算法進行求解。Singh等[5]針對具有模糊需求的CVRP，建立了該問題的數(shù)學(xué)模型，并以最小化成本為優(yōu)化目標，提出一種分支定界算法(Branch and Bound, B&B)進行求解。Lysgaard等[6]針對累積容量車輛路徑問題(Cumulative Capacitated Vehicle Routing Problem,CCVRP)，建立了問題的集劃分模型，以最小化路徑總成本為優(yōu)化目標，采用(Branch-and-Cut-and-Price,BCP)進行求解。Majid等[7]針對隨機需求的車輛路徑問題(Vehicle Routing Problem with Stochastic Demands, VRPSD)，建立了該問題的數(shù)學(xué)模型，以最小化總費用為優(yōu)化目標，采用整數(shù)L形(L-shaped)算法在分支切割框架下求解VRPSD，實驗結(jié)果驗證了所提算法的有效性。Christiansen等[8]針對VRPSD，以最小化成本為求解目標，提出了一種分支定價算法(Branch-and-price, B&P)進行求解。

當今社會持續(xù)推進節(jié)能減排，堅持綠色發(fā)展，考慮燃油消耗和碳排放等因素的綠色車輛路徑問題(Green Vehicle Routing Problem, GVRP)受到了重視。Andelmin等[9]針對求解目標為最小化行駛總距離的GVRP，并將該問題建模為集合劃分問題，提出了一種結(jié)合K-路(K-path)切割的列生成(Column Generation, CG)算法進行求解。Dabia等[10]研究污染路徑問題(Pollution-routing Problem, PRP)，以最小化運送成本為求解目標，運用CG算法對主問題進行了求解，并提出一種定制的標記算法解決了定價問題，實驗結(jié)果驗證了所提算法的有效性。?a?r?等[12]提出了GCVRP的混合整數(shù)規(guī)劃模型(Mixed Integer Programming, MIP)，并以最小化行駛距離為求解目標，同時提出了一種基于模擬退火(Simulated Annealing, SA)的精確解方法進行求解，作者運用改進B&P算法改進下界，并引用了一種基于SA的啟發(fā)式算法來獲得上界，實驗結(jié)果驗證了所提算法的有效性。Qiu等[13]研究碳排放限額交易下的污染生產(chǎn)路徑問題(Pollution Production-routing Problem, PPRP)，建立了該問題的數(shù)學(xué)模型并以最小化運營成本為求解目標，提出一種分支定價啟發(fā)式算法進行求解，仿真實驗證明了該模型具有降低二氧化碳排放水平和運營成本的潛力。由文獻調(diào)研可知，目前求解GCVRP的精確算法或與精確算法相結(jié)合的算法相對較少。

拉格朗日松弛算法(Lagrange Relaxation, LR)是一種求解復(fù)雜組合優(yōu)化問題的有效方法，LR主要思想是將問題的難約束引入目標函數(shù)得到松弛問題，再通過不斷更新迭代乘子，逐步求得松弛問題的解。運用傳統(tǒng)LR求解GCVRP，往往只能求得問題的下界，該下界一般為不可行解。因此，本文綜合LR和啟發(fā)式算法優(yōu)勢，采用兩者相結(jié)合的算法求解GCVRP，首先運用LR求得問題下界，再根據(jù)問題性質(zhì)，設(shè)計相應(yīng)啟發(fā)式算法對下界進行修復(fù)及優(yōu)化，得到問題的上界，通過上界和下界的間隙(Gap)，可以衡量算法性能，是一種求解GCVRP的有效方法。

綜上，本文針對實際生活中廣泛存在的GCVRP，建立以最小化運輸成本為優(yōu)化目標的混合整數(shù)規(guī)劃模型，并提出一種改進拉格朗日算法進行求解，具體為運用拉格朗日松弛算法得到GCVRP的松弛問題，進而得到原問題的對偶模型，并采用一種次梯度法求解該對偶模型，得到原問題的一個下界；然后，采用修復(fù)算法及鄰域搜索算法對下界進行修復(fù)和優(yōu)化，得到原問題的一個上界。最后，通過上下界的間隙大小以及對比Gurobi求解器所求解來衡量所提算法的有效性。

1 GCVRP問題描述與數(shù)學(xué)模型

1.1 GCVRP問題描述

GCVRP是在CVRP的基礎(chǔ)上進一步考慮了綠色因素，該問題建立在已知倉庫和客戶點坐標(坐標以km為單位)、客戶需求及車輛載重條件下，合理調(diào)度倉庫中的車輛為客戶送貨，使運輸費用最小。圖1為問題示意圖。

圖1 GCVRP問題示例Fig.1 GCVRP problem example

問題相關(guān)假設(shè)如下：

(1) 倉庫內(nèi)所有車輛均參與配送；

(2) 每輛車僅執(zhí)行一次配送任務(wù)；

(3) 車輛從倉庫出發(fā)，完成服務(wù)后返回原倉庫；

(4) 每個客戶只能被一輛車服務(wù)一次；

(5) 每輛車配送的貨物總量不能超過其載重；

(6) 車輛在配送過程中保持勻速行駛，忽略道路等因素對車輛的影響。

GCVRP相關(guān)符號定義見表1。

表1 符號及定義Table 1 Symbols and definitions

1.2 GCVRP混合整數(shù)規(guī)劃模型

其中：式(1)為目標函數(shù)；式(2)要求從倉庫出入車輛數(shù)目相等；式(3)要求倉庫內(nèi)所有車均參與配送，且配送完返回原倉庫；式(4)～(6)要求每個客戶只能被一輛車服務(wù)一次；式(7)為流平衡約束；式(8)要求倉庫的每輛車的載貨量不能超過其額定載重；式(9)為MTZ子環(huán)消除約束；式(10)表示車輛從倉庫出發(fā)時的載貨量等于所服務(wù)客戶的總需求量；式(11)表示車輛完成服務(wù)任務(wù)，返回倉庫的載貨量為0；式(12)表示某輛車服務(wù)完當前客戶點后的貨物量等于從該客戶點離開的貨物量；式(13)和(14)為決策變量。

1.3 碳排放量計算方法

車輛行駛過程中的碳排放與燃油的消耗直接相關(guān)，影響燃油消耗的因素較多，計算也相對復(fù)雜。本文在計算碳排放量上，采用文獻[14-16]中的計算方法，相關(guān)參數(shù)定義見表2。計算公式為

表2 參數(shù)定義與取值Table 2 Parameter definition and value

2 改進拉格朗日松弛算法

本節(jié)針對上節(jié)所提問題模型，提出一種結(jié)合拉格朗日松弛算法與啟發(fā)式算法的改進拉格朗日松弛算法進行求解。GCVRP求解流程如圖2所示。

圖2 GCVRP求解流程圖Fig.2 GCVRP solution flow chart

2.1 松弛載重約束

由于上節(jié)所提GCVRP具有NP難屬性，直接求解其最優(yōu)解非常困難，但求取該問題的下界較為簡單。LR是一種求解問題下界的有效方法，通過松弛問題中的難約束，可以大大減小問題的求解難度。通過對GCVRP進行分析，車輛的載重約束(8)影響了車輛的行駛距離，會影響目標值的求取，松弛這條約束會使原問題容易求解。本節(jié)將約束條件(8)松弛到目標函數(shù)中，得到拉格朗日松弛模型為

2.2 構(gòu)造可行解

求解對偶問題所得解一般是不可行解，因為可能違反車輛的載重約束條件，因此需要先對不可行解進行修復(fù)，得到原問題的可行解，再運用局部搜索操作對可行解進行優(yōu)化，得到原問題高質(zhì)量解。本小節(jié)采用修復(fù)算法和鄰域搜索算法對不可行解進行修復(fù)及優(yōu)化。

2.2.1 修復(fù)算法

2) 修復(fù)不可行解。

針對上述獲得的不可行解，運用以下算法進行修復(fù)，設(shè)車輛數(shù)為M，車輛編號為m，算法步驟如下：

(1) 將二進制編碼解轉(zhuǎn)換成排序解。

(2) 依次判斷車輛m(m=1,2,···,M)是否滿足載重約束，若滿足，將車輛m放入可行車輛集合X；否則，將車輛m中超過其服務(wù)能力的客戶點按服務(wù)先后順序摘出，放入剩余客戶點集合Y中，操作后，車輛m滿足載重約束，將車輛m放入可行車輛集合X中。直至判斷完所有車輛。

(3) 以最小化行駛距離為目標，并以車輛載重為約束條件，將剩余客戶點插入所有車輛(M輛)的所有位置，找出最佳插入點，將客戶插入，對剩余客戶按順序逐個執(zhí)行步驟3。若剩余客戶均分配完，輸出修復(fù)后的可行解，執(zhí)行步驟7；否則，執(zhí)行步驟4。

(4) 選出M輛車中剩余容量最大的車m1，將剩余客戶插入車m1的 最后一個位置，同時，從車m1中選出某一客戶點插入車m2(m2=1,2,···,M,m1≠m2)的最后一個位置，使兩車均滿足容量約束。對剩余客戶按順序逐個執(zhí)行步驟4。

(5) 若剩余客戶點分配完，輸出修復(fù)后的可行解，執(zhí)行步驟(7)；否則，執(zhí)行步驟(6)。

(6) 將剩余客戶點插入車輛m(m=1,2,···,M)的最后一個位置，同時，從車輛m選出某一客戶插入除車輛m的其他車輛的最后一個位置，使兩車均滿足容量約束。如果滿足，則插入；否則，直至判斷完所有車輛。對剩余客戶點按順序逐個執(zhí)行步驟(6)。

(7) 程序結(jié)束。

2.2.2 鄰域搜索算法

本小節(jié)針對上述獲得的可行解，設(shè)計了車輛間和車輛內(nèi)兩類共4種鄰域搜索操作，以進一步提高可行解的質(zhì)量。算法步驟如下：

(1) 設(shè)可行解為 Π，令當前最優(yōu)解Π?=Π，局部搜索次數(shù)為m axs， 車輛間搜索次數(shù)為vout，車輛內(nèi)搜索次數(shù)為vin， 車輛數(shù)為M，循環(huán)參數(shù)：l oop=1、Nout=1、Nin=1、Nv=1， 標志位 S I=1。

(2) 在 Π中隨機選取兩輛車，在兩輛車中分別隨機選取一個客戶，得到客戶r和t，將t客戶插到r客戶之前，若不滿足載重約束條件，令 Π =Π?，執(zhí)行步驟(3)；反之，若獲得更優(yōu)解 Π′， 則令Π?=Π′，Π =Π′，執(zhí)行步驟(4)，若沒獲得更優(yōu)解，執(zhí)行步驟(3)。

(3) 在 Π中隨機選取兩輛車，在兩輛車中分別隨機選取一個客戶進行交換，若不滿足載重約束條件，令Π =Π?，執(zhí)行步驟(4)；反之，若獲得更優(yōu)解Π′，則令Π?=Π′， Π =Π′，執(zhí)行步驟(4)。

2.3 次梯度算法

針對上述所提的拉格朗日對偶問題，本文采用次梯度算法進行求解。其求解基本思想是：根據(jù)已知條件，求解對偶問題，獲得問題的下界和決策變量值，進而計算更新次梯度以及拉格朗日乘子，通過不斷更新迭代求解對偶問題，直到找到滿足條件的對偶值。用次梯度算法求解對偶問題時，拉格朗日乘子λm(m∈C)更新公式為

2.4 改進拉格朗日松弛算法步驟

整個求解算法步驟如下：

(1) 初始化拉格朗日乘子 λ0=0 ，迭代次數(shù)n=0，設(shè)初始上界U B為原問題較大的可行目標值。

(2) 求解對偶問題(17)，獲得決策變量xn及下界LBn。

(3) 運用2.2.1中的修復(fù)算法對下界進行修復(fù)，若能修復(fù)得可行解 Π，執(zhí)行步驟(4)；否則上界U B保持不變，執(zhí)行步驟(5)。

(4) 運用2.2.2中的二階段算法對可行解 Π進行優(yōu)化，獲得優(yōu)化后的解 Π?，如果優(yōu)化后的目標值Z(Π?)

(5) 運用公式(18)、(19)和(20)更新次梯度、步長和拉格朗日乘子。

(6) 判斷是否滿足終止條件，若滿足，輸出最佳上界 UB 和下界L Bn，程序結(jié)束；否則繼續(xù)執(zhí)行步驟(2)～(5)。

3 實驗結(jié)果與分析

3.1 實驗環(huán)境及數(shù)據(jù)

3.2 結(jié)果分析

本文提出一種改進拉格朗日松弛算法求解GCVRP，與商業(yè)求解器Gurobi所求解進行了對比，求解結(jié)果見表3，從表可以看出：針對編號1～7的小規(guī)模算例，Gurobi求解器在1 800 s內(nèi)能求得4個算例的最優(yōu)解。小規(guī)模算例所求解平均G ap為5.78%，求解平均時間為825.277 s。改進拉格朗日松弛算法在合理時間內(nèi)求得了問題的滿意解，求解平均G ap為7.17%，求解平均時間為219.896 s。這也體現(xiàn)Gurobi求解器在求解小規(guī)模問題的優(yōu)勢，能大概率求解問題的最優(yōu)解，但求解時間較長。改進拉格朗日松弛算法求解效率較高，在較短時間內(nèi)能求得問題的上下界，能求得問題的滿意解。

表3 不同規(guī)模問題算法對比1)Table 3 Comparison of algorithms for different scale problems

針對編號為8～19的中大規(guī)模算例，隨著問題規(guī)模的變大，解空間更加復(fù)雜，Gurobi求解器弊端凸顯，12個算例的平均G ap為21.12%，平均求解時間為5 100 s，求解質(zhì)量和求解效率不高。而改進拉格朗日松弛算法通過松弛難約束，使問題變得容易求解，求解效率更高，求解質(zhì)量令人滿意，所求12個算例的平均G ap為7.87%，平均求解時間為2 422.016 s，求解質(zhì)量遠優(yōu)于Gurobi求解器。

綜合來看，針對19個不同規(guī)模的算例，改進拉格朗日松弛算法能較好地求解GCVRP，在合理時間內(nèi)，所有算例的平均G ap為7.61%。其在求得問題上界的同時，也能為原問題提供一個下界，可以衡量所求解的質(zhì)量，綜合了傳統(tǒng)拉格朗日松弛算法與啟發(fā)式算法的優(yōu)勢，是求解GCVRP的有效算法。

4 結(jié)論

本文在VRP的基礎(chǔ)上，進一步考慮車輛載重和綠色因素，建立了GCVRP的混合整數(shù)規(guī)劃模型，并以總費用最小為求解目標，提出一種改進拉格朗日松弛算法進行求解。首先，采用次梯度法求解對偶問題，得到原問題的下界；其次利用修復(fù)算法和鄰域操作對下界進行修復(fù)及優(yōu)化，得到原問題高質(zhì)量下界，實現(xiàn)上界下界并行求解。實驗結(jié)果表明，本文算法的綜合求解性能優(yōu)于Gurobi，是一種求解GCVRP的有效方法。后續(xù)將進一步考慮多車場多車型因素，并設(shè)計有效求解算法。