張圣忠,陳婷婷,孫榮庭,白 雪,高 超
(長安大學(xué) 經(jīng)濟與管理學(xué)院,陜西 西安 710064)
隨著工業(yè)化進程的加速,危險品運輸量在逐年上升[1]。據(jù)交通部統(tǒng)計,2017年國內(nèi)危險品運量達到11億t,且在以每年10%的速度增長,其中道路運輸是危險品運輸?shù)闹饕侄危?5%以上。不同于普通物品,危險品在運輸過程中容易出現(xiàn)側(cè)翻、碰撞、泄漏、燃燒、爆炸、中毒等現(xiàn)象,且一旦發(fā)生事故就可能造成人身傷亡、財產(chǎn)損毀和環(huán)境污染等嚴重后果[2]。這不僅會造成了運輸損失,還會對沿途人民的生命和周邊環(huán)境造成嚴重威脅。因此,在保障安全的前提下對危險品運輸車輛路徑進行合理的優(yōu)化就顯得十分重要。
危險品運輸車輛路徑優(yōu)化問題(vehicle routing problem for hazardous, HVRP)一般是將危險品運輸特性和VRP問題結(jié)合起來,考慮容量限制、總運輸時間限制和時間窗口限制等約束,以總運輸成本最小化和總風險最小化為目標。KARA等[3]首次提出一種考慮風險和成本的雙目標優(yōu)化模型, 使用庫恩-塔克條件和互補松弛條件將其轉(zhuǎn)換為單目標模型,并對模型進行求解。ZOGRAFOS等[4]提出了危險品運輸VRP的雙目標模型,使用線性加權(quán)方法將所提出的雙目標模型轉(zhuǎn)換為單目標模型,并對其進行求解。MENG等[5]提出了在有限的運行時間段、服務(wù)時間和節(jié)點等待時間窗約束下的車輛路徑模型,并提出了求解的動態(tài)規(guī)劃方法。PRADHANANGA等[6]提出了一種具有時間窗約束的雙目標VRP,并提出了一種元啟發(fā)式算法,以提供一組帕累托較優(yōu)解。ERKUT等[7]將危險品網(wǎng)絡(luò)設(shè)計問題歸結(jié)為樹選擇問題,采用樹狀網(wǎng)對雙層規(guī)劃模型進行分析,并用啟發(fā)式算法求解模型。辛春林等[8]針對HVRP問題給出了詳細的綜述,并提出了求解單目標和多目標問題的智能算法。
危險貨物運輸風險的度量也是HVRP問題中的重要因素。在以往的研究中,學(xué)者們一般利用“傳統(tǒng)風險模型”[9]來衡量危險品運輸過程中的風險,該模型通常以隨機分布來表示事故發(fā)生概率,以距事故現(xiàn)場特定距離內(nèi)受影響的居民人數(shù)[10]來表示事故后果。后來也有學(xué)者考慮到了運輸車輛的載貨量對風險帶來的影響,如PRADHANANGA等[11]使用危險品運輸車輛的額定載貨量和人口暴露數(shù)來衡量危險品的運輸風險。不過大多數(shù)研究者仍采用人口暴露數(shù)來衡量運輸風險,目前僅有少量研究考慮到了載貨量動態(tài)變化對風險的影響,如袁文燕等[12]提出了隨載貨量線性變化的風險測度方法,其所建模型以最小運輸距離和最小風險為目標,且未考慮服務(wù)時間窗的約束。綜上所述,目前針對危險品運輸路徑規(guī)劃問題的研究還不夠充分,與實際危險品運輸情景仍存在較大差距。為了進一步貼近現(xiàn)實情況,筆者考慮了HVRP問題中載貨量動態(tài)變化對運輸風險的影響,同時考慮配送時間窗約束,以最小化危險品運輸?shù)能囕v數(shù)、運輸總距離和運輸風險為目標函數(shù),建立了考慮載貨量的危險品運輸車輛路徑問題模型,決策危險品運輸車輛數(shù)和車輛服務(wù)路徑。
在求解HVRP優(yōu)化問題方面,常用傳統(tǒng)的精確算法往往由于問題復(fù)雜而造成計算量巨大,為此有學(xué)者提出分布估計算法求解多目標問題[13-14],可以提高模型的魯棒性,該算法可以有效應(yīng)用于車輛路徑優(yōu)化問題的求解[15-16]。因此,筆者設(shè)計一種基于概率模型的多目標進化算法對提出的HVRP問題進行求解,先采用NSGA2算法的非支配排序方法來獲取進化群體,再通過概率模型進行隨機采樣計算每個可行解的目標值和對應(yīng)的路線。針對考慮載貨量的危險品運輸車輛路徑優(yōu)化問題,通過概率模型來描述實際載貨量和風險之間的線性關(guān)系,并在選擇機制中融入NSGA2算法的概念,基于群體的宏觀進化方式,使其可以利用HVRP問題的結(jié)構(gòu)信息來產(chǎn)生更好的個體,擴大多目標進化算法的搜索范圍,保證算法解的可行性以及求得完整的Pareto前沿。
基于以上考慮,在采用實際載貨量建立新風險度量模型的基礎(chǔ)上,筆者研究了一種以最小化運輸距離、車輛數(shù)和運輸風險作為優(yōu)化目標且考慮服務(wù)時間窗約束的多目標路徑優(yōu)化問題,建立數(shù)學(xué)模型,提出一種基于概率模型的多目標進化算法對問題進行求解。并通過對3種按道路網(wǎng)絡(luò)節(jié)點形成的不同分布進行算例仿真,以驗證模型與算法的有效性。
研究問題可以描述為:在危險品道路運輸網(wǎng)絡(luò)中有一個配送中心、多個客戶需求點和若干危險品運輸車輛。在不超載的情況下,配送中心根據(jù)每個客戶節(jié)點的需求量在配送時間要求內(nèi)進行分配,直至滿足每個客戶需求,一輛車可為多個客戶需求點進行配送并確定配送車輛的數(shù)目。從配送中心同時派發(fā)多個危險品運輸車輛對各自分配的需求點進行服務(wù),服務(wù)完成后車輛必須返回到配送中心。主要假設(shè)有:①配送中心擁有足夠的庫存滿足轄區(qū)內(nèi)各個需求點需求;②各需求節(jié)點的需求量和時間窗要求可提前獲知且需求不可拆分;③配送中心及需求點之間的行駛距離、事故概率、暴露人口數(shù)已知;④所有配送車輛額定載貨量相同。
已知危險品道路運輸網(wǎng)絡(luò)系統(tǒng)G=(N,A),其中N={0,1,2,…,n}代表危險品配送中心0和該配送中心所有客戶需求節(jié)點i的集合,A代表連接需求節(jié)點及配送中心的所有道路的集合。各節(jié)點之間的距離為dij,每個需求節(jié)點i的需求量及配送時間窗要求分別為qi和[ETi,LTi],k輛具有容量W的配送車輛從配送中心0出發(fā),按照配送路線依次完成相應(yīng)需求節(jié)點的配送,最終返回到配送中心。
在HVRP問題中,運輸車輛從倉庫出發(fā)依次訪問若干個需求點后返回倉庫,配送車輛上危險品載重量的不斷變化與車輛的訪問路徑密切相關(guān)。
傳統(tǒng)的風險度量模型以Rij=pijφij來度量邊(i,j)上的道路風險,其中pij,φij分別為車輛在路徑(i,j)上發(fā)生事故的概率及事故導(dǎo)致的后果(通常用人口暴露數(shù)來表示),不考慮載貨量與事故影響后果的關(guān)系,因此無法準確衡量路徑風險?;谝陨峡紤],筆者假定風險影響后果與危險品的載貨量存在線性關(guān)系,于是有:
(1)
車輛k的第m段配送路徑上的載貨量為:
(2)
則車輛k的第m段配送路徑上的風險為:
(3)
那么車輛k配送路徑上的風險為:
(4)
該問題的優(yōu)化的目標為:
(5)
(6)
(7)
該問題的可行解必須滿足的約束條件有:
(8)
(9)
(10)
(11)
(12)
(13)
Tjk=Tik+tijxijk,ETj≤Tjk≤LTj,
?i∈N,?j∈C,?k∈K
(14)
xijk∈{0,1},?i,j∈N,?k∈K
(15)
yik∈{0,1},?i∈C,?k∈K
(16)
式中:目標函數(shù)(5)表示最小化運輸總距離;目標函數(shù)(6)表示最小化車輛數(shù);目標函數(shù)(7)表示最小化風險;約束條件(8)表示每個需求點只有一輛車為其服務(wù);約束條件(9)表示每輛車出發(fā)時最多只能訪問一個客戶;約束條件(10)表示每輛車到達時最多只能配送一個客戶;約束條件(11)表示車輛到達某客戶節(jié)點完成配送后必須從此節(jié)點出發(fā);約束條件(12)表示確保每個客戶被訪問一次,并返回到配送中心;約束條件(13)表示車輛載貨量約束;約束條件(14)表示配送時間窗約束,要求車輛必須在時間窗內(nèi)開始配送;約束條件(15)和(16)表示決策變量的取值范圍。
相關(guān)研究解決HVRP問題通常用遺傳算法進行求解[17],KALYANMOY等[18]提出基于擁擠距離的分布性方法的快速非支配排序遺傳算法NSGA2,但是NSGA2存在交叉、變異等遺傳操作對問題結(jié)構(gòu)破壞的設(shè)計缺陷,且無法考慮到變量之間的相關(guān)依賴關(guān)系,從而通過引入概率模型捕獲問題中的結(jié)構(gòu)關(guān)系,降低破壞程度,并在求解效率和解的質(zhì)量上有較大的優(yōu)勢?;诖?,筆者提出一種基于概率模型的多目標進化算法對HVRP問題進行求解,該算法通過NSGA2的非支配排序方法進行種群選擇和構(gòu)造非支配解集,將參與進化的個體進行構(gòu)建概率模型,并選取輪盤賭策略對概率矩陣進行更新、歸一化、隨機采樣等操作來求解多目標優(yōu)化模型。
對模型采用自然數(shù)定長編碼的方式進行求解,解碼得到的路徑由危險品車輛運輸?shù)呐渌椭行暮托枨簏c編號組成。例如1-12-8-10-21-1-20-11-2-1-13-4-25-1-3-16-7-1-24-5-14-1-22-23-1-6-9-18-1-15-17-19-1,該路徑的數(shù)字1為配送中心,其他數(shù)字表示車輛對需求點的配送順序,其中1-12-8-10-21-1表示同一輛車對不同需求點進行配送并返回到配送中心。通過編碼產(chǎn)生初始解時,首先采用NSGA2算法獲取種群的優(yōu)勢個體,既可得到運輸方案的可行解,又提高了算法的魯棒性;其次在其個體領(lǐng)域內(nèi)獲取部分解;最后與產(chǎn)生的其他解一起隨機生成初始解。
通過非支配排序的方法構(gòu)造問題的非支配解集,根據(jù)解之間的支配關(guān)系構(gòu)造,從配送中心出發(fā),為了保持種群在解空間的分布性和多樣性,將個體的擁擠距離作為種群中個體之間的比較準則,使得準Pareto域中的非支配解能均勻擴展到整個域,對可行解進行快速非支配排序并計算擁擠距離來進行比較操作,再根據(jù)相鄰解的擁擠距離大小選擇精英個體直接進入下一代種群參與進化。個體聚集密度的計算式為[19]:
(17)
式中:P[i]distance表示個體i的聚集密度距離;P[i]·m表示個體i在子目標m上的函數(shù)值;fk表示有k個子目標;r表示種群子目標的個數(shù)。
根據(jù)參與進化的個體進行構(gòu)建概率模型,在危險品運輸車輛的裝載不為空或者當前時刻不在配送時間窗范圍內(nèi)的情況下,采用輪盤賭策略對選取的精英個體進行采樣,選取第m段路徑配送首節(jié)點的車輛,對配送車輛k+1進行迭加,再計算節(jié)點ci為某條路徑首節(jié)點的概率qr,從中選取適應(yīng)度較高的個體,并統(tǒng)計相同車輛配送的節(jié)點ci和cj在同一路徑上相鄰的概率pr,計算方法如式(18)和式(19)所示。
(18)
(19)
算法具體流程:①參數(shù)初始化。構(gòu)造初始解并初始化非支配解集,令初始迭代次數(shù)gen=1,最大遺傳代數(shù)maxgen=100。②對目標進行解的構(gòu)造,構(gòu)造N個解;根據(jù)所有客戶的需求量,通過概率模型采用輪盤賭策略隨機選擇下一個被訪問的需求點j,直至配送車輛載貨量為0并記錄該配送路徑,并返回倉庫;進行下一路徑的選取。③進行非支配解的篩選,計算最優(yōu)的非支配解集第m段路徑的目標函數(shù)值并更新進化,在多次求得的非支配解集中再次求非支配解集,以達到較優(yōu)解集。④重復(fù)步驟②,直至gen=maxgen。基于概率模型的多目標進化算法流程如圖1所示。
圖1 基于概率模型的多目標進化算法流程圖
筆者通過數(shù)值仿真實驗來檢驗?zāi)繕四P偷男Ч?。實驗設(shè)置3組相同的數(shù)據(jù)規(guī)模和按需求點的地理位置或時間窗的3組不同的分布數(shù)據(jù),均考慮24個顧客需求點和1個倉庫,SOLOMON[20]算例為VRP領(lǐng)域公認的經(jīng)典算例,因此筆者也選取該算例進行研究。道路系統(tǒng)中的事故概率和人口密度信息隨機生成。算法的參數(shù)設(shè)置個體數(shù)為200個,最大迭代次數(shù)為100次,優(yōu)勢種群數(shù)為40個,采用Matlab編程實現(xiàn)該算法。
采用SOLOMON算例按網(wǎng)絡(luò)節(jié)點的地理位置或時間窗形成聚類分布、均勻分布和介于兩者之間3組不同分布的算例,分別為R101、RC101和C101數(shù)據(jù)。R101的車輛集中分布在運輸距離680~780 km之間,而C101的車輛數(shù)在各個運輸距離段呈現(xiàn)平均分布,且車輛數(shù)較少,RC101的車輛數(shù)在各個運輸距離段也呈現(xiàn)平均分布,但在每一個時間段內(nèi)的車輛數(shù)比C101集中。選取3組不同分布的顧客坐標、需求、時間窗和服務(wù)時間的信息作為基本算例,按車輛數(shù)進行分類后各自分布的Pareto解如圖2所示。
圖2 按網(wǎng)絡(luò)節(jié)點位置分布的Pareto解
不同分布生成的Pareto解的信息如表1所示。結(jié)合圖2和表1可知:①當顧客點相同時,考慮載貨量的聚類分布生成30個Pareto解,均勻分布生成11個Pareto解,在聚類與均勻之間的分布生成23個解,不考慮載貨量的Pareto解只有2個,這表明基于概率模型的多目標進化算法在考慮載貨量的情況下可以求得更多的Pareto解,對運輸過程中的動態(tài)變化較為敏感。②按網(wǎng)絡(luò)節(jié)點的地理位置形成的均勻分布生成的解空間相較而言更為分散,這表明基于概率模型的多目標進化算法可以對不同分布的個體聚類距離進行計算,根據(jù)顧客地理位置分布的分散程度,其對解的搜索范圍加大,導(dǎo)致運行時間增加。
表1 不同分布生成的Pareto解的信息
以25個網(wǎng)絡(luò)節(jié)點位置形成的聚類分布為例,對比分析兩種不同風險度量方式(即考慮載貨量和不考慮載貨量)下3個目標的Pareto解集風險值,如表2所示。
表2 不同風險度量方式下的Pareto解
對考慮載貨量的風險度量方式進行細節(jié)描述,通過R101算例進行仿真計算獲得一組Pareto解,記錄了對應(yīng)的3個目標值和路線,如表3所示。
表3 R101算例的Pareto解
從表3可以看出,隨著方案中車輛數(shù)目的增加,運輸風險不斷降低,但是運輸距離呈上升趨勢。這是因為在整體需求不變的情況下,車輛數(shù)增加會導(dǎo)致每輛車的配送節(jié)點減少,分配到每輛車上的危險品數(shù)量減少,完成配送任務(wù)至空車的速度加快,因此總風險下降;但是車輛數(shù)增加會降低單車配送效率,因此總距離上升。R101數(shù)據(jù)中顧客需求點位置坐標所形成70×80區(qū)域的運輸軌跡路線如圖3所示,如果決策者偏向于成本最低,可以選擇方案1的軌跡路線;如果決策者偏向于風險最低,可以選擇方案30的軌跡路線。
圖3 非支配解的軌跡路線
(1)考慮載貨量和服務(wù)時間窗的危險品車輛運輸路徑問題,建立以最小化車輛數(shù)、運輸總距離及運輸風險為目標的HVRP多目標優(yōu)化模型,并設(shè)計了一種基于概率模型的多目標進化算法來求解問題。模型引入了配送時間窗的約束,并將車輛數(shù)作為路徑優(yōu)化中考慮運輸成本的目標之一,從而決策者可以按車輛數(shù)目從Pareto解集中選擇滿意的危險品運輸路徑,且能夠為決策者提供更多的路徑選擇。
(2)通過數(shù)值分析發(fā)現(xiàn),概率模型對變量之間的關(guān)系進行估計的進化算法能夠獲得多目標問題的有效Pareto解。通過對顧客的地理位置或時間窗不同分布的算例求解,驗證了模型的有效性,計算結(jié)果和運行時間也表明了基于概率模型的多目標進化算法對解空間有著較強的搜索能力。
(3)考慮載貨量時,3種分布情況下得到的平均Pareto解集個數(shù)為21個;但不考慮載貨量時,3種分布的Pareto解僅有2個。通過對不同風險度量模型進行計算發(fā)現(xiàn),在不考慮載貨量的風險度量方式下,Pareto解集的風險值分布較單一,且同比平均風險值,其風險值波動范圍不超過2%;在考慮載貨量的風險度量方式下,Pareto解集的風險值分布更加多樣,其風險值在8.5%~12.9%范圍內(nèi)波動。因此,考慮載貨量的風險度量方式可以提供風險多樣化的決策方案,以滿足不同決策者的風險偏好。
(4)筆者在設(shè)計風險模型時,考慮到對模型復(fù)雜程度的影響,并沒有考慮危險品屬性、路況、天氣等多種因素對于危險品運輸風險的相關(guān)性,可以在未來進一步深入研究。