莊存波， 熊輝， 劉檢華， 唐承統(tǒng)

(北京理工大學 機械與車輛學院， 北京 100081)

0 引言

復雜產(chǎn)品是指客戶需求復雜、產(chǎn)品組成復雜、產(chǎn)品技術復雜、制造過程復雜、項目管理復雜的一類產(chǎn)品，如導彈、衛(wèi)星、火箭、飛機等[1]。裝配是復雜產(chǎn)品生產(chǎn)的最后環(huán)節(jié)，也是最為重要的環(huán)節(jié)之一，其直接關系到產(chǎn)品的質(zhì)量、壽命、性能和可靠性[2]。生產(chǎn)調(diào)度是指在一定約束條件下，把有限資源在時間上分配給若干個任務，以滿足或優(yōu)化一個或多個性能指標的過程[3]。生產(chǎn)調(diào)度是產(chǎn)品裝配過程的關鍵環(huán)節(jié)，也是裝配過程管理與控制的核心問題之一。復雜產(chǎn)品裝配多以流程為主組織生產(chǎn)，其裝配過程是典型的工作流。在航天器等復雜產(chǎn)品的裝配過程中，通常包括產(chǎn)品總裝和組件部裝兩個階段，每個階段均需繪制總裝工藝流程圖和部裝工藝流程圖。總裝裝配和部裝裝配多數(shù)是分布在不同的裝配車間，且二者的裝配過程相似，故本文僅研究單個階段的裝配調(diào)度問題。復雜產(chǎn)品裝配調(diào)度問題(總裝或部裝)可簡化為一個以最大完工時間最小為調(diào)度目標、針對關鍵路徑進行排產(chǎn)的多階段混合流水車間調(diào)度問題(HFSP)。HFSP為流水車間調(diào)度問題(FSP)與并行機調(diào)度問題的集合，也被稱為柔性FSP(FFSP)。與一般FSP相比， HFSP加大了機器的可選擇性，擴大了可行解的搜索范圍，是更加復雜的NP-hard問題。傳統(tǒng)HFSP按并行機類型可分為相同并行機HFSP[4]、均勻并行機HFSP[5]及不相關并行機HFSP[6]共3種，按階段數(shù)量可分為2階段HFSP、3階段HFSP及多階段HFSP[7].

自1954年Johnson[8]研究具有兩臺機器、工件加工順序已知的HFSP開始，國內(nèi)外學者對HFSP進行了研究，并取得了顯著成果。目前，HFSP求解算法主要包括精確算法、啟發(fā)式算法和元啟發(fā)式算法。精確算法中分支定界(B&B)算法使用最為廣泛，如Brah等[9]利用B&B算法求解HFSP時提出了一種對工件排列和機器分配進行枚舉的分支策略，Neron等[10]為了增強B&B算法搜索效率提出了一種基于推理的全局操作方法。精確算法在理論上能夠獲取問題的最優(yōu)解，但面對大規(guī)模調(diào)度問題時，其求解時間過長，故只適用于求解小規(guī)模調(diào)度問題。Gupta[11]提出了一種簡單的啟發(fā)式算法來求解第1階段含1臺機器、第2階段含2臺相同機器的2階段HFSP. Riane等[12]提出了2種新的簡單啟發(fā)式算法來求解3階段相同并行機調(diào)度問題。對于多階段HFSP，Ruiz等[13]和Pan等[14]在種群初始化階段利用Nawaz-Enscore-Ham(NEH)啟發(fā)式算法[15]生成少數(shù)優(yōu)質(zhì)個體以提高種群質(zhì)量，從而優(yōu)化算法尋優(yōu)性能。啟發(fā)式算法雖然能在較短時間內(nèi)獲得可行解，但難以保證解的質(zhì)量。為了獲得更優(yōu)質(zhì)的解，國內(nèi)外許多學者在近些年更多采用元啟發(fā)式算法來求解HFSP，如遺傳算法(GA)[16-17]、禁忌搜索(TS)[18-19]、模擬退火(SA)[20]、人工免疫系統(tǒng)(AIS)[21-22]、迭代局部搜索(ILS)[23]、粒子群優(yōu)化(PSO)[24-25]、人工蜂群(ABC)算法[26-27]和蟻群優(yōu)化(ACO)[28-29]等群體智能優(yōu)化算法(SIOA)，這些算法多用于求解多階段相同并行機HFSP. 2012年，Gandomi等[30]通過模擬磷蝦種群覓食行為提出了一種新穎的仿生SIOA，即磷蝦群(KH)算法。在函數(shù)優(yōu)化領域，KH算法及其改進算法在許多標準測試問題上被證明比一些熟知的元啟發(fā)式算法更好或非常具有競爭性[31]。KH算法具有較強局部探索能力，且控制參數(shù)少，易于實現(xiàn)，非常適用于并行計算。然而，KH算法無法一直很好地實現(xiàn)全局搜索，易陷入局部最優(yōu)，且KH算法更多依賴于隨機運動，易導致收斂速度緩慢[32-33]。目前，KH算法已在特征識別、人工神經(jīng)網(wǎng)絡訓練、比例—積分—微分(PID)控制、逆幾何設計、基于地圖的網(wǎng)絡路徑優(yōu)化和結(jié)構優(yōu)化等問題上得到廣泛研究與應用[34]。如Kowalski等[35]將KH算法應用于訓練人工神經(jīng)網(wǎng)絡，并通過仿真結(jié)果和算法比較得出KH算法在精度和時間上比其他元啟發(fā)式算法和傳統(tǒng)方法更有效。Yaghoobi等[36]提出了一種改進的混沌KH(ICKH)算法，將其應用于PID控制器參數(shù)的確定，通過對比發(fā)現(xiàn)ICKH算法比其他優(yōu)化算法性能更優(yōu)。KH算法在復雜產(chǎn)品裝配車間調(diào)度問題和HFSP上的應用研究尚無文獻可查。

本文以復雜產(chǎn)品裝配車間為研究對象，采用基于排列的編碼方法和基于啟發(fā)式規(guī)則的解碼方法，提出一種改進的離散KH(IDKH)算法，通過實驗設計方法分析不同參數(shù)設置對算法性能的影響，最后通過實例和算法比較驗證所提算法的高效性和魯棒性。

1 復雜產(chǎn)品裝配調(diào)度模型

1.1 問題描述

1.2 數(shù)學模型構建

基于上述分析，構建復雜產(chǎn)品裝配調(diào)度問題的數(shù)學模型，具體描述如下：

minCmax，

(1)

(2)

k∈{1,2,…,s-1},

(3)

(4)

zghk+zhgk≤1,g,h∈J,g≠h,

(5)

(6)

xjkm,zghk∈{0,1},

(7)

式中：Cmax為最大完工時間；g、h為產(chǎn)品編號；m∈AT，AT為裝配班組集合，AT={1,2,…,L}，L為裝配班組總數(shù)；pjkm為產(chǎn)品j在工序k、由班組m進行裝配操作所需要的時間；tjk為產(chǎn)品j的第k道工序開工時間；U為足夠大的正數(shù)；xjkm=1為產(chǎn)品j的第k道工序被分派到裝配班組m上，否則xjkm=0；zghk=1為在第k道工序上產(chǎn)品g先于產(chǎn)品h被裝配，否則zghk=0. (1)式表示調(diào)度性能指標，即工期最?。?2)式表示最后一個產(chǎn)品在最后一道工序s的完工時間為Cmax；(3)式表示同一個產(chǎn)品只有在前一道工序裝配完畢后才能開始下一道工序；(4)式表示每個產(chǎn)品在各個工序階段只能由一個裝配班組進行裝配操作；(5)式表示產(chǎn)品g和產(chǎn)品h間存在先于、后于、同時3種順序關系；(6)式表示同一個裝配班組在一個時刻只能裝配一個產(chǎn)品，即前一個產(chǎn)品裝配完成后才能開始下一個產(chǎn)品的裝配；(7)式表示變量的取值為0或1.

2 標準KH算法介紹

KH算法在優(yōu)化過程中的目標是磷蝦與食物之間距離最短及磷蝦種群密度最大[30]。磷蝦個體i(i∈{1,…,NP}，NP為種群數(shù)量)在t時刻的位置向量Xi由3個因素決定：1)受其他磷蝦個體誘導所產(chǎn)生的運動向量Ni；2)覓食運動向量Fi；3)隨機的物理擴散運動向量Di. 用拉格朗日模型來描述磷蝦個體的移動位置向量：

(8)

1)Ni由局部種群密度(即局部影響)、目標種群密度(即目標影響)和排斥種群密度(即排斥影響)3個因素決定。磷蝦個體i受其他磷蝦個體誘導所產(chǎn)生的運動向量可表示為

(9)

2)Fi由食物位置和磷蝦個體i的過往覓食經(jīng)驗2個因素決定。磷蝦個體i的覓食運動向量可表示為

(10)

3)Di為一個隨機過程，對于磷蝦個體i，有

Di=Dmaxδ,

(11)

式中：Dmax為最大擴散速度，取值范圍為[0.002,0.010]，一般取0.005；δ為隨機方向向量，每項都是位于[-1,1]區(qū)間的隨機數(shù)。

4)種群位置更新[30]：對于磷蝦個體i，第gen+1次迭代后的位置向量可表示為

(12)

式中：Δt為速度向量的比例因子，是非常重要的常數(shù)之一，需根據(jù)優(yōu)化問題仔細設定。Δt完全由搜索空間決定，可簡單表示為

(13)

式中：c為一個常數(shù)，取值范圍為[0，2];UBv和LBv為變量v的上界和下界。

5)遺傳算子：受進化算法的啟發(fā)，交叉算子和變異算子兩種遺傳繁殖機制被加入到KH算法當中，經(jīng)過仿真實驗發(fā)現(xiàn)，只將交叉算子應用到KH算法中綜合性能最優(yōu)。關于標準KH算法的詳細介紹見參考文獻[30].

3 IDKH設計

3.1 編碼與解碼

目前，求解HFSP的算法大多采用基于排列的方式對群體中個體進行編碼，即取所有工件序號排列作為一個個體。群體中的每一個個體都對應一個可行調(diào)度解，個體長度為工件數(shù)量n，工件號在排列中位置代表其在第一個階段加工順序。工件排序和后續(xù)機器分配由解碼規(guī)則決定[38-39]。

機器分配部分一般采用最先空閑機器規(guī)則[40]，但主要面向相同并行機調(diào)度問題。由于復雜產(chǎn)品裝配調(diào)度問題屬于不相關并行機調(diào)度問題，因此，本文基于該規(guī)則，提出如下改進策略：

1)對于每個工件j，選擇完工時間最短的機器作為工件j的加工機器，完工時間為am+pjkm，其中，am為工件j最早允許加工時間且am=max (rm,Cj,k-1)，rm為機器m釋放時間，Cj,k-1為工件j在上一階段的完工時間；

2)對于完工時間相同的機器，基于效率最大化原則，優(yōu)先選擇加工時間pjkm最小的機器；

3)對于完工時間和加工時間都相同的機器，為了使調(diào)度結(jié)果更加緊湊，選擇am-rm最小的機器。

工件排序部分采用如下改進策略：

1)第一階段按照編碼方式確定工件加工順序，后續(xù)階段基于先到先加工方式進行排序；

2)對于前一階段完成時間相同的工件，選擇當前階段加工時間短的工件優(yōu)先加工。對于復雜產(chǎn)品裝配調(diào)度問題，可用產(chǎn)品序號代替工件序號，用裝配代替加工，用裝配班組代替機器。

3.2 種群初始化

為了保證種群的多樣性以獲得全局最優(yōu)解，采用隨機方式生成種群個體的位置向量X=(X1,X2,…,XNP)，個體i位置向量為Xi，變量數(shù)為工件數(shù)n，上下界范圍可自定義。然而，標準KH算法只適用于連續(xù)優(yōu)化問題的求解，并不適用于離散優(yōu)化問題的求解，因此需將個體位置向量Xi轉(zhuǎn)化為工件序號的排列πi后才能夠被解碼。

為此，本文提出轉(zhuǎn)化策略，對于每一個磷蝦個體i，其對應的位置向量為

Xi=(Xi1,Xi2,…,Xin)，

(14)

3.3 局部搜索

針對標準KH算法收斂速度較為緩慢的問題，本文采用局部搜索策略來加快算法收斂速度，并增強算法局部開采能力。目前，常用的局部搜索策略包括：1)交換策略：從排列編碼中隨機選擇2個位置的元素進行交換；2)插入策略：從排列編碼中隨機選擇2個位置，將第1個位置中的元素放到第2個位置上。

為了不影響算法運算速度，在每一次迭代過程中，只對目前出現(xiàn)的最優(yōu)解執(zhí)行局部搜索策略。本文采用一種修改的變鄰域局部搜索策略，具體實現(xiàn)步驟如下：

1)對目前已知的最優(yōu)解Xbest執(zhí)行插入操作得到Xbest′.

2)對Xbest′執(zhí)行混合鄰域局部搜索策略(即在每一次局部搜索過程中，從交換策略和插入策略中隨機選擇一種，選擇概率均為50%)得到其鄰域X，如果鄰域X適應度值比Xbest′適應度值低，則用X替代Xbest′進行下一次混合鄰域搜索。該步驟迭代次數(shù)為n(n+1).

3)比較Xbest′與Xbest的適應度值，即比較f(Xbest′)與f(Xbest)，若f(Xbest′)

3.4 重啟操作

針對KH算法全局搜索性能有限而導致的容易陷入局部最優(yōu)問題，本文提出了一種重啟策略。具體運行流程為：1)當種群最優(yōu)解不更新代數(shù)Age達到設定常數(shù)Limit時，算法將執(zhí)行重啟操作；2)保留百分比為η的最優(yōu)個體，剩下(1-η)·NP個個體則被隨機生成的種群個體替代。重啟操作偽代碼如圖2所示。

3.5 IDKH算法流程

基于上述分析，構建如圖3所示的求解復雜產(chǎn)品裝配調(diào)度問題的IDKH算法運算流程。

4 仿真結(jié)果與算法比較

4.1 參數(shù)討論

采用MATLAB 2013平臺對IDKH算法進行編程，并在性能為Inter(R) Core(TM) i7-3770 CPU 和3.40 GHz RAM的計算機上運行IDKH算法，操作系

SetLimit,η;

IfAge>=Limit

Keep the current bestXbest

Randomly generate the left (1-η)*NPindividuals

While not termination

Repeat

Evaluate the fitness of the krill

Perform the three motions

Implement the crossover operator

Update the krill position

一切都和過往的經(jīng)驗不一樣了。在國家隊時，鄒習慣謹慎地談論目標。第一場職業(yè)賽前新聞發(fā)布會，一個泰國記者看到他笑瞇瞇的，感到吃驚。而現(xiàn)在，他會訓練自己不怒自威的樣子，擺一張臭臉，不光給媒體看，更多給對手看。

Evaluate the krill based on its new position

Until a number ofNPreplications

Sort all the krill and find the current bestXbest′

Iff(Xbest′)

Xbest=Xbest′

END If

END while

END If

END

圖2 重啟操作的偽代碼

Fig.2 Pseudo code of restart operation

統(tǒng)為Windows 7.

IDKH算法中，種群個數(shù)NP、控制時間間隔C、重啟操作常數(shù)Limit和最優(yōu)種群保留比例η等4個參數(shù)會對其性能產(chǎn)生影響。各參數(shù)水平取值如表1所示。

為了測試IDKH算法的性能，本文隨機生成一系列計算實例，產(chǎn)品(工件)數(shù)量J∈{10,15,20}，

表1 IDKH算法參數(shù)水平表Tab.1 Combinations of parameter values of IDKH

工序(階段)數(shù)S=5，裝配(加工)工時為取值范圍[3, 40]的離散隨機整數(shù)分布，每階段的并行機數(shù)均為3. 每種組合隨機生成3個實例，共計NI=3×3=9個實例。為每個實例選擇規(guī)模為L16(44)的正交實驗，IDKH算法在每種參數(shù)組合下均獨立運行20次，設定終止條件為最大迭代次數(shù)200，同時限制每次運行時間不超過20 s. 計算得到每個實例正交結(jié)果如表2所示。

表2 IDKH算法正交表Tab.2 Orthogonal array and response values

由表2可知，組合[3 2 4 1]的不同參數(shù)組合的平均值最小，是較為穩(wěn)定的、獲得最優(yōu)參數(shù)值的組合，此時NP=80，C=1.0，Limit=100，η=0.1. 根據(jù)表2，種群個數(shù)NP、控制時間間隔C、重啟操作常數(shù)Limit和最優(yōu)種群保留比例η的極差和重要程度如表3所示。

表3 IDKH算法各參數(shù)極差表Tab.3 Statistics of response values

由表3可知，Limit和η極差最大，其次為NP，C最小，這表明重啟操作對IDKH算法的性能影響最大。從數(shù)值上看，4個參數(shù)極差均不超過0.50%，說明IDKH算法魯棒性較強。不同參數(shù)在不同水平下對算法性能的影響趨勢圖如圖4所示。

由圖4(a)可知，NP越大，算法性能越優(yōu)，原因在于NP決定了群體在搜索空間的覆蓋范圍。當NP較小時，算法全局探索性能較差，易早熟收斂而陷入局部最優(yōu)；當NP較大時，群體能夠覆蓋更多搜索空間，更易獲得高質(zhì)量的解。NP較大時，計算量大、耗時長，同時變優(yōu)幅度隨著NP增加不斷變小(水平3和水平4差距不大，僅為0.03)，故NP不宜過大，NP=80為最佳。C決定了算法尋優(yōu)過程中在整個搜索空間的移動幅度，C過大易導致收斂速度過快而陷入局部最優(yōu)，C過小易導致收斂速度緩慢從而影響IDKH算法性能，由圖4(b)可知，C=1.0時IDKH算法整體性能最優(yōu)。對于重啟操作，Limit決定了重啟時機。當Limit較小時，IDKH算法還未達到局部最優(yōu)就被迫重啟，不僅重啟效果一般，且會影響IDKH算法性能。η決定了重啟時保留最佳個體比例，當保留比例較大時，IDKH算法易再次陷入局部最優(yōu)，導致重啟效果不明顯。由圖4(c)和圖4(d)可知，當Limit=100，η=0.1時，重啟效果最佳。綜上所述，組合[3 2 4 1]為最佳參數(shù)組合，即NP=80，C=1.0，Limit=100，η=0.1.

4.2 初始實驗

與離散KH(DKH)算法相比，IDKH算法主要是在局部搜索和重啟操作兩個方面進行了改進。

本文實驗中，考慮對比IDKH算法及其3種變型：DKH算法、包含重啟操作的DKH(DKH-RS)算法、包含局部搜索的DKH(DKH-LS)算法。為了驗證DKH算法、DKH-RS算法、DKH-LS算法和IDKH算法求解HFSP的性能，采用4.1節(jié)中的9個實例對算法性能進行比較，選擇最優(yōu)參數(shù)組合NP=80，C=1.0，Limit=100，η=0.1，每個實例均運行20次，每次運行時間不超過20 s. 如圖5所示為IDKH算法不同變型下的平均值、最優(yōu)值和標準差變化趨勢圖。

通過DKH算法和DKH-LS算法對比以及DKH-RS算法和IDKH算法對比發(fā)現(xiàn)，加入局部搜索操作后，算法的平均值、最小值和標準差都會得到優(yōu)化，說明局部搜索策略能夠全面提高算法尋優(yōu)性能和魯棒性。通過DKH算法和DKH-RS算法對比以及DKH-LS算法和IDKH算法對比發(fā)現(xiàn)，加入重啟操作后，算法的平均值和標準差得到了優(yōu)化，即算法平均尋優(yōu)性能和魯棒性得到了增強 但最小值略有下降，原因在于算法還未達到局部最優(yōu)就被迫重啟了，然而大多數(shù)情況下并不會出現(xiàn)這種情況，如DKH-LS算法僅在一個實例上比IDKH算法獲得的最小值小，其他均一致。綜合來看，IDKH算法性能最優(yōu)。

4.3 計算結(jié)果

4.3.1 文獻實例驗證

目前，求解HFSP比較常用的智能算法為GA[13,41]，近幾年涌現(xiàn)出了一些新穎算法，如萬有引力搜索算法(GSA)[42]、差分進化(DE)算法[43]、分布式估計算法(EDA)[38,44-45]。本文利用參考文獻[38,46]所提3個HFSP實例L1、L2、L3，對比IDKH算法和DKH算法與GA、GSA、DE算法、EDA在求解實例時的性能差異。其中：DE算法數(shù)據(jù)來源于參考文獻[43]，評價次數(shù)為10 000；GA、GSA、EDA數(shù)據(jù)均來源于參考文獻[45]，采用運行時間10 s為終止條件。對于IDKH算法，同樣設定終止條件為10 s，參數(shù)組合為NP=80，C=1.0，Limit=100，η=0.1，運算結(jié)果如表4所示。另外，由于智能算法均存在一些隨機因子， 運行10次的結(jié)果有時并不能準確反映出算法性能(如算法穩(wěn)定性)。因此，為了能夠更加準確地驗證所提算法性能，本文對DKH算法和IDKH算法運行70次，根據(jù)70次運行結(jié)果及其平均值得出10次運行結(jié)果。對于GA、GSA、DE算法、EDA(包括EDA_W算法[38]和EDA_I算法[45])，則采用參考文獻[45]所提供的10次運算結(jié)果，得出其平均值。

表4 3個實例在不同智能算法中的結(jié)果對比

由表4可知， IDKH算法和EDA_I算法結(jié)果要明顯好于其他算法，且都很穩(wěn)定。二者求解性能相差不大，但IDKH算法略優(yōu)于EDA_I算法，基本上每次都能獲得目前已知的最優(yōu)解，如：對于L2和L3，每次都能獲得目前已知的最優(yōu)解297和13；對于L2，在運行70次過程中，僅有一次沒有獲得已知的最優(yōu)解21，運行70次得到的平均值為21.014 3. 其次為GSA和DKH算法，再次為GA和EDA_W算法，最后為DE算法，從而可知IDKH算法在求解HFSP時是高效且穩(wěn)定的。

4.3.2 工程實例驗證

以某航天器裝配為例，在調(diào)度開始時刻，需要完成10個產(chǎn)品的裝配，其裝配工藝流程如圖6所示。裝配總工序數(shù)為7，關鍵路徑上的工序數(shù)為4，分別標記為S1、S2、S3、S4. 關鍵路徑上的每道工序均有3個裝配班組,分別標記為M1～M12，及5個裝配工位(滿足裝配工位數(shù)量大于等于裝配班組數(shù)量的前提)可供分配，不同裝配班組完成不同裝配工序所需要的裝配時間如表5所示。

基于該實例，分別用IDKH算法、DKH算法、EDA_W算法進行求解，每種算法均運行70次，每次最大運行時間仍為10 s.求解發(fā)現(xiàn)只有IDKH算法每次都能獲得目前已知的最優(yōu)解Cmax=195，其調(diào)度甘特圖如圖7所示。圖7中，矩形框的數(shù)字表示產(chǎn)品序號j.

表5 某航天器的裝配時間表Tab.5 Assembly time of a spacecraft h

5 結(jié)論

針對復雜產(chǎn)品裝配調(diào)度問題，引入了一種高效的群體智能優(yōu)化算法—KH算法，并采用局部搜索和重啟策略增強了KH算法的局部尋優(yōu)能力和全局搜索能力。針對已有參考文獻[38,46]中3個不同實例和1個工程實例，IDKH算法求解結(jié)果明顯比GSA、DKH算法、GA、EDA_W算法、DE算法等算法好，略優(yōu)于EDA_I算法，且基本上每次都能獲得目前已知的最優(yōu)解，反映出IDKH算法在求解HFSP和復雜產(chǎn)品裝配調(diào)度問題時是高效且穩(wěn)定的。