王智廣,王文亮,張同舉,魯 強(qiáng),劉偉峰(中國(guó)石油大學(xué)

(北京) 地球物理與信息工程學(xué)院，北京 102249; 2中國(guó)石油化工股份有限公司石油勘探開(kāi)發(fā)研究院，北京 100083; 3中國(guó)石油化工集團(tuán)公司 多波地震技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京 100083)

人工蜂群算法[1](Artificial Bee Colony，ABC)是一種模擬自然界蜜蜂群體覓食行為而提出的、相對(duì)比較新的群體智能算法.該算法模擬了真實(shí)世界蜜蜂群體采蜜行為，蜜蜂根據(jù)各自的分工進(jìn)行不同的采蜜活動(dòng).在采蜜過(guò)程中，由于獨(dú)有的角色轉(zhuǎn)換機(jī)制以及在自組織模式下實(shí)現(xiàn)信息的交流和共享，從而使蜜蜂尋找到最優(yōu)蜜源.該算法可解決一系列優(yōu)化問(wèn)題，包括函數(shù)優(yōu)化問(wèn)題、組合優(yōu)化問(wèn)題等[2,3].

在文獻(xiàn)[4]中，作者已經(jīng)對(duì)于選擇OpenCL并行編程模型來(lái)研究人工蜂群算法做了許多研究，并且文獻(xiàn)[4]對(duì)于人工蜂群算法的基本原理進(jìn)行了詳細(xì)的論述.本文在文獻(xiàn)[4]中提出的并行人工蜂群算法的基礎(chǔ)上，對(duì)于先前的并行人工蜂群算法進(jìn)行了很多改進(jìn)，采用OpenCL的本地內(nèi)存(local memory)、并行規(guī)約(reduction)等技術(shù)，提出了一種基于GPU改進(jìn)的并行人工蜂群算法.通過(guò)實(shí)驗(yàn)結(jié)果的分析，改進(jìn)的算法不僅提高了算法的執(zhí)行效率，得到了很好的加速比，而且其收斂性能得到明顯的改善，說(shuō)明本文提出的基于GPU改進(jìn)的并行人工蜂群算法是有效的.

1 人工蜂群算法性能分析

根據(jù)文獻(xiàn)[4]中對(duì)于人工蜂群算法基本原理的研究，在采蜜蜂和跟隨蜂階段，需要對(duì)鄰域進(jìn)行搜索，并且需要計(jì)算大量的適應(yīng)度值.隨著問(wèn)題維度的增加，需要耗費(fèi)更長(zhǎng)的時(shí)間來(lái)完成適應(yīng)度值的計(jì)算和搜索.

如表1所示，通過(guò)對(duì)于ABC算法主要功能函數(shù)在三個(gè)不同函數(shù)優(yōu)化問(wèn)題(Griewank函數(shù)、Rastrigin函數(shù)和Rosenbrock函數(shù))的測(cè)試，得到ABC算法主要功能函數(shù)的運(yùn)行時(shí)間(單位為s).

表1 ABC算法分析

通過(guò)表1的結(jié)果和對(duì)人工蜂群算法的研究，得到最耗時(shí)的部分SendOnlookerBees函數(shù)和SendEmployedBees函數(shù)，下面基于GPU改進(jìn)的并行人工蜂群算法主要針對(duì)采蜜蜂和跟隨蜂階段進(jìn)行并行化研究.

3 并行人工蜂群算法

通過(guò)文獻(xiàn)[4]對(duì)于人工蜂群算法基本原理的研究，結(jié)合GPU自身架構(gòu)特點(diǎn)，采用本地內(nèi)存、并行規(guī)約等技術(shù)來(lái)改進(jìn)并行人工蜂群算法的執(zhí)行效率.具體來(lái)說(shuō)，基于GPU改進(jìn)的人工蜂群算法并行化設(shè)計(jì)需要解決以下問(wèn)題：GPU中數(shù)據(jù)存儲(chǔ)、隨機(jī)數(shù)的GPU處理、采蜜蜂搜索過(guò)程的GPU設(shè)計(jì)、計(jì)算最優(yōu)適應(yīng)度值、跟隨蜂搜索過(guò)程的GPU設(shè)計(jì).

3.1 GPU中的數(shù)據(jù)存儲(chǔ)

設(shè)采蜜蜂和跟隨蜂種群大小都為N，用D表示搜索的維度大小，用N×D的二維數(shù)組存儲(chǔ)N個(gè)蜜源的位置信息(foods[N][D])，用N維度的列向量存儲(chǔ)適應(yīng)度值(fitness[N]).在文獻(xiàn)[4]中，把蜜源和適應(yīng)度值存儲(chǔ)在GPU的全局內(nèi)存(global memory)，但這樣讀取數(shù)據(jù)耗費(fèi)很多時(shí)間，計(jì)算結(jié)果很不理想.當(dāng)考慮到GPU的計(jì)算效率，我們使用本地內(nèi)存(local memory)來(lái)存儲(chǔ)部分?jǐn)?shù)據(jù)，如圖1所示，具體過(guò)程如下.

步驟一：初始時(shí)，把全局內(nèi)存中的數(shù)據(jù)(foods[N][D]和fitness[N])分配到本地內(nèi)存中來(lái)存儲(chǔ)，即每個(gè)工作組(work group)分配本地內(nèi)存，用來(lái)共享存儲(chǔ)數(shù)據(jù)——蜜源和適應(yīng)度值(s_foods和s_fitness).

步驟二：工作組中每個(gè)工作項(xiàng)(work item)分配私有內(nèi)存(private memory)，并且每個(gè)工作項(xiàng)都有自己的編號(hào)tid，用solution來(lái)記錄每個(gè)維度的蜜源信息(對(duì)應(yīng)于每個(gè)tid編號(hào)).之后，通過(guò)solution來(lái)計(jì)算出相應(yīng)tid所對(duì)應(yīng)新的適應(yīng)度值.

步驟三：每個(gè)工作項(xiàng)通過(guò)solution計(jì)算新的適應(yīng)度值進(jìn)行貪婪選擇，即新的適應(yīng)度值與之舊的適應(yīng)度值(s_fitness)進(jìn)行比較，如果新的適應(yīng)度值大于舊的，則在全局內(nèi)存中更新對(duì)應(yīng)的蜜源的位置信息和適應(yīng)度值.

圖1 數(shù)據(jù)存儲(chǔ)

3.2 隨機(jī)數(shù)的GPU處理

在人工蜂群算法采蜜蜂和跟隨蜂的搜索過(guò)程都需要使用大量的隨機(jī)數(shù)，而GPU本身不提供隨機(jī)數(shù)生成函數(shù)，若利用CPU生成隨機(jī)數(shù)，在傳給GPU會(huì)進(jìn)行大量的通信而影響效率.本文采用梅森旋轉(zhuǎn)算法[5](Mersenne Twister,MT)在GPU中生成大量隨機(jī)數(shù).定義一張1024×1024大小的隨機(jī)數(shù)數(shù)組rand，在CPU上初始化32個(gè)種子，并且傳遞到GPU，在GPU中生成隨機(jī)數(shù)數(shù)組，之后在采蜜蜂和跟隨蜂的搜索過(guò)程中進(jìn)行使用.整個(gè)過(guò)程只調(diào)用一次梅森旋轉(zhuǎn)算法來(lái)產(chǎn)生所需的隨機(jī)數(shù).

3.3 采蜜蜂搜索過(guò)程的GPU設(shè)計(jì)

在采蜜蜂尋找蜜源的過(guò)程中，每只采蜜蜂都在各自區(qū)域及其附近鄰域的進(jìn)行搜索，采蜜蜂之間依賴程度低，所以天然地存在并行性.在文獻(xiàn)[4]中，每只采蜜蜂對(duì)應(yīng)一個(gè)工作組(work group)，這樣計(jì)算效率很低.本文則采用每只采蜜蜂對(duì)應(yīng)一個(gè)工作項(xiàng)(work item).如圖2所示的OpenCL并行模型，一共有4個(gè)工作組，每個(gè)工作組有4個(gè)工作項(xiàng)，具體的過(guò)程如下.

步驟一：初始時(shí)，每個(gè)蜜源只有一只采蜜蜂，把蜜源的位置信息、適應(yīng)度值等從全局內(nèi)存拷貝到本地內(nèi)存中.

步驟二：每只采蜜蜂對(duì)應(yīng)一個(gè)工作項(xiàng)，并且獨(dú)立地進(jìn)行搜索工作，即通過(guò)本地內(nèi)存中蜜源位置信息隨機(jī)搜索，并計(jì)算出新的適應(yīng)度值.

圖2 人工蜂群算法GPU并行化模型

步驟三：等到搜索結(jié)束，對(duì)于每個(gè)蜜源通過(guò)適應(yīng)度值貪婪選擇，將新的適應(yīng)度值與舊的進(jìn)行比較，保留適應(yīng)度值較高的蜜源信息，并將全局存儲(chǔ)器中對(duì)應(yīng)蜜源的位置信息、適應(yīng)度值進(jìn)行更新.

3.4 計(jì)算最優(yōu)適應(yīng)度值

在CPU中計(jì)算最優(yōu)適應(yīng)度值很簡(jiǎn)單，只需將全部適應(yīng)度值進(jìn)行比較，求出最大值.但將這樣的算法移植到GPU則需要很長(zhǎng)的時(shí)間交換數(shù)據(jù)，所以傳統(tǒng)方法在GPU上無(wú)法直接計(jì)算.本文采用GPU并行規(guī)約(reduction)算法[6]來(lái)計(jì)算最優(yōu)適應(yīng)度值.并行規(guī)約算法可以將計(jì)算N個(gè)數(shù)的時(shí)間復(fù)雜度相對(duì)于原來(lái)串行算法O(N)降到O(logN).當(dāng)N的規(guī)模很大時(shí)，算法效率提升顯著.

3.5 跟隨蜂搜索過(guò)程的GPU設(shè)計(jì)

在跟隨蜂階段，由于其選擇過(guò)程是全局概率選擇，是一種全局輪盤賭選擇方式，相互之間有依賴，需要在整個(gè)種群搜集信息，不適合本文采用的GPU硬件進(jìn)行細(xì)粒度并行算法.本文采用局部選擇方式進(jìn)行選擇操作，相對(duì)于文獻(xiàn)[4]的左右鄰域選擇幾率相等的設(shè)計(jì)，本文采用優(yōu)先右鄰域的策略進(jìn)行局部選擇.將待處理的每只跟隨蜂i對(duì)應(yīng)搜索蜜源i，并選擇其左右鄰域.如圖3所示，具體過(guò)程如下.

步驟一：初始化數(shù)據(jù)，將全局內(nèi)存中蜜源位置信息和適應(yīng)度值等信息拷貝到本地內(nèi)存中.跟隨蜂0對(duì)應(yīng)蜜源位置0，其左右鄰域?yàn)镹-1、1，依此類推.

圖3 跟隨蜂局部選擇

步驟二：根據(jù)采蜜蜂搜索結(jié)果，計(jì)算出對(duì)應(yīng)的選擇概率prob值.之后，對(duì)于這3個(gè)蜜源分別產(chǎn)生3個(gè)隨機(jī)數(shù)(隨機(jī)數(shù)范圍[0,1])，如果“蜜源0”產(chǎn)生的隨機(jī)數(shù)大于選擇prob值，則把該蜜源標(biāo)記成1，不在進(jìn)行其他蜜源比較；否則選擇右鄰域“蜜源1”，同樣比較隨機(jī)數(shù)與prob值大小，如果滿足條件，則把該蜜源標(biāo)記成1.否則選擇左鄰域“蜜源N-1”進(jìn)行同樣操作，如果還不滿足條件，則終止搜索，并把trail值加1.

步驟三：根據(jù)步驟二的結(jié)果，對(duì)于已經(jīng)標(biāo)記的蜜源及其附近進(jìn)行搜索，即對(duì)于已經(jīng)標(biāo)記的蜜源計(jì)算相應(yīng)的適應(yīng)度值.等到搜索結(jié)束，通過(guò)適應(yīng)度值貪婪選擇，保留適應(yīng)度值高的蜜源，并將全局存儲(chǔ)器中對(duì)應(yīng)蜜源的位置信息、適應(yīng)度值等進(jìn)行更新.由于跟隨蜂的種群大小可能改變，可以根據(jù)實(shí)際情況設(shè)置鄰域大小.

3.6 人工蜂群算法的GPU實(shí)現(xiàn)

在解決上述問(wèn)題的基礎(chǔ)上，根據(jù)OpenCL并行編程模型，具體實(shí)現(xiàn)過(guò)程如下.

步驟一：初始化參數(shù).在CPU中，對(duì)foods[N][D]、fitness[N]、prob[N]、trial[N]以及梅森旋轉(zhuǎn)算法種子等參數(shù)進(jìn)行初始化.

步驟二：將數(shù)據(jù)和參數(shù)傳輸?shù)紾PU全局內(nèi)存中.

步驟三：產(chǎn)生隨機(jī)數(shù).采用3.2節(jié)提到的梅森旋轉(zhuǎn)算法一次性產(chǎn)生1024×1024大小的隨機(jī)數(shù)數(shù)組rand，存儲(chǔ)在GPU中.

步驟四：采蜜蜂階段.根據(jù)3.3節(jié)采蜜蜂搜索過(guò)程并行化方法，從全局內(nèi)存將蜜源的位置信息、適應(yīng)度值拷貝到本地內(nèi)存中，讓每只采蜜蜂在對(duì)應(yīng)的蜜源及其附近進(jìn)行搜索，并計(jì)算適應(yīng)度值，按照適應(yīng)度值大小進(jìn)行貪婪選擇，對(duì)于滿足條件的蜜源進(jìn)行更新蜜源的位置信息及其適應(yīng)度值，把數(shù)據(jù)拷貝回GPU全局內(nèi)存中.

步驟五：跟隨蜂階段.根據(jù)3.5節(jié)跟隨蜂搜索過(guò)程并行化方法，從全局內(nèi)存將蜜源的位置信息、適應(yīng)度值拷貝到本地內(nèi)存中，并按照可能性概率進(jìn)行右鄰域優(yōu)先的局部選擇，然后每只采蜜蜂在對(duì)應(yīng)的蜜源及其附近進(jìn)行搜索，并計(jì)算適應(yīng)度值，按照適應(yīng)度值大小進(jìn)行貪婪選擇，對(duì)于滿足條件的蜜源進(jìn)行更新蜜源的位置信息及其適應(yīng)度值，把數(shù)據(jù)拷貝回GPU全局內(nèi)存中.

步驟六：偵查蜂階段.若搜尋次數(shù)超過(guò)一定限制limit，仍沒(méi)找到具有更高適應(yīng)度值的蜜源，則放棄該蜜源，并且在GPU全局內(nèi)存中隨機(jī)產(chǎn)生一個(gè)新的蜜源.

步驟七：記錄到目前的全局最優(yōu)解，并跳轉(zhuǎn)至步驟四，直到滿足算法結(jié)束條件，把最終的全局最優(yōu)解等信息從GPU傳回CPU，進(jìn)行輸出顯示.

4 實(shí)驗(yàn)及其結(jié)果

4.1 測(cè)試的基準(zhǔn)函數(shù)

為了評(píng)價(jià)算法性能，我們利用本文的GPU人工蜂群算法執(zhí)行D. Karaboge和B. Basturk[7]中介紹的3個(gè)經(jīng)典基準(zhǔn)函數(shù)，分別為Greiwank函數(shù)、Rastrigin函數(shù)和Rosenbrock函數(shù)，具體函數(shù)表示以及取值范圍等可以參考文獻(xiàn)[4]實(shí)驗(yàn)及其結(jié)果部分.

4.2 實(shí)驗(yàn)環(huán)境及參數(shù)設(shè)定

本實(shí)驗(yàn)CPU采用Intel Xeon E5504處理器，8GB內(nèi)存.GPU為NVIDIA Tesla C1060處理器.操作系統(tǒng)為Red Hat Enterprise Linux 5.4，OpenCL 1.0 CUDA 3.2.串行和并行人工蜂群算法都設(shè)置最大循環(huán)次數(shù)為2000次.

4.3 實(shí)驗(yàn)結(jié)果及其分析

(1) 問(wèn)題維度的增加.實(shí)驗(yàn)中，假定蜜蜂種群的數(shù)量一定時(shí)(測(cè)試的時(shí)候選擇種群規(guī)模為6000)，問(wèn)題維度從30變化到90.如表2所列數(shù)據(jù)，f1、f2、f3分別表示Griewank函數(shù)、Rastrigin函數(shù)和Rosenbrock函數(shù)，表中的數(shù)據(jù)是3個(gè)函數(shù)在串行算法(CPU計(jì)算)和并行算法(CPU-GPU異構(gòu)計(jì)算)運(yùn)行時(shí)間(單位為s).通過(guò)表中的數(shù)據(jù)可以看出，在不同維度下，隨著問(wèn)題維度的增加，串行算法運(yùn)行時(shí)間也相應(yīng)的增加，在維度為90時(shí)，Griewank函數(shù)和Rastrigin函數(shù)運(yùn)行時(shí)間可以達(dá)到100 s以上；而并行算法運(yùn)行時(shí)間只在10 s左右變化，當(dāng)維度大于60時(shí)，運(yùn)行時(shí)間還有下降的趨勢(shì)，說(shuō)明GPU硬件適合進(jìn)行大規(guī)模密集型計(jì)算.

表2 維度與運(yùn)行時(shí)間關(guān)系表

如圖4所示，3個(gè)基準(zhǔn)函數(shù)隨著問(wèn)題維度的增加，加速比也相應(yīng)的增加，加速比變化范圍從0.5到16.從圖4中很直觀地看出，Griewank函數(shù)和Rastrigin函數(shù)加速比的效果比較顯著，Rosenbrock函數(shù)相對(duì)低一些，說(shuō)明本文的并行算法是有效的.同時(shí)，從圖4上可以看出，當(dāng)問(wèn)題的維度大于90時(shí)，加速比仍有增長(zhǎng)的趨勢(shì).

圖4 維度與加速比關(guān)系圖

(2) 種群規(guī)模的擴(kuò)大.實(shí)驗(yàn)中，假定問(wèn)題的維度一定時(shí)(測(cè)試的時(shí)候選擇問(wèn)題維度為60)，種群規(guī)模從600變化到6000，串行算法(CPU計(jì)算)比并行算法(CPU-GPU異構(gòu)計(jì)算)需要更多的時(shí)間.如圖5所示，f1、f2、f3分別代表Griewank函數(shù)、Rastrigin函數(shù)和Rosenbrock函數(shù).當(dāng)問(wèn)題規(guī)模增大時(shí)，由于算法中計(jì)算適應(yīng)度值等操作花費(fèi)大量的時(shí)間，串行算法Griewank函數(shù)和Rastrigin函數(shù)的運(yùn)行時(shí)間顯著上升，在種群規(guī)模為6000時(shí)運(yùn)行時(shí)間達(dá)到70～90 s之間；而并行算法3個(gè)函數(shù)的運(yùn)行時(shí)間增長(zhǎng)很緩慢，幾乎只在10 s左右變化.這說(shuō)明并行算法在種群規(guī)模增加的情況下，運(yùn)行時(shí)間不會(huì)受到很到影響，同時(shí)也說(shuō)明并行算法有很好的執(zhí)行效率.

圖5 種群規(guī)模與運(yùn)行時(shí)間關(guān)系圖

(3)算法的收斂情況.在精度相同的情形下，并行算法的收斂速度比串行算法快.本實(shí)驗(yàn)以Griewank函數(shù)、Rastrigin函數(shù)和Rosenbrock函數(shù)為例，設(shè)定維度為60、種群規(guī)模為6000情況下，如圖6所示，三個(gè)函數(shù)的并行算法基本在200次迭代就能到達(dá)比較好的效果，而串行算法要用至少500次迭代才能到達(dá)比較好的效果.事實(shí)上，并行算法比串行算法原本運(yùn)行時(shí)間就短，加之收斂速度快，非常有助于提高算法的執(zhí)行效率，減少算法的迭代次數(shù)和運(yùn)行時(shí)間.

但3個(gè)函數(shù)并行算法過(guò)程中，有時(shí)會(huì)出現(xiàn)陷入波動(dòng)情況，尤其在圖6(b)最為顯著，并行算法在迭代次數(shù)為500次左右，出現(xiàn)了波動(dòng)，有可能是在算法偵察蜂過(guò)程中達(dá)到limit限制，放棄原來(lái)蜜源，隨機(jī)初始化新蜜源產(chǎn)生的情況，對(duì)于并行算法有待進(jìn)一步研究和實(shí)驗(yàn).

總體來(lái)說(shuō)，并行算法的收斂性能方面好于串行算法，運(yùn)行時(shí)間有大幅縮減.

圖6 函數(shù)性能比較

4.4 本文算法實(shí)驗(yàn)結(jié)果與文獻(xiàn)[4]實(shí)驗(yàn)結(jié)果的對(duì)比分析

本文實(shí)驗(yàn)環(huán)境采用GPU為NVIDIA Tesla C1060處理器，文獻(xiàn)[4]中采用GPU為ATI Mobility Radeon HD 5470處理器，由于GPU架構(gòu)等方面的不同，實(shí)驗(yàn)數(shù)據(jù)沒(méi)有完全在相同的實(shí)驗(yàn)環(huán)境下進(jìn)行比較，此處只給出在加速性能結(jié)果方面的比較，如圖7所示.圖中f1、f2、f3分別代表本文中Griewank函數(shù)、Rastrigin函數(shù)和Rosenbrock函數(shù)，f1[4]、f2[4]、f3[4]表示文獻(xiàn)[4]中的實(shí)驗(yàn)結(jié)果.從圖中可以看出，本文中Griewank函數(shù)、Rastrigin函數(shù)的實(shí)驗(yàn)效果明顯好于文獻(xiàn)[4]，本文Rosenbrock函數(shù)的實(shí)驗(yàn)效果略低于文獻(xiàn)[4]的結(jié)果，但隨著維度的增加，Rosenbrock函數(shù)會(huì)有所增長(zhǎng)，而文獻(xiàn)[4]中3個(gè)函數(shù)的變化趨于平緩.總體講，本文提出的改進(jìn)算法與文獻(xiàn)[4]相比有較大提升.

圖7 加速比較

5 結(jié)語(yǔ)

本文在作者已發(fā)表的文獻(xiàn)[4]的基礎(chǔ)上，對(duì)于并行人工蜂群算法進(jìn)行了改進(jìn)，提出了一種基于GPU改進(jìn)的并行人工蜂群算法.該算法提高了算法的執(zhí)行效率，得到了很好的加速比，并且解決了之前并行版本收斂速度變慢的問(wèn)題.

對(duì)于本文并行人工蜂群算法的進(jìn)一步研究，可以從以下幾方面考慮:(1)雖然提高了算法的收斂速度，但收斂過(guò)程中，產(chǎn)生的波動(dòng)情況有待進(jìn)一步的深入的研究和實(shí)驗(yàn)；(2)進(jìn)一步優(yōu)化算法，提高算法的執(zhí)行效率；(3)本文圖7中闡述的Rosenbrock函數(shù)的實(shí)驗(yàn)效果略低于文獻(xiàn)[4]結(jié)果的問(wèn)題將在后續(xù)的研究中加以完善.

[1] Karaboga D,Basturk B.A powerful and efficient algorithm for numerical function optimization: artificial bee colony (ABC) algorithm[J]. Journal of Global Optimization,2007,39(3): 459-471.

[2] Karaboga D,Akay B. A survey: algorithms simulating bee swarm intelligence[J]. Artificial Intelligence Review,2009(31): 61-85.

[3] Karaboga D,Akay B. A comparative study of artificial bee colony algorithm[J]. Applied Mathematics and Computation,2009,214(1):108-132.

[4] 王文亮,王智廣,劉偉峰. 基于GPU加速的細(xì)粒度并行人工蜂群算法[J]. 微電子學(xué)與計(jì)算機(jī),2013,30(3):110-114.

[5] Matsumoto M,Nishimura T.Dynamic creation of pseudo-random number generators[J]. Monte Carlo and Quasi-Monte Carlo Methods,1998(1):56-69.

[6] 吳恩華,柳有權(quán). 基于圖形處理器(GPU)的通用計(jì)算[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào),2004,16(5): 601-612.

[7] Karaboga D,Basturk B. On the performance of artificial bee colony (ABC) algorithm[J]. Applied Soft Computing,2008,8(1): 687-697.