考慮碳排放的兩階段選址-路徑問題及其算法

湯希峰 ，何 杰 ，張 浩

（1.河海大學(xué)土木與交通學(xué)院，江蘇 南京 210098；2.東南大學(xué)交通學(xué)院，江蘇 南京 211189）

城市是生產(chǎn)、生活的主要集中地，也是物流活動最為集中的地方.兩層的物流網(wǎng)絡(luò)作為一種可行的城市物流解決方案正被國內(nèi)外許多大中城市所采用[1-4]，其核心問題之一是兩階段選址-路徑問題 （2ELRP）.經(jīng)典的2E-LRP 可以描述為：重型或中型貨車負責(zé)將貨物從一個位于城郊的物流園區(qū)運送至若干位于城市中心區(qū)周邊的配送中心，貨物在這些配送中心經(jīng)過重新配裝后再由輕型或微型貨車送達終端客戶，其目標(biāo)是尋求最優(yōu)的配送中心選址以及第一階段和第二階段貨車的最優(yōu)路徑以實現(xiàn)總的物流成本最小化[5].

眾所周知，貨車相較于客車容易導(dǎo)致交通事故和交通擁堵，更是污染物和碳排放的主要貢獻者.以我國為例，2019 年全國貨車保有量雖然只占汽車總量的10.71%，貨車CO、HC、NOx、PM 排放量卻占到了汽車排放總量的29.7%、26.3%、83.5%、90.1%[6].隨著生態(tài)環(huán)境保護問題日漸受到重視，如何優(yōu)化2 層網(wǎng)絡(luò)結(jié)構(gòu)的物流運作以減少貨車的碳排放引起了國內(nèi)外學(xué)者越來越多的研究興趣.Govindan 等[7]圍繞生鮮配送，建立了以物流成本最小化和環(huán)境影響最小化為目標(biāo)的2E-LRP 多目標(biāo)優(yōu)化模型，并設(shè)計了多目標(biāo)粒子群優(yōu)化算法（MOPSO）和改進的多目標(biāo)可變鄰域搜索算法（AMOVNS）2 種求解算法；胡大偉等[8]針對燃油車和電動車不同的行駛與排放特性，分別建立了燃油車和電動車的2 階段開放式選址-路徑問題模型，并提出了一種改進的模擬退火算法；劉成清等[9]以經(jīng)典的2 階段開放式選址-路徑問題模型為基礎(chǔ)，通過考慮車輛的燃油消耗和CO2排放以及車輛路徑選擇的靈活性，建立了基于路徑靈活性的兩階段開放式低碳選址-路徑問題模型，并運用Cplex進行了求解；Pitakaso 等[10]以燃油消耗最小化為目標(biāo)，建立了考慮碳排放的綠色兩階段選址路徑問題（G2ELRP）優(yōu)化模型，并設(shè)計了相應(yīng)的可變鄰域策略的自適應(yīng)搜索（VANSAS）算法.

一方面，鑒于車輛的碳排放受到車型、車速、載重、行駛距離以及駕駛行為等諸多因素的影響[11-12]，現(xiàn)有文獻大多基于車輛的燃油消耗進行碳排放的計算，所需參數(shù)較多，實用性不強；另一方面，現(xiàn)有的用于求解各種2E-LRP 的算法多適用于求解小規(guī)模問題 （客戶數(shù)量n≤100, 候選配送中心數(shù)量m<10），求解大規(guī)模問題 （n>100,m≥10） 時計算耗時較長，很難實際應(yīng)用.為此，本文基于一種簡單實用的以排放因子為主要參數(shù)的碳排放計算方法，建立了考慮碳排放的2E-LRP 優(yōu)化模型，并針對所研究問題的計算復(fù)雜性，設(shè)計了一種可快速求解大規(guī)模問題的兩階段混合算法.

1 模型構(gòu)建

1.1 問題描述

基于圖論，考慮碳排放的2E-LRP 可以描述為G={V,E}，V=V0∪Vs∪Vc，V0為物流園區(qū)，Vs= （1,2,…,m）為配送中心候選集，Vc= （m+1,m+2,…,m+n）為客戶集；E=E1∪E2，E1= {（i,j）:i,j∈V0∪Vs;i≠j}為物流園區(qū)與配送中心之間的路線集，E2= {（i,j）:i,j∈Vs∪Vc;i,j?Vs×Vs;i≠j} 為配送中心與客戶以及客戶與客戶之間的路線集；H為可供物流園區(qū)使用的重型或中型貨車的集合，每臺重型或中型貨車具有相同的載重能力，記為QH，滿載和空載時的碳排放系數(shù)分別為εvfH、εveH（單位： kg/km, 下同）；K為可供配送中心共用的輕型或微型貨車的集合，每臺輕型或微型貨車也具有相同的載重能力，記為QK（QK

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

1.2 數(shù)學(xué)模型

若用qijh和qijk分別為貨車h（h∈H）在路線 （i,j）∈E1、貨車k（k∈K）在路線（i,j）∈E2上行駛時的載重量；決策變量xijh=1 表示貨車h（h∈H） 經(jīng)過路線 （i,j）∈E1，否則xijh=0；yijk=1 表示貨車k（k∈K）經(jīng)過路線 （i,j）∈E2，否則yijk=0；zi=1 表示配送中心i（i∈Vs） 被選用，否則zi=0；uij=1 表示客戶j由配送中心i提供服務(wù)，否則uij=0；則以碳排放最小化為目標(biāo)的2E-LRP 數(shù)學(xué)模型可表示為

式 （1） 為目標(biāo)函數(shù)，表示從物流園區(qū)到配送中心以及配送中心到客戶2 個階段的貨車碳排放量之和最小，與經(jīng)典2E-LRP 模型的目標(biāo)函數(shù) （關(guān)于車輛通行成本或車輛行駛距離的函數(shù)） 不同，是一個關(guān)于車輛載重和車輛行駛距離的函數(shù)；式（2）、（3） 表示2 個階段中某一貨車要么不被使用，要么至多被使用一次；式（4）、（5） 保證2 個階段中貨車路線的連續(xù)性，同時保證貨車完成任務(wù)后返回原先出發(fā)的物流園區(qū)或配送中心；式（6） 表示客戶只能由一個配送中心提供服務(wù)；式（7） 表示客戶若由某一配送中心提供服務(wù)，則必被從該配送中心出發(fā)的貨車訪問一次；式（8） 表示某一配送中心要么不被選用，若被選用，所服務(wù)的客戶總需求量不能超過其作業(yè)能力；式（9）、（10） 表示2 個階段中貨車服務(wù)的客戶需求均不能超過其載重能力；式（11）、（12） 為2 個階段中路線的完整性約束；式（13）、（14）表示2 個階段中貨車載重量在相鄰路段上的數(shù)量關(guān)系；式（15） ～（18）為決策變量的取值約束.

2 算法設(shè)計

由于設(shè)施選址問題和車輛路徑問題都屬于NPhard 問題，所以本文研究的考慮碳排放的2E-LRP是更為復(fù)雜的NP-hard 問題，求解難度很大，尤其是隨著問題規(guī)模的增大，計算耗時呈指數(shù)增加.鑒于此，本文設(shè)計了一種可快速求解大規(guī)模問題的2 階段混合算法（TSHA ），算法流程如圖2 所示.第1 階段，通過將2E-LRP 轉(zhuǎn)化成不考慮車輛路徑的2 階段設(shè)施選址問題 （2E-FLP），直接求解配送中心選址方案和客戶分配方案；第2 階段，基于得到的配送中心選址和客戶分配方案，物流園區(qū)到被選用的配送中心以及每個配送中心到所分配客戶的車輛路徑問題就變成了獨立且較易求解的VRP 問題，進而加以解決.

圖2 兩階段混合算法的算法流程Fig.2 Flowchart of the TSHA

2.1 第1 階段算法

由于迭代算法在求解過程中非常低效，為提高所建模型的求解效率，TSHA 在第1 階段先將考慮碳排放的2E-LRP 轉(zhuǎn)化成不考慮車輛路徑的2EFLP，模型如下：

式 （19） 為目標(biāo)函數(shù)，表示客戶到配送中心以及配送中心到物流園區(qū)的以需求量為權(quán)重的距離之和最小，約束條件對應(yīng)原問題模型中的式 （6）、（8）、式（15）、（16）.不難看出，上述2E-FLP 模型為線性規(guī)劃模型，可直接運用Cplex 求解器進行求解，從而得到配送中心選址和客戶分配方案.

2.2 第2 階段算法

基于得到的配送中心選址和客戶分配方案，物流園區(qū)到被選用的配送中心以及每個被選用配送中心到所分配客戶的車輛路徑問題都變成了獨立的VRP 問題.為求解以碳排放最小化為目標(biāo)的VRP 問題，TSHA 在第2 階段采用了一種改進的蟻群算法，如圖3 所示.

圖3 改進蟻群算法的偽代碼Fig.3 Pseudocode of the improved ant colony algorithm

與經(jīng)典的蟻群算法相同[13]，改進的蟻群算法也包含2 個關(guān)鍵步驟：1） 主要負責(zé)路線的構(gòu)建，即基于當(dāng)前節(jié)點逐步選擇下一個將要訪問的節(jié)點；2） 主要負責(zé)信息素濃度的更新，具體包括局部信息素更新和全局信息素更新.

從式 （1） 碳排放的計算式可以看出，貨車的碳排放與其載重量和行駛距離密切相關(guān).為此，改進的蟻群算法采用了一種新的路線構(gòu)建規(guī)則，如式 （20） 、式（21） 所示.

式中：Ω為尚未被訪問的節(jié)點集，Ω∈Vs或Vc；τij為當(dāng)前節(jié)點i與未到訪節(jié)點j之間路線上的信息素濃度；ηij為節(jié)點i與節(jié)點j之間距離的倒數(shù)；α、β為重要度系數(shù)；q0為給定的常數(shù)；q為 [0, 1] 上的隨機數(shù).

當(dāng)隨機產(chǎn)生的q≤q0時，選擇使取值最大的j作為下一個將要訪問的節(jié)點；當(dāng)q>q0時，則根據(jù)的值所占的比例選擇下一個將要訪問的節(jié)點j，在TSHA 算法中采用輪盤賭的方式進行選擇.

每只螞蟻都按照上述的偽比例選擇規(guī)則構(gòu)建路線，當(dāng)螞蟻構(gòu)建出一條路線（即達到貨車的載重限制時，返回原先出發(fā)的物流園區(qū)或配送中心），然后重新出發(fā)構(gòu)建下一條路線，如此反復(fù)，直至遍訪所有客戶.當(dāng)螞蟻構(gòu)建出一個完整的路線方案后，對該路線方案中的每條路線進行局部信息素更新，第t+ 1 次信息素濃度為

式中：ρ（0 ≤ρ≤1） 為信息素的揮發(fā)系數(shù)；τ0為初始的信息素濃度，取值為1/（mEnns），Enns為基于最近鄰搜索算法得到的貨車碳排放量.

當(dāng)蟻群中每只螞蟻均構(gòu)建出一個完整的路線方案后，從中選擇一個碳排放量最小的路線方案作為當(dāng)次最優(yōu)的路線方案，為在當(dāng)次最優(yōu)路線方案和當(dāng)前最優(yōu)路線之間取得平衡，亦即鼓勵隨后的螞蟻在當(dāng)次最優(yōu)路線和當(dāng)前最優(yōu)路線附近的解空間繼續(xù)搜索，改進的蟻群算法采用Ting 等[14]推薦的策略對當(dāng)次最優(yōu)和當(dāng)前最優(yōu)路線方案中的每條路線進行全局信息素更新，如式（23）、（24）所示.

式中：EIW為當(dāng)次最差路線方案產(chǎn)生的碳排放.

如此，通過若干次循環(huán)后即可得到所有以碳排放最小化為目標(biāo)的VRP 的最優(yōu)路線方案.考慮到由物流園區(qū)提供服務(wù)的配送中心數(shù)量較少而由每個配送中心提供服務(wù)的客戶數(shù)量相對較多，對得到的配送中心到客戶的VRP 最優(yōu)路線方案，本文運用3 種局部搜索算法對其再進行改進：cross_relocation 算法將從某個配送中心出發(fā)路線上的客戶逐一插入從其他配送中心出發(fā)的路線中，若得到的路線可行且碳排放量小于原路線，則用改進后的路線替代原路線；cross_swap 算法將從不同配送中心出發(fā)的兩條路線上的客戶兩兩交換并重新生成新的路線，若新路線可行且碳排放量較小，則替代原路線；最后，再對每一條路線運用3-opt 算法進行改進.至此，改進后的路線方案即可作為最終的最優(yōu)路線方案.

3 算例測試

由于目前還沒有可用于測試2E-LRP 碳排放的標(biāo)準(zhǔn)算例集，所以本文選取了Prodhon 標(biāo)準(zhǔn)算例集[15]中全部6 個最大規(guī)模的算例稍加修改，用以測試本文所提出的模型和算法.Prodhon 算例集中6 個最大規(guī)模的算例均包括200 個客戶和10 個候選配送中心.算例的修改包括去掉原算例中的配送中心選址費用、增加不同類型貨車的碳排放系數(shù)，其中重型或中型貨車滿載和空載時的碳排放系數(shù)分別為2.392、1.638 kg/km，輕型或微型貨車滿載和空載時的碳排放系數(shù)分別為1.096、0.772 kg/km，其他數(shù)據(jù)保持不變.

設(shè)計的TSHA 算法采用MATLAB 2018b 編程并嵌入Cplex V12.7 求解器，改進的蟻群算法控制參數(shù)包括：求解物流園區(qū)到被選用的配送中心VRP 問題時，循環(huán)次數(shù)設(shè)定為選用的配送中心個數(shù)的1/2，蟻群中螞蟻的數(shù)量與選用配送中心個數(shù)相同；求解配送中心到被分配客戶的VRP 問題時循環(huán)次數(shù)均設(shè)定為被分配客戶數(shù)量的1/2，蟻群中螞蟻的數(shù)量均等于被分配的客戶的數(shù)量；重要度系數(shù)α= 2，β= 1；給定常數(shù)q0= 0.5；揮發(fā)系數(shù)ρ= 0.2.針對每個算例，獨立運行程序20 次，并分別記下計算結(jié)果和時間.

3.1 對TSHA 算法思想的驗證

為驗證TSHA 的算法思想，本文先編寫了基于同樣算法思想的TSHA-Ⅱ算法用以求解Prodhon 算例集中6 個最大規(guī)模的原算例.TSHA-Ⅱ與TSHA略微不同的是：1） 用物流成本代替碳排放作為目標(biāo)函數(shù)；2） 改進的蟻群算法中，路線構(gòu)建不考慮客戶的貨運需求量，即式 （20） 、（21） 中去掉啟發(fā)因子dj.TSHA-Ⅱ算法的求解結(jié)果與目前文獻中兩種最新算法的求解結(jié)果對比如表1 所示.表中：B為目前文獻所給出的最優(yōu)解；Cmin為基于某種算法得到的最優(yōu)解；Tmin為某種算法得到最優(yōu)解所需的時間；G為某個算法求得的最優(yōu)解與B之間的差異，即G=(Cmin-B)/B.

表1 TSHA-II 算法與LNS-2E、Two-stage 算法求解2E-LRP 原算例的結(jié)果對比Tab.1 Comparison of the results obtained by LNS-2E, two-stage algorithm, and TSHA-II in solving original 2E-LRP instances

通過對比可以看出：TSHA-Ⅱ算法求解大規(guī)模的2E-LRP 算例時，求解質(zhì)量稍遜于Breunig 等[16]提出的LNS-2E 算法，但明顯優(yōu)于Dai 等[17]提出的Twostage 算法；求解效率雖稍遜于后者，但大大優(yōu)于前者.這表明，TSHA-Ⅱ算法能夠在求解質(zhì)量與求解效率之間取得較好的平衡，可以在較短的時間內(nèi)得到較高質(zhì)量的解，TSHA 的算法思想可行有效.

3.2 對考慮碳排放的2E-LRP 算例的求解

算法思想通過驗證后，TSHA 再被用于求解考慮碳排放的2E-LRP 大規(guī)模算例，計算結(jié)果如表2所示.表中：Eavg為20 次運算得到的平均碳排放量；Emax、Emin分別為20 次運算得到的最大和最小碳排放量；Tavg為20 次運算的平均耗時.

表2 TSHA 算法求解考慮碳排放的2E-LRP 算例的結(jié)果Tab.2 Results obtained by TSHA in solving 2E-LRP instances considering carbon emissions

從表2 可以看出，TSHA 在求解考慮碳排放的2E-LRP 大規(guī)模算例時，求解結(jié)果的最大值、最小值和平均值以及計算時長的最小值和平均值之間相差很小，表現(xiàn)穩(wěn)定.可見，TSHA 可以作為求解考慮碳排放的2E-LRP 大規(guī)模算例的一種有效算法.

4 結(jié) 論

1） 相較于經(jīng)典的以物流成本最小化為目標(biāo)的2E-LRP 模型，本文基于一種簡單實用的以排放因子為主要參數(shù)的碳排放計算方法，建立了以碳排放最小化為目標(biāo)的2E-LRP 優(yōu)化模型.

2） 對Prodhon 標(biāo)準(zhǔn)算例集中全部6 個最大規(guī)模算例的測試，表明TSHA-Ⅱ算法比現(xiàn)有算法更具兼顧求解結(jié)果和求解效率的優(yōu)勢，也證明了TSHA算法思想的可行性和有效性.

3） TSHA 算法在求解考慮碳排放的2E-LRP 大規(guī)模算例時表現(xiàn)得高效而且穩(wěn)定，能夠滿足實際應(yīng)用場景的計算要求.

4） 本文提出的模型和算法有助于物流企業(yè)通過優(yōu)化城市物流運作實現(xiàn)車輛的減排目標(biāo)，服務(wù)國家“雙碳”戰(zhàn)略.