李家碧,韓曙光
(1.浙江理工大學 經濟管理學院,浙江 杭州 310018;2.浙江理工大學 理學院,浙江 杭州 310018)
隨著“構建智慧城市”理念的不斷深入,“無人化”和“去快遞員化”的配送模式脫穎而出.物流業的發展加重交通擁堵,道路交通狀況影響配送車輛行駛速度,繼而影響物流運輸效率,而配送的末端服務質量,會直接關聯客戶體驗.《浙江省基本公共服務體系建設2020 年工作要點》中明確,要推進智能收投終端和末端公共服務平臺建設的方針政策[1],給“無接觸配送”帶來發展機遇[2].
車輛路徑問題(vehicle routing problem,VRP)是由Dantzig 等[3]于1959 年提出的.隨著物流的不斷發展,VRP 得到大量關注,眾多學者在VRP 的基礎上進行大量的延申研究,引出VRP 的不同變體[4-5].針對帶時間窗的車輛路徑問題(vehicle routing problem with time windows,VRPTW),Solomon[6]于1987 年將時間窗的約束加入到VRP 中.Kang 等[7]將時間窗和懲罰函數相結合,對VRP 進行了更全面的研究.homberger 等[8]基于遺傳算法提出分別在變異和交叉操作上改進的2 種策略來求解VRPTW.Fan 等[9]針對具有模糊需求的時效綠色車輛路由問題,提出自適應鄰域搜索時間的策略和劣質解接受機制.針對電動車輛配送問題(electric vehicle routing problem,EVRP),Conrad 等[10]在VRP 上提出1 個電動車能夠在特定客戶點進行充電的EVRP 模型.更多學者從不同方面對EVRP 分析結合更多現實因素展開進一步研究[11].郭放等[12]研究了貨物需求的差異及電動車換電的策略.肖建華等[13]以道路限行為背景,研究電動車和燃油車的混合配送路徑優化問題.Desaulniers 等[14]考慮帶時間窗的電動車輛路徑問題的4 種變體.Ye 等[15]介紹電動車輛路徑問題模型的4 種類型,指出理論方法的發展趨勢.在考慮客戶等級問題中,Alexiou 等[16]考慮配送成本與不利情況的因素相結合,不同客戶等級設置相應權重,利用支樹結構搜索算法確定最優路線.Calvet 等[17]提出考慮市場細分策略,將統計學習技術與元啟發式框架相結合的混合方法.馬向國等[18]研究在配送需求隨機,當客戶不同等重要性下,通過自適應遺傳算法求解總成本最小化的數學模型.楊培穎等[19]研究減少成本的同時,并提高客戶滿意度的多目標VRPTW模型.王力鋒等[20]在考慮客戶等級劃分的基礎上,設計新的多目標冷鏈物流配送優化方法.
針對時變路況方面,Malandraki 等[21]研究時變速度的VRP,提出行駛時間分段函數.劉長石等[22]設計考慮交通擁堵指數的交通擁堵規避方法.葛顯龍等[23]考慮車輛配送過程中在2 個客戶點間路段等待的情況,以極小化碳排放量和車輛行駛時間為雙目標,建立優化模型.Keskin 等[24]考慮部分充電下的電動車輛路徑優化問題,并利用變鄰域搜索算法求解.Roberti 等[25]提出可以在短計算時間內解決20 個客戶實例的混合整數線性公式,以及基于通用變量鄰域搜索和動態規劃的三階段啟發式算法.張曉楠等[26]結合時變旅行時間和基于關鍵點構建的道路網路,建立以總旅行時間最小為目標的優化模型.
目前在大多數關于物流路徑優化的研究中,對于帶有時變路況的車輛路徑問題,尚未出現同時考慮到電動車輛和客戶等級的研究成果.不同等級的客戶能夠給物流公司帶來的收益不同,物流公司開展階梯式服務更有利于長久發展.在考慮客戶等級和時變路況2 個要素的基礎上,提出電動車輛路徑優化問題,建立以配送總成本極小化為目標函數的路徑優化模型,通過混合遺傳-模擬退火算法(genetic simulated annealing algorithm,GA-SA)對Solomon 數據集改造后的9 組數據進行數值實驗,并與遺傳算法(GA)、模擬退火(SA)進行對比分析,驗證GA-SA 的有效性等.
某物流企業有一個配送中心,負責給若干客戶點配送貨物,配送區域內有若干充電站點,并有一定數量的裝載容量相同和電池容量相同的無人電動車.已知每個客戶的等級分類、送貨需求、各節點之間的距離、道路狀況以及行駛過程中不同時間段對應的車速.從倉庫點出發,電動無人車對每個客戶點進行訪問,并回到倉庫點.無人車的實時載重量要滿足裝載容量的約束,剩余電量要足以支撐其完成服務,目標為極小化所有無人車行駛總成本,如圖1 所示,D為倉庫或配送中心.

圖1 考慮客戶等級和時變路況的城市無人物流配送示意圖Fig.1 Schematic diagram of urban unmanned logistics distribution considering customer grade and time-varying road conditions
1.2.1 客戶等級的劃分及等級評價指標體系建立
參考已有的客戶細分理論相關研究,構建客戶價值評價指標體系[27]如圖2 所示.

圖2 配送客戶的價值評價指標體系Fig.2 Value evaluation index system for delivery customer
1.2.2 客戶等級評價模型 將客戶等級的評分集假定為3 個級別V={V1,V2,V3} ,分別為VIP 客戶、主要客戶和普通客戶,并運用云模型的期望、熵和超熵3 個數字特性相應進行表示.運用云模型建立存在雙邊約束 [Cmin,Cmax] 的評分集,Cmin、Cmax分別為該評分集可取值的最小邊界值和最大邊界值,評價云模型的云參數計算式[28]為
式中:Ex為期望,En為熵,He為超熵.
根據公式可以計算如表1 所示.參考二八理論,企業利潤的80% 是由企業20% 的客戶所創造,企業創造利潤價值的客戶數量僅為企業客戶群體中的小部分.若是客戶評價值為區間邊界,則該客戶應取較低等級.假設客戶等級的指標集為U, 每個指標對應的權重為vi,vi≥0且v1+v2+···+vm=1,相應的評價矩陣為W.通過每位專家對權重V進行打分Xi,運用統計的方法,利用逆向云發生器計算出統計樣本對應的云參數Ex、En和He,得到權重矩陣V′.對m個客戶等級影響指標進行打分,利用逆向云發生器計算得到相應云的數字特征Ex、En和He,得到綜合評價矩陣W′.利用模糊合成算子計算綜合評價結果,得到綜合評價云模型,具體的云運算原理參考文獻[29].通過Matlab 將評價集云模型和等級綜合評價云模型分別仿真顯示出來,評價集云模型與綜合評價云模型距離最近的,也就是最終的客戶等級.

表1 客戶等級評價云模型相關數字特征值Tab.1 Customer rating cloud model related digital eigenvalue
模型基本假設如下:1)僅對客戶提供送貨服務,每個客戶服務時長相同,單個客戶送貨量不超過車輛最大裝載量;2)車輛均為同類型、裝載量及電池容量相同;3)只有1 個配送中心,并且所有車輛由此出發,為所有規定的交付目的地提供服務后返回起點;4)車輛出發時為滿電,途中可以多次去充電站充電,充電后電池均為滿電;5)客戶點位置信息和服務時間窗均為已知;6)充電站位置和數量已知,不考慮排隊充電;7)客戶點之間的車輛行駛時間和能源消耗與交付距離成正比;8)車輛在服務客戶的過程中無電量消耗.
1.4.1 目標函數 優化目標為極小化物流配送總成本:
(1)車輛的固定成本
式中:K為車輛集合;N為配送中心、客戶點集合與充電站集合的并集;c1為固定成本;=0~1,如果無人車k從客戶點i通過到客戶點j時為1.0,否則為0.
(2)車輛的運輸成本,該成本和行駛距離成正比.
式中:c2為單位行駛成本,di,j為客戶點i、j之間的距離.
(3)車輛的充電成本,由車輛到達充電站充滿電.
式中:c3為單位充電成本; α 為無人車的充電系數;R為電池最大容量;為到客戶點i時無人車k的電量;=0~1,如果無人車k通過充電站i充電為1.0,否則為0.
(4)懲罰成本.在物流配送的過程中,如果未在客戶的期望服務的時間窗內配送,而是在客戶可接受的服務時間窗內,將受到不同程度的懲罰;在客戶可接受服務的時間外才能送到,客戶可以選擇終止服務,懲罰成本將為無窮大.懲罰成本計算式為
式中:Ci為配送客戶i時產生的懲罰成本, θ1為早到的懲罰系數, θ2為晚到的懲罰系數,ai,k為無人車k到客戶點i的時間, E Ti、LTi為客戶i期望服務的時間窗, EETi、LLTi為客戶i可以接受服務的時間窗.故懲罰成本為
1.4.2 約束條件 客戶滿意度計算,與客戶時間敏感度系數相關.在給定 θ 為客戶整體滿意度的下限,使得配送路徑方案能夠達到一定的客戶滿意度,具體計算式為
式中:Ui為配送客戶i的滿意度, σ 為客戶時間敏感系數.
(5)路徑約束:
約束(10)為車輛在任一時間段h內在i和j這2 個客戶點之間行駛的距離不超過i、j兩點之間的距離.式中:為在時間段h內無人車k在客戶點(i,j)之間行駛的距離;H為所有時間段的集合;=0~1,如果無人車k在時間段h從客戶點i通過到客戶點j時為1.0,否則為0.
約束(11)為無人車在任何時間段在i和j之間行駛的總距離等于i,j兩點之間的距離.
約束(12)為消除路線中可能存在的所有子回路.
約束(13)為每個客戶點只被無人車服務一次,且只能由一輛無人車服務.
約束(14)為進入和離開的車流量相等.
(6)行程時間計算:
式中: [bh,eh] 為時間段h的開始時間和結束時間.約束(16)為無人車k在時間段h內的總行駛時間不大于時間段長度.
式中:li,k為無人車k從客戶點i離開的時間.約束(17)為時間段h的結束時刻減去無人車在該時間段內在客戶點 (i,j) 之間的行駛時間,晚于無人車離開客戶點i的時間.
約束(18)為時間段h的起始時間加上該時間段內在客戶點 (i,j) 之間的行駛時間小于或等于無人車到達客戶點j的時間.
約束(19)為從客戶點i出發后到達客戶點j的時間轉換關系.
式中:si為客戶i的服務時間,C為客戶點集合.約束(20)確保無人車在客戶點的出發時間必須大于或等于到達時間和服務時間之和.
(7)車輛載重約束:
式中:Q為車輛最大載重.約束(22)為無人車k不能攜帶超過其最大容量的貨物.
(8)電量約束:
式中:D為配送中心,aD,k為車輛從配送中心出發的時間,tmax為車輛返回配送中心的最晚時間.約束(23)確保在配送中心工作時間內提供所有服務.
約束(24)規定車輛離開充電站的時間小于或者等于無人車到達時間和從開始充電至完全充電的時間之和.
約束(25)為無人車從客戶i點出發到客戶點j的電量狀態提供邊界條件.
約束(27) 規定,車輛離開配送中心開始送貨時,電池已充滿電.
約束(28)為和 之間的限制關系.
約束(29)明確了決策變量的具體定義域.
遺傳算法和模擬退火算法都是基于概率分布的啟發式搜索優化算法但是擁有不同特性.SA由Metropolis 等[30]于1953 年提出,由于固體退火原理, 優化過程中以一定的概率接受較差解,有效防止陷入局部極值并逐漸趨于全局最優.GA是模擬生物進化過程的方法,使用達爾文的自然選擇理論和遺傳機制來尋找最優解,是密歇根大學的Holland[31]在1975 年開發的全局搜索優化算法.混合遺傳-模擬退火算法將2 種方法相結合,可以一定程度上彌補2 種算法的缺陷,在全局或者局部情況下都提高算法的搜索能力和效率.
GA-SA 的主要實現過程為1)對種群進行模擬退火操作;2)產生的新種群完成遺傳算法的交叉操作和變異操作,產生下一代種群;3)循環反復以上步驟直到滿足迭代次數,終止計算輸出最優解.GA-SA 的具體求解步驟如下.
1)初始化模型和算法參數.確定種群大小N,交叉概率Pc,變異概率Pm,初始溫度T0以及算法完成前的最大進化代數等.
2)初始化種群.初始種群由隨機產生的N個可行解Xi組成,進行編碼并計算種群適應度,在每一代中獲得最優個體和最高適應度.適應度函數為
式中:F(x) 為適應度函數,f(x) 為目標函數.
3) 溫度逐漸下降,進行隨機擾動.隨機置換2 個不同城市的坐標,計算新狀態適應度值F(x′).如果F(x′)>F(x) ,則接受新狀態;反之,則以一定概率選擇接受新狀態,如果接受則替換舊狀態,否則舊狀態保持不變.
4)交叉操作.設定概率p, 在染色體中任意選擇基因點,產生隨機數與概率p相比較.若小于p,則在同一個體中選取另一個基因點G2,將2 個基因點之前的部分倒置;反之,則在種群中選擇新個體,找到G1在該個體中前一個位置的點G3,再回到原來的個體中找到G3,將G1和G3之間的部分倒置.
5)變異操作.每個基因都有一定的變異概率與任意基因進行互換.
6)精英保留策略.在父代和子代中出現的精英個體,不進行交叉配對直接復制到下一代,保證最優個體不會被破壞.
7)更新當前最優解,判斷是否達到最大進化代數.如果滿足,則終止計算,解碼后輸出最優解,如果不滿足則返回步驟3),直到滿足算法終止條件.GA-SA 優化流程圖如圖3 所示.

圖3 遺傳-模擬退火算法的優化流程圖Fig.3 Optimization flowchart of hybrid genetic simulated annealing algorithm
GA-SA 存在循環嵌套,代入時間復雜度計算方法,循環的時間復雜度為循環部分的復雜度中循環運行次數的乘積,分析得到內循環的時間復雜度為O(n).對于外層的循環,為內部時間復雜度為O(n)的語句,再進行n次循環,則這段算法的時間復雜度為O(n2).當外循環的循環次數更改為m,時間復雜度為O(m?n).
從Solomon[20]標準算例(C1、R1、RC1)中分別選取25、50 和100 個客戶點的數據,模擬形成9 組測試算例.C 型數據是通過聚類操作之后生成的,客戶的分布區域較為集中,R 型數據為隨機生成的客戶分布,RC 型數據客戶的分布則是混合R 型數據和C 型數據分布的.應用Matlab R2021a,在Windows 10 操作系統的環境,Intel(R) Core(TM)i7-10710U CPU @ 1.10 GHz、內存16 GB 的筆記本電腦上進行仿真求解.
設定模型的相關參數如下:車輛載重為2 kg/輛,電池容為40 kwh/輛,耗電系數為4 km/kwh,充電系數為0.60 kwh/min,單位充電成本為0.75 元/kwh,車輛發車成本為20 元/輛,單位行駛成本為2 元/km,單位等待成本為0.40 元/min,早到和晚到的懲罰成本為2 元/min、3 元/min.
采取車輛在兩客戶點之間的配送路上不停靠策略, 即車輛從客戶點i行駛到客戶點j的過程中,始終處于行駛狀態,不在2 個客戶點中間路段停留.在時變的道路條件下,車輛速度的波動是已知的并且隨時間變化[32].參考杭州市交通擁堵大數據在線監測平臺提供的數據[33],模擬不同時間段車輛行駛速度如圖4 所示.圖中V為車輛在不同時間段的行駛速度,t為相對應的時間段.

圖4 車輛在不同時間段的行駛速度Fig.4 Vehicle speed in different time periods
設置種群規模N=100 最大進化代數G=3 000 ,交叉概率Pc=0.60 ,變異概率Pm=0.01 ,馬爾科夫鏈長度L=1.00 ,溫度衰減系k=0.99 ,用GA-SA、GA 和SA 各求解了9 組測試算例10 次,并使用當前的最優值和平均值進行比較和分析,求解結果如表2 所示.在B1、B2和B3列中,3 個算法求解10 次后的當前最優解;在Bˉ1、Bˉ2和Bˉ3列中,10 次運行結果的平均值; GAP1為當前最優解差值比較的百分比 GAP1=(B1-B2)/B2×100% , G AP2為平均值差值比較的百分比GAP3和 G AP4的運算參考 G AP1和 G AP2.從表2 的求解結果可以看出,無論是25 個、50 個或者100 個客戶點的規模,還是隨機分布、聚類分布和混合分布的客戶點類型,得到的當前最優解和平均解,GA-SA 算法均優于GA 和SA.與GA 和SA 對比,混合后的GA-SA 在不同算例中配送過程中產生的總成本均減少在-1.65%~-11.66%和-5.46%~-42.81%.由此可以得到混合后的GA-SA 是可行的、有效且優化效果顯著.
在U1、U2和U3列中,GA-SA、GA 和SA 分別運行10 次后的當前最優客戶滿意度.在和列中,10 次運行結果的平均客戶滿意度.其余運算仍參考 GAP1和 G AP2,求解結果如表3 所示.可以看出,GA-SA 優化后的車輛路徑優化方案達到的客戶滿意度明顯優于GA 和SA 求出的配送方案,在C 型客戶分布類型中較為明顯.結合表2、3 的結果可以發現,GA-SA 在優化的過程中和降低配送總成本的同時,也提高了配送客戶的整體滿意度.車輛數在配送過程中對每項成本有著不同程度的影響,使用上述的算例結果得到車輛數與成本的關系圖如圖5 所示.圖中Co為車輛在行駛過程中產生的不同成本,k為相對應的車輛數.從圖5 可以看出,當車輛數為6 時,配送方案的總成本最小.當車輛數減少時,雖然發車成本隨之減少,但存在配送無法滿足客戶的時間窗要求的情況,從而產生更多的懲罰成本.隨著車輛數的增加,車輛的固定成本、充電成本以及等待成本都會相應增加,因此配送總成本會出現明顯增高.合適的車輛數可以在配送路徑優化過程中有效地降低總配送成本.不同的道路狀況會改變車輛的行駛速度,進而影響到達客戶點的配送時間,規劃出不同的配送方案.通過將道路狀況設置為暢通(20 km/h)、時變和擁堵(50 km/h)3 種情況進行運行,得到的相應結果如圖6 所示.圖6 中U為配送完成的整體客戶滿意度.從圖6 可以看出,道路狀況由暢通逐漸變為擁堵時,配送方案的總成本在逐步增加,客戶滿意度卻在相應地減小.在車輛數相同的情況下,道路擁堵行駛時間和成本增加,且無法在客戶期望的時間窗內完成服務的情況增多,導致懲罰成本相應增加,最終造成配送總成本增加,客戶滿意度下降.

表3 不同方法考慮客戶滿意度的優化結果Tab.3 Optimization results of different methods considering customer satisfaction

圖5 車輛數對成本優化結果的影響Fig.5 Influence of vehicle number on cost optimization results

圖6 不同路況對成本與滿意度的優化影響Fig.6 Optimising impact of different road conditions on cost and satisfaction
客戶滿意度的閾值限制在配送方案中客戶滿意度的可接受下限.當設置不同的滿意度閾值時,即代表接受更低或者更高的整體客戶滿意度配送方案,配送的目標也會收到滿意度閾值的限制而產生不同結果.設置滿意度閾值為40%~100%時2 種算法得到的成本結果對比圖,如圖7 所示.可以發現,GA-SA 與GA 相比,在相同客戶滿意度的條件下求得的總成本較低,表明GA-SA 能夠有限求解問題且優于GA.當滿意度閾值低于80%時,總成本基本穩定,是因為滿意度閾值為更易達到的低值,在路徑規劃中對于總成本的波動影響不大.當滿意度閾值由80%增長至90%的過程中,GA-SA 和GA 得到的總成本分別增長了8.65%和7.29%;當滿意度閾值超過90%時,增長率分別為16.99%和11.39%,總成本增長較迅猛;當滿意度閾值增長為100%時,總成本最高.由于100%的客戶滿意度難以達到,可以看出80%~90%的客戶滿意度在制定配送方案時是更合適的選擇.

圖7 客戶滿意度對成本優化結果的影響Fig.7 Influence of customer satisfaction on cost optimization results
圖8 為劃分與不劃分客戶等級的配送路徑變化.通過圖8 可以看出,在劃分客戶等級進行配送后,全體客戶的整體滿意度由88.9%降至80.2%.由于高級客戶的時間敏感程度更高,存在優先配送的情況,高級客戶群體的整體滿意度增加.車輛的行駛成本、充電成本和懲罰成本都相應的增大,說明在路徑規劃中存在滿足高等級客戶的期望時間窗要求,盡管當前路況較為擁堵,仍優先配送高等級客戶而放棄近距離的其他等級客戶的情況.當配送資源不能滿足所有客戶時,路線規劃時可以考慮客戶分類,盡可能滿足高層次客戶的服務要求,并根據客戶的重要性合理分配配送資源.從全局的角度來說,能夠為企業制定更持久穩健的發展戰略,增強顧客黏性,帶來更穩定的收入來源.

圖8 電動車輛的路徑結果Fig.8 Path results of electric vehicles
在智慧物流背景下, 對客戶不同等級、時變速度和電動車輛等因素進行研究,針對道路狀況建立車輛行駛速度的分段函數,通過建立云等級評價模型對客戶進行分級,建立以物流配送總成本極小化為目標的模型,采用GA-SA 和GA、SA 對比求解.論文從標準的Solomon 算例中選取9 組不同類型和規模的時變電動車輛路徑規劃測試算例,測試了所建模型和算法的有效性,數值結果表明如下.1)在時變的道路行駛速度下,合理的出發時間,分等級地選擇服務客戶的順序,能明顯減少車輛的配送總成本,提高客戶的整體滿意度.2)根據GA-SA 與GA、SA 的對比實驗結果表明,GA-SA 具有更好的優化效果.考慮諸如車輛碳排放等環境因素,在經濟效益基礎上如何進一步優化和平衡,設計更優的快速算法,是后續值得研究的方向.