李曉輝,高 鐸,楊 晰,劉元東,趙 毅,董 媛
1(長安大學 電子與控制工程學院,西安 710054)
2(陜西航天時代導航設備有限公司,寶雞 721099)
由于機器人制造單元對制造行業的生產效率和生產質量有很大提高,近30年來,生產調度問題[1]已成為國內外自動控制領域、計算機科學等領域的研究熱點.現代智能制造市場的需求越來越個性化、多樣化,使傳統的量產模式受到了極大的挑戰.機器人制造單元在實現量產的同時,能夠保持小批量定制生產的靈活性,能夠更好地適應當前智能制造發展的需求,而合理的優化調度是機器人制造單元整體控制的一個重要方面.制造企業為了在激烈的競爭環境中生存和發展,需要提供能夠滿足客戶各種需求的資源調度方案.單一目標的尋優已經滿足不了人們日常生活的需要,工業上對于優化目標個數增多的要求使得求解問題的難度顯著增加.目前,調度問題的理論和實踐方法研究,仍不能滿足現實生產需求.傳統的優化算法在解決復雜多目標問題中[2–5]常常效率低下,難以達到預期.智能優化算法[6–9]通常是基于非劣解概念基礎上的現代啟發式算法,具有更好的處理高維多目標問題的能力.因此,在經典車間調整理論的基礎上,本文要針對Pareto 支配關系在高維多目標優化中的支配能力不足,將另外兩種支配關系與NSGA-III 算法相結合去對生產調度問題進行求解去解決機器人制造單元優化上的多目標協同調度優化問題.
本文主要研究的是有關job-shop 類型的機器人的智能制造車間的搬運問題.Job-shop 類型的機器人制造單元主要由加工站、搬運機器人和控制系統等組成.根據車間調度的工業路線特點,可將job-shop 類型的機器人制造單元調度過程描述為:有n個工件J={J1,J2,···,Jn} 將在m臺加工機器M={M1,M2,···,Mm}上進行加工,其中M0為裝載站兼卸載站.每個工件Ji有p道加工工序Oij={Oi1,Oi2,···,Oip},且每個工件都有預先設定好的工序加工順序,所有的工序都嚴格按照預定的加工工序順序在m臺機器上進行.現已知每道工序的加工時間、每道工序對應的加工機器和任意兩臺機器之間的搬運時間.如圖1所示:調度基本過程可描述為,首先,n個 工件在M0上等待加工,當加工任務開始時,機器手按照工序順序,當工件處于裝載區或工件上一道工序已完成時,從該工件所在位置上取走工件并搬運到該工件下一道工序對應加工機器上,根據當前加工機器狀態決定是否開始加工.當最后一道加工工序完成后并由機器手將該工件搬運至卸載站上時,則為加工任務結束.

圖1 以搬運機械人為中心的制造單元
近年來,在生產過程中人們越來越關注能夠使多項決策目標都能得到滿足的生產調度方案.通常考慮多個目標共同優化的方案.同時,由于對環境問題的重視,本文又引入了加工總能耗,生產總成本等一系列的綠色調度目標.以下是本文要優化的5 個目標:
1)最大完工時間:

2)加工總能耗:

3)交貨期提前量:

4)交貨期延遲量:

5)生產總成本:

其中,最大完工時間指所有工件完成加工后的總加工時間,N為工序總數;Ci為第i個工件的完工時間.加工總能耗指調度過程中機器加工總功率與機器手搬運總功率之和.M為機器總個數;Ni為第i個機器上的加工工序數;tij為第i個機器上第j道工序的加工時間;Pm,i為第i個機器的加工功率;ti,idle為第i臺機器空載時間;Pi,idle為第i臺機器空載功率;Nb為機器手搬運次數;tb,i為第i次 搬運所花費的時間;Pb為機器手拌搬運功率;tb,idle為機器手空載時間;Pb,idle為機器手空載功率.交貨期為制造企業交付加工產品的時間.若工件在理想交貨時間前完成生產則為提前,在理想交貨時間后完成則為延遲.N為工件個數;DTj為第j工件的理想交貨時間;RTj為第j工件的實際交貨時間.生產總成本為加工所需電費與機床使用成本之和.coste為電價(單位能耗成本);Etotal為加工總能耗;Cu,i為第i臺機器單位時間的使用成本.
NSGA-III[10]是Deb和Jain 在2014年提出來的用于處理多個優化目標的進化優化算法.其主要思路是在NSGA-II的基礎上引入了參考點機制,對于那些非支配并且接近參考點的種群個體進行保留,相比NSGAII 擁有優秀的處理更高維度目標的能力.
本文以NSGA-III的算法基本框架,全局搜索以遺傳算法為基礎,使用了一種符合調度問題特點的交叉變異策略,并提出了兩種更改支配關系的改進措施.改進后的NSGA-III 將具有更好的處理高維數據的能力,用于解決本文所提出的機器人制造單元的高維多目標調度問題.NSGA-III的算法流程如下.
Step 1.初始化,隨機產生大小為N的父代種群Pt.更新操作,父代通過交叉、變異產生大小為N的子代種群Qt,合并為大小為2N的種群Rt=Pt+Qt.
Step 2.計算目標值,通過非支配排序,將Rt劃分為不同的非支配層(F1,F2,···).
Step 3.選擇優先級高的非支配子集保留到下一代,用S t來保留.設Fl為臨界層,臨界層為加上該層解數目大于等于N時的非支配層.若數目等于N則S t=S t∪Fi,i=1,2,···,l,否則S t=S t∪Fi,i=1,2,···,l?1. 令Pt+1=S t.用K=N?|Pt+1|來 記錄需要從Fl中補充的解的個數.
Step 4.每一維度最大值作為極值點,最小值作為理想點,對每個解的目標值做歸一化處理.
Step 5.以極值點構造超平面,產生均布的C(M+H?1,H)個參考點,M為維度數,H為每維分割份數.
Step 6.對S t中和Fl層中每個解找出距離最近的參考點分別進行關聯.
Step 7.對于S t中關聯少的參考點,找出Fl層對應參考點關聯的一個最近的解進行補充.每添加一個解,就將該解從Fl中刪除,并將S t中關聯解的個數更新排序.重復Step 7 直至補充K個解.
Step 8.返回Step 2,直至迭代次數達到預設值,對當前解集進行非支配排序找出第一非支配層的解作為最終解集.
隨著目標數量的增加,算法對解的選擇壓力變小,基于Pareto 支配,解之間很難存在支配關系,這導致解集中很大部分的解都是非支配的.而算法本身就強調保留種群中的非支配解,所以在每一代種群中,就沒有太多的空間來保留新解.這樣就減緩了算法的搜索過程,使算法效率變得低下.使用改變支配關系的方法,將解的目標值進行調整,可以擴大每個解的支配空間,從而提高解之間的區分度.本文分別將Lorenz 支配和CDAS 支配與NSGA-III 結合來進行優化.第2.2.1 節和第2.2.2 節分別介紹了 Lorenz 支配和CDAS 支配.
2.2.1 Lorenz 支配

Lorenz 支配(式(6))是將標準化后的每一個解的目標值,按照從小到大進行排列,之后在一個新列表里,第i個位置存放前i個目標值的和(i=1,2,···,m),m為維度數.所有解的目標值f(x)都用Lorenz 修改過后,用此新的目標值f′(x)來進行非支配排序,修改過支配關系后,每個解的支配空間將會變大,每個解之間就有更好的區分度.
如圖2所示是一個雙目標的例子,對于解S,基于Pareto 支配區域為(a),經Lorenz 支配S的支配空間增加到(a)、(b)、(c).圖2中,X是一個解序列所求出對應的目標值的坐標.

圖2 Lorenz 支配關系
2.2.2 CDAS 支配

CDAS 支配(式(7))如圖3所示,式中,r為該解的目標值點與原點之間的距離;ωi為該點與原點連線和第i維的夾角;φi為該點與第i維上一個可控的夾角,其大小由系數S控制.當S=0.5時,CDAS 支配與Pareto支配效果相同;當S<0.5時,CDAS 支配將使原目標值點的支配區域擴大;當S>0.5時,CDAS 支配將使原目標值點的支配區域變小.圖3中,X是一個解序列所求出對應的目標值的坐標.

圖3 解x 修改前后第i 維目標值
通過試驗發現,當S取0.25 時,求解多目標優化問題的效果最佳.此時根據CDAS 支配將原目標值f(x)修改為f′(x),擴大該解的支配區域.
以雙目標優化問題為例,如圖4所示,當S=0.5,即Pareto 支配關系下解A、B、C 是互不支配的.取S=0.25,增大擴張角,從圖中虛線部分可以看到每個解的支配區域都有所擴大,解之間出現了更細微的序值分布,此時,3 個解從互不支配變為解B 支配解C.

圖4 CDAS 修改后解支配關系
為了驗證本文提出的改進算法的有效性,通過對本文算法與傳統的NSGA-Ⅲ算法進行對比.并從以下3 個方面來測試算法性能.
1)收斂性:反映解集與真實Pareto 前沿趨近程度.
2)多樣性:反映是否所得的Pareto 最優解集中的解都分布在 Pareto 前沿上.
3)均勻性:反映所得的Pareto 最優解集中的解在目標空間中分布的均勻程度.
本次實驗采用世代距離GD (generational distance)、反向世代距離IGD (inverted generational distance)和間距Spacing 作為評價算法性能的指標.其中,GD 用于評價解集收斂性,IGD 用于評價解集收斂性和多樣性的綜合性能,Spacing 評價解集均勻程度.
本文提出的改進優化算法以 PyCharm Python 3.8為開發工具編程實現,其中相關參數設置如下:實驗設計環境和配置為Intel(R)Core(TM) i5-7300HQ CPU@2.50、8.00 GB RAM、Windows 10 64 位操作系統.交叉概率設為1.0,變異概率設為0.5,迭代次數設為50,種群大小設為40.表1中,Ind表示案例編號,J表示該案例工件個數,M表示加工機器臺數,S表示加工工序總數,上標1為C-NSGA-III 算法測試結果,上標2為L-NSGA-III 算法測試結果,上標3為NSGA-III 算法測試結果,下標AVG為均值,SD為標準差,MIN為最小值.本文挑選了由小到中到大型的案例共10 個進行了實驗,用來測試本文所提出的算法的性能,各案例均按照一定規則隨機生成.每個案例運行10 次,取均值作為測試結果.
為了更好的描述最優方案的調度過程,可以通過甘特圖來分析工件的加工順序與機器人的搬運順序,如圖5為案例CS-3-16-5的甘特圖.該案例包含3 個工件,5 臺加工機器,16 道工序.圖中,橫坐標為時間軸,縱坐標為加工機器序號,M0為裝載卸載站.純色矩形代表加工過程,虛線矩形代表工件在當前機器上的等待時間,實線代表機器人帶載過程,虛線代表空載過程.“J00”代表工件0的第0 道工序,之后以此類推.則案例CS-3-16-5的加工過程可描述為:機器人將工件2 從M0 搬出至M3 進行工序J20的加工,搬運后機器人返回至M0 開始搬運工件0 至M1 進行工序J00的加工,隨后返回M0 搬運工件1 至M5 加工工序J10.隨后機器人轉至M3 取走在M3 上完成加工的且已在輸出緩沖區上等待的工件2,搬運至M2 加工工序J21.以此類推,直到所有工件加工結束.

圖5 案例CS-3-16-5 甘特圖
表1、表2、表3分別為C-NSGA-III、L-NSGAIII和NSGA-III 算法在GD、IGD、Spacing 上的對比結果.通過對比可以得出,C-NSGA-III和L-NSGAIII 都在GD 上表現良好,在絕大部分案例上數值都小于NSGA-III,且標準差不大,說明改進算法較為穩定.在GD 最小值上,C-NSGA-III 始終為0,總體可得兩種改進算法在收斂性上都要優于NSGA-III.但在IGD 值上,NSGA-III 幾乎始終占優,說明C-NSGA-III和LNSGA-III 在多樣性上會變差,分析原因是兩算法均加大了支配空間,它們的解相當于Pareto 關系下的解的子集,所以解的多樣性會變差.比較Spacing 值可知,兩改進算法在解集均勻程度上,較優于NSGA-III.總體來說,C-NSGA-III和L-NSGA-III 算法有效地提升了解集的收斂性和均勻性,在多樣性上表現一般.

表1 C-NSGA-III、L-NSGA-III 與NSGA-III 算法GD 指標對比表

表2 C-NSGA-III、L-NSGA-III 與NSGA-III 算法IGD 指標對比表

表3 C-NSGA-III、L-NSGA-III 與NSGA-III 算法Spacing 指標對比表
本文針對研究機器人制造單元的多目標協同調度問題,同時優化5 個目標:最大完工時間,加工總能耗,交貨期提前量,交貨期延遲量及生產總成本.針對job-shop的特點,分別將Lorenz 支配和CDAS 支配與NSGAIII 結合來進行優化并將結果與NSGA-III 進行比較,證明該算法有效.