基于多目標螢火蟲膜算法的學習效應生產調度問題

侯豐龍，葉春明，耿秀麗

（上海理工大學 管理學院，上海 200093）

隨著工業技術的發展和大數據環境的沖擊，生產調度作為生產環節的關鍵問題被再次廣泛地研究。生產調度問題主要包括單機調度問題（Single Machine Scheduling，SMS）、流水車間調度問題（Flow-shop Scheduling Problem，FSP）、作業車間調度問題（Job-shop Scheduling Problem，JSP）以及開放車間調度問題（Open-shop Scheduling Problem，OSP）等。在實際生產活動中，置換流水車間調度問題（Permutation Flow-shop Scheduling Problem，PFSP）是FSP問題中重要的一類，主要應用于以單件大批量生產為主的制造企業中，鑒于其實際應用，多目標置換流水車間（MPFSP）具有重要的研究意義。Lee等[1]采用遺傳算法，以最小化加權延遲時間和最大完工時間為目標研究流水車間調度問題。Lemesre等[2]提出精確并行方法求解雙目標置換流水車間調度問題，雙目標為最小化最大完工時間和總延遲時間。Arroyo等[3]以最小化最大完工時間、最大延遲時間和總流程時間為目標，研究置換流水車間調度問題，并提出了基于貪婪隨機自適應元啟發式算法。

當前學者對于MPFSP問題的研究主要集中在求解方法和多目標設定上，很少見到將操作者的不同學習效應應用到MPFSP研究中。在實際制造過程中，由于工人（機器）長時間加工相同或類似的工件，加工效率不斷提高，使得繼續加工的工件的加工時間逐漸減少，這種現象就稱為學習效應[4]。Biskup[5]首次把學習效應的概念引入生產調度領域，通過對單機問題的學習效應的研究認為，工件的實際加工時間與其在機器上的加工位置有關，并證明了引入學習效應下，極小化最大完工時間的單機問題多項式可解。Mosheiv等[6-7]研究了具有普遍性的工件學習曲線調度模型，證明了以極小化最大完工時間和總的加工時間為目標函數的單機調度問題，以及以到期時間安排和極小化總的加工時間為目標函數的異速并行機問題是多項式可解。Chen等[8]研究了具有學習效應的雙目標Flow-shop調度問題。

考慮到實際生產活動中，加工過程受環境的影響，本文把學習效應的因素引入到MPFSP問題的研究中，而當機器數量m＞2時，PFSP問題已是NP-難問題[9]。因此，具有學習效應的MPFSP問題也是NP-難問題。由于精確方法和啟發式算法的局限性很難有效處理此類問題，故本文提出一種新型的螢火蟲膜算法（FMA）。膜計算是Pǎun[10]根據生物細胞的結構抽象出的一種仿生計算模型，又被稱為P系統。Nishida[11-12]提出了具有嵌套結構的膜算法，并將其用來求解NP-完全問題，取得了比模擬退火算法更好的效果。本文把新型智能螢火蟲算法（FA）嵌入膜計算中，嘗試改進局部搜索策略。設計了求解MPFSP的新型螢火蟲膜算法，并與螢火蟲算法和粒子群算法進行對比，以驗證螢火蟲膜算法的性能。

1 基于學習效應的MPFSP模型建立

1.1 學習曲線描述

自從Wright[4]提出學習效應以來，一些經典的學習曲線模型相繼應用在生產調度領域，例如，對數-線性模型、S型模型、高原模型以及Dejong模型等，學習曲線模型如圖1所示。對數-線性學習模型是最常見也是應用最廣泛的模型，模型中的數據相關性可以通過統計獲得，它暗示無論最終的加工數量是多少，單元加工時間都會持續減少；高原模型認為，經過特定的工件加工數量后，學習效應會達到一定水平，并保持不變，想要改變這種水平必須經過大量的加工處理過程；S型學習模型近似于波音“駝峰”模型，模型假定學習經驗可以從一批（工件）傳遞到下一批（工件），不同之處在于，波音“駝峰”模型假定經過一定的初始加工量之后，模型曲線會恢復到線性或高原模型，而S型模型認為學習效應會以遞減的學習比率沿曲線運動，并逐步達到穩定值；Dejong模型是在高原模型和對數-線性模型之間提出一個折中的模型，它認為學習效應持續存在，只是在持續的加工學習過程中以一定的下降比率逐步達到穩定值。

圖1 不同學習效應模型曲線示意圖

Biskup[5]把對數-線性模型應用于排序問題，提出了基于位置的學習效應模型，即

式中：α≤0是學習效應因子；P j、Pjr分別為工件的基本加工時間和實際加工時間；r為工件的加工位置。

1.2 基于學習效應的MPFSP模型建立

設有n個工件在m臺機器上加工；工件的待加工序列為π（r1，r2，…，rn），其中，ri為排序中第i個位置上待加工工件對應的工件編號，工件ri在機器j上的正常加工時間為P（ri，j），i=1，2，…，n；j=1，2，…，m，若每個機器有相同的學習率，則根據學習曲線模型[5]，其實際加工時間為

式中：a為學習效應因子，且a≤0。

最大完工時間Cmax表示n個工件中最后加工的工件在最后一臺機器上完成加工的時間，最小的最大完工時間通常意味著機器的利用率高和產品的生產周期短。總流經時間是指所有工件的完成時間之和，記為

其中：Ci，m為工件i在機器m上的完工時間；ri為工件i的第1道工序可執行的時刻，稱為釋放時間或準備時間。對于MPFSP，假設ri=0，i=1，2，…，n，其總流經時間可簡化為

本文選取最大完工時間和總流經時間兩個指標作為多目標優化的目標函數，用三參數法可表示為Fm|prum|Cmax，Ftotal，則基于學習效應的置換流水車間調度問題的數學模型為：

2 基本算法描述

2.1 膜計算描述

膜計算是由Pǎun[10]受到生物細胞啟發提出的一種仿生計算模型。膜系統如圖2所示，膜1稱為表層膜，膜2、4、5為基本膜。膜計算優化方法基本思想是通過采用轉運型P系統，在各個基本膜內分別單獨采用算法進行尋優；然后，利用P系統的交流規則將各個基本膜中的最優個體輸送至表層膜中以引導算法進化，并將尋優結果返回各個基本膜內進行下一次尋優。

圖2 膜系統示意圖

膜計算模型用數學形式表示為

式中：O為非空有限集；μ為膜結構，由m個膜組成；w x為膜區域x中的多重集；Rx為區域x中擁有的規則集；i0=1，2，…，m表示輸出結果的膜區域，當i0=0時，表示結果輸出到外界環境中去。

膜計算方法能夠有效地避免陷入局部最優，提高優化算法搜索最優解的精度。

2.2 基本螢火蟲算法描述

自然界中多數種類的螢火蟲都有發光的特性，基本螢火蟲算法即是通過模擬自然界中螢火蟲的發光行為構建的智能優化算法，其假設：

（1）算法中的螢火蟲沒有性別之分，每一只螢火蟲都可被其他螢火蟲個體吸引。

（2）螢火蟲自身的光亮程度與對其他螢火蟲的吸引程度成正比例關系，并且，隨著距離增大其吸引程度減少。

（3）螢火蟲亮度與優化目標函數的性質相關[13]。在基本螢火蟲算法中，搜索和優化過程模擬的是螢火蟲個體的吸引和移動過程，用優化目標函數值的大小來衡量螢火蟲個體位置的優劣，個體位置的優勝劣汰過程也即求解算法過程中可行解的巡游過程。

利用螢火蟲算法求解過程中的關鍵參數定義[14-16]：

定義1螢火蟲的相對亮度計算公式：

式中：rij為螢火蟲i、j間的笛卡爾距離；γ為光強吸收系數，用來反映亮度在傳播過程中受到外界干擾而逐漸減弱的現象，在此設為常數；I0為r=0時的最大熒光亮度，與目標函數值成正比。

定義2螢火蟲個體的吸引度計算公式：

式中：β0為r=0時的最大吸引度，也即螢火蟲自身光源處的吸引度；γ、rij同式（10）。

定義3低亮度的螢火蟲i被高亮度的螢火蟲j吸引移動后的位置更新公式為

式中，x i、x j分別為螢火蟲個體i、j所處的空間位置；α（0≤α≤1）為步長因子；rand是在[0，1]上服從均勻分布的隨機數；α（rand-1/2）的設置是為了避免過早陷入局部最優而對螢火蟲個體搜索進行隨機擾動。

3 螢火蟲膜算法描述

3.1 編碼規則

針對MPFSP問題，利用基于ROV規則的隨機鍵編碼 方 式[17]，用 向 量x i=（x i，1，x i，2，…，x i，n）表示螢火蟲個體的位置，初始位置隨機產生。螢火蟲位置向量的每一個分量表示一個工件，對向量中的分量進行升序排序，通過對每個分量的位置與下標進行映射，將螢火蟲的位置向量和機器上各工件的加工序列π（r1，r2，…，rn）聯系起來。

3.2 多目標Pareto最優解集的構建

本文采用排除法構造目標函數Pareto最優解集，將螢火蟲進化群體Pop中的每個螢火蟲個體X依次與非支配集NDSet（初始狀態時為空）中的螢火蟲個體Y比較（當NDSet為空時，直接將個體X存入到NDSet中），如果X支配Y，說明Y是支配個體，將Y從NDSet中刪除（即排除掉）。如果X不被NDSet中任何一個螢火蟲個體所支配，則將X并入到NDSet中。

3.3 聚集距離

聚集距離的概念最早由Deb等[18]提出。聚集距離大的個體即其周圍解空間的密度小，該個體在選擇過程中具有更大的被選概率。一個個體的聚集距離即為該個體周圍與之相鄰的個體在每個子目標函數值上的距離差之和。如圖3所示，設有兩個子目標f1和f2，個體i的聚集距離是圖3中四邊形（虛線）的長、寬之和。個體i的聚集距離記作P[i]distance，個體i對應的子目標m的函數值為P[i].m，則圖3中個體i的聚集距離

在優化問題中，當有r個子目標時，個體i的聚集距離

邊界點的聚集距離設為無限大。

圖3 個體之間的聚集距離

3.4 局部搜索算法的改進

為了增加算法的尋優精度，在更新螢火蟲位置過程中增加互換和變異的局部搜索操作，并借鑒模擬退火算法[19]的思想，結合VNS鄰域搜索策略中的Insert和Inverse兩種鄰域結構，兩種鄰域結構如圖4所示，提出一種新型的局部搜索策略。

圖4 Insert和Inverse鄰域結構圖

具體過程：

（2）選定一個工件排序π，通過對工件排序進行互換、插入和逆序的鄰域搜索操作，得到新解π*。

（3）如果（Fπ*-Fπ）＜0，則用π*替換π，并更新螢火蟲位置；如果exp[（Fπ*-Fπ）/T k]＞rand（0，1），T k為第k次迭代的溫度，則用π*替換π，否則保留π不變。

（4）退溫操作，Tk+1=λ·T k，k=k+1。

（5）如果滿足T k+1＜Tend，則算法結束；否則，轉（2）。

3.5 膜內多目標螢火蟲算法的設計

基于膜計算優化的螢火蟲算法即螢火蟲膜算法是一種利用P系統設計的智能優化算法，計算過程如圖5所示，利用P系統在各個基本膜內分別單獨進行螢火蟲算法（FA）獲得局部最優；然后，利用交流規則將結果輸出表層膜中獲得當前全局最優，并把當前全局最優解分別返回各個基本膜內再次進化尋優，以此迭代，直到滿足終止條件。

圖5 螢火蟲膜算法圖解

對于多目標尋優問題，對膜內螢火蟲算法的設計：

（1）初始化螢火蟲種群。

（2）初始化非支配解集NDSet，從螢火蟲的初始位置中找出所有的非受支配解添加到NDSet中。

（3）循環迭代。

3.6 多目標螢火蟲膜算法的總流程

多目標螢火蟲膜算法求解過程：

（1）初始化一個由1個表層膜和m個基本膜組成的膜結構，即[[]1，[]2，…，[]m]0；

（2）初始化一個螢火蟲種群n，并將這n個螢火蟲平均分配到m個基本膜中；

（3）初始化各個基本膜內螢火蟲算法的各項參數；

（4）在m個基本膜中分別獨立運行螢火蟲算法，通過局部搜索策略尋求更優解，并將得到的非劣解集通過交流規則輸出到表層膜內；

（5）對表層膜中的所有解個體進行篩選，用排除法獲得新的非劣解集，若非劣解集中解的個數達到設定值，則利用聚集距離對其進行削減；

（6）當滿足運行次數時停止運行，輸出非劣解集；否則，更新非劣解集，返回執行（4）。

整體框架如圖6所示。

4 仿真實驗與分析

4.1 仿真測試

為了測試螢火蟲膜算法的性能，本文選取Car類基準測試問題[20]進行優化分析，并將測試結果與變鄰域搜索改進后的螢火蟲算法（FA）[21]以及粒子群算法（PSO）[22]的測試結果進行對比。仿真實驗環境為Windows7操作系統，CPU為Intel（R）Core（TM）i7-3770，內 存 為3.40 GB，編 程 環 境 為MATLAB R2 010 a。螢火蟲膜算法參數設置：基本膜個數m=4，每個基本膜中螢火蟲數目k=10，最大吸引度β0=1，光強吸收強度γ=1，步長因子α=0.3。螢火蟲算法中參數設置：螢火蟲數目k=40，最大吸引度β0=1，光強吸收強度γ=1，步長因子α=0.3。粒子群算法中粒子數為40，學習因子C1=C2=1.494 45，采用Wstart=0.9，Wend=0.4的線性遞減慣性權重[23]，3種算法最大迭代次數maxgen=300，各獨立運行20次。

圖6 多目標螢火蟲膜算法的整體框架

測試結果如表1所示，其中：C*為理論最優解；SR為尋優成功率；BRE為最優相對誤差；ARE為平均相對誤差；WRE為最差相對誤差。RE為所得結果C相對于C*的百分比誤差

由表1可見，螢火蟲算法無論在尋優成功率還是在相對誤差方面性能均明顯優于粒子群算法，而螢火蟲膜算法除了對Car6和Car8的尋優稍差于螢火蟲算法外，對其余測試問題的尋優好于螢火蟲算法。而且由于膜計算的并行性，使得時間上比螢火蟲算法的用時短很多。圖7所示為3種算法對Car5測試結果對比圖。由圖7可以明確地看出，螢火蟲膜算法的尋優性能優于螢火蟲算法，更明顯優于粒子群算法。因此，螢火蟲膜算法對PFSP有很好的尋優效果。

4.2 基于學習效應的MPFSP調度問題求解與分析

4.2.1 學習因子a=0時的MPFSP問題求解 測試問題、算法運行環境和參數設置同上，得到的Pareto最優解集如表2所示，Car1、Car3問題的Pareto最優解集及Pareto邊界如圖8所示。

表1 測試結果對比圖

圖7 Car5測試結果對比圖

表2 Pareto最優解集

4.2.2 基于不同學習效應因子的MPFSP問題求解與分析 利用多目標螢火蟲膜算法求解基于學習效應的兩目標PFSP問題，當學習率l分別為90%、80%、70%、60%和50%時，學習效應因子a分別 為：-0.152，-0.322，-0.515，-0.737，-1。算法對Car1～Car8測試問題進行測試，得到的部分結果集如表3所示。為了能夠明顯地看出學習效應對最大完工時間和總流程時間的影響情況，將Car3的測試結果繪成圖，如圖9所示。由圖9可見，隨著學習率的降低，兩目標函數值都減小的十分明顯，并且Pareto解集之間的直線距離有減小的趨勢。

圖8 Car1、Car3測試問題的Pareto邊界

圖9 不同學習率下Car3測試結果圖

表3 不同學習率下測試問題的部分Pareto結果集

利用敏感性分析方法對不同學習率下最大完工時間和總流程時間進行分析。敏感度系數為

式中：SFA為目標函數值F對于學習因素A的敏感度；ΔA/A為學習效應的變化率；ΔF/F為學習效應變化ΔA時，目標函數值對應的變化率。在此，對同一學習率下同一個目標所有解取平均值，然后用平均值作為敏感性分析中的評價指標進行計算。最大完工時間和總流程時間的敏感度系數如表4、5所示，敏感因素變化率即學習因素變化率，表4、5中，學習效應從100%分別上升到90%、80%…50%。

由表4可知，對于F1而言，除了在敏感因素變化率為-50%時，Car7對應的敏感度系數0.97小于1外，其余所有Car問題在不同學習率下的敏感性系數均大于1，其中部分敏感度系數大于2。這說明，學習率與最大完工時間同方向變化；學習率每變化1%，最大完工時間將變化1%以上；學習率越低，學習效果越明顯。

由表5可知，除了Car7和Car8的個別敏感度系數小于1外，剩余的Car問題在不同學習率下的敏感度系數均大于1。這說明，學習率對目標函數值F2的影響也很大。將表5中數據與表4中數據進行對比，可以看出，在不同的敏感因素變化率下，所有Car問題F2的敏感度系數均小于F1的敏感度系數。可見，學習效應對最大完工時間的影響稍大于對總流經時間的影響，最大完工時間對學習效應的存在更為敏感。

此外，由表4、5還可以看出，隨著學習率的降低，敏感度系數呈逐漸減小的趨勢，這與從不同學習率下測試問題的測試結果圖（見圖9）中觀察到的Pareto解集之間的直線距離有減小的趨勢這一現象相對應。這說明，學習效應對生產指標的影響會隨著學習率的降低而逐漸減小。

表4 最大完工時間的敏感度系數表

表5 總流程時間的敏感度系數表

5 結語

本文將學習效應模型加入到置換流水車間調度問題中，探索應用嵌套有螢火蟲算法的膜算法求解基于學習效應的多目標置換流水車間調度問題。對螢火蟲膜算法在搜索策略上進行一些改進，提高了算法全局尋優的性能，與螢火蟲算法和粒子群算法進行對比，驗證了螢火蟲膜算法在求解多目標置換流水車間調度問題中的有效性和可行性。探討學習效果對加工指標的影響程度，具有重要的理論價值和應用價值。研究表明，學習效應的存在對PFSP的時間指標影響顯著。在實際生產活動中，環境復雜，除了學習效應還有遺忘效應、不耐煩行為等均對生產調度產生重要的影響，綜合考慮各種效應影響是下一步研究的重點。此外，在調度優化問題方面，設計性能更好，更加合理、高效的智能算法，也是未來研究的重點方向。