Voronoi圖k階鄰近并行矩陣迭代算法

余婧，曹菡，靳朋飛

陜西師范大學(xué)計(jì)算機(jī)科學(xué)學(xué)院，西安 710062

Voronoi圖k階鄰近并行矩陣迭代算法

余婧，曹菡，靳朋飛

陜西師范大學(xué)計(jì)算機(jī)科學(xué)學(xué)院，西安 710062

針對(duì)Voronoi圖k階鄰近矢量法構(gòu)建復(fù)雜發(fā)生元困難，柵格法耗時(shí)長、精度受限等問題，提出了一種基于矩陣迭代的并行計(jì)算方法。以刀片機(jī)作為并行計(jì)算的硬件平臺(tái)，采用Arcgis軟件將MapInfo格式矢量數(shù)據(jù)轉(zhuǎn)換為柵格數(shù)據(jù)，實(shí)現(xiàn)了MPI并行環(huán)境中Voronoi圖k階鄰近的柵格計(jì)算新方法。實(shí)驗(yàn)結(jié)果表明，改進(jìn)后的Voronoi圖k階鄰近柵格并行算法明顯地提高了計(jì)算效率，且在柵格Voronoi圖精度較高時(shí)，運(yùn)行時(shí)間的拐點(diǎn)后移，加速比提高。

k階鄰近；Voronoi圖；矩陣迭代；并行計(jì)算；消息傳遞接口（MPI）

1 引言

空間鄰近問題是空間信息科學(xué)中一個(gè)重要的研究內(nèi)容，Gold[1-2]指出空間鄰近是與空間相關(guān)的諸多領(lǐng)域（GIS、空間推理等）內(nèi)一個(gè)經(jīng)常使用的概念。GIS領(lǐng)域中，Voronoi圖被應(yīng)用于空間數(shù)據(jù)模型、空間分析、空間查詢、地形建模與可視化、地形特征線提取等[3-6]方面。Voronoi拓?fù)浜涂臻g近鄰近等概念在空間內(nèi)插、空間覆蓋、數(shù)字化過程中的斷點(diǎn)捕捉、多邊形構(gòu)造、Buffer分析[4，7]等領(lǐng)域有著廣泛的應(yīng)用。

Gold將具有公共Voronoi邊的兩個(gè)空間目標(biāo)定義為具有相鄰關(guān)系[1]。Kainz等發(fā)現(xiàn)空間目標(biāo)之間的包含關(guān)系可以形成一種良好數(shù)學(xué)性質(zhì)的序關(guān)系，為Voronoi圖k階鄰近概念的提出奠定了理論基礎(chǔ)[8]。Zhao R等用Voronoi圖定義了空間目標(biāo)之間的k階鄰近關(guān)系并提出了三種生成算法：波浪法、對(duì)向法和穿越法[9]。波浪法適用于低階鄰近的計(jì)算；對(duì)向法在波浪法基礎(chǔ)上有所改進(jìn)，適合高階鄰近的計(jì)算。但此兩種算法，數(shù)據(jù)結(jié)構(gòu)復(fù)雜，計(jì)算難度大。穿越法計(jì)算效率高、數(shù)據(jù)結(jié)構(gòu)簡單，但所求的階數(shù)是近似解，只有特定空間模型中可得最優(yōu)解。李佳田[10]提出了基于Voronoi圖幾何鄰近關(guān)系的局部計(jì)算方法，但沒有考慮計(jì)算效率的問題。閆超德[11]歸納了Voronoi圖的首最鄰近遞歸收斂特性，但采用鏈?zhǔn)洁徑颖泶鎯?chǔ)空間目標(biāo)的鄰近關(guān)系，數(shù)據(jù)結(jié)構(gòu)較為復(fù)雜。本文首先介紹了Zhao R提出的Voronoi圖k階鄰近算法——波浪法和對(duì)向法，在此基礎(chǔ)上提出一種矩陣迭代的算法，采用鄰接矩陣存儲(chǔ)鄰近關(guān)系，適合大規(guī)模數(shù)據(jù)集的并行化處理，進(jìn)而提高計(jì)算效率。

2 Voronoi圖k階鄰近算法

Voronoi圖的生成方法分為矢量法和柵格法兩大類。矢量法的優(yōu)勢(shì)是生成的Voronoi圖精度高；問題和不足是：矢量法主要以點(diǎn)元生成有效的Voronoi圖，對(duì)線元、面元的處理破壞了空間目標(biāo)的完整性，計(jì)算量急劇增加；數(shù)據(jù)結(jié)構(gòu)復(fù)雜；存儲(chǔ)量要增加9～10倍[5，7]。柵格法較容易處理各種復(fù)雜空間目標(biāo)，但生成的Voronoi圖精度取決于柵格大小。鑒于矢量法構(gòu)造Voronoi圖的困難性，以下Voronoi圖k階鄰近生成算法均為柵格模式的討論。

設(shè)任意兩個(gè)空間目標(biāo)Pi，Pj之間的Voronoi區(qū)域的最少個(gè)數(shù)k為其Voronoi距離，記為vd(Pi，Pj)[12]。一般地vd(Pi，Pj)≥0。當(dāng)Pi=Pj時(shí)，vd(Pi，Pj)=0，即為兩目標(biāo)等價(jià)。當(dāng)vd(Pi，Pj)≥1時(shí)，兩目標(biāo)相鄰或者相離，值越大說明鄰近程度越弱?？筛鶕?jù)Voronoi距離對(duì)k階鄰近定義如下[9]：

定義1設(shè)Pi，Pj是空間目標(biāo)集合P中的任意兩個(gè)目標(biāo)，如果其Voronoi區(qū)域V(Pi)，V(Pj)存在，且vd(Pi，Pj)為k，則稱Pi與Pj之間存在k階鄰近關(guān)系。

特別地，當(dāng)Voronoi距離為1時(shí)，即兩目標(biāo)為1階鄰近關(guān)系時(shí)，結(jié)合點(diǎn)集拓?fù)浞椒ǎ挚煞譃橄嘟余徑椭苯余徑?/p>

定義2相接鄰近（Touch Adjacency），即兩個(gè)空間目標(biāo)的邊界具有公共部分，且Voronoi距離為1。

定義3直接鄰近（Immediate Adjacency），兩個(gè)空間目標(biāo)間不具有公共部分，但Voronoi距離為1。

Zhao R提出Voronoi圖k階鄰近生成算法——波浪法、對(duì)向法。根據(jù)公式（1）～（3），可知給定的目標(biāo)S的(n+1)階、n階與(n-1)階鄰近目標(biāo)間有下列關(guān)系式成立：

公式（4）準(zhǔn)確計(jì)算出Voronoi圖空間目標(biāo)間的k階鄰近關(guān)系，可作為基于Voronoi圖空間鋪蓋的k階鄰近的計(jì)算模型。但此兩種算法，采用鏈?zhǔn)洁徑颖泶鎯?chǔ)鄰近關(guān)系，數(shù)學(xué)模型復(fù)雜，不適合大規(guī)模數(shù)據(jù)集的并行化處理。

3 Voronoi圖k階鄰近矩陣迭代算法

空間關(guān)系計(jì)算方法的優(yōu)劣決定了Voronoi圖k階鄰近查詢效率的高低以及完備性與否。根據(jù)公式（1）～（3），提出了一種基于鄰接矩陣的迭代計(jì)算方法，核心思想如下：

以柵格形式存儲(chǔ)Voronoi圖k階鄰近關(guān)系分兩步：每個(gè)空間目標(biāo)賦予唯一的標(biāo)識(shí)，其對(duì)應(yīng)的Voronoi區(qū)域形成等值圖，特別地，公共邊界記為-1；空間目標(biāo)間的鄰接關(guān)系以二維表存儲(chǔ)。

步驟1設(shè)空間有N個(gè)目標(biāo)，形成N個(gè)Voronoi空間剖分區(qū)域。Voronoi圖k階鄰近矩陣中，第一行表示1號(hào)目標(biāo)，第二行表示2號(hào)目標(biāo)，依次類推，第N行表示N號(hào)目標(biāo)；同理可知列數(shù)。那么，定義第m行、第n列的數(shù)值為m、n目標(biāo)間的階數(shù)。當(dāng)兩目標(biāo)等價(jià)時(shí)，vd(Pi，Pj)=0。另外，任意兩個(gè)目標(biāo)的vd(Pi，Pj)均小于N。由此初始化矩陣對(duì)角線元素為0，其余元素全部賦為空間目標(biāo)數(shù)N。顯然，該矩陣為對(duì)稱矩陣?？紤]迭代時(shí)相互利用數(shù)據(jù)，且空間目標(biāo)數(shù)目遠(yuǎn)小于柵格總數(shù)，因而不實(shí)行壓縮存儲(chǔ)。

步驟2每個(gè)空間目標(biāo)對(duì)應(yīng)唯一的標(biāo)識(shí)，經(jīng)等速主動(dòng)擴(kuò)張后，形成Voronoi圖。擴(kuò)張過程中，兩個(gè)或兩個(gè)以上目標(biāo)共同占有同一柵格時(shí)，柵格標(biāo)記為-1，記錄相接鄰近(vd(Pi，Pj)=1)至k鄰近矩陣。擴(kuò)張完成后，全局掃描各個(gè)空間目標(biāo)的區(qū)域邊界，記錄直接鄰近(vd(Pi，Pj)=1)至k階鄰近矩陣。

步驟3高階迭代過程中，找目標(biāo)S的k階鄰近目標(biāo)，先到矩陣讀取S目標(biāo)的所有(k-1)階鄰近。由每個(gè)(k-1)階目標(biāo)出發(fā)，找出其一階鄰近對(duì)應(yīng)的目標(biāo)。S目標(biāo)與這些目標(biāo)在矩陣中的行列位置，比較矩陣中的數(shù)值是否比k大。如果是，找到一個(gè)k階鄰近目標(biāo)；如果不是不進(jìn)行k階更新操作。經(jīng)過多次迭代生成各個(gè)目標(biāo)的高階鄰近關(guān)系。當(dāng)k階鄰近矩陣中所有元素都小于目標(biāo)數(shù)目N時(shí)，即完成了k階鄰近矩陣的生成工作。

新算法一次性生成空間所有目標(biāo)的k鄰近關(guān)系，不必重復(fù)計(jì)算、收斂性好。生成k階鄰近矩陣后，只需進(jìn)行簡單的讀矩陣操作，即可查詢?nèi)我饽繕?biāo)的k階鄰近或者0到k階鄰近，也可查詢?nèi)魏文繕?biāo)對(duì)之間的k鄰近關(guān)系。此外，新算法適合大規(guī)模數(shù)據(jù)集的并行化處理，以下敘述的基于MPI實(shí)現(xiàn)方法，利用矩陣存儲(chǔ)k鄰近關(guān)系，有利于并行廣播、收集及規(guī)約等操作的實(shí)現(xiàn)。

4 Voronoi圖k階鄰近矩陣迭代的并行實(shí)現(xiàn)

柵格法能較好地處理復(fù)雜空間目標(biāo)，但生成的Voronoi圖精度取決于柵格大小、耗時(shí)長。通過調(diào)節(jié)柵格大小，細(xì)分格網(wǎng)提高Voronoi圖的精度。隨著柵格逐步細(xì)化，數(shù)據(jù)量進(jìn)一步增大，效率更低。針對(duì)以上問題，Voronoi圖k階鄰近算法引入并行處理技術(shù)，較好地保證了Voronoi圖的高精度和計(jì)算高效率。MPI（Message Passing Interface）是一種消息傳遞的編程模型，是消息傳遞并行程序設(shè)計(jì)的標(biāo)準(zhǔn)之一。基于MPI并行平臺(tái)，多進(jìn)程實(shí)現(xiàn)Voronoi圖k階鄰近并行矩陣迭代算法，算法流程如圖1所示。

圖1 Voronoi圖k階鄰近并行算法流程圖

首先，對(duì)柵格數(shù)據(jù)進(jìn)行數(shù)據(jù)分割，分配至各個(gè)子進(jìn)程。擴(kuò)張過程中，各子塊柵格邊界數(shù)據(jù)的計(jì)算需要用到其相鄰子塊柵格邊界數(shù)據(jù)，通信獲得其相鄰子塊柵格邊界數(shù)據(jù)。

其次，各進(jìn)程對(duì)分得的子塊柵格數(shù)據(jù)并行迭代擴(kuò)張。擴(kuò)張過程中，記錄下相接鄰近至k階鄰近矩陣。每個(gè)進(jìn)程維護(hù)各自的最鄰近矩陣。由于各進(jìn)程處理的目標(biāo)數(shù)和目標(biāo)大小均不等，各進(jìn)程處理時(shí)間不一定相同，只要有一個(gè)進(jìn)程未完成本次擴(kuò)張，其他進(jìn)程必須等待，以保證擴(kuò)張的同步性。各進(jìn)程在下一次擴(kuò)張迭代之前，子塊柵格數(shù)據(jù)需交換傳遞相鄰的邊界數(shù)據(jù)，保證生成Voronoi圖的正確性。各進(jìn)程循環(huán)進(jìn)行以上工作，直到各進(jìn)程擴(kuò)張結(jié)束。

最后，各進(jìn)程的子塊柵格數(shù)據(jù)合并得到Voronoi圖。全局掃描各Voronoi區(qū)域的邊界，記錄下直接鄰近至k階鄰近矩陣。對(duì)各進(jìn)程的最鄰近矩陣進(jìn)行規(guī)約處理，進(jìn)而得到全局最鄰近矩陣。依據(jù)上述矩陣迭代算法，得到高階鄰近矩陣，完成Voronoi圖k階鄰近生成算法。

5 實(shí)驗(yàn)及性能分析

實(shí)驗(yàn)平臺(tái)配置：惠普刀片機(jī)（型號(hào)：HP BL460c G6 E5540 6G 1P Svr），CPU為Intel Xeon E5540處理器（2.53 GHz，80 W，8 MB三級(jí)緩存），內(nèi)存為6 GB（3×2 GB）DDR3-1333Register。操作系統(tǒng)為Windows XP，集成開發(fā)環(huán)境為Visual Studio 2008，MPI為mpich2-1.2-win-ia32，協(xié)議為TCP/IP，開發(fā)語言為C++，采用Inter C++編譯器。實(shí)驗(yàn)數(shù)據(jù)采用MapInfo格式電子地圖數(shù)據(jù)，使用Arcgis軟件（Arcgis Desktop 9.2）將矢量數(shù)據(jù)轉(zhuǎn)換為柵格數(shù)據(jù)。

地圖數(shù)據(jù)樣本共有6個(gè)，分別是6、18、30、60、90、120個(gè)空間目標(biāo)的MapInfo格式電子地圖。采用第三方軟件Arcgis均將其轉(zhuǎn)換為柵格大小為0.000 5 m、0.000 3 m、0.000 1 m的3組數(shù)據(jù)，后依次簡稱精度1、精度2、精度3。

第1組實(shí)驗(yàn)如圖2所示，取精度1的6個(gè)地圖數(shù)據(jù)，對(duì)Zhao R提出的波浪法、對(duì)向法與并行矩陣迭代算法下1、2、4個(gè)進(jìn)程的運(yùn)行時(shí)間對(duì)比。

圖2 波浪法、對(duì)向法與本文算法運(yùn)行時(shí)間比較

由圖2可得出以下結(jié)論：

（1）隨著空間目標(biāo)數(shù)的增加，空白柵格數(shù)目減少，任何一種算法運(yùn)行時(shí)間總體呈下降趨勢(shì)。

（2）空間目標(biāo)數(shù)小于30時(shí)，波浪法較對(duì)向法運(yùn)行時(shí)間短；隨著空間目標(biāo)數(shù)的增加，空間目標(biāo)數(shù)大于30時(shí)，波浪法較對(duì)向法運(yùn)行時(shí)間長。這是由于波浪法適合較小規(guī)?？臻g目標(biāo)數(shù)的計(jì)算，對(duì)向法適合較大規(guī)模空間目標(biāo)數(shù)的計(jì)算。

（3）波浪法、對(duì)向法均比本算法下1個(gè)進(jìn)程運(yùn)行時(shí)間長；另外，隨著進(jìn)程數(shù)的增加，本算法取得了更短的運(yùn)行時(shí)間。

第2組實(shí)驗(yàn)如圖3～圖8所示，相同精度、不同空間目標(biāo)數(shù)、不同的進(jìn)程數(shù)之間的運(yùn)行時(shí)間和加速比比較。

圖3 精度1運(yùn)行時(shí)間比較

圖4 精度2運(yùn)行時(shí)間比較

圖5 精度3運(yùn)行時(shí)間比較

圖6 精度1加速比比較

圖7 精度2加速比比較

圖8 精度3加速比比較

由圖3～圖8可得出以下結(jié)論：

（1）空間目標(biāo)數(shù)增加，MPI并行計(jì)算時(shí)間減少。這是由于隨著空間目標(biāo)數(shù)的增多，待擴(kuò)張的空白區(qū)域相應(yīng)減少，運(yùn)行時(shí)間有所降低。

（2）隨著柵格的細(xì)化，Voronoi圖的精度顯著提高，計(jì)算規(guī)模成倍增加，多進(jìn)程的優(yōu)勢(shì)凸顯，運(yùn)行時(shí)間的拐點(diǎn)逐步后移，加速比逐步提高。精度1運(yùn)行時(shí)間的拐點(diǎn)在4個(gè)進(jìn)程處，加速比不超過2倍；精度2運(yùn)行時(shí)間的拐點(diǎn)在8個(gè)進(jìn)程處，加速比不超過3倍；精度3運(yùn)行時(shí)間的拐點(diǎn)在16個(gè)進(jìn)程處，加速比在4～8倍之間。

（3）受刀片機(jī)性能、并行同步等待及邊界通信等限制，在現(xiàn)有實(shí)驗(yàn)數(shù)據(jù)規(guī)模的基礎(chǔ)上，仍是16個(gè)進(jìn)程效率最高。

本文算法采用了鄰接矩陣這種較為簡單的數(shù)據(jù)結(jié)構(gòu)對(duì)Zhao R提出的波浪法、對(duì)向法進(jìn)行優(yōu)化，降低了計(jì)算難度，并且使用MPI對(duì)算法并行化。從第1組實(shí)驗(yàn)結(jié)果可知，本文算法較Zhao R提出的波浪法、對(duì)向法在計(jì)算效率上有所提高。第2組實(shí)驗(yàn)結(jié)果可知，隨著Voronoi圖柵格的細(xì)化，提高了Voronoi圖的精度，同時(shí)也增加了耗時(shí)，提高了多核處理器的利用率。由于較大的問題規(guī)模可提供較高的并發(fā)度；額外開銷（并行同步等待、邊界通信）的增加慢于有效計(jì)算（Voronoi圖k階鄰近的生成）的增加；Voronoi圖k階鄰近生成算法的串行分量比例隨著問題規(guī)模的增大而減小。

6 結(jié)束語

本文提出了一種基于空間目標(biāo)鄰接矩陣的迭代計(jì)算方法，采用并行柵格方法構(gòu)建Voronoi圖k階鄰近查詢。以鄰接矩陣為數(shù)據(jù)結(jié)構(gòu)，改進(jìn)了Zhao R提出的波浪法、對(duì)向法，新算法適合大規(guī)模數(shù)據(jù)集的計(jì)算。MPI并行計(jì)算模式顯著提高了Voronoi圖k階鄰近柵格生成算法的計(jì)算效率。隨著Voronoi圖的精度顯著提高，多進(jìn)程的優(yōu)勢(shì)凸顯，運(yùn)行時(shí)間的拐點(diǎn)逐步后移，加速比有所提高。下一步工作是利用并行計(jì)算的高效性及矩陣迭代算法的優(yōu)越性，做Voronoi圖k階鄰近的相關(guān)應(yīng)用。

[1]Gold C M.Dynamic spatial data structures—the Voronoi approach[C]//Canadian Conference on GIS，Ottawa，1992：245-255.

[2]Gold C M.Review：spatial tesselations-concepts and applicationsof Voronoi diagrams[J].International Journal ofGeographical Information Science，1994，8（2）：237-238.

[3]Mohammad K，Cyrus S.Voronoi-based k-nearest neighbor search for spatial network databases[C]//The 30th VLDB Conference，Toronto，Canada，2004：840-851.

[4]劉金義，劉爽.Voronoi圖應(yīng)用綜述[J].工程圖學(xué)學(xué)報(bào)，2004（2）：125-132.

[5]Aurenhammer F.Voronoi diagrams—a survey of a fundamental geometric data structure[J].ACM Computing Surveys，1991，23：345-405.

[6]劉萬增.平面離散點(diǎn)集拓?fù)溧徑€(wěn)定區(qū)域計(jì)算模型[J].測(cè)繪學(xué)報(bào)，2012，41（1）：127-132.

[7]Kekeng X，David T.Voronoi-based range and continuous range query processing in mobile databases[J].Journal of Computer and System Sciences，2011，77（4）：637-651.

[8]Kainz W.Spatial relationships topology versus order[C]// Proc of the 4th Int Symposium on Spatial Data Handling，1990：423-432.

[9]Zhao R，Chen J，Li Z.Voronoi-based k-order neighbour relations for spatial analysis[J].ISPRS Journal of Photogrammetry and Remote Sensing，2004，59（1）：60-72.

[10]李佳田.幾何鄰近空間關(guān)系形式化描述與計(jì)算的Voronoi方法[C]//Proceedings of the 4th China Association for Geographic Information System Congress，2007：407-411.

[11]閆超德，趙仁亮，陳軍，等.Voronoi圖的首最鄰近遞歸收斂特性及其應(yīng)用[J].武漢大學(xué)學(xué)報(bào)：信息科學(xué)版，2008，33（11）：1194-1197.

[12]陳軍.Voronoi動(dòng)態(tài)空間數(shù)據(jù)模型[M].北京：測(cè)繪出版社，2002：24-53.

YU Jing,CAO Han,JIN Pengfei

School of Computer Science,Shaanxi Normal University,Xi’an 710062,China

In view of the difficulty of vector method in building the Voronoi diagram k-order neighborhood with complex occurring elements,the problem of time-consuming and restricted accuracy with the raster method,this paper presents an iteration calculation based on the spatial objects adjacency matrix.The hardware is blade computer,MapInfo format vector data conversion for raster data by Arcgis software,and the new method implements the raster-based Voronoi diagram of k-order neighborhood in MPI parallel computing.Experiments show that MPI model significantly improves the calculation efficiency of the raster-based Voronoi diagram of k-order neighborhood.Experiments move knee point of running time back and get higher speed-up ratio,when the accuracy of raster-based Voronoi diagram is higher.

k-order neighbors;Voronoi diagram;iteration matrix;parallel computing;Message Passing Interface（MPI）

A

TP311

10.3778/j.issn.1002-8331.1205-0150

YU Jing,CAO Han,JIN Pengfei.Matrix iteration based parallel algorithm of k-order Voronoi diagram.Computer Engineering and Applications,2014,50（6）：102-105.

國家自然科學(xué)基金面上項(xiàng)目（No.41271387，No.40971213，No.41171310）；西安市科技計(jì)劃項(xiàng)目（社會(huì)發(fā)展引導(dǎo)計(jì)劃-軟科學(xué)研究項(xiàng)目）（No.SF1228-3）。

余婧（1987—），女，碩士研究生，主要研究方向?yàn)楦咝阅蹽IS計(jì)算；曹菡（1963—），通訊作者，女，教授，博導(dǎo)，主要研究方向?yàn)楦咝阅蹽IS計(jì)算、空間數(shù)據(jù)建模。E-mail：caohan@snnu.edu.cn

2012-05-18

2012-07-31

1002-8331（2014）06-0102-04

CNKI網(wǎng)絡(luò)優(yōu)先出版：2012-09-07，http://www.cnki.net/kcms/detail/11.2127.TP.20120907.0859.006.html