分布式多柔性裝配作業車間調度問題研究

魏光艷 葉春明

上海理工大學管理學院,上海,200093

0 引言

在經濟全球化的背景下,為應對快速變化的市場環境,大量制造企業由單一工廠生產轉變為多工廠生產,因此, 分布式工廠調度問題逐漸成為工業界和學術界的關注焦點。這種對多個制造單元進行調度決策的需求,催生了分布式調度的相關研究[1]。借助于智能生產系統,分布式制造可以實現對資源的統一配置,促進分散制造資源的組合優化與共享[2],因此,研究分布式車間調度問題具有重要的現實意義和應用價值。

分布式作業車間調度問題是作業車間調度問題的擴展,通常包括兩個子調度問題:工件在加工工廠之間的分配,以及對工序的操作順序進行排序[3-5]。由于柔性制造系統在分布式制造工廠中的普遍應用,分布式柔性作業車間調度問題(distributed and flexible job-shop scheduling problem, DFJSP)引起了廣泛關注。MENG等[6]提出了四種混合規劃模型和一種約束規劃模型用于求解DFJSP,以最小化最大完工時間為優化目標,在不同規模的DFJSP上進行了有效性驗證。吳秀麗等[7]構建了以總成本和提前/延期為優化目標的雙目標DFJSP模型,并利用改進的差分進化算法進行求解。LI等[8]以最小化最大完工時間、最大工作負載、總負載和提前/延誤為優化目標,構建了多目標DFJSP模型,并設計了混合禁忌搜索算法用于求解模型。其他一些學者在DFJSP的基礎上,考慮了一些更加復雜的約束,以使模型更適合實際生產環境,例如王凌等[9]、DU等[10]、張洪亮等[11]研究了考慮運輸時間的分布式綠色柔性作業車間的調度問題,并分別利用協同群智能算法、混合分布估計算法、改進的遺傳算法對模型進行求解;XU等[12]考慮部分工件的操作需要外包的情況,構建了多目標DFJSP模型,并將遺傳算法和禁忌搜索算法相結合對模型進行求解;LUO等[13]同時考慮了運輸時間和允許工件操作外包的情況,構建了多目標DFJSP模型,并提出了一種模因算法用于求解模型;郭晨等[14]以最小化生產成本和拖期時間為目標,建立了分布式裝配柔性模糊車間調度模型,并利用混合分布估計算法進行求解。還有一些學者結合實際生產案例對DFJSP進行了研究,例如TANG等[15]以砂型鑄造作業車間為例,構建了以最小化最大完工時間為目標的DFJSP模型,通過試驗證明了所提出的混合教學優化算法的有效性;SANG等[16]基于分布式智能工廠生產系統,構建了多目標DFJSP模型,并提出了一種模因算法用于求解模型,最后利用機床部件生產廠的生產信息進行了試驗。

以上學者針對DFJSP以及帶運輸時間或允許工序外包的DFJSP模型進行了研究,提高了模型的實用性。文獻[15-16]則分別基于砂型鑄造作業車間和機床部件生產廠的調度需求,對具體案例進行了建模和求解。然而,該類研究僅考慮了調度過程中的機器柔性,而在集裝配和加工一體化的分布式多柔性裝配作業車間調度問題(distributed and multi-flexible assembly job-shop scheduling problem, DMFAJSP)中,還需要考慮工人安排柔性和工序順序柔性[17-19]。

目前已有一些研究,在調度過程中考慮了工序順序柔性[20]和工人安排柔性。例如,胡瑞淇等[21]基于表達式樹結構建立了考慮工序順序柔性的作業車間調度問題(flexible job-shop scheduling problem, FJSP)模型,設計了一種工序順序的隨機生成方式,并以遺傳算法為例設計了求解流程;ZHANG等[22]研究了考慮工藝順序柔性的FJSP,構建了帶裝配的FJSP模型,并利用分布式蟻群算法進行求解;董海等[23]構建了具有機器、工人和工序多柔性的FJSP模型,并設計了離散回溯搜索算法進行求解;VITAL-SOTO等[24]利用有向無環圖表示工序的操作優先級,建立了具有排序柔性的雙資源約束FJSP模型,提出了帶有創新算子的精英非支配排序遺傳算法用于求解模型。

這些研究提供了工序順序柔性的表示方式,為構建DMFAJSP模型提供了參考。然而,上述考慮工序順序柔性和工人安排柔性的調度研究都限定在單個柔性作業車間范圍內,缺乏對分布式多柔性制造系統的調度研究。因此,本文將對考慮工序順序柔性、工人安排柔性和機器選擇柔性的DMFAJSP進行研究。

DMFAJSP需要對四個子調度問題進行決策,包括將工件分配給工廠、為工序執行順序排序、為工序安排操作工人和為工序選擇機器。由于DMFAJSP的調度過程和可行解的結構較為復雜,本文將采用模因算法(memetic algorithm, MA)進行求解,該算法可以同時進行全局搜索和局部搜索。

目前已有一些學者基于遺傳算法如NSGA-Ⅱ[13,25-26]、NSGA-Ⅲ[16]等,結合各類局部搜索算子設計了模因算法,并通過大量試驗證明了模因算法在求解DFJSP問題上具有良好的性能。與遺傳算法通過交叉變異等操作產生種群的方法不同,分布估計算法(estimation of distribution algorithm, EDA)通過空間采樣和統計學習,預測未來的搜索空間[27-28],從而利用概率模型解決多維調度問題。EDA算法具有顯著的全局搜索優勢,但同時存在局部搜索能力較弱和易收斂早熟的缺點[14]。而模因算法通過將基于種群的全局搜索和基于個體的局部搜索相結合,可以有效提高算法的尋優能力[29]。

因此,本文嘗試根據DMFAJSP的特點,提出一種以分布估計算法作為全局搜索組件,以一些鄰域搜索算子作為局部搜索組件的多維模因算法(multi-dimensional memetic algorithm, MDMA)。通過與其他算法進行對比試驗,驗證所提算法在求解DMFAJSP模型上的有效性。

1 DMFAJSP問題及模型

1.1 具有工序順序柔性的作業車間調度問題

工件的工序順序柔性分為圖1所示的5種柔性等級[17],可以看出:當圖1c中每個優先級中只有一道工序時,與圖1a情況相同;當圖1c中只存在一種優先級時,與圖1b情況相同;而圖1d和圖1e則可以看作圖1c的組合與變形。因此本文以圖1c工序順序非完全確定性工件為代表,研究考慮該類工序順序柔性的DMFAJSP。

(a)工序順序完全確定

1.2 分布式多柔性裝配作業車間調度問題

在分布式多柔性裝配作業車間的環境中,裝配和加工過程被統一作為工序進行調度。DMFAJSP具有地理分布式、多柔性以及集加工和裝配一體化的特點,除了上述提到的工序順序柔性外,本文還考慮了工人安排柔性和機器選擇柔性。該問題中有F個柔性工廠U,每個工廠Uf(f=1,2,…,F)中有Wf名工人,每名工人均會操作多臺機器,Uf中的機器數量為Mf,每一臺機器均具有多種加工功能。假設該問題中有N個工件J需要處理,工件Jn的工序有Rn個優先級,工件Jn中優先級為r的工序O有Knr個。其中,同一工件的工序若具有同樣的優先級,則工序加工無次序要求;否則,低優先級的工序需在所有高優先級的工序加工完成后才可以進行加工。在DMFAJSP中,工序由不同的工人和機器組合進行加工時,加工時間可能不同。

1.3 相關假設

本文以最小化最大完工時間和總能耗為優化目標,相關的假設如下:

(1)每個工件只能由一個工廠進行加工;

(2)屬于同一工件的工序,加工順序受優先級約束,不同工件間的工序無優先級約束;

(3)每個工人均會操作多臺機器,每臺機器均可以加工多道工序;

(4)一個工人同一時刻只能操作一臺機器,一臺機器同一時刻只能加工一道工序;

(5)工序加工不可中斷,直至該工序加工結束。

1.4 符號定義

DMFAJSP模型的數學描述中使用的符號及其含義如表1所示。

表1 符號表

1.5 數學模型

目標函數:

minMK=max(Tn) ?n

(1)

minTE=EP+EI

(2)

其中,

(3)

?f,w,m,n

(4)

?w,n,r,k

(5)

?w,n,r,k

式(1)和式(2)為目標函數,式(1)表示最小化最大完工時間,其中Tn的計算方式如式(3)所示;式(2)表示最小化總能耗,為所有機器的加工能耗和空轉能耗的總和,所有機器的總加工能耗的計算如式(4)所示,總空轉能耗的計算如式(5)所示。

約束條件:

(6)

(7)

(8)

(9)

(10)

Efwmnrk=Sfwmnrk+Pfwmnrk?f,w,m,n,r,k

(11)

Sfwmnrk≥max(Efwmn(r-1)k|k=1,2,…,knr)

(12)

?f,w,m,n,r∈[2,Rn]

Lwm′-Lwm≥V·(Gwmm′-1) ?w,m,m′

(13)

Efwmn′r′k′-Efwmnrk≥Pfwmn′r′k′+V·(Dnrkn′r′k′-1)

(14)

?f,w,m,n,r,k

(Sfwmnrk-Lwm)·Cfwmnrk≥0

(15)

(16)

(17)

約束條件中式(6)和式(7)表示每個工件只能分配給一個工廠,每個工廠可以分配到0至N個工件;式(8)表示工人為多能工,均會操作多臺機器;式(9)表示每臺機器都可以處理多道工序;式(10)保證每道工序只由一名工人操作一臺機器加工;式(11)保證工序加工不中斷,直到該工序結束;式(12)表示屬于同一工件的工序,需在所有高優先級的工序加工完成后,才可以進行低優先級工序的加工;式(13)和式(14)表示同一時刻,每名工人只能操作一臺機器,每臺機器只能加工一道工序;式(15)表示每道工序只有在約定工人空閑時才可以開始加工;式(16)表示每名工人只能操作其會操作的機器且機器只能加工其可以加工的工序;式(17)保證所有工件的工序總數為所有工件中各優先級所包含的工序數量的總和。

1.6 問題案例

為了更好地理解DMFAJSP,以2個柔性工廠的生產環境為例{U1:[[w1,w2],[m1,m2,m3]],U2:[[w1,w2],[m1,m2]]},每個工廠各有2名操作工人以及分別有3臺和2臺機器,有3個工件{J1:[[2],[2],[1]],J2:[[1],[3],[1]],J3:[[1],[2],[2],[1]]}待分配,其中工件J1的工序具有3個優先級,每個優先級中分別有2道、2道和1道工序,其他工件的優先級和工序數量的表示方法與J1類似。工廠U1中工人W1所對應的Q11矩陣如表2所示,其他工人對應的Q矩陣與Q11類似。該案例在U1和U2中的一種調度安排結果如圖2所示。

表2 案例問題中U1的Q表

(a)U1

2 應用MDMA算法求解DMFAJSP

本文針對DMFAJSP可行解的多維特點構建了MDMA算法,該算法綜合運用全局搜索和局部搜索進行尋優。其中全局搜索組件以分布估計算法為核心,局部搜索組件包括3種鄰域搜索算子。

2.1 編碼與解碼

由于DMFAJSP問題可以根據需求將其分解為四個子問題:為工件分配工廠,為工序加工順序排序,為工序分配加工工人以及為工序選擇加工機器,因此,本文設計了多維編碼方案來表示DMFAJSP模型的可行解。如圖3中的編碼圖所示,每個解包含4個變量,各變量的釋義如下:

圖3 編碼與解碼

UJ表示工件-工廠的分配方案,該變量中的數字代表工廠編號,位置索引表示工件編號;

WS表示工人序列,該變量中的數字代表工人編號,工人編號的排列順序與OS對應,即WS某位置上的工人將被安排加工OS對應位置上的工序;

MS表示機器序列,該變量中的數字代表機器編號,機器編號的排列順序與OS和WS對應,即MS某位置上的機器被安排由WS中對應位置上的工人操作,一起加工OS中對應位置的工序;

圖3中的甘特圖展示了編碼示意圖中各變量彩色代碼的解碼結果,即工件J1的各工序在U1中被工廠內的工人和機器加工的調度安排。

2.2 概率模型

EDA算法中的概率模型用于描述候選解的分布信息,在迭代過程中通過對精英解和前代概率矩陣的統計學習完成概率矩陣的更新,從而實現解的進化。本文根據DMFAJSP的多維編碼方案,構建了概率模型UJM、OSM和WMM對應的概率矩陣MUJM、MOSM、MWMM用于抽樣生成種群。

2.2.1概率模型的構建

MUJM是工件-工廠分配概率矩陣,其元素unf表示第n個工件Jn分配到第f個工廠Uf進行加工的概率。MUJM在第0代的初始狀態為

(18)

第g次迭代過程中的更新方式為

(19)

(20)

?n,f,z∈[1,μ·P]

其中,α為學習率;μ為每代種群的精英選擇率,本文設置為0.2;P為種群。

MOSM是工序序列概率矩陣,其元素sni表示工件Jn的工序出現在第i個位置上的概率,其中i∈[1,K]。MOSM的初始狀態為

(21)

在迭代過程中的更新方式為

(22)

?n,i

(23)

?n,i,z∈[1,P]

其中,β為學習率。

(24)

?f,w,m,n,k

在迭代過程中的更新方式為

(25)

?f,w,m,n,k

yfwmnk(g)=(1-γ)·yfwmnk(g-1)+γ·ρ

(26)

?f,w,m,i

(27)

(28)

其中,γ為學習率。

2.2.2概率模型的通用性證明

為保證概率模型的通用性,需確保以下條件成立:

本文利用數學歸納法給出與條件(1)、(2)和(3)相關的模型通用性證明,如引理1、引理2和引理3所示。

證明:當g= 0 時,有

(29)

(30)

證明成立。

證明:當g=0時,有

(31)

(32)

證明:當g= 0時,有

(33)

(34)

(35)

2.3 種群

2.3.1初始化種群

初始化種群采用全局生成方式、局部生成方式和隨機生成方式,按照6∶3∶1的比例[9,13]生成。

全局生成方式:將工件分配給加工其所有工序總用時最少的工廠。

局部生成方式:為工件隨機分配一個工廠,然后在該工廠為每道工序選擇可用的且用時最少的工人和機器組合進行加工。

隨機生成方式:將工件隨機分配給一個工廠,并為每一道工序隨機安排加工機器和操作工人。

2.3.2抽樣生成種群

在迭代過程中,通過采樣2.2節中的概率模型生成UJ、OS、WS和MS變量組成新的個體。為方便操作,分別對UJM、OSM和WMM模型的概率矩陣MUJM、MOSM和MWMM中的元素進行累加,形成輔助矩陣MUJSM、MOSSM和MWMSM。第g代輔助矩陣中的元素u′nf(g)、s′ni(g)和y′fwmnk(g)如下所示:

(36)

(37)

(38)

在一次迭代過程中,利用輔助矩陣生成變量UJ、OS、WS和MS的過程如算法1、算法2和算法3所示。

算法1抽樣生成UJ

輸入:MUJSM

輸出:UJ

01: forn←1 toNdo

02:θ1←random(0,1)

03: forf←1 toFdo

04: ifθ1

05:UJ[n]←f

06: end if

07: end for

08: end for

算法2抽樣生成OS

輸入:MOSSM

輸出:OS,OK

01: forn←1 toNdo

03:θ2←random(0,1)

04: fori←1 toKdo

05: ifθ2

06:OS[i]←n

07:OK[i]←k

08: forj←1 toNdo

09:ossni← 0

10: end for

11: end if

12: end for

13: end for

14: end for

算法3抽樣生成WM

輸入:UJ,OS,OK,MWMSM

輸出:WM

01: fori←1 toKdo

02:n←OS[i]

03:f←UJ[n]

04:k←OK[i]

05:θ3←random(0,1)

06: forw←1 toWfdo

07: form←1 toMfdo

08: ifθ3

09:MS[i]←m

10:WS[i]←w

11: end if

12: break

13: end for

14: end for

15: end for

2.4 局部搜索算子

本文設計了3種鄰域結構NS1、NS2、NS3,以提高種群多樣性和收斂速度。鄰域結構的性能取決于鄰域結構的連通性和規模[30],只有對關鍵路徑上的工序進行操作才有可能減少最大完工時間[31]。因此,本文針對UJ變量設計了鄰域結構Ⅰ,并基于關鍵路徑構建了針對OS、WS和MS變量的鄰域結構Ⅱ和Ⅲ。具體操作如下所示:

(1)鄰域結構Ⅰ。選擇完工時間最晚工廠,將其中最后完工的工件分配給最早完工的工廠。工序序列保持不變,將原工廠中為該工件工序分配的工人和機器編號剔除,并在新工廠中隨機為該工件的工序安排工人和機器。

(2)鄰域結構Ⅱ。如圖4所示,首先,對原個體中的變量對應關系進行存檔。然后,在關鍵路徑(OS變量中填充陰影的工序)中隨機選擇一個關鍵工序,將其移動到頭關鍵工序之前(或尾關鍵工序之后),并為移動后的工序選擇加工時間最小的工人和機器。其他工序的加工工人和機器依照之前存檔的對應信息保持不變。

圖4 鄰域結構Ⅱ示意圖

(3)鄰域結構Ⅲ。首先,對原個體中的變量對應關系存檔。然后,在關鍵路徑上隨機選擇一個與頭關鍵工序(或尾關鍵工序)不屬于同一工件的工序,并將該工序與頭關鍵工序(或尾關鍵工序)進行位置交換。所有工序對應的加工工人和機器與存檔中的對應信息保持一致。

2.5 算法框架

本文設計的MDMA框架如圖5所示。首先,根據2.3.1節中的初始化種群策略生成初始種群。接著,對初始化種群進行局部搜索,根據貪心策略保留初始個體或其鄰域解。然后進行全局搜索,按照2.2.2節中描述的方法通過快速非支配排序和擁擠距離篩選精英解。根據2.3.2節的方法,利用新的概率模型抽樣生成新一代種群。再對新種群循環進行局部搜索和全局搜索,直到達到迭代停止條件,輸出最終的帕累托前沿解。

圖5 MDMA框架

3 試驗與分析

本文試驗采用python3編程語言進行編寫,運行環境為Intel Core i7, 3.6 GHz,4×8 GB RAM。為驗證MDMA在求解DMFAJSP方面的有效性,將MDMA與EDA、NSGA-Ⅱ[32]和HEDA(hybrid estimation of distribution algorithm)[14]進行對比試驗。本文算例為基于BRANDIMARTE[33]提出的MK01-10算例,擴展形成的20個試驗算例,如表3所示。

表3 試驗算例

3.1 參數設置

MDMA需要確定的參數包括:種群規模popsize、精英選擇率μ和3個概率模型的學習率α、β、γ。通過田口試驗[34]確定參數的最優組合,其中每個因素設置4個水平,并采用L16(45)正交表進行試驗設計。測試算例為表3中的算例10,迭代停止時間設置為100 s,每組參數組合獨立運行10次。參數試驗的評價指標選用多目標偏差和(objectives deviation sum, ODS)[13],其計算方法如下所示:

(39)

i=1,2,…,9n=1,2,…,in

其中,in為算法基于參數組合i得到的前沿解總數量;bjn為第n個前沿解第j個目標的值。參數組合及10次運行的平均ODS如表4所示,參數水平趨勢如圖6a所示,根據參數水平趨勢確定MDMA的最優參數組合為:α=0.5,β=0.3,γ=0.4,μ=0.2,popsize=150。

表4 MDMA算法參數組合及平均ODS

(a)MDMA

此外,為保證對比算法的公平性,對比算法的最優參數組合同樣采用田口試驗設計方法進行設置。對比算法的參數及水平如表5所示,其中NSGA-Ⅱ的參數cr為交換率,mu為變異率;HEDA中的f為縮放因子。基于對比算法的參數數量,分別采用L16(45)、L9(33)、L16(44)正交表進行參數組合試驗。各算法的參數水平趨勢如圖6所示,EDA的最優參數組合為:α=0.5,β=0.2,γ=0.4,μ=0.2,popsize=50;NSGA-Ⅱ的最優參數組合為:popsize=50,cr=0.6,mu=0.1;HEDA的最優參數組合為:α=0.4,f=0.2,cr=0.2,popsize=150。

表5 對比算法參數水平

3.2 對比試驗及結果分析

為測試所提算法MDMA對求解DMFAJSP模型的性能,本文將MDMA與EDA、NSGA-Ⅱ和HEDA進行對比。其中EDA算法為不包含局部搜索算子的MDMA,因此通過對比MDMA與EDA,可以驗證MDMA中局部搜索算子的有效性。所有算法的初始種群生成方式一致。所有算法在每個算例上獨立運行10次。由于NSGA-Ⅱ是解決車間調度類問題最為經典的算法之一,為保證試驗公平,對于每個算例,所有算法的運行時間統一設定為NSGA-Ⅱ迭代100次所需的時間。

對比試驗的評價指標選擇常見的世代距離(generational distance, GD)、反轉世代距離(inverted generational distance, IGD)和解集覆蓋率(coverage,C)指標,各指標的計算公式如下:

(40)

(41)

(42)

試驗結果如表6～表8和圖7所示。其中,表6展示了對比試驗中每個算法在每個算例上運行10次得到的最小MK、TE的非支配解。對比試驗的平均GD、IGD指標結果和每個算例的運行時間如表7所示,C指標如表8所示。此外,利用箱形圖展示了對比算法C指標的分布情況,如圖7所示。

表6 對比試驗最小單目標值的非支配解

表7 對比試驗的平均GD和平均IGD

表8 對比試驗的C指標

圖7 C指標箱線圖

表6中第1列表示算例編號,第2列表示MDMA算法在每個算例上10次運行結果中得到的最小MK對應的非支配解,第3列表示其得到的最小TE對應的非支配解,后續4～9列分別表示EDA、NSGA-Ⅱ和HEDA獲得的最小MK和TE對應的非支配解。從表6中可以看出,MDMA在75%的算例上取得的MK最小值優于其他算法,且在85%的算例上取得的TE最小值優于其他算法。由表7可以看出,除了算例11,MDMA的GD和 IGD指標均優于EDA,說明MDMA求得的解在收斂性和多樣性方面均優于EDA。結合表8和圖7中C(MDMA,EDA)、C(EDA,MDMA)指標及分散情況,MDMA獲得的非支配解占比均高于EDA,可以驗證MDMA中局部搜索算子可以顯著提高MDMA的解的質量。另外,從表7、表8和圖7中可以看出MDMA的GD、IGD和C指標均優于NSGA-Ⅱ和HEDA,說明MDMA算法所求得的解在收斂性、多樣性以及非支配解占比方面均優于NSGA-Ⅱ和HEDA求得的解。驗證了MDMA算法在求解DMFAJSP問題上具有明顯優勢。

4 結語

本文針對分布式多柔性裝配作業車間的調度問題,構建了以最小化最大完工時間和最小化總能耗為優化目標的數學模型。該模型考慮了調度過程中存在的機器選擇柔性、工人安排柔性和工序順序柔性。為了求解該模型,本文提出了一種多維模因算法(MDMA)。首先,根據調度模型的特點,設計了多維編碼方案,用于表示問題模型的可行解。隨后,在MDMA算法的全局搜索階段,根據編碼方案設計了3種概率模型,用于抽樣生成種群。并且,為了驗證概率模型的通用性,采用數學歸納法進行了推理和證明。在MDMA的局部搜索階段,設計了3種鄰域結構,其中有2種鄰域結構基于關鍵路徑進行設計,以有效提高算法的局部搜索效率。最后,在20個分布式多柔性裝配作業車間調度問題(DMFAJSP)算例上,將所提算法與分布估計算法、遺傳算法和混合分布估計算法進行了對比。結果驗證了所提MDMA算法對于求解DMFAJSP問題具有顯著優勢。

此外,本文在研究分布式多柔性裝配作業車間的調度問題時,僅考慮了靜態的調度過程。后續將對更復雜的動態DMFAJSP模型進行研究。