求解旅行商問題的混合粒子群優化算法

沈繼紅，王侃

(1.哈爾濱工程大學理學院，黑龍江 哈爾濱 150001;2.哈爾濱工程大學自動化學院，黑龍江 哈爾濱 150001)

優化問題可以自然分為2類:一類是連續變量的優化問題;另一類是離散變量的優化問題，即所謂的組合優化問題.旅行商問題(travel salesman problem，TSP)是組合優化問題中的一個著名NP難題，TSP因其典型性已經成為許多啟發式搜索、優化算法的間接比較標準.同時TSP也是一個具有廣泛的應用背景與重要理論價值的組合優化難題，對求解該問題高效的全局優化算法的研究，一直被科學界和工程界所高度重視.

TSP問題的求解方法歸納起來可以分為得到最優解的精確算法和找到近似解的近似算法.完全枚舉法、動態規劃法和全局搜索算法屬于精確算法.TSP問題精確算法的運行時間是指數級復雜度，難以適應大規模的實例，隨著對TSP問題的認識加深，精確算法的研究越來越少.近年來受到自然界的啟發，人們提出了各種各樣的計算智能方法，如人工神經網絡、遺傳算法、蟻群優化算法、粒子群優化算法和人工免疫系統等.智能優化算法為解決TSP問題提供了新的思路，它們被廣泛地應用于各種NP難題的優化問題求解，雖然不能保證獲取最優解，但在問題規模較大時也可以在可行時間內找到滿意的解.

粒子群優化算法(particle swarm optimization，PSO)是一種群智能優化方法，它是由美國社會心理學家Kennedy和電氣工程師R.Eberhart在1995年提出的，它利用了生物群體中信息共享的思想，其概念簡單、易于實現，同時又有深刻的智能背景，既適合科學研究，又適合工程應用.因此，PSO一經提出，就引起了眾多學者的關注，得到了非常廣泛的應用.為解決組合優化問題，Kennedy等［1］首先提出了PSO算法的離散二進制版;Clerc［2］提出了求解TSP問題的離散粒子群優化算法，對TSP問題的求解重新定義粒子的位置、速度和相關運算，但其性能與其他算法相比仍有不小的差距;高尚等［3］在粒子群算法中加入遺傳算法思想，構造了混合算法;Hendlass等［4］通過對離散PSO的每個粒子增加記憶功能，成功解決了一個小規模的TSP;王康平等［5］通過引入“交換子”和“交換序”的概念，給出了另一種解決TSP的PSO方法，為求解TSP問題提供了新的思路.但在算法的收斂速度方面以及在城市規模較大的情況下，現有文獻中的粒子群算法都存在著一定的缺陷，本文試圖通過與其他智能優化算法的結合來解決這一問題.光學尋優算法［6］是2007年沈繼紅教授提出的模擬光自然屬性的智能優化算法，利用在可行域中填充正方形介質，模擬光的折射以及反射現象，通過最基本光學定律找到最優值，算法迭代機理簡單、收斂速度快，具有很強的并行計算能力.2010年，李焱等［7］給出了基于正六邊形網絡的光學尋優算法，在高維迭代中具有良好的仿真效果，李加蓮等［8］給出了光學尋優算法的基本理論證明，并與其他算法進行了比較，完善了算法的理論體系.本文通過正方形網絡光學尋優算法的搜索機制形成初始粒子群，加入混沌優化的思想，并引入“交換子”概念，利用離散粒子群算法求解TSP問題，提出了一種新的解決TSP問題的光學混沌粒子群算法.

1 混合粒子群優化算法的構建

1.1 旅行商問題以及標準粒子群算法

TSP的描述十分簡單，即尋找一條最短的遍歷N個城市的路徑，其數學描述如下:

設有N個城市的集合C={c1，c2，…，cN}，每2個城市之間的距離為d(ci，cj)∈R+，其中 ci，cj∈C(1≤i，j≤N)，求使目標函數

達到最小的城市序列(c∏(1)，c∏(2)，…，c∏(N))，其中，∏(1)，∏(2)，…，∏(N)是1，2，…，N的全排列.

1998年，Shi等［9］給出了標準的 PSO 算法的數學描述:設搜索空間為D維空間，粒子數為n，第i個粒子的位置用 xi=(xi1，xi2，…，xiD)表示;第i個粒子的速度變化率用vi=(vi1，vi2，…，viD)表示;第i個粒子迄今為止搜索的得最好位置為pi=(pi1，pi2，…，piD)，記為pbest，整個粒子群迄今為止搜索到的最好位置為 pg=(pg1，pg2，…，pgD)，記為 gbest，對于每一次迭代，第i個粒子在第D維運動的表達式如下:

式中:c1、c2為正常數，稱為加速因子;rand()為［0，1］之間的隨機數;w稱為慣性因子.第d維的位置和速度的變化范圍為［－xdmax，xdmax］和［－vdmax，vdmax］，如果在某一維中迭代的xid、viid超過了取值邊界則按照邊界取值.

1.2 光學尋優算法

光學尋優算法借鑒了費馬定理，利用光在傳播過程中自動尋優的機制，將光的折射與反射原理與最優化的尋優過程聯系起來，給出一種新的最優化搜索算法.這種算法將坐標空間設想為填充了具有不同折射率的介質的空間，如圖1所示，將搜索路徑設想為光的傳播路徑，通過光的折射原理，使搜索方向趨向于目標函數值減小的方向，通過光的反射原理，改變搜索方向，使得搜索在折射無法進行時繼續下去.

圖1 在可行域中填充介質Fig.1 Filling medium in feasible region

針對以下問題進行研究:

式中:f(X)是正函數，即?(x，y)∈M，f(x，y)＞0.X是可行解，M是f(X)的可行域，R×R是二維實數空間.P(i)為第i次的搜索方向，h、τ為網格步長.對于第i次迭代，令搜索以P(i)方向在矩形分塊Di中搜索到點X(i)=(xi，yi)，并在X(i)點改變搜索方向，到達X(i+1)點，方向的迭代關系滿足:

式中:αi為Di中的入射角，αi+1為Di+1中的折射角，vi為光在Di中的傳播速度，可設定為Di中X(i+1)=(xi+1，yi+1)點的函數值.vi+1為光在Di+1中的傳播速度，可設定為Di+1中X(i+1)=(xi+1，yi+1)點的函數值.

式中:αi為Di中的入射角，αi+1為Di+1中的反射角.圖3～4是加入了多個臨界面的水平以及豎直方向的搜索路徑.光學尋優算法能快速找到最優搜索路徑.

圖2 基于反射和折射搜索方向的更新Fig.2 Updating searching direction based on reflection and refraction

圖3 增加多個臨界面光路的更新Fig.3 Updating light line after adding

圖4 增加水平豎直臨界面Fig.4 Adding horizontal and many critical surfaces vertical critical surfaces

1.3 TSP問題中的光學尋優思想

光學尋優算法遵循費馬原理，費馬原理表明，光在介質中從一點向另一點傳播時，總是沿著時間最少的路徑，光的傳播在不同介質中速度不同，運用光學尋優的思想可以很容易地得到粒子群的一組性能良好的初始粒子，大大減少算法的迭代次數.

定義1 當前城市密度.從當前城市到其他未被遍歷過的城市的最短路徑.

定義2 前沿城市密度.去除已經遍歷的城市，剩余城市路徑的平均值.

定義3 TSP問題中的反射.從當前城市開始，隨機選擇一條與未被遍歷城市的路徑.

定義4 TSP問題中的折射.從當前城市開始，選擇與未被遍歷城市之間最短的路徑.

光學尋優初始化初始值的過程如下:

1)隨機選擇一個城市作為出發點，選擇與當前城市之間的最短路徑城市作為下一個遍歷點.

2)計算當前城市密度與前沿城市密度.

3)如果當前城市密度小于前沿城市密度，則發生折射，選擇與未被遍歷城市之間最短路徑城市作為下一個遍歷點;如果當前城市密度大于前沿城市密度，隨機選擇一個未被遍歷的城市作為下一個遍歷城市.

4)重復3)的過程直到所有的城市都被遍歷.

5)將每一個城市都作為起點，依據光學尋優的思想形成N個城市序列，作為混沌粒子群算法的初值.

2 光學混沌粒子群算法

2.1 TSP問題中的混沌粒子群算法

2.1.1 TSP問題中混沌粒子群算法的相關定義

定義5 交換子.2個城市序列［10］Xi=［xi1xi2… xim］與Xj=［xj1xj2… xjm］，如果2個序列在相同的位置，數值不相同，即 xia≠xja，稱(xia，xja)為城市序列的交換子，即為Vij(xia，xja).

定義6 交換序列.由交換子組成的序列V=［V1V2… Vn］，其中n為2個城市對應序列相同，但數值不同的位置個數.

定義7 粒子的位置.粒子的位置是由城市序列X=［X1X2… Xm］表示，m為城市的個數;

粒子的速度.粒子的速度V=［V1aV2b… Vmn］，其中Vmn表示交換子，速度為交換序列.

2.1.2 交換子與交換序列的運算法則

1)位置與交換子的加法.

位置與速度的加法形成新的城市序列:設X=［X1X2… Xm］為城市序列，Vij(Xi，Xj)為交換，則

X=［X1X2… XjXiXm］為新形成的城市序列.

例1:

2)位置與位置的減法.

位置與位置的減法形成交換序列即生成新速度:Vij=Xi－ Xj，其中 Xi、Xj為城市序號.先找到與第1個城市序列中第1個元素相同的第2個城市序列位置，形成交換子v(1，i)，然后將此交換子作用在第1個序列上得到新的第1個序列，再找到新的第1個城市序列與第2個城市序列數值相同的第1個位置，形成交換子v(2，i)，依次進行下去，得到2個城市序列的交換序列.

例2:

3)交換子的數乘.

速度的數乘具有概率意義，例如Via=c·Vjb，其中c∈［0，1］是一個常數，在計算Via時，對Vjb中的每一維速度Vjn生成一個(0，1)的隨機數.

4)混沌思想.

通常一類非常簡單卻又廣泛應用的混沌系統是Logistic 映，其定義如下［11］:

式中:zk為實值序列，u為參數，研究表明當3.571 448≤u≤4時，該混沌映射處于混沌狀態，把3.571 448≤u≤4稱為混沌區域.Logistic混沌映射所生成的序列具有如下混沌特性:①非周期的序列;②該混沌序列不收斂;③zk可以遍歷整個(0，1)區域.

本文取u=4，設城市數量為n，按照順序隨機地生成(0，1)的n個隨機數，構成向量序列Z(k)=［z1(k)z2(k)… zn(k)］，將 z1，z2，…，zn按照從小到大排列1，2，…，n的下標號構成城市序列Xi=［Xi1Xi2… Xin］，按照式(1)生成向量，然后對Z(k+1)中元素從小到大進行排列下標號1，2，…，n構成新的城市序列 X(i+1)=［X(i+1)1X(i+1)2… X(i+1)n］.混沌遍歷的引入，有效地增加了城市序列的多樣性.

例3:

生成的城市序列為Xi=［4 1 2 3 4 6］，運用式(3)，

生成新的城市序列為Xi+1=［6 4 3 5 1 2］.

針對TSP問題，結合混沌思想的粒子群的更新公式變為:

式中:?1、?2為(0，1)的數，f(x)為 xi(k)向量從小到大排列下標形成的向量函數，Xi(k)為當前城市序列，Xi(k+1)為新生成的城市序列，Pbesti為當前個體最好序列，Gbesti為全局最好序列.Pbesti－ Xi(k)=Vp(k)，Gbesti－Xi(k)=VG(k)表示為交換序列，?1(Pbesti－Xi(k))表示當概率小于?1時發生交換操作，當概率大于?1不進行操作，城市序列保持不變.

2.2 光學混沌粒子群算法步驟

1)利用1.3中光學尋優的思想初始化粒子群，從每一個城市出發，得到N個初始值良好的城市序列;隨機生成混沌序列Z(1)并運用函數f(x)得到城市序列Xf(1).

圖5 混合粒子群算法的流程Fig.5 A flow sheet of mixed particle swarm algorithm

2)如果滿足最優條件，最短城市距離不再變——或者達到最大迭代次數轉到5)，如不滿足條件對初始種群進行更新，找到個體最好位置Pbesti以及全局最好位置Gbesti.

3)根據式(4)更新城市序列;

①根據式(1)和(3)計算混沌生成序列Z(k)并運用函數f(x)得到新的城市序列Xf(k);

②運用例2的思想計算交換序列，Pbesti－ Xi(k)=Vp(k)，Gbesti－ Xi(k)=VG(k);

③根據式(2)計算 φ1(Pbesti－Xi(k))+φ2(Gbesti－Xi(k));φ1、φ2為執行操作的控制概率.

④根據以上3個結果結合式(4)得到新的城市序列Xi(k+1).

4)迭代次數增加，更新粒子的位置轉到2).

5)輸出最優城市序列，并輸出最短距離.

2.3 混合算法時間復雜度分析

算法的時間復雜度是對算法運行時間的度量，用來表示算法的計算效率的高低.算法的時間復雜度的大小在一定程度上反映了算法性能的優劣.忽略硬件及環境因素，假設每次執行時硬件條件和環境條件是完全一致的.設城市數量為n，運行迭代次數為m，則光學混沌粒子群算法的時間復雜度量級為O(m(2n2+1)).

以下來說明該算法的時間復雜度數量級為O(m(2n2+1)).根據混合算法的特點首先應用光學尋優算法初始化城市序列，并應用混沌方法生成新的城市序列，然后應用粒子群算法的核心思想不斷地更新城市序列直到滿足條件.

1)以一次迭代為例，城市數量為n，隨機選取城市作為城市序列初始點，并應用光學尋優算法進行初始化，需要n(n－1)次運算，所以時間復雜度為O(n2－n).

2)應用混沌思想更新城市序列，對N個序列都進行更新，時間復雜度為O(n).

3)n個序列用來計算適應度函數的時間復雜度為O(n)，在此基礎上進一步確定個體極值以及群體極值，計算交換子序列從而對每一個城市序列進行更新，此步驟一共運行的時間復雜度為O(n2).

4)在一次迭代完成之后判斷是否達到終止條件，操作的時間復雜度為O(1).

通過以上的分析可以得出1)～4)的時間復雜度為O(2n2+1)，假設整體算法的迭代次數為m，則混合算法的整體運行時間為O(m(2n2+1)).因此整體算法的時間復雜度與城市的序列以及算法的運行代數有關，如何在保證精確度的前提下減少迭代次數也就成了控制算法運行時間的關鍵因素.

2.4 光學混沌粒子群算法收斂性分析

根據混沌序列的更新以及例3可得f(xi(k))是線性映射，可以寫成城市序列與交換子的線性組合的形式:f(xi(k))=X(k)+βv(k).β∈［0，1］，線性系統變成如下形式:

為了方便計算與分析，首先將空間簡化為一維空間，將模型轉化為式(5):

令 ?=?1+?2，記 Pbesti=pi(t)，Gbesti=pg(t)，并對模型進行簡化可得表達式:

整理記為式(6):

通過標準粒子群算法的數學描述可以將系統變為離散線性系統，這里假定在t次迭代之后粒子找到最優位置時，pbest、gbest將保持不變，式(6)中pi(t)、pg(t)不隨著時間變化.

定理1當w＜1+β，β＜?1+?2＜2w+2+β時，線性定常離散系統(6)漸近穩定，并且系統收斂.

證明通過系統(6)以及模型(5)的表達式可將xi(t)消去得到如下的差分表達式:

對式(7)求特征方程可得

為了對式(7)進行穩定性分析，用雙曲線性變換將離散系統轉化為線性系統，令帶入式(7)得:

根據勞斯－赫爾維茨判據［12］，很容易得到表達式:

可得，當w＜1+β，β＜?1+?2＜2w+2+β時，線性定常離散系統(6)漸進穩定.

當系統(6)穩定，

當參數滿足條件(8)，|A|＜1，所以，

3 數值仿真

首先運用30個城市的標準TSP測試數據Oliver30對算法性能進行評估，運用蟻群算法、遺傳算法、模擬退火方法、禁忌搜索方法［13］、混沌粒子群算法以及光學混沌粒子群算法，在最大迭代次數Mt=300的限定下分別對測試城市進行計算，根據定理1設定 ?1=0.5，?2=0.5，w采用線性遞減原則，t為當前的迭代次數，w=0.6 － (t/Mt)*0.5，混沌系數u=4.每種算法運行20次，對比如圖6所示(X、Y的散點坐標圖，取計算最好效果).各種算法計算的最好結果、最差結果以及平均迭代次數如表1所示.

表1 對于Oliver30的算法效果對比Table 1 Algorithm contrast effects of Oliver30

從表1結果中可以看出，這4種經典智能算法精度有所差異，在有限次迭代的情況下得到的最優解效果一般，如果設置較高的迭代次數，精度能達到滿意的要求，不過迭代次數較高，往往容易收斂到局部最優點.本文構建的混沌粒子群算法以及加入光學原理的混沌粒子群算法具有良好的精度，對比于基本的粒子群優化算法改進效果顯著，在最大迭代次數上限為300的情況下，可以找到精確的結果，其中光學混合粒子群算法性能更為突出，每次運行都能找到最優值423.740 6.在迭代次數上，蟻群算法具有較快的收斂速度，但容易陷入局部最優點，影響算法的精度.同比與其他的5種算法，光學混沌粒子群算法在收斂速度上有很大的優勢.從圖7中可以看出光學混沌粒子群算法具有很強的收斂速度.

圖6 7種算法Oliver30效果對比Fig.6 Seven contrast figures of Oliver30

圖7 Oliver30 4種算法進化曲線Fig.7 Four evolutionary curves of Oliver30

本文提出的光學混沌粒子群算法的GUI界面如圖8所示，得到的最優城市序列為(28，27，26，25，24，15，14，8，7，11，10，21，20，19，18，9，3，2，1，6，5，4，13，12，30，23，22，17，16，29)，平均迭代 113 次迭代找到最優解，混沌粒子群算法得到的最好結果為424.691 8，平均迭代次數為272，遠遠大于光學混沌算法，無論是迭代時間、算法精度都要劣于加入光學尋優思想的粒子群算法.光學尋優算法初始化的本質就是要減少迭代次數，提高算法精度，原因就在于用光學尋優思想形成的初始解與最優序列之間部分區域序列是相似的，甚至是完全相同的，這樣在運用粒子群算法迭代時就減少了迭代次數.混沌方法具備的全局遍歷性在保證算法的精度前提下，加強了算法跳出局部最優點的能力，增強了算法的適用性.

圖8 光學混合算法GUI效果圖Fig.8 GUI interface

表2是對測試問題eil51，參數設置與Oliver30實驗相同，最大迭代次數300，每種算法計算20次所得的數據.

表2 對于51個城市問題eil51的算法效果對比Table 2 Algorithm contrast effects of eil51

在TSP城市規模增大的情況下，光學混沌粒子群算法的優勢得到明顯體現，無論是算法精度還是算法的收斂速度都要好于其他算法.從以上結果可以看出本文提出的基于光學原理的混合粒子群算法相比于其他優化算法具有良好的效果，能成功解決中大型TSP路徑優化問題并保持較高的精度，加入的光學尋優思想能大大節約迭代次數，保證算法在規定的迭代次數內找到最優解.

表3是對測試問題CH130的對比數據，CH130是相對復雜的130個城市的TSP問題，最大迭代次數設定為500，每種算法運行20次，誤差計算是對比于CH130給出的精確解6 110，應用最優路徑計算得出，所有數據如表3所示.

表3 對于130個城市問題CH130的算法效果對比Table 3 Algorithm contrast effects of CH130

從結果可以看出，對于較大規模的TSP問題光學混合算法效果顯著，誤差在所有同類算法中最低，無論是從迭代時間搜索精度，還是誤差大小同比于其他算法都有很大的優勢.與標準粒子群算法以及混沌粒子群算法相比，改進效果顯著.圖9為最優搜索形成的TSP效果圖.

圖9 應用光學混合算法求解CH130最短路徑6 210.422 0Fig.9 The shortest path 6 210.422 0 of CH130 using the algorithm in this paper

4 結束語

本文提出了一種結合光學尋優算法、混沌思想的混合粒子群算法，通過光學尋優思想形成最優初值，利用加入混沌的粒子群算法成功解決了TSP問題，該算法迭代次數少、收斂速度快，對比于其他智能優化算法具有明顯的優勢，并用實驗表明加入光學尋優思想的搜索方式大大減少了算法的迭代次數，并在一定程度上提高了算法的精度，為高效率解決大規模TSP問題提供了新的思路.

［1］KENNEDY J，EBERHART R.A discrete binary version of the particle swarm algorithm［C］//Proceedings of the World Multiconference on Systemic，Cybernetics and Informatics.Piscataway，USA:IEEE Service Center，1997:4104-4109.

［2］CLERC M.Discrete particle swarm optimization［C］//New Optimization Techniques in Engineering.Berlin:Spinger-Verlag，2004:204-219.

［3］高尚，韓斌，吳小俊，等.求解旅行商問題的混合粒子群優化算法［J］.控制與決策，2004，19(11):1286-1289.GAO Shang，HAN Bin，WU Xiaojun，et al.Solving traveling salesman problem by hybrid particle swarm optimization algorithm［J］.Control and Decision，2004，19(11):1286-1289.

［4］HENDLASS T.Preserving diversity in particle swarm optimization［J］.Lecture Notes in Artificial Intelligence，2003(2718):155-199.

［5］XIE Shenli，TANG Min，DONG Jinxiang.An improved genetic algorithm for TSP problem［J］.Computer Engineering and Application，2002，38(8):58-60.

［6］SHEN Jihong，LI Yan.Light ray optimization and its parameter analysis［C］//Proceedings of the 2009 International Joint Conference on Computational Science and Optimization.Kunming，China，2007:918-922.

［7］沈繼紅，李焱.基于正六邊形網格的光線尋優算法［C］//中國運籌學會第十屆學術交流會論文集.南京，中國，2010:89-94.SHEN Jihong，LI Yan.Light ray optimization on hexagonal grid［C］//Proceedings of the 10th ORSC.Nanjing，China，2010:89-94.

［8］SHEN Jihong，LI Jialian.The principle analysis of light ray optimization［C］//2010 Second International Conference on Computational Intelligence and Natural Computing.Wuhan，China，2010:154-157.

［9］SHI Y，EBERHART R C.A modified particle swarm optimizer［C］//Proceedings of the Congress on Evolutionary Computation.Anchorage，USA，1998:69-73.

［10］ZHANG Guoping，WANG Zhengou，YUAN Guolin.A chaotic search method for a class of combinatorial optimization problems［J］.Systems Engineering Theory ＆ Practice，2001，21(5):102-105.

［11］梁艷春，吳春國.群智能優化算法理論與應用［M］.北京:科學出版社，2009:17-21.

［12］鄭大中.線性系統理論［M］.北京:清華大學出版社，2009:213-251.

［13］ANDRIES P E.計算智能導論［M］.北京:清華大學出版社，2010:111-123.