基于GPU的并行Voronoi圖柵格生成算法

Voronoi圖是一種空間分割算法。其是對空間中的n個(gè)離散點(diǎn)而言的，它將平面分割為n個(gè)區域，每個(gè)區域包括一個(gè)點(diǎn)，此區域是到該點(diǎn)距離最近的點(diǎn)的集合。由于Voronoi圖具有最鄰近性，鄰接性等眾多性質(zhì)和完善的理論體系，其被廣泛的應用在地理學(xué)、氣象學(xué)、結晶學(xué)、航天、機器人等領(lǐng)域。

Voronoi圖的生成主要有矢量方法和柵格方法。矢量法中，典型的方法有增量法、分治法和間接法。分治法是一種遞歸方法，算法思路簡(jiǎn)單，但是很難在應用過(guò)程中實(shí)現動(dòng)態(tài)更新。間接法則是根據其對偶圖Delaunay三角網(wǎng)來(lái)構造Voronoi圖，因此其性能的高低由所采用的Delaunay三角網(wǎng)的構造算法所決定。增量法通過(guò)不斷向已生成的Voronoi圖中增加點(diǎn)來(lái)動(dòng)態(tài)構建Voronoi圖。相對于前兩種方法，增量法構造簡(jiǎn)單并且容易實(shí)現動(dòng)態(tài)化，所以被廣泛應用。矢量方法的優(yōu)勢是生成Voronoi圖精度高，但是存在存儲復雜，生長(cháng)元只能是點(diǎn)和線(xiàn)，以及難以向三維及高維空間擴展等問(wèn)題。因此本文重點(diǎn)研究了Voronoi圖的柵格生成方法，首先比較了常見(jiàn)的柵格方法生成Voronoi圖的優(yōu)缺點(diǎn)，然后結合CUDA的出現，提出一種基于GPU的 Voronoi圖并行柵格生成算法。

1 柵格法簡(jiǎn)介

柵格方法生成Voronoi 圖主要是將二值圖像轉化為柵格圖像，然后確定各個(gè)空白柵格歸屬。主要方法有兩類(lèi)，一類(lèi)以空白柵格為中心，計算每個(gè)空白柵格到生長(cháng)目標的距離，以確定其歸屬，常見(jiàn)的方法有代數距離變換法，逐個(gè)空白柵格確定法等;另一類(lèi)以生長(cháng)目標為中心，不斷擴張生長(cháng)目標的距離半徑，填充其中的空白柵格，直到將整個(gè)圖像填充完成，主要有圓擴張法，數學(xué)形態(tài)學(xué)距離變換法等。代數距離變換法對距離圖像進(jìn)行上行掃描(從上到下，從左到右)和下行掃描(從下向上，從右到左)兩次掃描，計算出每個(gè)空白柵格最鄰近的生長(cháng)目標，以此生長(cháng)目標作為其歸屬。此方法中柵格距離的定義直接影響了空白柵格的歸屬和Voronoi圖的生成精度，通常使用的柵格距離定義有街區距離、八角形距離、棋盤(pán)距離等。距離變換的柵格生成方法精度低、耗時(shí)長(cháng)，所需要花費的時(shí)間和柵格的數量成正比，當柵格為n×n大小時(shí)，其時(shí)間復雜度為O(n×n)。圓檢測法以生長(cháng)目標為圓心，以一定的步長(cháng)為初始半徑，所有生長(cháng)目標同時(shí)對其構成的圓內的空白柵格進(jìn)行覆蓋。通過(guò)不斷擴大生長(cháng)目標的半徑，將會(huì )有越來(lái)越多的空白柵格被各個(gè)圓所覆蓋，直到最終覆蓋完整個(gè)圖像。數學(xué)形態(tài)學(xué)距離變換法與圓檢測法類(lèi)似，其思想來(lái)源于數學(xué)形態(tài)學(xué)中膨脹操作，膨脹操作起到了擴大圖像的效果，通過(guò)不斷的對生長(cháng)目標進(jìn)行膨脹操作，最終擴張到所有的空白柵格。這兩種方法有個(gè)共同的缺點(diǎn)，在每次擴張后，都需要判斷整個(gè)柵格圖像是否已完成擴張，而這需要遍歷柵格圖像，十分耗時(shí)。

2 GPU下的柵格生成方法

2.1 CUDA編程模型與GPU

CUDA是一個(gè)并行編程模型和一個(gè)軟件編程環(huán)境，其采用了C語(yǔ)言作為編程語(yǔ)言，提供了大量的高性能計算指令開(kāi)發(fā)能力，使開(kāi)發(fā)者能夠在GPU的強大計算能力上建立起一種更加高效的密集數據計算解決方案。

CUDA將CPU作為主機端，GPU作為設備端，一個(gè)主機端可以有多個(gè)設備端。其采用CPU和GPU協(xié)同工作的方式，CPU主要負責程序中的串行計算的部分，GPU主要負責程序中的并行計算的部分。GPU上運行的代碼被稱(chēng)為內核函數，其能夠被GPU上內置的多個(gè)線(xiàn)程并行執行。一個(gè)完整的任務(wù)處理程序由 CPU端串行處理代碼和GPU端并行內核函數共同構成。當CPU中執行到GPU代碼時(shí)，其首先將相關(guān)數據復制到GPU中，然后調用GPU的內核函數，GPU中多個(gè)線(xiàn)程并行執行此內核函數，當完成計算后，GPU端再把計算的結果返回給CPU，程序繼續執行。通過(guò)將程序中耗時(shí)的且便于并行處理的計算轉移到GPU中使用GPU并行處理，以提高整個(gè)程序的運行速度。CUDA是以線(xiàn)程網(wǎng)格(Grid)，線(xiàn)程塊(Block)，線(xiàn)程(Thread)為三層的組織架構，每一個(gè)網(wǎng)格由多個(gè)線(xiàn)程塊構成，而一個(gè)線(xiàn)程塊又由多個(gè)線(xiàn)程構成，如圖1所示。在GPU中，線(xiàn)程是并行運行的最小單元，由此可見(jiàn)，當存在大量的線(xiàn)程時(shí)，程序的并行程度將會(huì )十分高。目前的GPU上一個(gè)網(wǎng)格最多包含65535×65535個(gè)線(xiàn)程塊，而一個(gè)線(xiàn)程塊通常有512個(gè)或1024個(gè)線(xiàn)程，所以理論上可以對65535×65535×512個(gè)柵格同時(shí)進(jìn)行計算。

2.2 并行Voronoi圖柵格生成算法

傳統的柵格生成算法中，不論是采用以空白柵格為中心確定其歸屬的方法，還是以生長(cháng)目標為中心通過(guò)不斷增長(cháng)生長(cháng)目標半徑對空白柵格進(jìn)行覆蓋的方法，他們在計算每個(gè)空白柵格距離時(shí)，只能通過(guò)遍歷柵格，逐一處理。而柵格處理過(guò)程中的一個(gè)重要特點(diǎn)是，各個(gè)柵格的計算并不依賴(lài)于其他柵格的計算結果。即各個(gè)柵格的計算是相互獨立的，而由于CPU的串行性，導致了各個(gè)柵格只能順序處理，降低了處理速度。

圖1GPU組織架構

由于GPU下的多個(gè)線(xiàn)程都是硬件實(shí)現的，各個(gè)線(xiàn)程的處理都是并行的，因此將柵格距離的計算分散到GPU端各個(gè)線(xiàn)程，必然能夠提高其生成速度。為了并行處理柵格化圖像，可以采用如下的想法，將每一個(gè)柵格點(diǎn)對應于一個(gè)線(xiàn)程，此線(xiàn)程計算此柵格到所有的生長(cháng)目標的距離，取最小距離的生長(cháng)目標作為其歸屬。即采用一個(gè)線(xiàn)程用來(lái)確定一個(gè)空白柵格歸屬的方法。

確定方法后，就需要對GPU端內核函數進(jìn)行設計，由于內核函數是并行處理的執行單元，其設計方式直接決定了GPU端的程序運行效率。因此如何設計良好的內核函數是提高并行速度的關(guān)鍵。本文采用如下方式進(jìn)行內核函數的設計，假設共分配了K個(gè)并行處理線(xiàn)程，柵格規模為M×N，設A[i]為第i個(gè)線(xiàn)程處理的柵格編號。當K

由于顯卡上的內存是動(dòng)態(tài)隨機存儲(DRAM)，因此最有效率的存取方式，是以連續的方式存取。當采用第一種方式時(shí)，看似是一種連續的存取方式，實(shí)際上恰好是非連續的，當第i個(gè)線(xiàn)程處理第i個(gè)柵格時(shí)，由于處理需要一定的時(shí)間，此時(shí)GPU自動(dòng)將下個(gè)一線(xiàn)程i+1需要的內存數據取出給其使用，此時(shí)下一個(gè)線(xiàn)程的內存數據卻是在i+C處，內存變成了間斷存取。而在使用第二種方式進(jìn)行處理時(shí)，恰好是一種連續的存取方式，由于第i個(gè)線(xiàn)程正在處理第i個(gè)柵格數據，此時(shí) GPU為第i+1個(gè)線(xiàn)程準備數據，而此時(shí)的數據正好為第i+1內存處。滿(mǎn)足了內存的連續存取特性。因此本文采用第二種方式，內核函數的設計偽代碼如下：