對具有二部圖特點的二次指派問題的研究

（蘭州交通大學 交通運輸學院，蘭州 730070）

0 引言

二次指派問題[1,2]（quadratic assignment problem，QAP）是組合優(yōu)化中最難的問題之一，現(xiàn)實生活中的許多問題都可以以它為形式化模型。它已被應用到校園規(guī)劃、醫(yī)院布局、車間調(diào)度等諸多生產(chǎn)服務領域，因此研究該問題具有很強的實際意義。

目前，對大規(guī)模的二次指派問題基本都依賴于智能算法求解。Izabela在文獻[3]中介紹了幾種求解QAP的方法，包括：遺傳算法、禁忌搜索、分支定界和模擬退火。Shigeyoshi[4]提出一種基于對稱多處理的并行蟻群優(yōu)化算法對二次指派問題進行求解。Guo[5]針對二次指派問題，將粒子群算法嵌入蟻群算法中，提出一種混合蟻群優(yōu)化算法。文獻[6]提出一種具有并行化的狡猾螞蟻系統(tǒng)來求解QAP，并將魯棒禁忌搜索融合到CAS中。李引珍等[7]通過分析模型特征，提出一種帶時間約束的分支定界算法。Dokeroglu等[8]在網(wǎng)格上提出一種高性能多起點超啟發(fā)式算法，用于解決QAP，它通過一種進化選擇機制來實現(xiàn)。文獻[9]基于足球的不同概念，提出了金球算法與模擬退火相結合的應用（GBSA）來求解QAP；由于不可接受的運動，模擬退火搜索在局部最優(yōu)中被阻擋，GBSA提出的策略指導模擬退火搜索逃離局部最優(yōu)并以有效的方式探索搜索空間。上述文獻表明，智能優(yōu)化算法對求解大規(guī)模QAP的優(yōu)勢更為明顯。

因蟻群算法具有并行化、正反饋、自組織、魯棒性強等先天優(yōu)越性，所以蟻群算法是用于求解QAP最有效的算法之一。文獻[10～13]的研究內(nèi)容涉及最優(yōu)迭代最大最小螞蟻算法[10]、目標引導蟻群優(yōu)化算法[11]、預處理快速螞蟻系統(tǒng)[12]以及改進蟻群算法[13]等多種蟻群算法。戢守峰等[14]提出一種混合分布估計算法（HEDA），針對HEDA概率模型，提出一種概率矩陣初始生成機制和擾動操作，用于提高算法的全局搜索能力。然而，在現(xiàn)有的文獻中為簡化模型減少參數(shù)，而將啟發(fā)式因子排除只考慮信息素，使得算法前期收斂速度過快，容易陷入局部最優(yōu)；或把距離因子與流量因子的乘積當作啟發(fā)式因子，雖然體現(xiàn)了貪婪的思想，但并未反映QAP的本質特點；在生成初始解及種群初始化時，現(xiàn)有的文獻采取隨機生成，因為不能產(chǎn)生質量較好的初始化種群會導致收斂速度過慢；另外，在現(xiàn)有的大多數(shù)文獻中，在算法的后期解的多樣性丟失，使其最后得到的最優(yōu)解喪失全局代表性。

基于以上分析，文章結合QAP的本質特點，提出一種基于熵收斂的改進蟻群算法（EC-IACA），用于求解具有二部圖特點的二次指派問題。該算法主要貢獻有4點：1）利用了QAP的距離和流量矩陣，設計了基于貪婪思想和設施交互次數(shù)相結合的方法生成初始種群，提高初始解質量，并在每輪迭代時，用上輪得到的最優(yōu)解來設置初始值，減少每輪搜索的盲目性；2）針對QAP的特點，使用偽隨機比例規(guī)則，并引入先驗概率，設計了具有QAP特點的狀態(tài)轉移策略，包括自適應信息素更新策略、構造動態(tài)啟發(fā)式信息以及對參數(shù)ɑ和β的動態(tài)調(diào)整；3）局部搜索策略采用2-opt的鄰域結構，并在算法后期引入信息素軌跡平滑機制來增強解的多樣性，用于提高算法的局部搜索能力、提高解的質量、并使最后得到的最優(yōu)解具有全局代表性；4）引入“信息熵”的概念，通過對信息素的熵收斂來分析整個算法的收斂性。這些均是已有的研究鮮少涉及到的。

1 二次指派問題的描述

本文對具有二部圖特點的QAP描述如下：V1={1，2，…，n}為n項任務或設施的集合，V2={n+1，n+2，…，2n}為n個機器或者位置的集合，E為V1中節(jié)點與V2中節(jié)點之間的連接構成的集合。設G=（V1，V2，E），這樣二部圖G中每一條邊對應著某項任務或設施的指派。對該問題的求解目標是考慮將一組設施指派至一組位置，使其各運行單元間的無論和運輸成本最??；D=[dij]n×n為距離矩陣，dij為位置i與為位置j之間的歐氏距離；F=[fij]n×n為設施或任務間的物流量矩陣，fij為設施i與設施j之間的流量。數(shù)學表達式如下：

其中：p={p（1），p（2），…，p（n）}是1～n的一個排列，p（i）和p（j）表示位置i和j所指派的設施或任務；Π是所有排列構成的集合。

2 改進蟻群算法設計

2.1 編碼

根據(jù)二次指派問題的特點，考慮到設施和位置的雙映射，編碼方式采用有序串實數(shù)編碼。如圖1所示。

圖1 問題編碼方式

2.2 具有QAP特點的初始解構架策略

在算法的搜索過程中，初始種群距最優(yōu)解的距離直接影響著算法總的運行時間。根據(jù)距離矩陣D=[dij]n×n和流量矩陣F=[fij]n×n的信息，有以下兩種策略：

策略1：提出一種基于假定應急物流中心的貪婪法規(guī)則。首先計算δj，即位置j到其他所有位置的距離之和。

其次，對于流量矩陣F也可以計算出設施i與所有設施之間的流量總和ωi，即：

近似的認為，δj值越小，位置j與假定應急物流中心越近，根據(jù)貪婪思想，我們總是想在靠近假定應急物流中心的位置指派流量大的任務或設施，通過將經(jīng)過升序排序之后的設施Ri對應指派給經(jīng)過降序排序之后的Lj可得到一組完整的初始解。

策略2：從未指派的設施中選取一個與最近已指派的設施交互次數(shù)最多的設施；然后，將它指派至能夠使得總成本最小增加的空余位置。重復此步驟，直到所有的設施被指派完畢，可獲得一組完整的初始解。

選用以上兩種策略中的任意一個，都可獲得初始解。策略1的優(yōu)點是運用貪婪的思想將設施或任務的的最大流量指派至離假定物流中心最近的位置；缺點是這種指派方法未考慮對迄今為止已指派的設施所產(chǎn)生的影響。策略2的優(yōu)點是選擇與最近已指派設施交互次數(shù)最多的設施，并對其進行指派，使得目標函數(shù)的增量最小；缺點是設施的選擇既沒有考慮到它與所有其他設施的相互作用，也沒有考慮到所選擇的位置對迄今尚未指派的設施可能具有的重要性。因此。文章考慮將這兩種策略相結合來產(chǎn)生一個完整的初始解，進而通過隨機互換操作來產(chǎn)生初始種群。初始解的產(chǎn)生見算法1，種群初始化偽代碼如圖2所示。

圖2 種群初始化偽代碼

2.3 具有QAP特點的初始位置的選擇

現(xiàn)有的求解二次指派問題的蟻群算法或其他啟發(fā)式算法在構造解的過程中，往往會隨機選擇一個任務或者設施將其指派至位置1，而這種方法明顯存在不足：如果所有的螞蟻都從同一初始位置1開始迭代，會使得搜索的效果明顯下降[10]；同時，如果在位置1指派了錯誤的設施或任務則會導致此次迭代是無效的。因此，本文根據(jù)QAP的特點結合蟻群算法并行性特點對初始位置構造進行設計。

從第二代開始，構造初始位置時，從當前最優(yōu)解p*中隨機選擇一個指派對作為本輪的起點，而不是隨機選擇設施或任務。采用此方法，螞蟻走過的路徑為可表示為：

Location：i→…→n→1→…→j→…→i-1；

Facility：pt+1（i）→…→pt+1（n）→pt+1（1）…→pt+1（j）→…→pt+1（i-1）。

改進的設計在構造解初始位置過程中，丟失全局最優(yōu)解的概率會明顯降低；同時由于對起點位置的正確選擇避免了螞蟻選擇初始任務的盲目性，加強了正反饋效果。

2.4 具有QAP特點的狀態(tài)轉移策略設計

本文在狀態(tài)轉移過程中使用偽隨機比例規(guī)則，即確定性選擇和隨機性選擇相結合的選擇策略。在構造可行解過程中，為建立當前路徑與探索新路徑之間的平衡，本文引入先驗概率（prior probability）q0。具體描述如下：

當處于第t代螞蟻k為設施i指派位置j時，先隨機生成q，如果q的值小于等于q0，則從設施i到所有可行的位置j中找出{[τis（t）]α[ηij（t）]β}最大的位置j為設施i指派，如式(4)所示；如果隨機數(shù)q大于q0，螞蟻k則按輪盤賭選擇策略按照概率逐一將各個設施指派至不同的位置，將設施i指派至位置j上的概率為，如式(5)所示。

其中：q是均勻分布在區(qū)間[0，1]中的一個隨機變量；q0為（0，1）間的常數(shù)；τij（t）表示第t輪迭代中將設施i指派至位置j上的信息素濃度；ηij（t）表示第t輪迭代中將設施i指派至位置j上的啟發(fā)式信息；α表示信息素權重，β表示啟發(fā)式因子相對權重；allowedk表示螞蟻k當前可選位置的集合；J是根據(jù)以上公式給出的概率分布產(chǎn)生出來的一個隨機變量；Jk（i）表示允許放置設施i的位置集。

2.4.1 啟發(fā)信息的改進

對蟻群的求解過程分析可知，在搜索過程中已得到的部分解和當前未遍歷完指派集的狀態(tài)以及啟發(fā)信息具有一定的聯(lián)系：所選路徑包含的有效指派對和部分解會影響解的適應度，而且未遍歷完的組合集以所選設施點為起點的有效指派對的多少也很大程度影響了后續(xù)搜索過程中所產(chǎn)生的新的有效指派對的機會。由于在搜索過程中部分解是動態(tài)變化的，所以文章提出用動態(tài)啟發(fā)式信息的方法來引導螞蟻求解，計算公式如式(6)所示。

式(6)表示在構造解的過程中，邊Eij啟發(fā)信息會被動態(tài)修改，它不但依賴于某時刻螞蟻k構建的部分解，而且還會涉及邊Eij在后續(xù)搜索過程中潛在的作用，因此，有兩部分組成。其中：||A||表示未遍歷結束的有效指派集中指派對的個數(shù)；Nk（Ri，j）表示螞蟻k在到達節(jié)點Ri之后，螞蟻所走過的邊與Eij組合所產(chǎn)生的有效指派對的個數(shù)，計算如式(7)所示；Uk（Ri，j）表示以邊Eij為入邊，節(jié)點Ri之后的各邊為出邊所組成的有效指派對個數(shù)，計算如式(8)所示。

其中：nu表示第u個設施點Ru出邊的個數(shù)；C（Eij，Euv）表示當指派對（Eij，Euv）屬于未遍歷完的有效指派||A||時，C（Eij，Euv）=1，否則為0（其中j∈ni，v∈nu）；μ和υ表示權重。

2.4.2 具有QAP特點的信息素更新策略

因為是將信息素布置于邊（i，j），即將任務i置于位置j的邊，所以信息素的分布空間也是一個二部圖。為使螞蟻在搜索過程中既不過分集中也不過分分散，本文運用一種新的信息量更新策略來根據(jù)解空間情況對各路徑上信息量進行自適應動態(tài)調(diào)整，進而避免了早熟和局部收斂，提高全局搜索能力。信息素更新公式為：

其中：m為改進解的個數(shù)，c為常數(shù)；τmax和τmin根據(jù)文獻[15]確定。

其中：Q表示信息素強度，為常數(shù)；Lgb表示當前最優(yōu)解對應得的目標值，表示當前最優(yōu)解對信息素的增量。

2.4.3 參數(shù)α、β的自適應調(diào)整

現(xiàn)有文獻中，為減少算法的復雜性，參數(shù)α、β取值會提前給定，并在搜索過程中始終保持不變；但是在同一問題不同階段α、β的不同取值會促進算法收斂，同時也會避免算法陷入局部最優(yōu)。為此本文提出對α、β的自適應調(diào)整，如式(11)所示。

螞蟻每搜索完一次，通過設定參數(shù)θ來調(diào)節(jié)參數(shù)α、β在路徑選擇中的比例。

2.5 具有QAP特點的局部搜索策略

為提高算法的局部搜索性能，本文使用最簡單的2-opt鄰域結構，因為2-opt是一種無記憶算法，在每次迭代中達到最優(yōu)解的概率是相同的，而改進之后的蟻群算法能夠從最優(yōu)解中獲取信息，并在個迭代中很好的利用這一信息，因此將二者結合可增強算法的整體性能。

現(xiàn)分析QAP的2-opt鄰域結構，如圖3所示。

由式(1)知，解pgb的目標函數(shù)值的變化僅僅取決于進行了互換操作的元素。因此將式(1)分解，可得到目標函數(shù)的另一種形式：

圖3 2-opt鄰域結構示意圖

因此，通過交換兩個設施位置獲得的目標函數(shù)值的變化量Δ（pgb，i，j）由式(14)計算。

可根據(jù)Δ（pgb，i，j）與0的關系來判斷進行元素互換操作是否有意義，如果Δ（pgb，i，j）<0，則認為此次交換有效，并得到相應的目標函數(shù)值和解的編碼。同時在計算變化量時，只需進行16n-16次乘加操作，即采用此方法的算法復雜度為O（n）；而如果用式(1)來計算pnew，1目標函數(shù)值則需進行2n2-1次乘加操作，即算法復雜度為O（n2）。因此，使用2-opt鄰域操作可有效降低計算目標函數(shù)的算法復雜度，提高局部搜索的效率。當采用2-opt鄰域結構時，評價所有個體的算法復雜度為O（n3），而如果采用式(1)來評價所有個體的算法復雜度為O（n4）。

綜上分析，所提出的2-opt鄰域結構有利于提高局部搜索的效率，也有利于對可能存在優(yōu)秀解的pgb鄰域進行精細的搜索。

2.6 多樣性--信息素初始值的設定與重置

如果不考慮當前最優(yōu)解對信息素矩陣的增強，當所有螞蟻都集中在同一個指派方案上時，信息素矩陣中的最大值將趨向于Q×n/（1—ρ）。因此，設定信息素的初始值為τ0=Q×n/（1—ρ）。

如果算法求得的最優(yōu)解在N次循環(huán)內(nèi)沒有明顯改進時，則表明此時算法的多樣性丟失，通過對信息素重設置來提升算法的多樣性和提高解的搜索空間，參考最大最小蟻群系統(tǒng)（MMAS）中的信息素軌跡平滑機制。更新公式如式(15)所示。

2.7 算法的收斂過程分析

本節(jié)通過計算信息素濃度的熵來論證算法的收斂性，即“熵收斂”。因為當算法的解空間分布趨于穩(wěn)定時，算法的系統(tǒng)信息熵處于收斂狀態(tài)，因此，通過對算法的熵收斂進行分析進而判斷出解收斂。在蟻群優(yōu)化算法中，因為反映解的分布和信息素概率分布特征的兩個統(tǒng)計指標期望和方差分別近似相等，因此兩者的收斂代數(shù)近似相等，又因為反映信息素概率分布趨于穩(wěn)定的依據(jù)為信息熵I（t）所以使得算法的最優(yōu)迭代次數(shù)與信息素概率分布穩(wěn)定的收斂代數(shù)近似相等。

對于二次指派問題，信息熵I（t）可定義為：

其中：pij（t）表示狀態(tài)發(fā)生的概率，0≤pij（t）≤1，定義為：

當所有的τij（t）值處于相同值時，即在信息素的初始化階段我們獲得的I（t）的上限I+；當所有的τij（t）值都收斂于τmin或τmax時，即在一種極端情況下我們獲得I（t）的下限I-。在分析過程中我們使用歸一后的信息熵IN（t）。

其中：r=τmax/τmin；IN（t）∈[0，1]。

3 仿真計算實驗與靈敏度分析

仿真計算實驗的算例從QAPLIB數(shù)據(jù)集[16]中選取，編程軟件為MATLAB，所有算法均獨立重復運行10次；基本參數(shù)設置：n為問題規(guī)模，q0=0.4，c=20，pop_size=n，ɑmax=8，βmax=10，初始階段ɑ=2、β=5，ρ=0.1，ε=0.5，μ取0.6，υ取0.4，N=100，pbest=0.005，Q=100，進化代數(shù)為1000代。

3.1 靈敏度分析及熵收斂性

以Bur26c為例，通過在算法搜索過程中ɑ、β值的動態(tài)變化繪制曲線圖，如圖4、圖5所示；并同時繪制該算例的收斂圖，如圖6所示。通過對比分析三個圖，可知在100、430代時由于ɑ、β波動幅值較大，可以看出最優(yōu)值有明顯變化；在600代之后隨著最優(yōu)值收斂，ɑ、β在一定范圍內(nèi)波動。因此，在文章設計的算法中ɑ位于[3，5]、β位于[3，4]較好。

QAPLIB數(shù)據(jù)集根據(jù)文獻[17]被分為4類，文章從每一類中隨機選取一個算例，根據(jù)式(13)～式(15)所示，計算I（t），并繪制熵收斂圖，如圖7所示。從圖7可知這4個算例的熵變化規(guī)律，Nug30和Tai30a在300代左右時，熵值開始趨于穩(wěn)定；Tai40b和Kra30b分別在500和800代左右時，熵值開始趨于穩(wěn)定。因此，此時算法的解空間也開始趨于穩(wěn)定即已收斂。

圖4 α值變化趨勢

圖5 β值變化趨勢

圖6 Bur26c收斂圖

圖7 熵收斂圖

3.2 具有QAP特點初始解構架的性能驗證

為驗證所提出的構造初始種群的方法比隨機初始種群更為有效，本文選取兩個指標：目標函數(shù)值的均值（avg）和最好值（best）。算例隨機選取，算法運行10次，結果如表1所示。

由表1知，使用所提出的具有QAP特點初始解的構架策略可明顯改善初始種群的質量，因此可有效提高算法的搜索效率。

表1 兩種初始種群的比較

3.3 EC-IACA算法的總體性能

為了觀察EC-IACA算法的總體性能，考慮將ECIACA與9組算例的已知最好解進行比較，比較結果如表2所示。為表示方便，表中誤差的單位為%。

表2 EC-IACA與已知最優(yōu)解比較結果

由表2知，在與已知最優(yōu)解比較中，除Wil100和Tho159外，其余算例的誤差和平均誤差都控制在1%之內(nèi)。所以用設計的EC-IACA求解QAP是有效的，解的質量是可接受的。

為與其他算法比較EC-IACA的性能，下面用ECIACA算法與文獻[10]、文獻[12]中所提算法進行比較，如表3所示。從表3的實驗結果可知，對于這四類算例，文章設計的EC-IACA算法的解的質量總體上要優(yōu)于文獻[10]和[12]，例外的只有Tai35a，Tai50a，Tai100b，Sko81，Ste36a。對比文獻[12]，文獻[12]對模擬現(xiàn)實問題和現(xiàn)實問題的求解有較大的優(yōu)勢，而EC-IACA算法在解的總體質量上也稍優(yōu)于文獻[12]；至于隨機產(chǎn)生無結構問題和隨機產(chǎn)生網(wǎng)格距離問題EC-IACA算法總體優(yōu)于文獻[12]。對比文獻[10]，EC-IACA算法全面優(yōu)于文獻[10]，例外的只有Kra30b。綜上，文章設計的EC-IACA算法在求解質量和穩(wěn)定性上具有明顯優(yōu)勢，是一種有效的優(yōu)化算法。

4 結語

通過詳細分析QAP的特點，文章設計提出了ECIACA算法，首先，為了提高解的初始種群質量，提出將基于應急物流中心的貪婪法規(guī)則與已指派設施的交互性相結合的策略用于生成初始種群，并利用得到的最優(yōu)解來設施初始值，避免了搜索的盲目性；其次，為避免算法過早的出現(xiàn)停滯、陷入局部最優(yōu)、收斂速度過慢以及提高解的質量和全局搜索能力，使用偽隨機比例規(guī)則，引入先驗概率，分別對對狀態(tài)轉移策略、信息素、啟發(fā)式信息的計算進行改進設計；然后，提出了一種2-opt鄰域結構的局部搜索策略，用于提高算法的局部開發(fā)能力；最后，引入信息素軌跡平滑機制，來增強解的多樣性與全局搜索的深度，并通過信息熵來分析算法的收斂性。仿真實驗和算法比較驗證了所提出的EC-IACA算法在求解的精度、魯棒性有效性方面與其他算法相比具有一定優(yōu)勢。下一步將繼續(xù)探蟻群算法參數(shù)的自適應性和設施間流量不確定下的二次指派問題。

表3 EC-IACA算法的測試結果