董 海，吳 瑤

(沈陽大學 1.應用技術學院；2.機械工程學院，遼寧 沈陽 110044)

隨著全球變暖加劇和工業(yè)全球化，企業(yè)對于環(huán)境、社會效益的重視日漸提高。近年來，在Listes等[1]最早研究了閉環(huán)供應鏈網絡設計并提出了隨機規(guī)劃模型之后，越來越多的學者在進行閉環(huán)供應鏈(closed loop supply chains，CLSC)網絡設計時開始考慮環(huán)境和社會效益等因素。Pati等[2]利用混合整數(shù)目標規(guī)劃模型設計一個多產品紙張回收網絡，用以降低總成本，并通過提高廢紙回收率實現(xiàn)環(huán)境效益。Dehghanian等[3]利用生命周期分析(life cycle analysis，LCA)和層次分析法(analytic hierarchy process，AHP)建立一個多目標數(shù)學規(guī)劃模型，使經濟效益和社會效益達到最大化。Devika等[4]也通過同時考慮CLSC中經濟、環(huán)境和社會因素，建立混合整數(shù)規(guī)劃模型。Chen等[5]從可持續(xù)性的角度考慮太陽能行業(yè)的供應鏈成本和環(huán)境問題，并研究集成CLSC網絡的設計問題。Sahebjamnia等[6]在保護環(huán)境和最大化社會效益的前提下，研究輪胎行業(yè)可持續(xù)CLSC網絡設計的混合整數(shù)線性規(guī)劃模型，用以優(yōu)化總成本。呂寶龍等[7]在考慮利潤最大化的同時，兼顧碳稅情況和產品綠色度，研究閉環(huán)供應鏈網絡Nash均衡模型。

現(xiàn)實生產中的CLSC網絡的設計，有很多的不確定性，而這些因素往往會造成企業(yè)風險較大，成本較高等問題。學者在CLSC網絡的設計上，大多采用隨機規(guī)劃、模糊規(guī)劃和魯棒控制等方法降低不確定性的影響。Ozceylan等[8]提出不確定性條件下CLSC網絡的多產品、多周期的混合整數(shù)模糊數(shù)學模型，并借助線性隸屬度函數(shù)實現(xiàn)參數(shù)和目標的確定性轉化。Ozkir等[9]開發(fā)一種強調客戶滿意度和凈利潤的多目標CLSC網絡模型，并采用模糊集理論進行求解。Ozceylan等[10]建立戰(zhàn)術和戰(zhàn)略決策場景下的CLSC網絡規(guī)劃模型，并利用模糊交互方法解決CLSC網絡設計問題。Tsao等[11]提出一個基于CLSC網絡的可持續(xù)多目標規(guī)劃模型，引入一個交互式的兩階段模糊概率法來處理各種不確定性。Fakhrzad等[12]提出一個基于模糊數(shù)的多目標混合線性整數(shù)規(guī)劃模型，達到成本最小，碳排放最小且供應商可靠性最強的目標。但以上文獻均采用悲觀值或樂觀值處理模糊問題，結果具有極端性。

CLSC網絡設計是一個NP-hard問題，解決此類問題多采用智能算法。Pishvaee等[13]通過基于特殊局部搜索的模擬退火算法解決了多級RL網絡的規(guī)劃模型。Wang等[14]將環(huán)境問題整合到一個整數(shù)CLSC模型中，并提出一種基于生成樹結構的遺傳算法求解所提出的模型。Fakhrzad等[15]在解決多目標CLSC網絡設計的NP-hard問題時，采用非支配排序遺傳算法進行求解。但傳統(tǒng)的智能算法步驟繁瑣且用時較高，2016年Mirjalili等[16]開發(fā)的鯨魚優(yōu)化算法(whale optimization algorithm, WOA)具有原理簡單，人工參數(shù)少和收斂速度優(yōu)等特點。

綜上，本文構建基于Me測度的CLSC網絡的模糊規(guī)劃模型，并提出一種基于帕累托(Pareto)最優(yōu)的改進鯨魚算法(Pareto improved whale optimization algorithm，PIWOA)求解綠色CLSC網絡優(yōu)化設計問題，最后通過實例驗證所提模型及算法的優(yōu)越性。

1 模型建立

1.1 問題描述

本文研究一個多層級、多周期的CLSC網絡，如圖1所示。其中包含7個節(jié)點：原材料供應商I、制造中心M、配送中心D、客戶集合C、分離中心S、處理中心O和回收中心R。正向流動包含原料供貨給制造中心生產，再由配送中心分發(fā)給客戶。逆向物流包括客戶中心退回的產品由分離中心挑揀，可回收再利用的貨物發(fā)到回收中心轉運給制造中心進行再制造，而無法再利用的產品發(fā)給處理中心進行統(tǒng)一銷毀或填埋等處理。本文主要解決供應鏈運營成本最低，碳排放最小和社會效益最大的閉環(huán)供應鏈網絡多目標優(yōu)化問題。

圖1 CLSC網絡結構圖Figure 1 CLSC network structure diagram

1.2 模型假設及符號定義

本文假設如下：1) 產品回收利用能夠節(jié)約成本；2) 每個節(jié)點設施均有容量限制；3) 制造與再制造的生產成本一致，再制造只省掉一部分原材料費用；4) 客戶集合和原材料供應商的位置已知且固定，對于各類產品需求也已知；5) 各節(jié)點間運輸只采用一種運輸方式?；诮Ｐ枰?，定義符號如表1～4。

表1 模型中集合定義Table 1 Collection definitions in the model

本文建立如下模型Z1、Z2、Z3。

目標函數(shù)1(Z1)表示最小化網絡成本，分別包含原材料成本、各個節(jié)點的建設費用、配送運輸費用、正逆向物流的節(jié)點處理費用以及回收節(jié)省的原料費用。目標函數(shù)2(Z2)表示最小化環(huán)境污染，分別包含建設時CO2排放量和運輸時CO2排放量。目標函數(shù)3(Z3)表示最大化社會效益。約束條件主要包含各節(jié)點流量約束和各節(jié)點能力約束，具體模型如式(1)～(3)，約束條件如式(4)～(15)。

表2 模型中參數(shù)定義Table 2 Parameters in the model

表3 0-1變量Table 3 Variables 0-1

表4 決策變量Table 4 Decision variables

式(4)～(9)和 式(13)中 ξi、 ξd、 ξr、 ξo、 ξs、 ξm、ξc?[0,1]，表示節(jié)點I、D、R、O、S、M、C的能力約束的置信水平。式(10)～(12)表示配送中心、再制造中心、分離中心的流量守恒，式(14)～(15)表示客戶中心的可用產品數(shù)量約束和決策變量取值約束。

2 模型處理

為減弱模型中的不確定性，本文采用模糊規(guī)劃方法對已建模型和約束進行模糊處理。以往文獻的模糊測度均采用極端值表示。例如，Pos表示可能性測度，極喜；Nec表示必要性測度，極悲，此類方法數(shù)據(jù)測定不夠準確。而本文的模糊測度采用Xu等[17]提出的Me測度法，此類方法能夠實現(xiàn)在2個極端值間靈活變化，減小結果的誤差性。具體定義如下[18]，其中 λ是事物發(fā)生的概率參數(shù)，A表示事件。

其中， λ∈(0,1)。下式中，引入Me測度水平值K， K ∈[0.5,1]。

對應上下界估計模型如下。

據(jù)此目標函數(shù)1的清晰等價模型Lower和Upper如Z4所示。

式(31)～(34)同式(10)、(11)、(12)、(14)。

Lower

s.t.

式(35)～(38)同式(31)～(34)。

3 求解方法

本文中所采用的解決方法是WOA算法。近幾年來，WOA算法因其操作模式簡單，人工參數(shù)量少，已應用于許多大規(guī)模求解問題。但在傳統(tǒng)的WOA算法中，其結果容易陷入局部最優(yōu)的狀態(tài)，使問題難以得到最優(yōu)值?；诖祟悊栴}，本文引入變異因子，并結合實數(shù)編碼和Pareto優(yōu)化方法，對WOA算法進行改進，得到改進后的鯨魚算法(PIWOA)，該算法力求在解決多目標問題的基礎上，提高全局搜索能力的同時又避開了局部最優(yōu)值。

3.1 原始的鯨魚算法

WOA的基本思想是通過模擬鯨魚獨特的捕食行為和社會行為來進行優(yōu)化。通過3個數(shù)學模型模擬了包圍獵物、氣泡網攻擊和尋找獵物的過程。

3.1.1 包圍獵物

該階段選擇一個搜索代理，將其定義為前一代中的最佳解決方案。其他搜索代理試圖避開最優(yōu)代理來實現(xiàn)全局搜索，具體描述如下。

式中， Π 代表當前迭代次數(shù)； D W是包圍步長；AW和 CW是系數(shù)向量； XW(Π)是當前鯨魚位置向量； XW′(Π)是 當前最優(yōu)鯨魚位置向量；XW(Π+1)是每一次迭代中需要更新的最佳位置向量。系數(shù)向量A W和C W的計算如下。

其中，?稱為收斂因子，其隨迭代次數(shù)增加從2線性遞減到0；rand為[0,1]之間的隨機數(shù)。?的求解公式如下。

Π為當前迭代次數(shù)，Πmax是最大迭代次數(shù)。

3.1.2 泡沫網攻擊

鯨魚在捕食的時候通常采用螺旋運動的方式先包圍獵物，再進行捕獵。

在這一階段討論了2種方法：收縮環(huán)和螺旋更新位置。收縮環(huán)的機制與全局搜索相似，并取值AW[?1,1]。另一種方法是螺旋更新位置,根據(jù)當前位置和最優(yōu)代理構造對數(shù)螺旋曲線，使得搜索代理逐步逼近最佳位置。具體描述如下。

式中， D W′指當前鯨魚與最佳位置(也就是與獵物位置)的距離；Λ為對數(shù)螺旋形狀的常數(shù)量；l是介于 [ ?1,1]之間的隨機數(shù)。

3.1.3 搜索獵物

搜索獵物采取魚群隨機性位置更新的方法，數(shù)學模型如下。

式中， Xrand是 隨機所得的位置向量，如果 A W超出 [?1,1]的范圍，就通過隨機尋找鯨魚個體，并根據(jù)個體的方位，排查其他鯨魚位置，面向全局搜索找到合適的獵物。

3.2 變異操作

為了更準確地找到全局最優(yōu)位置，避免陷入局部最優(yōu)并保證在搜索后期具有較快的收斂速度，本文對原有式中收斂因子 ?進行深層考慮，該因子可以在搜索前期使最小值與迭代次數(shù)成正比，同時又能將 ?在增大到較大值時與迭代次數(shù)成反比。考慮以上內容，本文將隨進化迭代次數(shù)而非線性變化的收斂因子更新，公式如下。

其中，?start與?final為 ? 的 初值與終值，? ＞0，稱為非線性調節(jié)系數(shù)[19]。

3.3 編碼與變量處理

模型中決策變量只能取0-1，利用rand函數(shù)和round函數(shù)隨機取值0/1。本文采用平均隨機數(shù)生成函數(shù)配合約束系數(shù)，得到決策變量。在計算中，各個變量均采用實數(shù)編碼方法進行編碼。將個體適應度子集F編碼為群體中的解 Xτ(τ =1,···,θ)，假定子集包含 θ個個體適應度的值，即F={F(ρ1) ，F(xiàn)(ρ2),···,F(ρθ) }，若F(ρθ+1) ＞F(ρθ) ，則F(ρθ+1)被選擇。

3.4 多目標優(yōu)化處理

解決多目標問題時，可將模型轉化成如下。

式中，ρ的集合ZM為問題的決策空間；Fτ(ρ),(τ=1,2,3,···,θ)為 目標函數(shù)；θ維函數(shù)稱為問題的目標空間。多目標優(yōu)化問題無法達到每項最優(yōu)，只能得到相對最優(yōu)解集，而最優(yōu)解的集合稱為非支配集。

經處理過后的PIWOA算法如圖2所示。圖中P是用于確定選擇搜索環(huán)繞還是螺旋更新位置的因子。假設收縮包圍機制和螺旋位置更新概率相同，均為0.5。

圖2 PIWOA算法流程Figure 2 PIWOA algorithm flow

4 算例分析

本文以某玻璃制造廠商為例，此供應鏈由原料供應商、制造中心、配送中心、客戶、分離中心、處理中心和回收中心7個節(jié)點組成。為評價本文采用的PIWOA算法的性能，本文將PIWOA與WOA和遺傳算法(GA)在Matlab軟件上運行并進行比較。Matlab的運行計算機環(huán)境為：處理器為Intel(R)Core(TM)i7-8550u CPU＠ 3.9 GHz，內存8 GB，操作系統(tǒng)為Window10的手提電腦。

本文針對目標函數(shù)共生成10個測試問題如表5，參數(shù)范圍如表6，并運用PIWOA算法、WOA算法和遺傳算法對每個測試問題運行20次，得到表7。在表7中權重ζ1、 ζ2、 ζ3的取值是根據(jù)不同測試問題中相關成本所占的比重而確定的。將種群規(guī)模都為200的PIWOA算法、WOA算法和GA算法設置相同數(shù)據(jù)進行400次迭代，目標函數(shù)1的Lower與Upper函數(shù)值與迭代次數(shù)曲線如圖3、4所示，函數(shù)平均值與計算時間曲線如圖5所示。

表5 測試問題Table 5 Test questions

表6 參數(shù)取值范圍Table 6 Parameters scope

表7 算法計算結果Table 7 Results of algorithm calculation

圖3 Z1的Lower最優(yōu)值與迭代次數(shù)圖Figure 3 The lower optimal value of Z1 and the graph of iteration times

圖4 Z1的Upper最優(yōu)值與迭代次數(shù)圖Figure 4 The upper optimal value of Z1 and the graph of iteration times

圖5 Z1的平均函數(shù)值與求解時間圖Figure 5 Average function value and solution time diagram of Z1

根據(jù)表7中3組數(shù)據(jù)對比可得，使用PIWOA算法時，總成本Z1整體數(shù)值最小，CO2排放量Z2整體數(shù)值最小，社會效益Z3整體數(shù)值最大，方差最小，穩(wěn)定性最強。而WOA算法相比于GA算法，總成本Z1整體數(shù)值偏大，CO2排放量Z2整體數(shù)值偏多，社會效益Z3整體數(shù)值趨同，但方差更小，所以WOA算法穩(wěn)定性更強。根據(jù)圖3、4可得，無論Lower函數(shù)還是Upper函數(shù)，采用PIWOA算法對比WOA算法和GA算法總是能更快地跳出局部最優(yōu)，精確尋找到最優(yōu)解且成本函數(shù)值總是相對更低。而WOA算法與GA算法相比，成本雖略高些，但WOA算法可以更快搜索到最優(yōu)解。同時，通過圖5可知PIWOA的算法時間最短，其次是GA和WOA。

5 結論

針對綠色閉環(huán)供應鏈網絡設計的模糊優(yōu)化問題，本文得出以下結論。

1) 建立以成本最小、CO2排放最小和社會效益最大的目標函數(shù)，并采用Me測度和三角模糊數(shù)對模型和約束進行模糊處理，通過實例驗證，采用模糊處理后的閉環(huán)供應鏈網絡優(yōu)化模型可以減弱模糊影響并取得較好的結果。

2) 采用PIWOA算法解決NP-hard問題，并將該方法與WOA算法和GA算法進行對比。結果表明，PIWOA算法計算時間較短，能夠較快跳出局部最優(yōu)且穩(wěn)定性強，在處理閉環(huán)供應鏈網絡設計問題時，具有較好的優(yōu)越性。

3) 實際生產中，供應鏈成員成分更為復雜，很多貨物要考慮時間約束，下一步的研究工作主要考慮貨物種類多樣性對于閉環(huán)供應鏈網絡設計時的影響，以保證閉環(huán)供應鏈網絡設計更加貼近現(xiàn)實生產。