利用歷史信息和限制算子求解MOFJSP

吉訓生，蔡益青

1.江南大學 物聯網工程學院，江蘇 無錫214122

2.物聯網技術應用教育部工程研究中心，江蘇 無錫214122

1 引言

當下，智能制造日益成為制造業發展的重要趨勢和核心內容。從機器層面上來講，智能制造要求機器具有柔性的生產制造能力，以實現“高效化”“綠色化”的生產目標[1]。柔性作業車間調度問題（Flexible Job-shop Scheduling Problem，FJSP）假設工序可在多臺機器上加工，突破機器唯一性限制，是典型的柔性生產制造問題[2]。相比傳統的作業車間調度問題（Job-shop Scheduling Problem，JSP），FJSP 是復雜度更高的NP-hard 難題，也是智能制造及組合優化領域的研究熱點之一[3]。

為了實現綠色、高效生產，建設綠色工廠，必須要求合理調度柔性作業車間內的生產過程。通常，認為以滿足消費者要求，降低生產成本、減小污染氣體的排放等多個目標的FJSP是一個多目標柔性作業車間調度問題（Multi-Objective Flexible Job-shop Scheduling Problem，

MOFJSP）。MOFJSP需要優化多個性能指標（最小化最大完工時間、最小化生產、最小碳排放量等），而多個性能指標之間往往存在沖突。對于一個解來說，無法使所有性能指標都為最優值。因此，求解MOFJSP的目的是為了得到各性能指標的均衡解，稱為Pareto 最優解[4]。多目標進化算法（Multi-Objective Evolutionary Algorithm，MOEA）因其在求解MOFJSP 時可獲得一致性的Pareto最優解的特點，已成為求解此類問題常用的方法[5]。DEB 等[6]利用解之間Pareto 支配關系和擁擠距離，提出一種經典的MOEA——帶精英策略的非支配排序遺傳算法（Non-dominated Sorted Genetic Algorithm-II，NSGA-II）；張超勇等[7]將NSGA-II 作為求解MOFJSP 的策略，驗證了NSGA-II在解決MOFJSP有良好的性能。Zhang等[8]提出一種基于問題分解的多目標進化算法（Multi-Objective Evolutionary Algorithm based on Decomposition，MOEA/D），此算法的求解速度都遠高于NSGA-II；隨后，Li等[9]將穩定配對思想引入MOEA/D，提出了基于穩定配對選擇策略和分解相結合的多目標進化算法（Multi-Objective Evolutionary Algorithm based on Decomposition and Stable Matching，MOEA/D-STM），成功解決了MOEA/D“丟失解”的缺點；Wang 等[10]在多目標遺傳算法（Multi-Objective Genetic Algorithm，MOGA）的基礎上引入了免疫原則和熵，以保證種群的分布性，從而解決了早熟問題，提高收斂速度。

參考文獻[10]的做法，本文在MOEA/D-STM 基礎上提出一種基于歷史信息和限制算子多目標進化算法（Multi-Objective Evolutionary Algorithm based on Decomposition with Historical information and Limited operators，MOEA/D-HL），并利用MOEA/D-HL 求解MOFJSP。此算法在保證收斂性的同時，提高了種群的多樣性，增大進化強度，增強全局搜索能力，避免早熟。實驗結果表明，本文算法在求解MOFJSP問題上比較可行且有效。

2 MOFJSP模型

MOFJSP 可用集合( M,J,O,T )來描述，其中M 代表機器，J 代表工件，O 代表工序，T 代表加工時間。通常設定有m 臺機器( Mk,k ∈{1 ,2,…,m} )和n 個工件(J={ J1,J2,…,Jn} )，每個工件有多道工序，且工序的先后順序是確定的，每道工序Oij(i=1,2,…,n,j=1,2,…,ni，ni表示工件i 所含的工序數）可由多臺機器加工，所用時間tijk由機器能力唯一確定[11]。一個包括3 個工件、4 臺機器的MOFJSP問題描述如表1所示。

表1 3×4 MOFJSP問題加工時間表

本文將最小化最大完工時間、最小化生產成本和最小化能耗作為優化目標，構建一個MOFJSP 模型[12-13]。求解MOFJSP 即在滿足電力資源約束和工序先后約束的前提下，根據目標最優要求，把工序分配到各機器上并確定工序的加工序列和開始時間[14]，達到高效、綠色生產的目的。

MOFJSP優化目標如下：

式中，bijk為工件i 的第j 道工序在機器k 上的開始加工時刻；tijk為工件i 的第j 道工序在機器k 上加工所用時間；Sijk為決策變量，若工序Oij在機器k 上加工則Sijk=1，否則，Sijk=0。

式中，γi為材料費；νijk為人工費，為工件的加工費用，表示使用第v 種供電方式在單位時間內為正常運行（額定負載）狀態下的機器k 供電所需的成本，其中，v=1,2,3 分別表示使用主電網、分布式能源和燃機三種供電方式為機器提供能源；為機器空載運行費用，tk為機器k 在整個加工過程中的空載運行時間表示使用第v 種供電方式在單位時間內為空載狀態下的機器k 供電所需的成本。

求解MOFJSP滿足如下約束條件：

（1）在同一時段內，一臺機器只能加工一個工件，即在某時刻h( )h ＞0 ，若1 則在i ≠p 或j ≠q 條件下，Spqk≠1 恒成立。

（2）同一個工件，只有在前工序完成后，才可加工下道工序，即：

其中，Cij表示工序Oij的完成時刻。

（3）任何一道工序一旦開始加工，就不可打斷，即：

除了滿足上述約束條件，還應做如下假設：

（1）某工件當前工序加工完成后，應立即運送至加工下一道工序所用機器前的緩沖器，若機器空閑應立即開始加工，否則在緩沖器中等待加工。

（2）工件在機器間的傳送時間、在機器上的加工時間、裝卸時間、調整時間統一視為加工時間。

（3）在加工開始前，所有機器處于正常待命運行狀態，且所有工件同優先級。

3 求解MOFJSP的MOEA/D-HL

利用MOEA/D-STM 求解MOFJSP 時，在迭代過程中使用了無限制穩定配對選擇策略，通過該策略選擇的種群收斂性優于多樣性，這將降低最終種群的性能[15]。另外，在交叉操作時只使用了父代信息，這會降低收斂速度。由于以上原因，本文選擇有限制穩定配對策略作為選擇策略，以平衡迭代過程中被選種群的收斂性與多樣性；另外設置一個外部存儲器用來存儲每個個體的歷史信息。歷史信息表示個體上一次迭代時的狀態。如圖1所示分別表示第g+2 次進化操作時，x2、x4表示當前狀態，則分別表示x2、x4的歷史信息；同理，分別表示x2、x4在g+1 次進化操作時的歷史信息。在交叉操作時，借助個體的歷史信息以達到加速收斂的目的。利用MOEA/D-HL 求解MOFJSP的過程如圖2所示。

圖1 個體歷史信息示意圖

3.1 染色體編碼與解碼

圖2 MOEA/D-HL算法流程圖

編碼是算法與MOFJSP 之間的橋梁[16]。本文選擇基于工件工序和機器的雙層整數編碼，第一層為工序編碼，第二層為機器編碼。若有m 個待加工工件，每個工件有ni道工序，則染色體的每層編碼長度為。具體步驟為：（1）在[0，1]之間，隨機初始化長度為∑ni的實數數組；（2）根據每個位置上實數值的排名，確定其代表的工序；（3）通過工序編碼的具體數值，確定所用的機器，并得到機器編碼。

如圖3所示，染色體表示3個待加工工件，每個工件有2～3 道工序，有3 臺機器可用。隨機初始化一組實數數組，通過比較每個基因位上的值確定排名，排名前3位的代表工件1 的三道工序，排名4、5 的基因位代表工件2 的兩道工序，排名6、7 的基因位代表工件3 的兩道工序，因此工序編碼為[1 3 2 1 2 3 1]，表示工序[O11,O31,O21,O12,O22,O32,O13]。在確定工序編碼后，從每道工序可用的加工庫中，選擇機器，并確定機器編碼[2 1 3 2 2 2 1]，表示加工每道工序所用的機器分別為[M2,M1,M3,M2,M2,M2,M1]。

圖3 染色體編碼

解碼的具體步驟為：首先根據工序編碼確定工件及工序，然后從機器編碼中找出對應的加工機器，并從時間矩陣中找到機器加工相應工序所需時間。通過上述操作，一條合法的染色體可通過解碼操作得到一個可行、完整的調度方案。

3.2 進化操作

進化操作包括交叉和變異操作。交叉操作決定算法的全局搜索能力，是算法的關鍵[17]。如圖1所述，本文交叉操作使用的種群歷史信息存放在外部存儲器中。算法開始時，初始化父代種群，同時將父代種群個體存放在外部存儲器中，以初始化個體的歷史信息。隨后進行交叉操作，步驟如下：首先，從父代種群中選擇三個染色體e1、e2、e3，并從外部存儲器中選擇染色體e1的歷史信息h1；其次，選擇從左向右依次在工序編碼的基因位處隨機產生r ∈(0,1]，若r ＜pc（pc為交叉概率），在父代染色體上選擇此基因位為交叉點；最后，將兩交叉點之間的基因進行交叉操作，并根據工序編碼調整機器編碼。如圖4所示，選擇a、b 為交叉點，對a、b 之間的基因位進行交叉操作：yr=e1+0.5×( e2-e3)+0.5( e1-h1)，生成子代染色體y 后，通過排序操作確定工序編碼，并根據父代和子代信息確定機器編碼。經上述交叉操作得到的子代染色體滿足約束條件，是可行解。

圖4 交叉操作示意圖

變異操作有利提高算法的局部優化能力和種群多樣性，防止算法陷入局部收斂。本文選用變異操作過程如圖5所示：在工序編碼里隨機選擇基因a、基因b 進行交換，此交換將改變工序編碼中每個基因所代表的工序，改變工序編碼中每個基因所代表的工序，如c 位置的基因1 在原染色體表示工序O13，在新染色體變為表示O12，d 位置的基因3在原染色體表示工序O31，在新染色體變為表示O32，同時機器編碼也隨之變化。經此變異操作得到的新染色體能夠滿足約束條件。

圖5 變異操作示意圖

當一次進化操作結束后，將父代種群個體按順序送入外部存儲器，更新外部存儲器所含歷史信息，為下一次交叉操作提供歷史信息。

3.3 選擇操作

MOEA/D-STM 因其在構造偏好矩陣時，只用了解與子問題之間的歐氏距離和切比雪夫距離這兩類信息，造成選擇到的解多樣性較差。如圖6（a）所示，解擁擠在p1、p2之間，在下一次交叉和變異操作時，由于多樣性的不足，使得生成的新染色體改變較小，一方面使收斂速度降低，另一方面容易使尋優過程陷入局部最優。為了提高父代種群的多樣性，加速收斂，本文將解的位置信息引入到子問題對解的偏好矩陣的構造過程中，以增加那些離Pareto Front 遠但分布均勻的解的被選擇概率。主要操作有：（1）處理目標空間與生成待選解集；（2）獲取解的位置信息，將位置信息代入自適應轉移函數中，以得到限制信息；（3）構造帶限制信息的子問題對解的偏好矩陣，以及解對子問題的偏好矩陣，并根據兩偏好矩陣的信息選擇父代種群。

（1）處理目標空間與生成待選解集

設目標空間為F(s)=[f1(s),f2(s),…,fd(s)]∈Rd,d 為目標空間的維數。使用一組均勻的權向量w={w1,w2,…,，將目標空間劃分為N 個子空間，其中wt=(wt,1,…,wt,j,…,wt,d)∈Rd,wt,j≥0，且1，同時得到對應的子問題P={ p1,…,pl,…,pN}。將待選種群S={s1,…,sN,…,s2N}中的所有染色體通過目標函數值的計算，映射到目標空間中，得到待選解集X={x1,…,xN,…,x2N}，其中，x=[ f1( s),f2( s),…,fd( s )]∈Rd。

（2）獲取限制信息

本文選用解相對于子問題的角度作為位置信息。將d 維目標空間F( s )=[ f1( s),f2( s),…,fd( s )]∈Rd轉化為個二維空間，其中，表示從d 個變量中選出2個所有的組合方法。對于二維空間Fc( s )=[ fu( s),fv( s )],c=1,2,…,,u,v ∈[1 ,2,…,d ]，解x 在此空間中的目標值分別為fu( s),fv( s )，子問題ph對應的權向量wh在此二維空間的分量為wh,uv=( wh,u,wh,v)，則夾角分量為：arctan(| fu( s )-ωh,u|| fv( s )-ωh,v|),θuv∈[0 ,π/2]。角度θ 定義為：子問題ph與解x 在個二維平面上的夾角分量之和。

構造一個自適應轉移函數，并引入位置信息θ，即：

（3）選擇操作

將上述限制信息引入到子問題pr,r=1,2,…,N 對候選解x,x ∈X 的偏好值計算式中，如式（7）。利用式（7）可計算出子問題對所有候選解的偏好值，按升序排列則構成子問題對解的偏好矩陣ψp(N×2N)維 中的一行。

候選解x,x ∈X 對子問題pr,r=1,2,…,N 的偏好值計算式如式（8）。利用式（8）可計算出解對所有子問題的偏好值，按升序排列則構成解對子問題的偏好矩陣ψX( 2N×N維 )中的一行。

其中，‖ · ‖為歐氏距離，F( si)是將F( xi)歸一化處理后的目標函數值，其第l 個目標函數值)/l=1,2,…,d ，其中，為最差點，可通過,l=1,2,…,d 計算。

利用兩個偏好矩陣的信息和延遲接受程序[18]，得到最終配對結果如圖6（c）所 示：} 。很明顯，此種方法選擇的解比不加限制信息的配對策略選擇的解多樣性更好。

4 實驗設計與結果分析

NSGA-II 因其快速排序、快速估算擁擠距離以及容易比較擁擠距離的三個特點，成為CIMS 領域最常用的多目標算法之一。因此，本文將與NSGA-II 比較。NSGA-II和MOEA/D-STM需要設置的參數參照文獻[19]。MOEA/D-HL 參數設置：種群染色體數目N=40，迭代次數K=400。為驗證交叉概率pc，變異概率pm，限制算子控制參數L，鄰域大小T 對生成調度方案的影響，本文進行了各相關參數的靈敏度分析。

4.1 靈敏度分析

采集某制桶公司實際生產車間的原始數據，取整處理后的數據如表2、3。車間設備名稱與編號列于表2，6種桶，每種桶加工所需的20工序、可用的加工機械及所用加工時間列于表3。

以制桶車間實際訂單為例，設置不同的實驗參數，以三個目標值為參考，所得實驗結果如表4所示。

從表5 可知，三種算法對pc、pm的值不敏感。另一方面，MOEA/D-HL中限制算子的控制參數對實驗結果有較大的影響。因此，為了得到最佳實驗結果，本文設置交叉概率pc=0.8，變異概率pm=0.6，限制算子控制參數L=1，鄰域大小T=10。

4.2 實例仿真及結果分析

用3種算法求解該調度問題，不考慮工件在機器間轉移的時間，實驗參數不變。收斂性能見表5，MOEA/D-HL 所需最大完工時間、生產成本和C 排放量最優解分別為4 105 s、4 223.6元和39.5 g。從加工實時性上考慮，MOEA/D-HL性能優于NSGA-II、MOEA/D-STM，分別提高了0.8%，1.4%；在生產成本控制上，MOEA/D-HL較NSGA-II、MOEA/D-STM 提高了0.8%和1.8%；碳排放量是衡量綠色生產車間的一個重要指標，通過MOEA/D-HL 得到的調度方案對環境更加友好，較NSGA-II、MOEA/D-STM減少了2.5%和4.8%的碳排放量。另外，如圖7 所示，MOEA/D-HL 收斂速度在整個進化過程中優于另外兩種算法，也即MOEA/D-HL能更快地得到調度方案。

圖6 MOEA/D-STM和MOEA/D-HL兩者對比

表2 設備及編號

表3 工序與加工設備

表4 參數靈敏度分析

表5 實例仿真結果

圖7 最大完工時間平均收斂曲線

5 結論

本文利用歷史信息和帶限制算子的選擇策略，在保證良好收斂性的同時，加快了收斂速度。在生成子代的過程中，利用個體的歷史信息，生成較大的進化量，在選擇過程中，本文利用解在目標空間中的位置信息得到自適應限制信息，并將其加入到偏好矩陣的構造過程中，以獲得分布性能好的父代種群，從而加快收斂速度。以某企業生產訂單為例，驗證MOEA/D-HL在解決實際生產車間調度問題時，收斂性和收斂速度優于MOEA/DSTM和NSGA-II，能夠快速地獲得收斂性能好的調度方案指導實際生產。