楊華龍,陸 婷,辛禹辰
1.大連海事大學 交通運輸工程學院,遼寧 大連116026
2.大連海事大學 物流研究院,遼寧 大連116026
庫存路徑問題(Inventory Routing Problem,IRP)是VMI(Vendor Managed Inventory)模式下,供貨商整合庫存管理、車輛路徑?jīng)Q策的綜合性優(yōu)化問題。供貨商為滿足位置較遠地區(qū)的多零售商(客戶)需求,節(jié)省配送成本,需要在靠近這些客戶位置租用配送中心,以便為客戶提供配送服務,這就出現(xiàn)了二級IRP[1]。此時,供貨商需要考慮三方面成本:一是配送中心庫存成本;二是客戶處的庫存成本;三是配送成本。由于配送和庫存之間存在效益背反規(guī)律,且在不同銷售區(qū)和不同時段,客戶需求常常會出現(xiàn)隨機波動等情況[2],因此隨機需求下的多客戶二級IRP 優(yōu)化已成為VMI 模式下供貨商經(jīng)營決策的核心與瓶頸問題。
國內外學者針對IRP展開了許多有益的研究。Soysal等[3]在隨機需求的情況下,針對易腐品,建立了多周期的IRP模型,研究將單一產品配送給多個顧客的輸出型配送網(wǎng)絡。Soysal等[4]進一步研究多個供貨商服務多個客戶的物流網(wǎng)絡結構,建立橫向整合策略下的易腐品IRP模型,研究考慮橫向整合戰(zhàn)略對整條供應鏈的影響。在此基礎上,Micheli等[5]研究擁有異質車隊的多對多網(wǎng)絡配送系統(tǒng),重點強調了碳排放策略,如碳限額、碳稅、限額交易等對IRP決策的影響。此外,Cheng等[6]研究將單一產品配送給多個顧客的輸出型配送網(wǎng)絡,建立了同時考慮環(huán)境問題和異質車隊的IRP模型,并證明了異質車隊配送可以節(jié)約成本。以上研究都是針對供貨商直接交付給顧客的情形,并沒有考慮二級IRP的情形。
二級IRP 是IRP 的擴展問題,考慮一個具有兩階段間接交付和路線決策的兩級系統(tǒng)。Zhao 等[7]提出了由一個供貨商、一個中央倉庫和一組客戶組成的三級分銷系統(tǒng),提出了用于運輸?shù)墓潭ǚ謪^(qū)策略(Fixed Partition Policy,F(xiàn)PP)和用于庫存補充的二次冪(Power-of-Two,POT)庫存策略。Rahim 等[8]研究同樣的物流系統(tǒng),提出兩階段啟發(fā)式解決方案,第一階段最大程度減少系統(tǒng)總庫存成本,第二階段對車輛路徑進行優(yōu)化。Guimar?es等[9]研究了二級多配送中心庫存路徑問題,強調不同的庫存策略對供貨商產品輸入和最終產品輸出的影響,并運用分支定界法和自適應大規(guī)模鄰域搜索算法對問題進行求解。Rohmer等[10]提出了關于易腐品的二級庫存路徑問題的模型,主要考慮易腐品的產品生命周期問題,提出了自適應大規(guī)模鄰域搜索算法進行求解。
在上述已有研究中,研究問題從單個供貨商到多個供貨商、從直接交付和聯(lián)合交付等方面不斷深化。但是以上關于二級IRP的研究都是假定同質車隊,當需求波動較大時,運用同質車隊并不能與配送需求很好地匹配,此時采取異質車隊配送可以節(jié)約成本[11]。因此,研究異質車隊的二級IRP更加貼合現(xiàn)實。
此外,IRP 是典型的NP 難問題[12],因此二級IRP 也是NP難問題。近年來,粒子群算法以其概念簡單清晰,參數(shù)較少,魯棒性好等優(yōu)點被廣泛應用于VRP 及相關領域[13]。鑒于此,本文結合客戶需求波動因素,創(chuàng)新性地開展以下研究:一是提出干線運輸和配送相結合的異質車隊二級IRP優(yōu)化問題,構建異質車隊二級IRP模型;二是設計改進的粒子群算法對模型進行求解,以期為解決VMI下的二級IRP提供決策參考。
考慮某供貨商通過一個配送中心為區(qū)域內的多個客戶提供VMI服務。在這個供應鏈系統(tǒng)中,一方面,供貨商為配送中心庫存補貨,由于供貨商至配送中心距離遠,供貨數(shù)量大,集貨過程(即干線運輸)中通常采用運量大的單一類型運輸工具,因此只要在銷售周期內零售商的總需求給定,則干線運輸總成本便不會發(fā)生變化。另一方面,供貨商通過其配送中心給零售商配送產品,配送車輛完成配送作業(yè)后,便返回配送中心。供貨商、配送中心與零售商組成的供應鏈網(wǎng)絡庫存與配送結構如圖1所示。
圖1 二級IRP供應鏈網(wǎng)絡
由圖1 可見,二級IRP 決策包含虛線圈內的干線運輸決策、配送中心庫存決策、各零售商庫存決策,以及配送中心至零售商的配送車輛路徑?jīng)Q策。
由于銷售期內不同時段零售商的需求量隨機波動,為了提高車輛的使用效率,供貨商需要以配送中心和各零售商的庫存與配送路徑總成本最小的原則,使用異質車隊進行配送。由此可見,二級IRP優(yōu)化的目標是為了最小化配送中心與各零售商兩級庫存成本,以及配送中心向各零售商的車輛配送成本。決策內容包括:(1)產品從供貨商運送給配送中心的時間和數(shù)量;(2)配送中心將產品配送給零售商的時間和數(shù)量;(3)基于不同車型,如何組成配送中心—零售商的配送車輛和路徑。
為了便于問題的解決,本文做以下假設:
(1)在銷售期內每個時段零售商需求量已知,但在不同時段會發(fā)生隨機變化;
(2)只考慮配送中心和零售商的庫存持有成本,不考慮訂貨成本;
(3)配送中心和零售商的初始庫存為0;
(4)配送中心和零售商不允許缺貨;
(5)配送中心有足夠的各類型車輛可供派遣。
(1)集合
ΔC={1,2,…,N}表示零售商的集合;
Δ={0,1,…,N}表示配送網(wǎng)絡節(jié)點的集合,0表示配送中心;
ΔT={1,2,…,T}表示時段的集合。
(2)參數(shù)
Lij:表示從配送網(wǎng)絡節(jié)點i 到節(jié)點j 的距離,i,j ∈Δ;
Dit:表示在t 時段零售商i 對商品的需求量;
Qk:表示k 型車輛的最大載重量,k ∈{1,2,…,K};
Hi:表示配送網(wǎng)絡節(jié)點i 處單位商品的存儲成本,i ∈Δ;
Ii:表示配送網(wǎng)絡節(jié)點i 處的最大庫存容量,i ∈Δ;
F:表示每升燃油的成本;
W :表示司機每小時的薪水;
Rk:表示k 型車輛單位租車成本;
Q:表示干線運輸車輛最大載重量;
C:表示每輛干線車每次運輸?shù)目偝杀尽?/p>
(3)決策變量
wt:表示t 時段供貨商給配送中心運送產品的重量;
ut:表示t 時段供貨商給配送中心運送產品車輛數(shù);
zit:表示t 時段末,節(jié)點i 處商品的庫存數(shù)量;
由于本文采用異質車隊進行配送,車型不同會影響車輛載重、車輛速度等各方面差異,為準確計算不同車型在配送過程中的燃油消耗成本,參見文獻[5],引入公式如下:
式(1)中,F(xiàn)Ck表示第k 種車型從節(jié)點i 行駛到節(jié)點j過程中的燃油消耗量,λ 為燃油熱值系數(shù),φk為第k 種車型的發(fā)動機系數(shù),S 為車輛行駛速度,γk為第k 種車型的車輛傳動效率系數(shù),βk為第k 種車型的空氣阻力系數(shù),τ 為車輛的阻力系數(shù),μk為第k 種車型的整備質量。
于是,可建立考慮異質車隊的庫存路徑優(yōu)化模型如下:
目標函數(shù):
約束條件:
目標函數(shù)式(2)表示總成本最小化,第一項為配送中心和零售商的庫存成本,第二項為干線運輸成本,第三項為配送中心至零售商的配送司機薪水,第四項為配送中心至零售商的租車成本,第五項為配送中心至零售商的車輛燃油消耗成本;式(3)表示配送中心和零售商的初始庫存為0;式(4)和式(5)表示配送中心和零售商當前時段末庫存數(shù)量與上一階段末庫存數(shù)量及需求的關系約束;式(6)表示配送中心和零售商的最大庫存水平約束;式(7)表示本時段配送中心的庫存數(shù)量與本時段配送給零售商的產品數(shù)量的關系約束;式(8)表示零售商最大庫存水平的補貨策略約束;式(9)表示如果沒有車型為k 的第v 輛車訪問零售商i,則車型為k 的第v 輛車向零售商i 交付的產品數(shù)量為0;式(10)和式(11)表示配送給配送中心和零售商的裝載量約束;式(12)表示配送中心至零售商每個時段每個零售商最多被訪問一次;式(13)表示配送中心至零售商的配送任務只能由一輛車來執(zhí)行;式(14)表示車輛平衡等式;式(15)表示連接兩個決策變量之間的關系;式(16)表示產品平衡等式以及消除所有子回路;式(17)表示車輛不會從配送中心直接回到配送中心;式(18)~(21)表示決策變量非負;式(22)和式(23)表示0-1變量。
本文構建的上述模型是混合整數(shù)規(guī)劃模型,對于小規(guī)模算例,可以通過整數(shù)線性規(guī)劃求解器來解決,而對于較大規(guī)模的實例,需要設計啟發(fā)式算法求解。為此,結合模型的特征,本文設計改進的粒子群算法對模型進行求解。
設在一個D 維的搜索空間中,由n 個粒子組成種群X=(X1,X2,…,Xn),每個粒子在D 維搜索空間中的位置代表問題的一個解,第i 個粒子的位置為Xi=(xi1,xi2,…,xiD),速度為Vi=(Vi1,Vi2,…,ViD),個體極值為Pi=(Pi1,Pi2,…,PiD) ,種群的群體極值為Pg=(Pg1,Pg2,…,PgD),粒子在m+1 代的搜索空間中更新自身的速度和位置可如下確定[14]:
式中,ω 為慣性權重;d ∈(1,2,…,D);i ∈(1,2,…,n);m為粒子群算法當前迭代次數(shù);c1和c2為加速度因子,是非負的常數(shù);r1和r2是分布在[0,1]之間的隨機數(shù)。
慣性權重的設置影響著算法收斂速度和結果,較大的慣性權重能增強算法的全局搜索能力,而較小的慣性權重則增強算法的局部搜索能力。為避免算法較早陷入局部極值,本文采用線性遞減的慣性權重,其取值如下:
式中,ω 為慣性權重,ωmax是慣性權重最大值,ωmin是慣性權重最小值,m 是粒子群算法當前迭代次數(shù),mmax是粒子群算法的最大迭代次數(shù)。
此外,粒子群算法存在容易早熟和迭代后期收斂較慢的缺陷[15],因此針對粒子群算法的收斂速度、收斂精度和早熟問題,本文基于帳篷映射[16]對粒子群算法進行混沌局部搜索改進,運用帳篷映射模型對種群進行初始化,使粒子均勻地分布在空間內,增大粒子群的多樣性,以增加局部搜索能力。帳篷映射模型表達式為:
式中,p 為混沌局部搜索的當前迭代次數(shù),p ∈(1,2,…,n);為第j 維變量上的混沌變量。
因此,改進的粒子群算法步驟如下:
步驟1 設置參數(shù)。根據(jù)問題規(guī)模,設置種群大小、粒子群算法的最大迭代次數(shù)、混沌搜索的最大迭代次數(shù)、學習因子等相關參數(shù)。
步驟2 種群初始化。隨機產生粒子位置和速度,再利用混沌局部搜索進一步初始化粒子群的速度和位置,進行以下混沌局部搜索:
式中,xmax,j和xmin,j分別為第j 維變量的搜索上下界。
(3)判斷當前混沌局部搜索迭代次數(shù)p 是否達到混沌局部搜索的最大搜索次數(shù),如果未達到,則返回(2);否則,將混沌搜索的新結果作為新的粒子。
步驟3 根據(jù)目標函數(shù)式(2)計算步驟2產生的粒子適應度值。
步驟4 根據(jù)步驟3 計算的粒子適應度值尋找個體極值和群體極值。
步驟5 更新粒子。根據(jù)式(24)和式(25)進行粒子的速度和位置更新,并根據(jù)更新后的粒子按照目標函數(shù)式(2)重新計算粒子的適應度值。
步驟6 根據(jù)更新后新種群中的粒子適應度值更新個體極值和群體極值。
步驟7 終止條件判斷。判斷粒子群算法的迭代次數(shù)m 是否達到粒子群算法的最大迭代次數(shù)mmax,如果未達到,返回步驟4;否則,終止算法,輸出最優(yōu)值。
設有一家供貨商,在6 周的決策期內,供貨商通過其配送中心D 為位于8 個不同銷售區(qū)域的零售商C1~C8 提供VMI 服務,供貨商和零售商位置數(shù)據(jù)取自文獻[8]。每個零售商每周的隨機需求量數(shù)據(jù)和配送車輛的燃油消耗相關系數(shù)數(shù)據(jù)取自文獻[5]。供貨商至配送中心的集貨過程統(tǒng)一采用30 500 kg 車型進行運輸,其每車每次的干線運輸成本為3 000元,配送中心給各零售商配送產品的過程有三種車型,載重量分別為4 000 kg(車型1)、12 500 kg(車型2)和17 236 kg(車型3),車輛的行駛速度均為80 km/h。假設配送中心以及零售商的最大庫存水平分別為100 000 kg和10 000 kg,配送中心庫存持有成本為0.03 元/(kg·周),零售商庫存持有成本由均勻分布U[0.01,0.05] 元/(kg·周)隨機生成。
本文在CPU 為Intel?CoreTMi7-7700 3.60 GHz,內存為8 GB的電腦上利用Matlab 2016a軟件進行數(shù)值算例分析,得到供貨商使用異質車隊時(方案1)的二級IRP配送車輛行駛路徑、平均裝載率和總行駛距離的計算結果如表1所示。同時,根據(jù)每個時段零售商只能配送一次的限制以及結合零售商各時段的需求量數(shù)據(jù)情況,可知供貨商不能使用車型1組成的同質車隊進行配送。為此,為了進行異質車隊與同質車隊配送效果的對比,本文分別計算供貨商使用車型2(方案2)和車型3(方案3)組成的同質車隊進行配送的結果,亦見表1所示。
由表1中配送車輛行駛路徑可以看出,使用異質車隊時,會根據(jù)配送量合理安排配送車輛,減少車輛的租用成本。根據(jù)其行駛距離的計算結果顯示,方案1車輛總行駛距離小于方案2車輛總行駛距離,但略高于方案3 車輛總行駛距離。當考慮其平均裝載率時,方案1 車輛平均裝載率明顯高于另兩種方案,從而可以有效減少車輛空間的浪費。究其原因,是由于采用本文提出的方案1,配送中心可以根據(jù)顧客的需求選擇合適的配送車輛進行產品的配送。三種方案下各時段末的客戶庫存數(shù)量計算結果如表2所示。
由表2可以看出,在方案1中,配送中心與零售商的總庫存累計為84 900 kg,其中配送中心的庫存量累計為70 200 kg,各零售商總庫存量累計為14 700 kg;在方案2中,配送中心與零售商的總庫存累計為84 900 kg,其中配送中心的庫存量累計為61 900 kg,各零售商總庫存量累計為23 000 kg;在方案3 中,配送中心與零售商的總庫存累計為84 900 kg,其中配送中心的庫存量累計為65 200 kg,各零售商總庫存量累計為19 700 kg。說明在三種方案中,配送中心與零售商的庫存總量一樣,但在方案1中,庫存主要集中在配送中心處。究其原因是在方案1中車型相對多樣,配送中心會根據(jù)客戶需求合理安排車型,避免不必要的配送,從而減少配送成本,使零售商處庫存降低,因此配送中心處的庫存較高。
表1 三種方案車輛行駛路徑和總距離
表2 使用異質車隊與同質車隊各時段末的客戶庫存數(shù)量 kg
計算三種方案配送中心庫存成本、零售商庫存成本、租車成本、燃油成本和總成本等指標,得到結果如表3所示。
表3 三種方案的各項成本 元
從表3的總成本數(shù)據(jù)可以看出,在方案1中,使用異質車隊來完成配送比同質車隊更節(jié)約成本,盡管在庫存持有成本(配送中心與零售商庫存持有成本之和)、干線運輸成本及司機成本方面,三種方案的結果相差不大,但在租車成本和燃油消耗成本方面,方案1都有明顯的降低。這是因為零售商的需求波動變化會導致每一期產品配送量會有差異,使用異質車隊進行配送,能夠合理地配置車輛,使車輛租用成本和燃油消耗成本普遍減少,從而使總成本減少。此外,從環(huán)保的角度來說,也會降低車輛行駛中的排放污染。
為了驗證本文改進粒子群算法的有效性,將本文改進粒子群算法(算法1)與傳統(tǒng)粒子群算法(算法2)、Cplex最佳可行解算法(算法3)和Cplex精確解算法(算法4)進行比較,結果如表4所示。
表4 算法結果對比
由表4可見,本文提出的改進粒子群算法比傳統(tǒng)粒子群算法和Cplex最佳可行解算法的求解效果更好,其不僅與精確解相差不大(與Cplex精確解算法的誤差率為0.62%),而且其求解時間優(yōu)勢明顯。
圖2和圖3為本文算法與傳統(tǒng)粒子群算法運行結果。
圖2 改進粒子群算法迭代變化
圖3 粒子群算法迭代變化
從圖2 和圖3 的算法迭代變化對比可以看出,改進粒子群算法能有效地引導算法跳出局部最優(yōu),故算法改進效果明顯。
為了分析當產品需求變化波動不同時,考慮異質車型對二級庫存路徑問題的影響,本文在保持產品的需求總量為72 000 kg 不變的情況下,令產品的需求標準差分別為200 kg至1 600 kg,共8種情況,對算例進行敏感性分析。本文針對每種情況,運用不同的車型運行10次,根據(jù)10 次的平均值,得到使用異質車隊和同質車隊時二級IRP庫存路徑總成本對比結果,如圖4所示。
圖4 不同需求波動下總成本對比
由圖4可以看出,無論產品的需求波動程度怎樣變化,使用異質車隊時二級IRP庫存路徑總成本都要低于使用同質車隊時的總成本。這是由于產品需求波動變化程度會直接影響產品送貨量的大小,間接影響產品的庫存水平,由此對成本產生影響。當使用異質車隊配送時,即從三種車型中選取合適的車型來滿足配送需求,不僅可使車輛運送的路線減少,降低車輛的燃油消耗成本以及駕駛成本,提高車輛的裝載率,還可以降低配送中心和零售商的庫存成本,從而降低二級IRP庫存路徑總成本。
VMI 模式下的二級IRP 無論在理論上還是實際應用中都具有非常重要的研究價值。針對這一問題,本文考慮了客戶需求頻繁波動的實際情況,建立了異質車隊的二級IRP 優(yōu)化模型,并設計了改進粒子群求解算法,算例分析驗證了模型和算法的適用性和有效性。研究表明:(1)使用異質車隊,不僅可以提高配送車輛的裝載率,還能夠減少零售商的庫存水平;(2)本文提出的異質車隊二級IRP 模型,可以有效地降低二級IRP 庫存配送系統(tǒng)總成本;(3)不論零售商的需求如何波動,相比于同質車隊,異質車隊都會降低供貨商二級IRP庫存路徑總成本。研究結論可為供貨商VMI決策提供有益的參考。
本文在研究二級IRP 時,假設產品種類是單一的。然而,在實際配送中,產品種類也存在多樣化的情形。因此,考慮多種產品的異質車隊二級IRP優(yōu)化將是下一步的研究方向。