考慮運(yùn)輸時間的多目標(biāo)柔性作業(yè)車間調(diào)度研究

陳 魁，畢 利

(寧夏大學(xué) 信息工程學(xué)院，銀川 750021)

1 引 言

柔性作業(yè)車間調(diào)度問題(Flexible job shop scheduling problem，F(xiàn)JSP)是經(jīng)典車間調(diào)度問題(Job shop scheduling problem，JSP)的延伸.已知FJSP是一個NP-hard問題.目前，F(xiàn)JSP的求解方法主要集中在群體智能優(yōu)化算法，如遺傳算法[1]和粒子群算法[2]等.

粒子群算法(Particle Swarm Optimization，PSO)相較于遺傳算法具有收斂速度快、容易實(shí)現(xiàn)等優(yōu)點(diǎn)，在求解FJSP中被廣泛應(yīng)用.丁舒陽[2]提出了一種離散粒子群優(yōu)化算法求解FJSP，將不同的交叉策略引入PSO使其離散化.Shao[3]通過PSO的全局搜索與模擬退火算法的局部搜索相結(jié)合求解FJSP的帕累托解，采用擁擠距離來識別粒子的適應(yīng)度.Song[4]采用多目標(biāo)粒子群算法中加入變鄰域搜索策略提高了PSO的開發(fā)能力，通過分配規(guī)則和調(diào)度規(guī)則提高初始粒子的質(zhì)量.Edilson[5]將PSO和隨機(jī)重啟爬山算法結(jié)合求解多目標(biāo)FJSP.Rim[6]在求解FJSP時，提出了一種兩級粒子群優(yōu)化算法，將FJSP的兩個子問題機(jī)器選擇和操作排序分成兩層處理.Huang[7]提出了一種基于教和學(xué)的混合遺傳粒子群優(yōu)化算法來解決多目標(biāo)FJSP.仲于江[8]提出了一種將小生境技術(shù)和PSO相結(jié)合求解多目標(biāo)FJSP的優(yōu)化方法，基于小生境技術(shù)計(jì)算粒子的刪除概率對非支配解的外部存檔進(jìn)行更新.Li[9]提出了一種求解多模態(tài)優(yōu)化的小生境粒子群算法，在該算法中加入平衡因子，使得小生境劃分更均勻.

工件移動是柔性制造動態(tài)性的最直接體現(xiàn)，考慮工件運(yùn)輸時間的FJSP更貼切實(shí)際生產(chǎn)[10].趙寧[11]建立了考慮運(yùn)輸時間的關(guān)鍵鏈優(yōu)化方式，將運(yùn)輸時間集成到經(jīng)典FJSP的析取圖模型中，結(jié)合鄰域搜索算法實(shí)現(xiàn)快速尋優(yōu)，獲得FJSP在考慮運(yùn)輸時間下的近優(yōu)解.田旻等[12]設(shè)計(jì)了粒子群遺傳混合算法求解考慮運(yùn)輸時間的FJSP，并在種群初始化階段設(shè)計(jì)了一種考慮運(yùn)輸時間的初始化方式.Zhang[13]提出一種左移插入解碼的改進(jìn)遺傳算法求解運(yùn)輸時間約束的FJSP.Dai[14]提出改進(jìn)遺傳算法求解在考慮運(yùn)輸時間下的最小能耗和完工時間的FJSP.Huang[15]結(jié)合免疫遺傳算法和變鄰域搜索算法，提出了一種混合算法求解考慮運(yùn)輸是時間和劣化效應(yīng)的FJSP.Huang[16]將遺傳算法與模擬退火算法相結(jié)合，解決了考慮運(yùn)輸時間的MOFJSP.

綜上，在大多數(shù)研究中，PSO求解多目標(biāo)FJSP時算法全局搜索能力不足，并且未有效解決PSO的早熟問題；考慮運(yùn)輸時間的FJSP求解僅在解碼時考慮運(yùn)輸時間，而忽略了運(yùn)輸時間對粒子初始解質(zhì)量的影響.本文考慮到多目標(biāo)FJSP的解空間分布，建立基于Pareto支配的多目標(biāo)優(yōu)化模型，設(shè)計(jì)了一種PSO和小生境結(jié)合的算法來求解具有運(yùn)輸時間的多目標(biāo)柔性作業(yè)車間調(diào)度問題.算法改進(jìn)主要體現(xiàn)在3個方面：首先，在初始化階段就考慮工件的運(yùn)輸時間，初始化選擇工件的加工機(jī)器時動態(tài)選擇完工時間最短的機(jī)器；其次，在粒子群算法的優(yōu)化過程中，引入鄰域搜索算法增強(qiáng)PSO的局部探索能力；再次，提出了一種無須共享半徑的小生境劃分技術(shù)，克服了利用小生境技術(shù)劃分粒子種群時共享半徑難以確定的問題，并通過該技術(shù)將粒子群均勻劃分成多個小生境，提高PSO的全局搜索能力，有效避免PSO的早熟問題.最后，通過在Kacem數(shù)據(jù)集上做考慮運(yùn)輸時間的多目標(biāo)FJSP仿真實(shí)驗(yàn)，證明本文方法在求解不同規(guī)模運(yùn)輸時間下的FJSP時的優(yōu)越性.

2 問題描述與建模

2.1 問題描述

考慮運(yùn)輸時間的FJSP可以描述為：有n個工件在m臺機(jī)器上進(jìn)行加工，每個工件至少有一道工序，每道工序可以在至少一臺進(jìn)行加工，同一工件在不同機(jī)器之間的運(yùn)輸時間是確定的.該問題一方面要受到工序順序和設(shè)備資源約束，另一方面受到工件在不同機(jī)器間傳輸時間的約束[12].因此，合理地安排工序的加工順序，并為工序選擇合適的加工機(jī)器，可以使考慮工件運(yùn)輸時間的FJSP的目標(biāo)達(dá)到最優(yōu).為簡化待研究的問題，本文采用如下假設(shè)：

1)工序的加工過程不能中斷；

2)同一時刻同一機(jī)器只能加工一道工序；

3)每道工序在不同機(jī)器上的加工時間是確定的；

4)如果同一工件的前后工序在不同機(jī)器上加工，需要考慮工件的移動時間，否則，不需要考慮移動時間；

5)每個工件的首道工序之前和末道工序之后不考慮運(yùn)輸時間.

表1為一個柔性作業(yè)車間調(diào)度問題實(shí)例的加工時間表；表2為該實(shí)例的工件在機(jī)器之間的運(yùn)輸時間.

表1 加工時間表

表2 運(yùn)輸時間表

2.2 數(shù)學(xué)模型

符號說明：

Oij表示工件i的第j道工序；

Oi(j-1)表示工序Oij的前一道工序；

Oi′j′表示與工序Oij在同一臺機(jī)器上加工的前一道工序；

Pijh表示工序Oij在機(jī)器h上的加工時間；

Thk表示機(jī)器h與機(jī)器k之間的運(yùn)輸時間；

Sijh表示工序Oij在機(jī)器h上的開始加工時間；

Cijh表示工序Oij在機(jī)器h上的加工完成時間；

Cj表示工件j的完工時間；

Cmax表示所有工件的最大完工時間；

本文綜合考慮FJSP的最大完工時間(Cmax)、機(jī)器最大負(fù)載(Wm)和機(jī)器總負(fù)載(Wt)，建立最小化目標(biāo)的數(shù)學(xué)模型，該模型的目標(biāo)函數(shù)為公式(1)～公式(3)：

(1)

(2)

(3)

約束條件：

Cijh=Sijh+Pijh·xijh

(4)

Sijh=max(Ci(j-1)k+Tkh，Ci'j'h)

(5)

公式(4)表示工件按照工藝順序進(jìn)行加工，工序的完工時間等于開始時間和加工時間之和；公式(5)表示在加工過程中工序受工藝和運(yùn)輸時間的約束，只有當(dāng)機(jī)器不被占用且工件到達(dá)機(jī)器時才能開始加工.

2.3 多目標(biāo)優(yōu)化

本文采用基于Pareto支配的多目標(biāo)優(yōu)化方法，其中Pareto支配的定義為：對于決策空間中的兩個決策變量a和b，對?i∈{1，2，…，G}有fi(a)≤fi(b)，且?i∈{1，2，…，G}使得fi(a)

為測試算法的性能指標(biāo)，令算法A得到的Pareto解集為PFA，所有算法得到的Pareto解集為PFT，解集PFT中解的個數(shù)為nT；算法的非支配性指標(biāo)QS和收斂指標(biāo)GD定義如下：

1)解的非支配性QS：令解集PFA和解集PFT的交集有mA個，則算法A得到Pareto解集的QSA指標(biāo)定義為：

(6)

公式(6)可知，若QSA的值越接近1，則算法A求得解的非支配性越高，在所有算法求得非支配解集中占的比例越大.

2)當(dāng)代距離指標(biāo)GD：令解集PFA中解的個數(shù)為nA，解集PFA距離真實(shí)非支配解前沿PFT的距離GDA定義為：

(7)

公式(7)中Di表示解集PFA中第i個解到PFT中最近點(diǎn)的歐式距離，可知GDA的值接近0說明算法A求得解集到真實(shí)Pareto前沿距離越小，解集PFA收斂效果越好.

3 小生境粒子群算法求解考慮運(yùn)輸時間的FJSP

3.1 編碼與解碼

本文采用仲于江[8]提出的整數(shù)編碼實(shí)現(xiàn)PSO在非連續(xù)域上求解問題.粒子編碼分為兩個向量，即機(jī)器選擇向量MS和工序排序向量OS，其中，OS向量中從左到右工件號p出現(xiàn)的次數(shù)q表示該工件的工序號，記為Opq；MS向量從左到右按工件的工序順序排列，元素值表示對應(yīng)工序所選擇的機(jī)器號.

以表1的3*3FJSP為例，對某個可行解進(jìn)行粒子編碼如圖1所示，圖1中OS向量表示加工順序，MS向量表示對應(yīng)工序所選擇的加工機(jī)器，例如MS向量的第1個元素“1”表示工序O11選擇在機(jī)器1上加工.

圖1 粒子編碼

在對粒子進(jìn)行解碼時，先根據(jù)OS向量得到待加工工序，再根據(jù)MS向量找到該工序選擇的加工機(jī)器，在機(jī)器上從左向右尋找空閑時間段，若該時間段可以滿足開始時間和加工時間的要求，將工序放在該時間段加工，如果沒有滿足要求的空閑時間段，則將工序放到機(jī)器末尾加工.

3.2 動態(tài)初始化

初始解的優(yōu)劣直接影響算法的求解效率和求解質(zhì)量[12].考慮到運(yùn)輸時間會直接影響初始解的優(yōu)劣，本文設(shè)計(jì)了一種動態(tài)初始化的方法，在初始化過程中考慮加工時間和運(yùn)輸時間，為每道工序選擇完工時間最短的機(jī)器進(jìn)行加工.具體方法為，隨機(jī)初始化OS向量，從左向右依次選擇OS向量中的工序，為當(dāng)前工序選擇完工時間最短的機(jī)器記錄到MS向量的對應(yīng)位置；其中機(jī)器選擇向量MS確定方法的偽代碼如圖2所示.

圖2 MS向量的偽代碼

3.3 離散粒子群算法

基于丁舒陽等[2]提出的離散粒子群算法(Discrete Particle Swarm Optimization，DPSO)，本文引入自適應(yīng)慣性權(quán)重，設(shè)計(jì)一種自適應(yīng)DPSO的粒子位置更新方式，如公式(8)所示.

(8)

F1操作算子表示粒子受當(dāng)前位置影響的慣性運(yùn)動，在DPSO中體現(xiàn)為粒子變異.具體操作為分別對OS向量和MS向量進(jìn)行調(diào)整，其中基于OS的變異為：隨機(jī)交換OS向量中兩個位置的元素，實(shí)現(xiàn)對工序加工順序的調(diào)整；基于MS的變異為：隨機(jī)選擇一道工序，在該工序的可選機(jī)器集中隨機(jī)選一個與當(dāng)前機(jī)器號不同的機(jī)器號替換MS向量中的對應(yīng)元素，如果該工序只能在一臺機(jī)器上加工，則不做任何改變.

F2表示粒子受個體的當(dāng)前最優(yōu)位置影響進(jìn)行的粒子運(yùn)動，F(xiàn)2操作算子在OS向量中進(jìn)行POX交叉策略，在MS向量中進(jìn)行RPX交叉策略，其中POX交叉策略和RPX交叉策略的具體操作見文獻(xiàn)[2].F3表示粒子受當(dāng)前全局粒子最佳位置影響進(jìn)行的粒子運(yùn)動，F(xiàn)3操作算子采用RPX交叉策略對MS向量進(jìn)行調(diào)整.

公式(8)中ω表示慣性權(quán)重，本文采用自適應(yīng)的ω變化策略，如公式(9)所示，ω的值隨進(jìn)化代數(shù)在[ωmin，ωmax]區(qū)間內(nèi)線性減小，公式(9)中t表示當(dāng)前進(jìn)化代數(shù)，G表示總代數(shù).

(9)

3.4 局部搜索算法

DPSO在粒子運(yùn)動過程中受自身位置的慣性影響導(dǎo)致算法搜索不穩(wěn)定，故本文引入鄰域搜索算法提高DPSO的穩(wěn)定性，具體實(shí)現(xiàn)為：每代粒子位置更新完成之后，對每個粒子進(jìn)行局部搜索，即每個粒子生成它的鄰域解，如果鄰域解比粒子當(dāng)前解優(yōu)，則用鄰域解代替當(dāng)前解，否則不做改變.對該鄰域搜索算法重復(fù)進(jìn)行Km次，Km即為局部搜索代數(shù).

本文在局部搜索算法中針對工序排序和機(jī)器選擇設(shè)計(jì)了不同的鄰域解的產(chǎn)生方式.工序排序部分通過倒置、插入和交換操作來移動工序順序調(diào)整OS向量產(chǎn)生鄰域.例如，對粒子[31321213]的第3個位置和第6個位置進(jìn)行倒置操作后的粒子為[31212313]，進(jìn)行插入操作后的粒子為[31321213]，進(jìn)行交換操作后的粒子為[31221313].機(jī)器選擇部分隨機(jī)選擇MS向量中一個位置，確定該位置對應(yīng)的工序，在其對應(yīng)工序的可加工機(jī)器集上隨機(jī)選擇一臺與當(dāng)前機(jī)器不同的機(jī)器替換，若該工序只有一臺可加工機(jī)器，則另隨機(jī)選擇MS向量的一個位置，同上進(jìn)行機(jī)器替換.

3.5 小生境粒子群算法

小生境技術(shù)的核心為將種群劃分成多個小生境，常見的劃分方式是以某一粒子為小生境中心，根據(jù)小生境半徑確定其他粒子是否在該小生境中.這種小生境劃分方式需要確定小生境半徑，而且不能均勻劃分小生境，在以往的研究中并沒有確定小生境共享半徑的有效方法.本文設(shè)計(jì)了一種不需要小生境半徑的種群劃分方法，具體操作為：先根據(jù)小生境數(shù)量確定每個小生境應(yīng)該包含的粒子數(shù)S，再求出種群的Pareto解GB，以GB為中心找到與GB歐式距離最近的S-1個粒子，將包含GB在內(nèi)的這S個粒子劃分到一個小生境中；同理，求出剩余種群的Pareto解，在剩余種群劃分出S個粒子的小生境；以此類推，直到所有粒子均劃分完畢.該小生境劃分方式的偽代碼如圖3所示.

圖3 小生境劃分的偽代碼

粒子群優(yōu)化初期粒子間分散比較均勻，可通過小生境技術(shù)將種群劃分為更多相同大小的小種群，保證算法的全局搜索；隨著進(jìn)化代數(shù)的增加，需要減少小生境數(shù)量使算法收斂到最優(yōu)解，還需滿足當(dāng)前種群的均勻劃分.設(shè)計(jì)了一種滿足上述要求的小生境數(shù)量隨進(jìn)化代數(shù)遞減的方法，假設(shè)種群規(guī)模P=100，進(jìn)化代數(shù)G=100，對該規(guī)模粒子群的小生境數(shù)量遞減方式的代碼如圖4所示，圖5為該實(shí)例的小生境數(shù)量變化趨勢圖.

圖4 小生境數(shù)遞減的偽代碼

圖5 小生境數(shù)遞減圖

3.6 算法整體步驟

小生境粒子群算法(NPSO)求解考慮運(yùn)輸時間的多目標(biāo)FJSP時的算法流程圖如圖6所示.

圖6 算法流程圖

4 仿真實(shí)驗(yàn)

本文采用Matlab R2018b編程，運(yùn)行環(huán)境為Intel(R)Core(TM)i7-7700 CPU @3.6GHz，RAM 8GB；本文在標(biāo)準(zhǔn)算例Kacem數(shù)據(jù)集的8×8、10×10和15×10三組實(shí)例上進(jìn)行仿真，分別驗(yàn)證動態(tài)初始化的有效性和NPSO求解考慮運(yùn)輸時間的多目標(biāo)FJSP的有效性和可行性.Kacem標(biāo)準(zhǔn)算例沒有運(yùn)輸時間，故隨機(jī)生成0-1、1-5之間隨機(jī)分布的兩類運(yùn)輸時間[16].為便于實(shí)驗(yàn)對比，本文采用文獻(xiàn)[16]的兩類運(yùn)輸時間數(shù)據(jù).仿真實(shí)驗(yàn)的主要參數(shù)設(shè)置為：種群規(guī)模P=100，進(jìn)化代數(shù)G=100，鄰域搜索代數(shù)Km=10，學(xué)習(xí)因子c1=0.5，c2=0.8.

4.1 驗(yàn)證動態(tài)初始化

驗(yàn)證考慮了運(yùn)輸時間的動態(tài)初始化對算法求解質(zhì)量的影響，具體的實(shí)驗(yàn)設(shè)計(jì)為：分別在0-1和1-5兩類運(yùn)輸時間下的3組Kacem實(shí)例中通過隨機(jī)初始化和動態(tài)初始化得到兩組初始粒子群，將隨機(jī)初始化得到的粒子群作為對照組.在0-1運(yùn)輸時間下的Kacem算例在兩種初始化方式得到的初始解和經(jīng)過DPSO優(yōu)化得到的最優(yōu)解的結(jié)果見表3；在1-5運(yùn)輸時間下的Kacem算例通過兩種初始化方法得到的初始解和最優(yōu)解的結(jié)果見表4.表3和表4中的值為10次實(shí)驗(yàn)的平均值.

表3 不同初始化的初始解和最優(yōu)解(0-1運(yùn)輸時間)

表3和表4可知，通過動態(tài)初始化得到的初始解均優(yōu)于隨機(jī)初始化得到的初始解，且通過動態(tài)初始化種群經(jīng)過DPSO優(yōu)化后的最優(yōu)解均優(yōu)于隨機(jī)初始化種群優(yōu)化后的最優(yōu)解.

為了直觀體現(xiàn)動態(tài)初始化對粒子群算法優(yōu)化過程的影響，以實(shí)例15×10為例，比較兩類運(yùn)輸時間下動態(tài)初始化和隨機(jī)初始化兩種方式產(chǎn)生的初始種群，在DPSO優(yōu)化過程中各粒子完工時間平均值的收斂過程.圖7表示在0-1運(yùn)輸時間下實(shí)例15×10完工時間均值的收斂圖，圖8表示在1-5運(yùn)輸時間下實(shí)例15×10完工時間均值的收斂圖；可見通過動態(tài)初始化產(chǎn)生的初始種群在DPSO優(yōu)化過程中的收斂效果優(yōu)于隨機(jī)初始化產(chǎn)生的種群.

圖7 實(shí)例15×10完工時間的收斂圖(0-1運(yùn)輸時間)

4.2 驗(yàn)證NPSO求解考慮運(yùn)輸時間的多目標(biāo)FJSP

驗(yàn)證小生境粒子群算法(NPSO)在求解考慮運(yùn)輸時間的多目標(biāo)FJSP時的性能，對Kacem的3組案例分別在0-1和1-5兩類運(yùn)輸時間下進(jìn)行仿真，將NPSO與DPSO及Huang[16]提出的GASA進(jìn)行比較.

表5和表6分別給出了Kacem測試集中8×8、10×10、15×10三組實(shí)例在兩類運(yùn)輸時間下各算法求得的非支配解集.

表5 各算法在Kacem數(shù)據(jù)集上的最優(yōu)解(0-1運(yùn)輸時間)

表6 各算法在Kacem數(shù)據(jù)集上的最優(yōu)解(1-5運(yùn)輸時間)

表7和表8為兩類運(yùn)輸時間下各算法在Kacem數(shù)據(jù)集上求解的性能指標(biāo)對比結(jié)果.

表7 各算法求解的性能指標(biāo)對比結(jié)果(0-1運(yùn)輸時間)

表8 各算法求解的性能指標(biāo)對比結(jié)果(1-5運(yùn)輸時間)

由表7可知，NPSO在0-1運(yùn)輸時間下的3組算例上求得的解集的非支配性指標(biāo)QS均大于DPSO和Huang[16]提出的GASA，說明NPSO算法得到的解集的非支配性較高，解的質(zhì)量較好.NPSO得到解集的收斂指標(biāo)GD在8×8和10×10算例中均接近0，在15×10算例中的收斂指標(biāo)GD為0.57，明顯小于DPSO和GASA的1.55和1，說明NPSO得到的非支配解集距真實(shí)Pareto前沿的距離最近，該解集相較于其他兩個算法得到的解集的收斂性較好.

由表8可知，NPSO在求解1-5運(yùn)輸時間下的Kacem算例時，得到解集的非支配性指標(biāo)QS均小于DPSO和GASA的QS，說明NPSO算法得到非支配解集的質(zhì)量最好.在收斂性指標(biāo)GD的比較中，NPSO和GASA算法在算例8×8上得到解集的GD都接近0，NPSO在算例15×10中得到解集的GD為0.32均小于GASA和DPSO的1.85和2.39，說明NPSO算法得到解集的收斂性較好.

綜上，NPSO在兩類運(yùn)輸時間下的FJSP算例中得到非支配解集的收斂性和分布性均優(yōu)于未改進(jìn)的DPSO和Huang[16]提出的GASA，尤其在求解高維度算例時NPSO得到解的質(zhì)量更高.

5 總 結(jié)

本文針對考慮運(yùn)輸時間的多目標(biāo)FJSP，提出一種將小生境技術(shù)和離散粒子群算法相結(jié)合的優(yōu)化算法，構(gòu)建了以最大完工時間、最大機(jī)器負(fù)載、總機(jī)器負(fù)載為最小優(yōu)化目標(biāo)的多目標(biāo)優(yōu)化模型.首先，在粒子群初始化階段，提出了一種動態(tài)初始化方法產(chǎn)生考慮工件加工時間和運(yùn)輸時間的初始種群，提高初始種群的質(zhì)量.其次，針對傳統(tǒng)粒子群算法的不穩(wěn)定性和易早熟現(xiàn)象，采用鄰域搜索算法提高粒子群算法在優(yōu)化過程中的穩(wěn)定性，并將小生境技術(shù)與粒子群算法結(jié)合，給出一種以Pareto解為中心多次迭代的種群劃分方法，降低全局最優(yōu)解對粒子位置更新的影響，擴(kuò)大粒子搜索空間，有效避免算法陷入局部最優(yōu)解.最后，在Kacem標(biāo)準(zhǔn)數(shù)據(jù)集上做不同規(guī)模運(yùn)輸下的仿真實(shí)驗(yàn)，結(jié)果表明，動態(tài)初始化可提高粒子群算法的優(yōu)化效率，小生境粒子群算法在求解考慮運(yùn)輸時間的多目標(biāo)FJSP時是可行有效的.