馬艷芳，張 文，李宗敏，閆 芳，郭凌云

（1.河北工業(yè)大學(xué) 經(jīng)濟(jì)管理學(xué)院，天津 300401；2.四川大學(xué) 商學(xué)院，成都 610064；3.重慶交通大學(xué) 經(jīng)濟(jì)管理學(xué)院，重慶 400074；4.重慶市環(huán)衛(wèi)集團(tuán)有限公司，重慶 401121）

0 引言

研究發(fā)現(xiàn)，建立臨時(shí)中轉(zhuǎn)站，構(gòu)建高效的醫(yī)院-中轉(zhuǎn)站運(yùn)輸模型是解決城市醫(yī)療垃圾運(yùn)輸問題的關(guān)鍵方向之一［1］。Liu 等［1］結(jié)合環(huán)境影響評價(jià)指南，采用免疫算法建立中轉(zhuǎn)站選址模型；然后采用蟻群-禁忌混合算法，建立醫(yī)院與中轉(zhuǎn)站之間的路徑優(yōu)化模型；最后用武漢市作實(shí)例取得了良好的驗(yàn)證效果。Babaee Tirkolaee 等［2］研究了具有時(shí)間窗口的選址-路徑問題（Location-Routing Problem，LRP），構(gòu)建以廢物收集車輛的總行駛時(shí)間最小、時(shí)間窗口違規(guī)數(shù)量和居住在處理設(shè)施地點(diǎn)周圍的人口數(shù)量最小化為目標(biāo)的多目標(biāo)模型，為高效地收集、運(yùn)輸和處理醫(yī)療廢物給各國家提供了參考。趙佳虹等［3］根據(jù)環(huán)境系統(tǒng)的動(dòng)態(tài)性，設(shè)計(jì)時(shí)變“環(huán)境-人口”風(fēng)險(xiǎn)度量模型，建立風(fēng)險(xiǎn)最小和成本最小的多目標(biāo)選址-分配模型，可根據(jù)決策者偏好提供多個(gè)有效優(yōu)化方案。這些研究為醫(yī)療廢物管理提供了新思路，但大部分僅考慮單級(jí)選址路徑。雖然廢棄物暫存點(diǎn)或中轉(zhuǎn)站的建立十分必要，但是最終端的處理地點(diǎn)（焚燒站）也需重視；因此本文考慮兩級(jí)選址-路徑問題（Two-Echelon LRP，2E-LRP）。

2E-LRP 是復(fù)雜的組合優(yōu)化問題，眾多學(xué)者致力于該類問題的研究。Jacobsen 等［4］是首先考慮選址-路徑問題中多階存在性的經(jīng)典代表，提出了三個(gè)建設(shè)性啟發(fā)式求解模型。雖然有不少學(xué)者提出運(yùn)用純精確算法解決兩級(jí)問題［5-7］，但由于2E-LRP 包含多個(gè)非確定性多項(xiàng)式（Non-deterministic Polynomial，NP）問題，常常會(huì)加入不止一種啟發(fā)式法，尤其在求解子問題選址-路徑問題（LRP）的過程中使用。馬艷芳等［8］在求解多個(gè)決策主體下的LRP 時(shí)設(shè)計(jì)了遺傳和粒子群混合算法；胡大偉等［9］在求解兩級(jí)選址-路徑模型時(shí)，采用改進(jìn)的兩階段模擬退火（Simulated Annealing，SA）算法，求解時(shí)間明顯少于商業(yè)求解軟件CPLEX。鯨魚優(yōu)化算法（Whale Optimization Algorithm，WOA）雖然廣泛應(yīng)用于各個(gè)領(lǐng)域［10-12］，但關(guān)于2E-LRP 的研究較少。本文將改進(jìn)WOA 與SA算法結(jié)合，用于求解所提兩級(jí)多目標(biāo)選址-路徑（Two-Echelon Multi-Objective LRP，2E-MOLRP）模型。

本文結(jié)合目前國內(nèi)具體情況，借鑒前人的研究思路，研究在物流網(wǎng)絡(luò)中同時(shí)考慮垃圾中轉(zhuǎn)站和焚燒站的兩級(jí)選址-路徑問題。除了考慮經(jīng)濟(jì)目標(biāo)之外，同時(shí)把焚燒站運(yùn)行時(shí)對周圍居民產(chǎn)生的影響量化為負(fù)效應(yīng)值作為目標(biāo)函數(shù)。為求解所提2E-MOLRP 模型，設(shè)計(jì)鯨魚-模擬退火多目標(biāo)算法WOA-SA；改進(jìn)初始化方式，采用非線性動(dòng)態(tài)慣性權(quán)重系數(shù)和非支配排序方法以得到更好的結(jié)果；最后，基于2E-LRP、LRP 兩組基準(zhǔn)案例進(jìn)行算法對比分析，驗(yàn)證WOA-SA 的尋優(yōu)能力、收斂性和穩(wěn)定性，再結(jié)合天津市實(shí)際情況驗(yàn)證兩級(jí)多目標(biāo)模型的適用性。

1 問題描述

1.1 兩級(jí)選址-路徑問題

兩級(jí)選址-路徑問題（2E-LRP）的兩級(jí)物流網(wǎng)絡(luò)由若干個(gè)備選焚燒站、若干備選中轉(zhuǎn)站以及一些日常垃圾產(chǎn)生點(diǎn)（以下簡稱垃圾產(chǎn)生點(diǎn)）組成。收集過程中，垃圾首先在產(chǎn)生點(diǎn)被車輛收集，然后被送至?xí)捍嬷修D(zhuǎn)站，被簡單處理后運(yùn)到焚燒站處置，其中：焚燒站與中轉(zhuǎn)站及相應(yīng)路徑構(gòu)成一級(jí)運(yùn)輸網(wǎng)絡(luò)；中轉(zhuǎn)站和垃圾產(chǎn)生點(diǎn)及相應(yīng)收集路徑構(gòu)成二級(jí)運(yùn)輸網(wǎng)絡(luò)。具體網(wǎng)絡(luò)結(jié)構(gòu)如圖1 所示。

圖1 兩級(jí)物流網(wǎng)絡(luò)結(jié)構(gòu)Fig.1 Two-echelon logistics network structure

1.2 負(fù)效應(yīng)

垃圾收集與處置是為城鎮(zhèn)居民創(chuàng)造更好的居住環(huán)境，但焚燒站或中轉(zhuǎn)站在運(yùn)行時(shí)會(huì)產(chǎn)生污染，某些特殊背景下，這種有害健康的設(shè)施更會(huì)帶來心理上的排斥；垃圾處理設(shè)施對周邊居民產(chǎn)生的負(fù)面影響應(yīng)當(dāng)被考慮。

負(fù)效應(yīng)是指設(shè)施運(yùn)行對人體健康、生命或生存造成的負(fù)面損害。大多數(shù)人將討厭效應(yīng)定義為與設(shè)備容量正相關(guān)、與兩個(gè)節(jié)點(diǎn)之間的距離負(fù)相關(guān)的函數(shù)，并假設(shè)負(fù)效應(yīng)在所有方向上一致傳播［13］。垃圾在焚燒站處理過程中會(huì)產(chǎn)生廢氣，氣體傳播速度與風(fēng)向和風(fēng)速有關(guān)，所以并不是所有的居民區(qū)都受到同樣的影響，因此，本文考慮所選節(jié)點(diǎn)的氣候，建立了與風(fēng)向和距離有關(guān)的負(fù)效應(yīng)函數(shù)。

定義焚燒站為點(diǎn)(xi，yi)(i∈B)，附近居民區(qū)節(jié)點(diǎn)為(xj，yj)(j∈P)，dij為焚燒站與居民區(qū)節(jié)點(diǎn)之間的歐氏距離，兩點(diǎn)間直線與當(dāng)前風(fēng)向的夾角為θij，ωij(dij)是與距離有關(guān)的負(fù)效應(yīng)參數(shù)。若兩點(diǎn)間距離一定，θij在0°到90°之間增長時(shí)，風(fēng)速在兩點(diǎn)間直線方向分量減小，則對居民區(qū)產(chǎn)生的負(fù)效應(yīng)值也會(huì)減小；若只有距離為變量時(shí)，居民區(qū)受到最大的負(fù)效應(yīng)影響值為Mi，在距離閾值[di，1，di，2]內(nèi)隨著速率mb逐漸減小，在di，1以內(nèi)負(fù)效應(yīng)值均為Mi，在di，2以外均為0。此外，負(fù)效應(yīng)值還與焚燒站的規(guī)模正相關(guān)，處理的廢棄物數(shù)量越多，對周邊造成的影響越大。

負(fù)效應(yīng)參數(shù)ωij(dij)的示意圖如圖2 所示，其公式為：

圖2 負(fù)效應(yīng)參數(shù)ωij(dij)的二維與三維示意圖Fig.2 2D and 3D schematic diagrams of obnoxious effect parameter ωij(dij)

2 2E-MOLRP模型構(gòu)建

2.1 模型假設(shè)

1）候選焚燒站、中轉(zhuǎn)站和垃圾產(chǎn)生點(diǎn)的數(shù)量、地理位置、固定成本已知，各級(jí)運(yùn)輸車輛固定成本已知，垃圾產(chǎn)生點(diǎn)的垃圾量已知。

2）各級(jí)物流設(shè)施容量已知，分配給各級(jí)物流設(shè)施的總需求量不能超過設(shè)施的容量限制。

3）各級(jí)配送車輛容量已知，各運(yùn)輸路徑上車輛載重不能超過容量限制。

4）各節(jié)點(diǎn)間距離按照歐氏距離計(jì)算，各級(jí)運(yùn)輸車輛單位長度運(yùn)輸成本已知。

5）每個(gè)垃圾產(chǎn)生點(diǎn)只能被一輛車服務(wù)一次，各級(jí)運(yùn)輸車輛從物流設(shè)施節(jié)點(diǎn)出發(fā)并返回相同物流設(shè)施節(jié)點(diǎn)。

2.2 符號(hào)說明

本文所使用的符號(hào)說明如表1 所示。

表1 參數(shù)定義Tab.1 Parameter definition

2.3 模型建立

式（1）～（2）是模型的目標(biāo)函數(shù)：式（1）是負(fù)效應(yīng)目標(biāo)函數(shù)，最小化垃圾焚燒站對周圍居民產(chǎn)生的環(huán)境影響，與垃圾焚燒站的容量、周圍居民數(shù)量、和負(fù)效應(yīng)參數(shù)ωij(dij)有關(guān)。式（2）為二級(jí)網(wǎng)絡(luò)相關(guān)的經(jīng)濟(jì)成本，由六部分組成，依次為：垃圾焚燒站的開放成本、中轉(zhuǎn)站的開放成本、第一級(jí)配送網(wǎng)絡(luò)中的車輛固定成本和運(yùn)輸成本、第二級(jí)配送網(wǎng)絡(luò)中的車輛固定成本和運(yùn)輸成本。

式（3）～（4）保證開放的垃圾焚燒站和中轉(zhuǎn)站數(shù)量合理：

式（5）表示回收設(shè)施在一、二級(jí)路徑上回收的垃圾數(shù)量不超過其自身容量：

式（6）保證每個(gè)被選擇的垃圾中轉(zhuǎn)站只能被一輛一級(jí)配送車輛服務(wù)一次，未被選擇的中轉(zhuǎn)站不能被服務(wù)：

式（13）～（15）保證一個(gè)垃圾產(chǎn)生點(diǎn)只分配給一個(gè)垃圾中轉(zhuǎn)站，且只有開放的垃圾中轉(zhuǎn)站才能有車輛進(jìn)出：

式（22）表示如果一輛車的當(dāng)前載重加上下一點(diǎn)需求若小于車輛容量，它就會(huì)服務(wù)下一點(diǎn)：

式（23）保證垃圾中轉(zhuǎn)站發(fā)出的車輛總數(shù)不能超出其擁有的數(shù)量：

式（24）～（25）表示垃圾產(chǎn)生點(diǎn)的垃圾量和弧（i，j）上的車輛載荷均非負(fù)：

式（26）保證開放的垃圾中轉(zhuǎn)站容量能夠滿足垃圾產(chǎn)生點(diǎn)的需求：

3 非支配WOA-SA設(shè)計(jì)

2E-MOLRP 是NP 難問題，因?yàn)榘藥讉€(gè)已知復(fù)雜問題：兩級(jí)設(shè)施選址問題、兩級(jí)車輛路徑問題、多目標(biāo)優(yōu)化和有容量約束選址-路徑問題，因此很難用精確方法求解。本文設(shè)計(jì)了非支配WOA-SA 求解模型，算法的流程如圖3 所示。

圖3 本文算法流程Fig.3 Flow of proposed algorithm

3.1 第一階段：WOA-SA

第一階段要求解的問題可視為單層選址-路徑問題（LRP），中轉(zhuǎn)站的選址極其重要；涉及大量垃圾點(diǎn)的分配，路徑的優(yōu)化可節(jié)約大量物流成本，合理的選址也利于一級(jí)網(wǎng)絡(luò)成本優(yōu)化。為更好地求解第二級(jí)物流網(wǎng)絡(luò)問題，采用混合啟發(fā)式算法WOA-SA，在改進(jìn)的WOA 之后運(yùn)用SA 算法以增強(qiáng)找到的最佳解決方案。

WOA 是澳大利亞學(xué)者M(jìn)irjalili 等［14］根據(jù)座頭鯨的狩獵方式提出的一種群智能優(yōu)化算法。WOA 和其他智能優(yōu)化算法一樣，在求解復(fù)雜優(yōu)化問題時(shí)，存在迭代后期種群個(gè)體多樣性減少、易陷入局部最優(yōu)的缺陷。

Kirkpatrick 等［15］提出基于單一解的模擬退火（SA）算法可以克服局部最優(yōu)陷于停滯的缺陷，該算法利用一定的概率來接受較差的解，從而加強(qiáng)跳出局部最優(yōu)的能力。Mafarja等［16］證明了WOA 和SA 的結(jié)合方式優(yōu)于WOA 的性能，并且大量實(shí)驗(yàn)結(jié)果表明在WOA 之后再使用SA 比在WOA 中嵌入SA 的效果更好。本文借鑒此思路加以改進(jìn)作為第一階段算法，求解中轉(zhuǎn)站選址問題及第二級(jí)運(yùn)輸路徑。

3.1.1 編碼與解碼

本階段研究垃圾產(chǎn)生點(diǎn)（客戶）的垃圾收集路徑和中轉(zhuǎn)站選址兩個(gè)子問題，車輛向量、順序向量和車輛所屬中轉(zhuǎn)站的中心向量表示一個(gè)完整的解。

X=｛X1，X2，…，XN｝表示鯨魚群的位置，N是客戶數(shù)，X1，X2，…，XN是區(qū)間（1，K+1）生成的隨機(jī)數(shù)，K是車輛個(gè)數(shù)。鯨魚位置編碼X的整數(shù)部分表示二級(jí)路徑上每位顧客對應(yīng)的被服務(wù)車輛編號(hào)，在解碼中由Xk表示，即每位顧客對應(yīng)的被服務(wù)車輛。鯨魚位置編碼X的小數(shù)部分表示顧客被服務(wù)的順序，在解碼中由XS表示，降序排列，即小數(shù)部分越小對應(yīng)顧客被服務(wù)的優(yōu)先級(jí)越高。R是中轉(zhuǎn)站個(gè)數(shù)，XR為隨機(jī)生成的K個(gè)（1，R）內(nèi)整數(shù)，表示每輛車所屬中轉(zhuǎn)站。該編碼方式能夠滿足模型中對于中轉(zhuǎn)站的一些約束，例如：保證被選中的中轉(zhuǎn)站數(shù)量至少為1 但不超過候選個(gè)數(shù)，每個(gè)客戶只能被一個(gè)開放的中轉(zhuǎn)站服務(wù)，滿足約束中的式（4）、式（13）～（17）、式（21）和式（23）。同理，焚燒站的數(shù)量和它與中轉(zhuǎn)站之間的關(guān)系也能足約束式（3）、式（6）～（8）。

設(shè)有4 個(gè)中轉(zhuǎn)站、10 個(gè)客戶和3 輛車，則編碼解碼的示例圖如圖4 所示。

圖4 編碼與解碼示例圖Fig.4 Schematic diagram of coding and decoding

從圖4 可以看出，在第二級(jí)物流網(wǎng)絡(luò)中選擇開放編號(hào)為1 和2 的中轉(zhuǎn)站，共派出3 輛車服務(wù)10 個(gè)客戶。例如，Xk中1對應(yīng)客戶2、3、5，對應(yīng)XS中數(shù)字8、6、4，對應(yīng)XR的第一個(gè)數(shù)字為1，則車輛1 從中轉(zhuǎn)站1 出發(fā)先后服務(wù)客戶5、3、2。此解的完整表示如圖5 所示。

圖5 2E-LRP的決策示例Fig.5 Decision example of 2E-LRP

3.1.2 初始化種群

WOA 是基于種群的算法，標(biāo)準(zhǔn)WOA 一般隨機(jī)生成初始化種群。本文希望得到一個(gè)良好的初始化種群，為后續(xù)提供一個(gè)良好的開端，避免計(jì)算資源的浪費(fèi)；也希望種群具有良好的多樣性，降低算法復(fù)雜度且不易陷入局部最優(yōu)。因此，本文采用不止一種初始化方法來解決當(dāng)前的問題。首先，采用聚類的方法來分配客戶給中轉(zhuǎn)站；然后，同時(shí)采用隨機(jī)方法與Clarke 和Wright（Clarke and Wright，CW）節(jié)約算法生成初始車輛路徑（其中：80%的初始種群是隨機(jī)的；20%由CW生成）。CW 中節(jié)約距離計(jì)算公式為：

聚類方法如下。

步驟1 計(jì)算所需最少回收中心個(gè)數(shù)。

步驟2 隨機(jī)從R中選取Num個(gè)中轉(zhuǎn)站，將選定的序號(hào)放入集合C中，每個(gè)集合都有一個(gè)對應(yīng)的空集合Si。

步驟3 計(jì)算客戶點(diǎn)I=(xi，yi)與每個(gè)選中中轉(zhuǎn)站R=(xj，yj)之間的距離并排序，把客戶分配給最近的中轉(zhuǎn)站R，即客戶I∈SR。計(jì)算兩點(diǎn)間的距離：

步驟4 計(jì)算每個(gè)集合SR中的客戶數(shù)量，并按降序排序，選擇排名最高且未被選中的中轉(zhuǎn)站開放。

步驟5 根據(jù)中轉(zhuǎn)站容量分配客戶。

3.1.3 WOA-SA更新流程

在WOA 的三個(gè)捕捉獵物階段游走位置不同，因此算子的更新公式也各不相同。當(dāng)|A| ＜1 時(shí)，下一代計(jì)算公式為：

其中：t為當(dāng)前迭代次數(shù)；X*、X分別是當(dāng)前最佳和目前個(gè)體的位置向量；A、C、B是系數(shù)；r為［0，1］的隨機(jī)向量；a從2 到0遞減；b是對數(shù)螺線形狀的常數(shù)；l在［-1，1］均勻分布；p為［0，1］的一個(gè)隨機(jī)數(shù)，表示每種包圍方式的概率。圖6 舉例說明式（30）的更新過程。

圖6 WOA的更新過程Fig.6 Update progress of WOA

為了執(zhí)行全局搜索，Xrand為隨機(jī)數(shù)，通過隨機(jī)更新而不是依賴迄今為止找到的全局最佳個(gè)體解來選擇下一代個(gè)體位置。此時(shí)A?(-1，1)，具體公式如下：

由此可以看出WOA 設(shè)置的參數(shù)較少，全局與局部搜索的平衡由參數(shù)A決定，而根據(jù)式（31）可知它跟a密切相關(guān)：a越小代表算法局部搜索能力越好；反之全局能力增強(qiáng)。在復(fù)雜的函數(shù)優(yōu)化中，無法保證a線性減小的規(guī)律易導(dǎo)致算法陷入局部最優(yōu)。針對上述不足，引入余弦非線性控制因子，前期a從2 緩慢減小利于充分的全局搜索，后期下降速度加快促進(jìn)算法收斂，表達(dá)式如下：

此外，文獻(xiàn)［17］指出權(quán)重因子也會(huì)影響算法性能：當(dāng)權(quán)重因子較大時(shí)，其搜索范圍較大，有利于全局搜索；當(dāng)權(quán)重因子較小時(shí)，可在較小范圍內(nèi)執(zhí)行高精度搜索。為了增強(qiáng)WOA 的尋優(yōu)能力，采用非線性動(dòng)態(tài)慣性權(quán)重系數(shù)，更新公式如下：

其中ω(t)i代表第t代第i個(gè)搜索代理的權(quán)重，公式如下：

其中：f(t)avg、f(t)min、f(t)max是適應(yīng)度值的平均、最小和最大值；ω1、ω2是初始值；Tmax為最大迭代次數(shù)。當(dāng)個(gè)體優(yōu)于平均水平時(shí)，權(quán)重減小保留該個(gè)體，加快算法收斂；反之促使它向最優(yōu)個(gè)體靠近，增強(qiáng)算法全局搜索能力。

在WOA 得到最優(yōu)解之后，利用SA 算法避免陷入局部最優(yōu)解，具體流程可參考文獻(xiàn)［16］，此處不再贅述。

3.2 第二階段：精確求解一級(jí)網(wǎng)絡(luò)

第二階段求解時(shí)已知中轉(zhuǎn)站信息，簡化第一級(jí)物流網(wǎng)絡(luò)問題為多目標(biāo)選址-路徑問題。先計(jì)算所有被選中的中轉(zhuǎn)站需求，然后精確求解目標(biāo)函數(shù)式（1）～（2），適應(yīng)度函數(shù)的設(shè)置在下文詳細(xì)介紹。兩個(gè)目標(biāo)會(huì)產(chǎn)生非劣解，因此使用快速非支配排序與擁擠距離計(jì)算方法［18］：建立一個(gè)帕累托解集的存檔E，把得到的解放入存檔；通過非支配排序，根據(jù)非支配關(guān)系將每組解決方案排列到不同的前沿，擁擠距離計(jì)算機(jī)制用于確保非支配解的分布盡可能均勻。此外，每當(dāng)帕累托解集存檔存儲(chǔ)已滿時(shí)，在存檔的刪除操作中再次使用擁擠距離計(jì)算。算法1 描述了擁擠距離的計(jì)算過程。最后給出一組完全非支配的帕累托解集供決策者參考。

3.2.1 適應(yīng)度函數(shù)

在WOA-SA 中，每個(gè)個(gè)體的適應(yīng)度值能夠反映目標(biāo)值，據(jù)此判斷個(gè)體優(yōu)劣。本文考慮設(shè)施和車輛的容量約束，對于違反約束的個(gè)體設(shè)置一個(gè)足夠大的懲罰值，避免生成不可行解。在第一階段，只考慮中轉(zhuǎn)站的選址-路徑問題，求解的目標(biāo)函數(shù)是式（2）中的一部分，因此該適應(yīng)度函數(shù)的設(shè)置也適用于第一階段算法。當(dāng)滿足容量約束時(shí)，以目標(biāo)函數(shù)式（1）～（2）作為適應(yīng)度函數(shù)，當(dāng)違反容量約束時(shí)，則采用帶有懲罰值的適應(yīng)度函數(shù)：

其中：P1、P2、P3、P4和P5是懲罰因子，P1對應(yīng)負(fù)效應(yīng)函數(shù)，P2和P3對應(yīng)違反焚燒站和中轉(zhuǎn)站容量的情況，P4和P5對應(yīng)違反一、二級(jí)車輛容量的情況。通過適應(yīng)度函數(shù)的設(shè)置，保證模型中約束的實(shí)現(xiàn)，如與車輛容量有關(guān)的約束為式（9）～（11）、式（18）～（20）、式（22）～（26），設(shè)施容量約束為式（5）和式（12）。

3.2.2 擁擠距離計(jì)算

擁擠距離計(jì)算的偽代碼見算法1。

4 基準(zhǔn)案例測試

首先，選用兩組經(jīng)典基準(zhǔn)進(jìn)行測試，驗(yàn)證本文算法的性能，包括一組2E-LRP 算例［19］和一組LRP 算例［20］。兩個(gè)算例集合都屬于國際通用案例，已經(jīng)被諸多學(xué)者引用測試［21-23］，均可在http：//prodhonc.free.fr 下載。其次，結(jié)合天津案例，驗(yàn)證考慮負(fù)效應(yīng)的兩級(jí)多目標(biāo)模型的合理性。本文構(gòu)建的數(shù)學(xué)模型與算法使用Matlab R2012a 編碼，在Dell XPS 13 筆記本運(yùn)行，Intel Core i5-7200U 處理器2.70 GHz，操作系統(tǒng)為Windows10（64 位）。

4.1 算例參數(shù)設(shè)置與描述

本文算法涉及的參數(shù)主要有：迭代次數(shù)、種群規(guī)模、慣性權(quán)重、初始溫度、終止溫度、降溫速度和較差解的接受概率等。通過實(shí)驗(yàn)測試并參考文獻(xiàn)［8，16］，設(shè)置迭代次數(shù)為200，種群規(guī)模為30，ω1=0.4，ω2=0.9，T0=1 000，T=1/1 000，降溫速度q=0.96，較差解的接受概率設(shè)定為p=，內(nèi)循環(huán)迭代次數(shù)設(shè)置為50。

本文選用Nguyen 等［19］提出的2E-LRP 算例中的一組來驗(yàn)證所提算法的有效性，該算例由Prins 等［24］提出的一級(jí)LRP 算例改進(jìn)而來。在LRP 基礎(chǔ)上添加坐標(biāo)（0，0）作為主倉庫，將原始倉庫作為中轉(zhuǎn)站，并且假設(shè)從主倉庫到中轉(zhuǎn)站的距離是兩倍的歐幾里得距離。選取的24 個(gè)案例具有以下特征：客戶數(shù)量n∈{20，50，100}；中轉(zhuǎn)設(shè)施數(shù)量m∈{5，10}；客戶聚類個(gè)數(shù)β∈{1，2，3}；二級(jí)車輛容量Q∈{70，150}；一級(jí)車輛的容量是最大中轉(zhuǎn)站容量的1.5 倍；一、二級(jí)車輛的運(yùn)輸成本分別通過距離乘以200 和100 并四舍五入計(jì)算得到。每個(gè)算例的命名格式為“n-m-β”或“n-m-βb”，后綴“b”表示算例中Q=150。

考慮到一級(jí)物流網(wǎng)絡(luò)中涉及焚燒站和中轉(zhuǎn)站的數(shù)量點(diǎn)較少，而單級(jí)選址-路徑問題更廣泛，所以選取LRP 基準(zhǔn)案例進(jìn)行測試。除了Prins 等［24］提出的LRP 算例，較為流行的還有Barreto 等［20］提出的LRP 算例，他總結(jié)了諸多前人研究的算例，并改進(jìn)部分車輛路徑問題算例。

本文選取11 組算例，數(shù)據(jù)與Prins 算例［24］的主要區(qū)別在于無車輛固定成本，客戶數(shù)量n∈[20，100]，設(shè)施個(gè)數(shù)m∈[5，10]，車輛容量各不相同，具體可參考文獻(xiàn)［20］，每個(gè)算例的命名格式為“人名-n-m”。

4.2 算法對比分析

為評估本文算法的性能，將本文算法求解基準(zhǔn)案例的結(jié)果分別與目前已知最優(yōu)解（Best Known Solution，BKS）和現(xiàn)有其他算法的求解結(jié)果進(jìn)行了對比。對比結(jié)果如表2 所示，其中：差值Gap為本文算法最優(yōu)計(jì)算結(jié)果Best與BKS的比較，計(jì)算公式為Gap=100%*(Best-BKS)/BKS。

選取的對比算法為：多起點(diǎn)局部搜索（Multi-Start Iterated Local Search，MS-ILS）算法［19］、帶學(xué)習(xí)過程的貪婪隨機(jī)自適應(yīng)搜索法（Greedy Randomized Adaptive Search Procedure × Learning Process，GRASP×LP）［21］、帶路徑重連的貪婪隨機(jī)自適應(yīng)搜索法（Greedy Randomized Adaptive Search Procedure × Path-Relinking，GRASP×PR）［21］、SA［22］和帶路徑重連的多起點(diǎn)混合啟發(fā)式（Multi-start Hybrid Heuristic with Path Relinking，MHH-PR）算法［23］?？梢钥闯?，MS-ILS 算法整體表現(xiàn)良好，在大規(guī)模算例中求得最優(yōu)解次數(shù)最多，其次是SA。與MS-ILS 算法相比，本文算法求解大規(guī)模算例的能力較弱，只有100-10-2b 和100-10-3 兩個(gè)案例求解結(jié)果優(yōu)于其他對比算法；但本文算法在中小規(guī)模算例中均能找到最優(yōu)解，且總體上平均值為188 718，僅次于表現(xiàn)最優(yōu)的MS-ILS（188 308），差值的平均值僅為0.37%。

表3 為不同算法對Barreto 等［20］提出的11 個(gè)LRP 基準(zhǔn)案例的求解結(jié)果對比。對比算法為多起點(diǎn)偏差隨機(jī)（Multi-Start Biased-Randomized Heuristic，MSBRH）算法［25］和有偏隨機(jī)變鄰域搜索（Biased-Randomized Variable Neighborhood Search，BR-VNS）算法［25］。表3 中記錄了各算法運(yùn)行10 次的平均值、最優(yōu)值以及本文算法與最優(yōu)值的差值，每個(gè)案例的最優(yōu)結(jié)果和平均值用粗體標(biāo)識(shí)。從表3 可以看出，除Ch69-100×10 外，本文算法能找到所有案例的最優(yōu)解，且多次運(yùn)行結(jié)果的平均值也都是3 個(gè)算法中最優(yōu)的。

表3 不同算法求解Barreto案例的結(jié)果比較Tab.3 Results comparison of different algorithms solving Barreto case

結(jié)合表2～3 可以看出，本文算法在求解顧客規(guī)模小于100 的時(shí)候均能找到最優(yōu)解；在大規(guī)模算例中雖不能全部找到最優(yōu)解，但是與其他算法差距不大，具有一定的競爭力。

表2 不同算法求解Prins案例的結(jié)果對比Tab.2 Results comparison of different algorithms solving Prins case

4.3 算法收斂性及穩(wěn)定性驗(yàn)證

為探究提出的WOA-SA 在求解不同規(guī)模算例的能力，記錄算法求解案例時(shí)每代運(yùn)行結(jié)果并繪制迭代圖。選取顧客數(shù)量為21、50、88 和100 的4 個(gè)案例，收斂曲線如圖7 所示。可以看出50 代以內(nèi)算法快速下降，100 代左右已收斂于最優(yōu)值附近，100 到150 代之內(nèi)下降趨勢減弱逐漸到達(dá)最優(yōu)結(jié)果。

圖7 本文算法求解Barreto案例的收斂曲線Fig.7 Convergence curves of the proposed algorithm on Barreto cases

為進(jìn)一步驗(yàn)證算法性能，檢驗(yàn)算法穩(wěn)定性。將4.2 節(jié)中35 個(gè)基準(zhǔn)案例均運(yùn)行10 次的結(jié)果與BKS 的差值繪制為箱型圖，以觀察不同算法求解案例差值的離散程度。從圖8（a）可以看出，求解2E-LRP 案例時(shí)，本文算法的異常值最少；雖然矩形的高度不是最短（僅次于MS-ILS），但是遠(yuǎn)小于MHHPR 算法，且總體高度明顯低于GRASP×LP 和GRASP×PR 算法。從圖8（b）可以看出，求解LRP 案例時(shí)，本文算法存在一個(gè)異常值（Ch69-100×10），但是矩形高度很低，這說明數(shù)據(jù)分布集中。

此外，繪制3 種算法對Barreto 案例中前10 個(gè)案例的求解結(jié)果均值的條形圖，如圖8（c）所示。本文算法結(jié)果在大多數(shù)情況下都低于其他兩種算法。

圖8 算法穩(wěn)定性分析結(jié)果Fig.8 Algorithm stability analysis results

綜上所述，WOA-SA 無論是求解2E-LRP 還是LRP 案例，在不同規(guī)模中都具有較強(qiáng)的穩(wěn)定性。

4.4 實(shí)例應(yīng)用分析

為驗(yàn)證所提模型與算法的適用性，以天津市垃圾回收問題為案例，解決垃圾暫存中轉(zhuǎn)點(diǎn)與焚燒站兩者的選址以及中間路徑規(guī)劃問題。

4.4.1 天津市案例描述

外賣打包和快遞便利了生活，但也導(dǎo)致了垃圾數(shù)量急劇增長，其回收與處理工作已引起政府高度重視。具體到各市區(qū)的社區(qū)地理位置已知且穩(wěn)定，選取垃圾點(diǎn)個(gè)數(shù)為50，候選中轉(zhuǎn)站個(gè)數(shù)為6，候選焚燒站個(gè)數(shù)為3。備選焚燒站（A）和中轉(zhuǎn)站（B）的經(jīng)緯度信息、開放成本和容量的具體信息見表4。它們均擁有足量的車輛以供調(diào)配，一級(jí)車輛能力為5 000 kg，固定成本400 元；二級(jí)車輛容量為1 000 kg，固定成本150 元。垃圾點(diǎn)的經(jīng)緯度信息和需求量信息見表5，為了表述方便，用編號(hào)代指真實(shí)地名，且只展示了前10 個(gè)地點(diǎn)。

表4 候選站點(diǎn)的信息Tab.4 Candidate stations information

表5 垃圾點(diǎn)的部分信息Tab.5 Part information of waste stations

相關(guān)的公路網(wǎng)、行政區(qū)域、人口信息數(shù)據(jù)由天津市交通廳和天津市區(qū)域地理數(shù)據(jù)庫獲取。除經(jīng)濟(jì)成本外，決策者還考慮焚燒站對周圍居民的負(fù)效應(yīng)影響，提出的負(fù)效應(yīng)函數(shù)與當(dāng)?shù)貧夂蛳嚓P(guān)，因此在中國天氣網(wǎng)等平臺(tái)獲取天津市風(fēng)向信息。天津市風(fēng)向具有明顯的季節(jié)變化特征：冬季，西北風(fēng)、北風(fēng)多；夏季有東南、南風(fēng)；春秋兩季，西南風(fēng)比較多。設(shè)t=｛1，2，…，9｝分別代表9 種常見的風(fēng)向，其中第9 種風(fēng)是靜風(fēng)，每年有45 天是靜風(fēng)，其余8 種風(fēng)向各40 天。此外，負(fù)效應(yīng)函數(shù)根據(jù)距離劃分三段式，選取兩個(gè)閾值分別為0.5 km和1.0 km［15］。

4.4.2 2E-LRP模型分析

運(yùn)用本文算法求解天津市案例，設(shè)置種群規(guī)模為30，迭代200 次，運(yùn)行10 次記錄結(jié)果。

第一階段算法求出的第二級(jí)物流網(wǎng)絡(luò)最優(yōu)選址決策為開放候選點(diǎn)序號(hào)為2，3，5 的暫存中轉(zhuǎn)站，經(jīng)濟(jì)成本為91 319.30 元。其他選擇可有1，2，5（91 449.12）或者2，4，5（91 398.33），都需派發(fā)8 輛車服務(wù)50 個(gè)顧客。雖然第一種的成本目前最小，但是還要累加一級(jí)物流網(wǎng)絡(luò)的經(jīng)濟(jì)成本，如中轉(zhuǎn)站4 距離焚燒站3 較近而離焚燒站2 比較遠(yuǎn)，而中轉(zhuǎn)站1 距離焚燒站3 最遠(yuǎn)而離1 和2 較近。

第一級(jí)物流網(wǎng)絡(luò)中開放焚燒站的成本是一樣的，因此經(jīng)濟(jì)成本的決定因素在于開放的中轉(zhuǎn)站與焚燒站之間的距離。從圖9（a）可以看出，焚燒站2 距離大多數(shù)中轉(zhuǎn)站較近，服務(wù)路程短，但是否開放編號(hào)為2 的焚燒站還應(yīng)考慮另一個(gè)負(fù)效應(yīng)目標(biāo)。

4.4.3 多目標(biāo)分析

一級(jí)物流網(wǎng)絡(luò)決策者只需開放一個(gè)候選焚燒站即可，因考慮到設(shè)施運(yùn)行對周邊居民產(chǎn)生的負(fù)面影響，需平衡經(jīng)濟(jì)成本和負(fù)效應(yīng)值這兩個(gè)目標(biāo)。運(yùn)行結(jié)果如圖9，表6 給出幾個(gè)推薦方案的數(shù)據(jù)。

圖9 天津市案例信息圖Fig.9 Information maps of Tianjin case

表6 天津市案例的帕累托解集Tab.6 Pareto solution set of Tianjin case

從表6 可以看出，焚燒站2 的經(jīng)濟(jì)成本明顯低于焚燒站1、3，但是焚燒站2 周圍居民區(qū)較1、3 更為密集，負(fù)效應(yīng)值最大；焚燒站1 和3 的經(jīng)濟(jì)成本和負(fù)效應(yīng)值無明顯差異。綜上所述，若決策者更注重社會(huì)影響則可開放焚燒站1 或3；若可接受焚燒站2 的負(fù)效應(yīng)影響則可節(jié)約許多經(jīng)濟(jì)成本。

負(fù)效應(yīng)值與設(shè)施附近居住人口數(shù)量有關(guān)，該參數(shù)穩(wěn)定無法通過算法優(yōu)化，因此選取帕累托解集中有最優(yōu)路徑的解展示，如圖9（c）所示。一級(jí)網(wǎng)絡(luò)用虛線表示路徑，派出至少2輛車服務(wù)；二級(jí)網(wǎng)絡(luò)路徑用實(shí)線表示，共8 輛車，分別用不同顏色展示。它們的裝載量分別為：900、930、980、970、910、950、900 和860 kg，均在容量范圍內(nèi)且利用率合理。

4.4.4 靈敏度分析

為驗(yàn)證負(fù)效應(yīng)函數(shù)的合理性，對相關(guān)參數(shù)設(shè)施容量（QB）和附近的人數(shù)（p）進(jìn)行了靈敏度分析，分析結(jié)果見表7，第一列是參數(shù)變化幅度，以最優(yōu)成本291 595 為基準(zhǔn)。

表7 模型參數(shù)靈敏度分析Tab.7 Sensitivity analysis of model parameters

總體而言，負(fù)效應(yīng)值與兩個(gè)參數(shù)呈正相關(guān)，但在不同的變化區(qū)間內(nèi)波動(dòng)幅度不同。當(dāng)參數(shù)降低到一定程度時(shí)，負(fù)效應(yīng)值趨于穩(wěn)定，因?yàn)樵O(shè)施一旦運(yùn)行就會(huì)對周圍環(huán)境產(chǎn)生影響。當(dāng)需求未超出設(shè)施的容量時(shí)，人口增長導(dǎo)致負(fù)面影響的緩慢增加。但增長到一定程度（115%）時(shí)，負(fù)效應(yīng)值激增，可能產(chǎn)生垃圾的數(shù)量超過現(xiàn)有設(shè)施的容量，需開放新設(shè)施。當(dāng)其他參數(shù)不變時(shí)，設(shè)施容量小范圍變化時(shí)目標(biāo)值變化幅度穩(wěn)定，若容量過大，即使處理垃圾不多波及半徑也會(huì)增大。

此外，還對二級(jí)車輛載重（QK）和垃圾量（q）進(jìn)行分析，以驗(yàn)證參數(shù)對經(jīng)濟(jì)目標(biāo)的影響。車輛載重QK越小，使用數(shù)目越多，經(jīng)濟(jì)成本越高。但是，當(dāng)QK增加到一定程度（115%）時(shí)，車輛的固定成本將非常高，造成經(jīng)濟(jì)浪費(fèi)。還有，參數(shù)q與經(jīng)濟(jì)成本呈正相關(guān)，但是當(dāng)q增加到110%時(shí)，經(jīng)濟(jì)成本顯著增加，這表明需求超過了目前設(shè)施容量，必須開設(shè)一個(gè)新的設(shè)施來服務(wù)客戶。

靈敏度分析的結(jié)果清楚地表明應(yīng)在人口稀少的區(qū)域建設(shè)適當(dāng)規(guī)模的垃圾處理設(shè)施，并應(yīng)選擇適當(dāng)大小的車輛。

5 結(jié)語

本文研究垃圾的收集、運(yùn)輸、處理問題，構(gòu)建兩級(jí)多目標(biāo)選址-路徑模型，通過優(yōu)化中轉(zhuǎn)站和焚燒站的選址及路徑問題，以達(dá)到負(fù)效應(yīng)最低和經(jīng)濟(jì)成本最小化的目標(biāo)。

為求解提出的模型，設(shè)計(jì)改進(jìn)的兩階段WOA-SA 非支配算法：1）初始化方式，先聚類分配客戶到中轉(zhuǎn)站，后用隨機(jī)（20%）和CW（80%）同時(shí)初始化路徑；2）采用非線性動(dòng)態(tài)慣性權(quán)重系數(shù)調(diào)整收斂速度以增強(qiáng)WOA 的尋優(yōu)能力；3）采用WOA 和SA 并行結(jié)構(gòu)加強(qiáng)算法跳出局部最優(yōu)的能力；4）運(yùn)用快速非支配排序方法確定帕累托解集以供決策者參考。

為驗(yàn)證算法求解2E-LRP 的能力，選取24 個(gè)2E-LRP 基準(zhǔn)案例，與已知5 個(gè)對比算法的結(jié)果進(jìn)行對比分析。在13 組數(shù)據(jù)中，本文算法均優(yōu)于其他對比算法，總體平均提升了0.37%。為驗(yàn)證算法求解LRP 的能力，選取11 個(gè)LRP 基準(zhǔn)案例，與其他對比算法相比，本文算法能找到除Ch69-100×10外其他所有案例的最優(yōu)解，且總體差值平均值（0.08%）是最優(yōu)的。測試了算法收斂性和穩(wěn)定性，結(jié)果表明算法在前期快速收斂于最值附近，且在不同規(guī)模的算例中表現(xiàn)穩(wěn)定。

此外，設(shè)計(jì)基于距離和風(fēng)向的負(fù)效應(yīng)分段函數(shù)，結(jié)合當(dāng)?shù)貧夂?，考慮設(shè)施與顧客點(diǎn)直線方向和風(fēng)向的夾角，使其更加貼合實(shí)際情況。為驗(yàn)證所提模型和算法的適用性，選取天津市為分析對象，給出三種不同選擇方案的負(fù)效應(yīng)值以及經(jīng)濟(jì)成本供決策者參考，并給出最優(yōu)經(jīng)濟(jì)成本時(shí)的路徑方案。

實(shí)驗(yàn)結(jié)果表明，本文模型與算法具有廣泛的適用性，為生活垃圾的收集處置給各地政府提供決策支持，也為其他類型的垃圾處置提供思路。但是本文也具有一定的局限性，比如醫(yī)療類廢棄物的處理不僅需要考慮設(shè)施的影響，還需要考慮運(yùn)輸風(fēng)險(xiǎn)，未來可就此進(jìn)一步研究。