武子科，潘 攀+，彭 誠，呂秀莎，梁子涵，張洪光

(1.中國空間技術(shù)研究院 北京東方計量測試研究所，北京 100093； 2.北京郵電大學 電子工程學院，北京 100876)

0 引 言

儀表檢測工作正由傳統(tǒng)的人工檢測轉(zhuǎn)變?yōu)樽詣踊瘷z測。自動化儀表檢測中，為節(jié)省批次任務時間同時滿足儀表檢測順序約束的資源分配問題，構(gòu)成了一個作業(yè)車間調(diào)度問題[1](job shop scheduling problem，JSSP)。由于作業(yè)車間調(diào)度問題的NP-Hard特性[2]和No-Free-Lunch特性對于不同應用場景沒有通用的求解方法[3]，目前已有較多在經(jīng)典啟發(fā)式算法上做出改進的算法，如并行迭代的雙種群混合遺傳算法[4]和混合并行布谷鳥搜索算法[5]，使用鄰域搜索的改進Jaya算法[6]和局部鄰域搜索算法[7]，以及結(jié)合幾種算法優(yōu)勢的ALPS-GA算法[8]和混合NSGA-II算法[9]等。鯨魚優(yōu)化算法由Mirjalili等[10]提出，是一種源于座頭鯨群狩獵的仿生學算法[11,12]。在作業(yè)車間調(diào)度問題中，已有使用量子旋轉(zhuǎn)門[13]、使用萊維飛行[14]的鯨魚優(yōu)化算法，但在具體場景中應用的研究還剛剛起步。

自動化儀表檢測工作中的儀表種類較多，且檢測差異較大，其調(diào)度問題對比經(jīng)典調(diào)度問題具有數(shù)據(jù)復雜(不同儀表檢測時長相差較大)、穩(wěn)定性要求高的特點。因此，需要有能夠同時在有效性和穩(wěn)定性方面表現(xiàn)更好的算法來提升自動化儀表檢測工作的效率。本文根據(jù)自動化儀表檢測中的作業(yè)車間調(diào)度問題特點，提出了一種基于生命力選擇的精英鯨魚優(yōu)化算法，達到了提高檢測效率、節(jié)省檢測成本的目的。首先，介紹了作業(yè)車間調(diào)度問題的相關(guān)定義和模型；然后，提出了基于生命力選擇的精英鯨魚優(yōu)化算法；最后，在標準數(shù)據(jù)集和北京東方計量測試研究所自動化儀表檢測的實際調(diào)度問題中進行了仿真實驗，驗證了算法的有效性和穩(wěn)定性。

1 作業(yè)車間調(diào)度模型

作業(yè)車間調(diào)度問題模型由工件集J={J1,J2,…,Jn} 和對其進行加工的機器集M={M1,M2,…,Mm} 組成，每個工件Ji均按照自身的加工順序在對應機器Mk上執(zhí)行相應的工序Oik。相關(guān)的參數(shù)定義見表1。對于該問題有如下假設：

(1)每個工件需按照自身加工順序在對應機器上完成加工；

(2)工件工序有固定的加工時間，且加工過程不能中斷；

(3)一臺機器同一時刻只能加工一個工件；

(4)一個工件同一時刻只能被一臺機器加工；

(5)0時刻所有機器和工件均可用；

(6)不考慮機器故障及工件到達延遲等問題。

本文使用整數(shù)規(guī)劃模型[5]來描述作業(yè)車間調(diào)度問題，并以最小化Cmax為優(yōu)化目標，如式(1)～式(4)所示。式(1)為目標函數(shù)；式(2)表示每個工件的加工順序?qū)谠摴ぜ付ǖ募庸ろ樞?；?3)表示每臺機器一次只能加工一個工件；式(4)表示每個工件在每臺機器上的完工時間為非負數(shù)

(1)

cik-tik+G(1-aihk)≥cih,i=1,2,…,n;h,k=1,2,…,m

(2)

cjk-cik+G(1-xijk)≥tjk,i,j=1,2,…,n;k=1,2,…,m

(3)

cik≥0,i=1,2,…,n;k=1,2,…,m

(4)

表1 參數(shù)定義

在自動化儀表檢測中，儀表在多個自動化檢測設備之間進行檢測，將儀表視為工件、自動化檢測設備視為機器，則其中的調(diào)度問題模型與作業(yè)車間調(diào)度問題模型一致。因此，對于作業(yè)車間調(diào)度模型的優(yōu)化方法同樣能夠應用于實際自動化儀表檢測工作中。

2 求解作業(yè)車間調(diào)度問題的精英鯨魚優(yōu)化算法

2.1 基本鯨魚優(yōu)化算法

鯨魚優(yōu)化算法(WOA)是一種群體智能優(yōu)化算法，該算法模擬了自然界中座頭鯨群的捕食過程，包括包圍獵物、氣泡網(wǎng)攻擊和搜尋獵物3種行為。鯨魚優(yōu)化算法使用座頭鯨代表所求可行域中的可行解，從一組隨機解開始，在每次迭代中，選取當前最優(yōu)解作為獵物位置，其余座頭鯨隨機選擇執(zhí)行其中一種行為來更新自身位置，實現(xiàn)向獵物的靠近，直至完成捕食過程。

2.1.1 包圍獵物

座頭鯨識別獵物的位置，并將其包圍的過程可以用以下數(shù)學公式進行表述

D=|C·X*(t)-X(t)|

(5)

X(t+1)=X*(t)-A·D

(6)

A=2a·r1-a

(7)

C=2·r2

(8)

其中，“·”表示點乘，t為當前迭代次數(shù)；式(5)中D表示座頭鯨與獵物之間的距離， “|·|” 表示取絕對值，X(t) 表示座頭鯨在第t次迭代中的位置，X*(t) 表示獵物的位置；式(6)表示座頭鯨個體位置更新的方式；A和C表示系數(shù)矢量，分別用式(7)和式(8)計算，其中，a在迭代過程中從2線性減少到0，r1和r2為[0,1]中的隨機向量。

2.1.2 氣泡網(wǎng)攻擊

在氣泡網(wǎng)攻擊過程中，座頭鯨首先會吐出氣泡將獵物包圍在氣泡中，然后通過收縮圍圈和螺旋上升同時進行的方式靠近獵物。收縮圍圈過程與包圍獵物過程相同。螺旋上升過程用以下數(shù)學公式進行描述

D′=|X*(t)-X(t)|

(9)

X(t+1)=D′·ebl·cos(2πl(wèi))+X*(t)

(10)

式(9)中D′表示座頭鯨與獵物之間的距離，式(10)中b是用于限定對數(shù)螺旋形狀的常數(shù)，l是在[-1,1]中的隨機數(shù)。

為了使收縮包圍和螺旋上升同步進行，座頭鯨在捕食過程中有50%的可能性通過收縮包圍或螺旋上升來進行更新位置。此行為用以下數(shù)學公式進行描述

(11)

式中：p是[0,1]中的隨機數(shù)。

2.1.3 搜尋獵物

為了保證鯨魚優(yōu)化算法的全局搜索能力，當 |A|≥1時，每個座頭鯨根據(jù)隨機選取的其它座頭鯨的位置來更新自身位置信息。此行為用以下數(shù)學公式進行描述

D=|C·Xrand-X(t)|

(12)

X(t+1)=Xrand-A·D

(13)

其中，Xrand是隨機選取座頭鯨的位置。

當 |A|<1時，選擇當前最優(yōu)解作為獵物位置來更新其余座頭鯨的位置，并且根據(jù)p的值進行收縮包圍或螺旋上升。

2.2 基于生命力選擇的精英鯨魚優(yōu)化算法

鯨魚優(yōu)化算法具有結(jié)構(gòu)簡單、求解精度高、穩(wěn)定性好的特點，但仍存在收斂速度慢、容易陷入局部最優(yōu)等問題。因此，本文提出了精英鯨魚優(yōu)化算法，使用生命力選擇方法來確定種群中連續(xù)表現(xiàn)較差的個體，并使用表現(xiàn)較好個體的變異體或隨機生成的新個體進行替換。該方法可以及時淘汰連續(xù)表現(xiàn)較差的個體，避免其浪費過多的計算成本；反之，對連續(xù)表現(xiàn)較好的個體，則進行深入挖掘，提升種群搜索效率。通過反復迭代，不斷更新鯨群位置，使得種群逐漸向最優(yōu)解收斂。

2.2.1 算法框架

精英鯨魚優(yōu)化算法框架如算法1所示。

算法1：精英鯨魚優(yōu)化算法

輸入：生命力上限Vmax，生命力下限Vmin，選擇壓力Sp

(1)初始化種群位置Xi(i=1,2,…N)

(2)計算每個個體的適應度

(3)令X*為當前最佳適應度

(4)while(t<最大迭代次數(shù))

(5) fori=1 toNdo:

(6) 生成[0,1]中的隨機數(shù)r1,r2,p

(7) 按式(7)、 式(8)計算A,C的值

(8) if1(|A|≥1)

(9) 根據(jù)式(13)更新第i個個體位置(搜尋獵物

(10) 過程)

(11) else if1(|A|<1)

(12) if2 (p<0.5)

(13) 根據(jù)式(6)更新第i個個體位置(包

(14) 圍獵物過程)

(15) else if2(p≥0.5)

(16) 根據(jù)式(10)更新第i個個體位置(螺

(17) 旋上升過程)

(18) end if2

(19) end if1

(20) end for

(21) 計算每個個體當前的適應度

(22) 根據(jù)個體適應度的變化，計算個體生命力

(23) 根據(jù)算法2，使用生命力選擇方法替換生命力較

(24) 差的個體

(25)end while

2.2.2 編碼方法

由于精英鯨魚優(yōu)化算法為連續(xù)函數(shù)的優(yōu)化算法，而作業(yè)車間調(diào)度問題為離散的組合優(yōu)化問題，因此，應對精英鯨魚優(yōu)化算法中輸出的種群進行離散化處理。本文首先將連續(xù)編碼進行升序排列，然后將連續(xù)編碼分配給按照序號升序排列的工件，最后將工件序號填充到連續(xù)編碼的原位置，得到對應的離散編碼。在該編碼方案下，所有生成的解均為可行解，每個離散編碼序列代表一個可行的調(diào)度方案，由工件編號組成，每個工件序號代表該工件工序的順序執(zhí)行。

考慮3個工件3臺機器的作業(yè)車間調(diào)度問題，每個工件都需在3臺機器上執(zhí)行一道對應的工序，工件的加工順序見表2。以算法中某次的連續(xù)編碼(0.15-0.26-0.36-0.69-0.34-0.59-0.06-0.77-0.82)為例，首先將其升序排列為(0.06-0.15-0.26-0.34-0.36-0.59-0.69-0.77-0.82)，然后將(0.06-0.15-0.26)分配給工件一，將(0.34-0.36-0.59)分配給工件二，將(0.69-0.77-0.82)分配給工件三，最后將工件序號填充到連續(xù)編碼的原位置，即可得到離散編碼(1-1-2-3-2-2-1-3-3)，由于每個工件都有其固定的工序執(zhí)行順序，則編碼中工件序號的第N次出現(xiàn)，可以表示該工件的第N個工序。則該離散編碼對應的執(zhí)行順序為(O12-O13-O21-O31-O22-O23-O11-O33-O32)，具體過程如圖1所示，甘特圖如圖2所示。

本文使用一組從-100到100隨機生成的個體作為初始種群。

表2 3×3問題加工順序

圖1 編碼的離散化和對應的調(diào)度方案

圖2 離散編碼(1-1-2-3-2-2-1-3-3)對應調(diào)度順序的甘特圖

2.2.3 生命力選擇方法

生命力選擇方法是一種廣泛適用于優(yōu)化問題的選擇算子[15]，其原理為，使用個體在迭代中取得的進步作為個體的評價標準，如果一個個體在一代中取得了進步，則對該個體進行保護，如果該個體在長期迭代中未能得到進步，則可能被放棄。生命力選擇方法通過識別和保護具有旺盛生命力的個體，在種群早熟時拋棄雖然具有較好適應度，但是難以繼續(xù)取得進步的個體，并提供產(chǎn)生新個體的機會，加強了全局搜索范圍、避免了算法提前收斂。

在搜索過程中，種群每個個體的初始生命力被設置為floor((Vmax+Vmin)/2)， 其中floor(·)為向下取整函數(shù)。在每次迭代中，首先將適應度最高的個體生命力設置為Vmax，其余個體根據(jù)當前適應度和上一代適應度計算生命力，當前適應度小于上一代適應度時，個體生命力計算公式如式(14)所示；當前適應度大于上一代適應度時，個體生命力計算公式如式(15)所示。計算完種群中所有個體的生命力后，將根據(jù)個體生命力做出調(diào)整，生命力較低的個體將被替換。具體算法過程如算法2所示

(14)

(15)

算法2：生命力選擇方法

輸入：生命力上限Vmax，生命力下限Vmin，選擇壓力Sp

(1)將種群中適應度較高的一半組成集合PopA， 適應

(2)度較低的一半組成集合PopB

(3)對于PopA中的每個個體:

(4) if1Vi==Vmin

(5) 生成一個隨機數(shù)r3∈(0,1)

(6) if2r3

(7) 對該個體執(zhí)行逆序變異操作

(8) else if2r3≥Sp

(9) 對種群中適應度最高的個體執(zhí)行逆序變

(10) 異操作，并替換掉當前個體

(11) end if2

(12) end if1

(13)對于PopB中的每個個體:

(14) if3Vi==Vmin

(15) 生成一個隨機數(shù)r4∈(0,1)

(16) if4r4

(17) 從PopA中隨機選擇一個個體執(zhí)行逆序

(18) 變異操作，并替換掉當前個體

(19) else if4r4≥Sp

(20) 隨 機生成一個新個體，并替換掉當前個體

(21) end if4

(22) end if3

其中，算法中的逆序變異操作指在個體中選擇兩個位點，將兩個位點間的編碼進行逆序重排。如圖3所示，隨機選取編碼中的位置3和位置7，則所選連續(xù)編碼為(0.36-0.69-0.34-0.59-0.06)，對應離散編碼為(2-3-2-2-1)，則轉(zhuǎn)換得到的調(diào)度方案為(O21-O31-O22-O23-O11)。對位置為3-7之間的元素進行逆序重排，則得到的連續(xù)編碼為(0.06-0.59-0.34-0.69-0.36)，由于離散編碼與連續(xù)編碼的大小對應，所以得到的離散編碼為(1-2-2-3-2)。因為調(diào)度方案與離散編碼并非是一一對應關(guān)系，而是需要根據(jù)工件的實際工序做出調(diào)整，所以得到的調(diào)度方案并非是原方案的逆序重排，而應是與離散編碼相對應的(O11-O21-O22-O31-O33)。

圖3 逆序變異過程

2.2.4 性能分析

圖4顯示了本文提出EWOA算法的搜索過程。使用二維平面表示問題的搜索解空間，灰色五角星代表真實最優(yōu)解的位置，圓圈代表種群中的個體。在算法初始階段，對種群進行初始化，使所有個體均勻分布于搜索空間(圖4(a))。在WOA搜索階段，將當前最優(yōu)個體標記為最優(yōu)解位置(圖4(b)中灰色圓)，所有個體隨機的對其進行螺旋包圍，或向其它個體靠近，進行全局的搜索。而在一次WOA搜索過程結(jié)束后，使用生命力選擇算法對種群進行進一步的優(yōu)化，淘汰掉長期無進步的個體(圖4(c)中帶網(wǎng)格的圓)，使用新生的個體(圖4(c)中黑色圓)對有搜索潛力的區(qū)域做進一步挖掘，避免提前收斂的同時提高局部的搜索性能。EWOA算法引入的生命力搜索雖然在每一代的搜索中增加了對生命力的計算工作，但避免了由無搜索潛力的個體帶來的計算成本，所以并不會造成過大的計算負荷，可以在有限的計算資源中提升算法的搜索性能。

圖4 EWOA算法搜索過程

因此，本文提出的EWOA算法從理論上提升了傳統(tǒng)WOA算法的性能，能夠更好地適用于優(yōu)化問題的求解。

3 實驗仿真與討論

3.1 實驗設計

為測試精英鯨魚優(yōu)化算法求解作業(yè)車間調(diào)度問題的性能，本文使用了作業(yè)車間調(diào)度問題標準數(shù)據(jù)集和北京東方計量測試研究所自動化儀表檢測數(shù)據(jù)集，進行了仿真實驗。標準數(shù)據(jù)集包括OR-Library(http://people.brunel.ac.uk/～mastjjb/jeb/orlib/files/jobshop1.txt)中的3個FT實例、4個LA實例和5個ABZ實例，共計12個實例。北京東方計量測試研究所自動化儀表檢測實驗室由10個機械臂自動檢測環(huán)節(jié)O1～O10組成，30個不同的待檢測件均需按不同的順序要求完成相應檢測工作。圖5顯示了實驗室布局情況，和其中兩個執(zhí)行不同檢測任務的自動檢測環(huán)節(jié)實物圖。本文選取該實驗室實際檢測中的10組數(shù)據(jù)作為10個實例，使用編號P1-P10進行區(qū)分(https://github.com/ZikeW/JSSP.git)。

圖5 自動化儀表檢測流水線布局

本文選取近年來求解車間調(diào)度問題的高性能算法，即灰狼優(yōu)化算法[16](GWO)、擴展遺傳算法[17](EGA)、改進Jaya算法[18](IJA)和鯨魚優(yōu)化算法[10](WOA)作為對比算法。實驗中，所有算法的種群規(guī)模設置為50；GWO算法的鄰域搜索參數(shù)設定為10，局部搜索參數(shù)設定為30；EGA算法的交叉參數(shù)設為0.7，變異參數(shù)設為0.1；IJA算法的鄰域搜索參數(shù)設為0.9；EWOA算法的選擇壓力設定為0.7，生命力上限設定為11，生命力下限設定為1。所有算法的終止條件均設置為種群進化800代。

針對每個問題，每種算法都運行50次，通過比較50次Cmax的平均值(Mean)，來測試算法的有效性。此外，使用誤差棒來表示50次Cmax的標準差，來比較算法的穩(wěn)定性。

3.2 實驗結(jié)果與分析

如圖6所示，在作業(yè)車間標準數(shù)據(jù)集測試問題中，GWO、EGA和IJA等3種求解車間調(diào)度問題的成熟算法和EWOA均表現(xiàn)良好，有著較好的Cmax平均值，而WOA算法未能取得較好的結(jié)果。在對實驗結(jié)果進一步分析后可以發(fā)現(xiàn)，除WOA外的4種算法在所有實例中的搜索結(jié)果的最小值與文獻中記錄的實驗下限相同或接近。以實例LA31為例，文獻中記錄的最小值為1784，4種算法搜索到的最小值均為1784，GWO算法搜索到的平均值為1920.3，標準差為39.2；EGA算法搜索到的平均值為1870.7，標準差為41.6；IJA算法搜索到的平均值為1831.1，標準差為38.6；EWOA算法搜索到的平均值為1836.0，標準差為36.5。該實例上的實驗結(jié)果表明，幾種算法均能搜索到實例的最小值，而EWOA算法有著較小的Cmax平均值和最低的標準差。在所有實例的實驗結(jié)果中，EWOA算法實驗結(jié)果的平均值與文獻中記錄的最小值(BKS)較為接近，且有著較小的標準差，雖然部分算法的搜索性能超過了EWOA算法，但EWOA算法的綜合表現(xiàn)仍是最優(yōu)的。以上實驗結(jié)果說明：①EWOA算法相對于傳統(tǒng)的WOA算法能夠在相同的迭代次數(shù)內(nèi)獲得更好的實驗結(jié)果，有著顯著的搜索性能提升；②盡管有很多算法在標準車間調(diào)度問題數(shù)據(jù)集中取得了優(yōu)秀的結(jié)果，EWOA算法仍是非常有競爭力的算法之一，有著相對更好的有效性和穩(wěn)定性。

圖6 標準實例實驗結(jié)果

圖7顯示了在北京東方計量測試研究所自動化儀表檢測數(shù)據(jù)集中的實驗情況。WOA算法未能在終止條件內(nèi)完成收斂，有著極大的平均值和標準差，不能適用于問題的求解。在標準數(shù)據(jù)集中有著優(yōu)秀表現(xiàn)的GWO、EGA和IJA等3種算法，表現(xiàn)都較為普通。而EWOA算法在所有問題中都取得了最好的Cmax平均值和標準差，相對于其它算法有著明顯較高的有效性和穩(wěn)定性優(yōu)勢。以上實驗結(jié)果說明：

EWOA算法相比于GWO、EGA、IJA和WOA算法，更適合于作業(yè)車間調(diào)度實際問題的求解，可以有效應用于自動化儀表檢測實際任務中，能夠提升工作效率，節(jié)省檢測成本。

圖7 自動化儀表檢測數(shù)據(jù)集實例實驗結(jié)果

兩組實驗結(jié)果表明，作業(yè)車間調(diào)度問題的復雜性對算法的有效性和穩(wěn)定性有一定影響，能夠適用于某種場景的算法未必能夠在另一種場景中有著同樣優(yōu)秀的表現(xiàn)。在自動化儀表檢測場景中，本文提出的EWOA算法能夠很好適應場景的復雜性，同時能夠滿足場景的有效性和穩(wěn)定性要求，可以應用于實際調(diào)度問題的求解。

4 結(jié)束語

鯨魚優(yōu)化算法是一種新興的智能優(yōu)化算法，多用于連續(xù)優(yōu)化問題，在離散的調(diào)度問題中應用較少。本文提出了一種精英鯨魚優(yōu)化算法，驗證了其在作業(yè)車間調(diào)度問題的可行性。在標準數(shù)據(jù)集和北京東方計量測試研究所的自動化儀表檢測數(shù)據(jù)集中進行了仿真實驗，驗證了提出算法的有效性和穩(wěn)定性，可以應用于自動化儀表檢測的實際工作中。

在未來自動化儀表檢測工作調(diào)度問題的研究中，儀表檢測任務臨時增減、工件運輸?shù)仁录绊憴z測工作的統(tǒng)籌調(diào)度，如何解決這些在線的自動化儀表檢測調(diào)度問題，是進一步的研究方向之一。