面向異構架構的傳遞閉包并行算法

肖 漢，郭寶云，李彩林，周清雷

（1.鄭州師范學院 信息科學與技術學院，鄭州 450044；2.山東理工大學 建筑工程學院，山東 淄博 255000；3.鄭州大學 信息工程學院，鄭州 450001）

0 概述

傳遞閉包運算在圖論、網絡、計算機形式語言、語法分析以及開關電路中的故障檢測和診斷領域都有著廣泛的應用價值［1-2］。根據定義，關系傳遞閉包的計算是通過多次進行集合復合運算完成，運算量很大。同時，假如二元關系在某種情況下發(fā)生了改變，其中的某些序偶增加或減少，需要按照原方法將變化的二元關系重新計算來得到新關系的傳遞閉包，運算量則進一步增大［3-5］。這樣容易造成大量數(shù)據無法實時處理，最終使整個應用系統(tǒng)處理的時間增加，因此如何快速有效地處理傳遞閉包問題成為了一個急需解決的問題［6-7］。

開展傳統(tǒng)的利用CPU 集群的高性能計算是解決大規(guī)?？茖W計算問題的常用方法，然而集群的并行計算性能對于CPU 的更新?lián)Q代的依賴性很大。由于CPU 芯片單位面積內的晶體管集成度越來越高，散熱和能耗問題凸顯，致使提升CPU 的速度放緩，發(fā)展陷入瓶頸［8-10］。為了更快地增強計算能力，計算機硬件設計的異構化的趨勢越發(fā)明顯［11-12］。由若干不同架構的CPU 處理器和協(xié)處理器共同工作，通用計算處理器與多個加速器設備互連構成的異構計算系統(tǒng)逐漸成為主流［13-14］。

開放式計算語言（Open Computing Language，OpenCL）是一個面向異構硬件平臺的、免費的、開放的行業(yè)標準。遵循OpenCL 規(guī)范的不同架構的硬件，提供需要的編譯和運行平臺，就能夠在OpenCL 平臺上開發(fā)普適的應用系統(tǒng)，為多核CPU、CPU+GPU、DSP 和多GPU 等異構計算提供良好的研發(fā)平臺［15-16］。

本文基于開放式計算語言平臺，提出一種基于CPU+GPU 的高效傳遞閉包并行算法，并采用具有可移植性的OpenCL 架構來實現(xiàn)該算法。對在不同數(shù)據集下和不同體系結構下的算法和加速比進行分析。

1 相關研究

近年來，很多學者對傳遞閉包運算進行了研究。文獻［17］用一階有界傳遞閉包模糊邏輯來刻畫模糊有窮自動機。文獻［18］研究了稠密圖條件下采用XHop 方法，對傳遞閉包進行高壓縮比存儲和有效查詢的算法。文獻［19］提出改進的傳遞閉包求解方法，并在傳遞閉包改進的求解方式基礎上，設計了傳遞閉包的增量式更新方法。文獻［20］證明只要函數(shù)在Jensen 的J層次結構的某個多項式級別中是統(tǒng)一可行的，則相對于函數(shù)參數(shù)的傳遞閉包，其在任意集上是安全遞歸的。文獻［21］提出改進的Floyd-Warshall算法，其中最耗時的部分（描述程序循環(huán)中自相關的傳遞閉包）是通過依賴距離向量計算，減少了傳遞閉包計算時間。文獻［22］基于程序依賴圖的傳遞閉包，提出一種在瓦片內生成具有任意順序循環(huán)的并行代碼的方法。文獻［23］使用循環(huán)嵌套依賴圖的傳遞閉包來執(zhí)行原始矩形瓦片的校正，生成并行無同步代碼。文獻［24］通過應用依賴圖的傳遞閉包，提出生成Nussinov RNA 折疊算法的并行代碼的加速因子。文獻［25］通過MPI 并行化實現(xiàn)了Warshall 方法，進而快速求取了關系R的傳遞閉包R+。

文獻［26］通過向量化方法將循環(huán)結構并行化，實現(xiàn)了傳遞閉包并行算法。文獻［27］在二叉樹并行計算模型上，實現(xiàn)了一種基于MPI 的傳遞閉包并行算法。文獻［28］通過實現(xiàn)傳遞閉包并行算法，提高了在圖形和高維數(shù)據中挖掘中心對象算法的收斂性。文獻［29］利用MPI 提出了基于VLSI 的傳遞閉包并行算法。文獻［30］通過合并Dijkstra 單源最短路徑方法中的貪婪技術的特征和傳遞閉包屬性來找到所有點對最短路徑，并在MapReduce 平臺上實現(xiàn)了ex-FTCD 算法。

綜上所述，目前大部分研究工作是通過優(yōu)化算法本身從而實現(xiàn)對傳遞閉包算法的快速計算，有些則利用傳統(tǒng)的向量化和CPU 集群的MPI 并行計算方式設計傳遞閉包算法。但是，性能加速效果在這些相關研究中表現(xiàn)的均不明顯。同時，算法研究和平臺設計局限于單一類型，對于多算法和多平臺的系統(tǒng)性能評估不多。本文將根據傳遞閉包算法特性和OpenCL 架構的特征，研究異構協(xié)同計算下的傳遞閉包并行算法，以及在多種計算平臺上算法的性能移植。

2 傳遞閉包算法

2.1 OpenCL 異構編程模型

OpenCL 是一種面向開放的、通用并行編程的、跨平臺的行業(yè)標準，軟件開發(fā)人員可以方便地將CPU、GPU和其他各類計算設備接入系統(tǒng)計算［31-32］。OpenCL 標準是編程語言與編程框架的集合體，人們可以基于硬件抽象層API和面向數(shù)據的異構編程環(huán)境進行OpenCL系統(tǒng)的開發(fā)和優(yōu)化應用。OpenCL 框架主要由OpenCL平臺層、OpenCL 運行時環(huán)境和OpenCL 編譯器3 個部分組成［33-35］。平臺層允許用戶收集可用的OpenCL 設備信息。開發(fā)者可以查詢特定設備的詳細資料，比如緩存大小、存儲器結構、核心數(shù)量等。OpenCL Runtime提供了管理設備存儲器、運行kernel、在設備與主機之間傳輸數(shù)據等一系列API［36-39］。OpenCL 編譯器創(chuàng)建包含OpenCL kernel 的可執(zhí)行程序，把kernel 編譯成設備能夠識別的代碼。

2.2 算法定義

圖的傳遞閉包可以采用布爾矩陣的平方法計算。首先假定A是一個m點有向圖的m×m的布爾鄰接矩陣，當且僅當有向圖中從頂點i到頂點j之間有一條邊時，矩陣元素aij為1。然后利用矩陣乘法對布爾矩陣的傳遞閉包A+求解［40-41］。設I是單位矩陣，大小為m×m的關系矩陣為B=A∪I。使矩陣B的第i行上的元素與矩陣B的第j列上的元素按順序分別相乘再相加，得到新的關系矩陣B的第i行第j列的元素（關系矩陣B在不斷更新），即B的定義如下：

得到新的關系矩陣B重復上一步進行循環(huán)，即依次計算，即執(zhí)行p＜logam次，得到布爾矩陣的傳遞閉包A+［42-44］。由此可知，算法的時間復雜性在最壞的情況下為O(m3logam)，當m非常大時，該算法運算將非常耗時［45］。

2.3 算法的并行特征分析

算法的可并行性高低與算法自身存在的數(shù)據依賴性有關。如果算法運算前后依賴性越強，則算法的可并行性就越低，反之，如果算法運算前后依賴性越弱，則算法的可并行性就越高，算法并行化后進行并行計算的性能就會越好。圖1 所示是一個有向圖的傳遞閉包算法的可并行性分析。

圖1 有向圖Fig.1 Directed graph

根據圖1 的5 個頂點的有向圖表示出布爾矩陣Aij，計算布爾矩陣Aij的閉包(Aij)÷過程如下：

在(Aij)÷的計算過程中可以發(fā)現(xiàn)，每一計算步驟中的任意一個元素的計算過程與其他元素計算互不影響，相互之間并沒有依賴性。因此，可以在計算某一個元素值時，同時對其他元素值進行運算。結合OpenCL 的計算模型，將每個元素的計算過程放入工作項中，每個工作項計算得出相應元素的結果。若矩陣中每個元素計算結束，則本次計算結束，如果需要繼續(xù)迭代，則再次重復以上過程。

3 傳遞閉包算法并行映射模型

3.1 并行算法

有向圖布爾矩陣A的傳遞閉包可以利用B=(A+I)的自乘logam次得到。設定工作空間中的工作組和工作項排成m×m的二維陣列，即其坐標為(tx，ty)。每個工作組用數(shù)組as和數(shù)組bs存儲矩陣B中相應子矩陣，Pvalue 保存的是每次子矩陣計算之后得到的值，數(shù)組C為每次計算完成之后最終數(shù)據。傳遞閉包并行算法描述如算法1 所示。

3.2 并行算法整體并行化思路

基于OpenCL 的傳遞閉包并行執(zhí)行流程如圖2所示。

圖2 傳遞閉包并行算法實現(xiàn)流程Fig.2 Implementation procedure of transitive closure parallel algorithm

傳遞閉包并行算法執(zhí)行過程如下：

1）在主機端根據對應的頂點數(shù)，初始化布爾矩陣A，并保存初始化后的布爾矩陣。

2）初始化OpenCL 平臺。

3）創(chuàng)建上下文，并在目標設備上創(chuàng)建命令對象。為了協(xié)調內核計算，在上下文和計算設備之間利用clCreateCommandQueue 命令建立一個邏輯鏈接。

4）讀入源程序文件，并創(chuàng)建和編譯程序對象。根據上下文中的設備特性，利用運行時編譯系統(tǒng)構建程序對象。

5）設置存儲器對象和數(shù)據傳輸。在全局內存中創(chuàng)建buffer 存儲器對象，然后將存儲器訪問任務加入到命令隊列，最后通過clCreateBuffer 將布爾矩陣A、B從CPU 端隱式地傳輸?shù)皆O備端的全局內存中。

6）建立內核對象。在指定的一個內核對象中將內核參數(shù)和內核函數(shù)通過clCreateKernel 封裝進來。

7）設置需要傳遞的內核對象參數(shù)。

8）創(chuàng)建kernel 函數(shù)，調度kernel 執(zhí)行。

9）循環(huán)調用kernel 函數(shù)對數(shù)據進行相應的處理，計算矩陣乘積大小。

10）將顯存端完成運算任務后的結果復制到主機端內存，并且釋放設備端顯存空間，將最終的計算結果保存到對應的文件中。

3.3 算法的并行方案設計

在設計傳遞閉包并行算法時，矩陣乘法的并行計算采用了工作項分塊的方法實現(xiàn)，計算原理如圖3 所示。工作組中的每個工作項讀取矩陣B中的一行和矩陣B中的一列，將行、列中對應元素相乘之后再相加，得到新矩陣B的對應位置元素值，即每個工作項對應計算新矩陣B中的一個元素。以上操作循環(huán)經過p＜logam次后得到有向圖的傳遞閉包A+的形式矩陣。

圖3 傳遞閉包算法中的矩陣相乘Fig.3 Matrix multiplication in the transitive closure algorithm

按照相互之間無重疊的劃分原則，整個矩陣B將被劃分成若干個計算區(qū)域。計算區(qū)域可作為一個基本處理單位，由工作組處理。文中采用二維工作空間進行設計，從數(shù)據層面上看，每個工作組在x，y方向上的維度均為BLOCK_SIZE。工作空間在x、y方向上共有個工作組，每個工作組中執(zhí)行了BLOCK_SIZE×BLOCK_SIZE 個工作項。

矩陣乘法中每一對元素間的乘-加計算由一個工作項負責。在內核函數(shù)中循環(huán)完成矩陣B第i行元素與矩陣B第j列元素的乘-加運算，并將乘-加的結果賦給Pvalue。在該kernel 函數(shù)中，矩陣B中的每一個元素從全局存儲器中讀取了m次，造成了時間上的大量延遲。

3.4 優(yōu)化設計

GPU 全局存儲器屬于片下存儲器，存儲空間較大，但具有較高的訪存延遲。而本地存儲器是GPU片上的高速存儲器，它的緩沖區(qū)駐留在物理GPU上。因此，本地存儲器的訪存延遲要遠遠低于全局存儲器，大量工作項的并行執(zhí)行能夠在一定程度上掩蓋全局存儲器操作的延遲。

將矩陣相乘后得到的新矩陣B分解成小矩陣塊，每一個工作組負責計算一個小矩陣塊。若矩陣B的大小是m×m，則新矩陣B=B×B。假設m=b×b，將新矩陣B分為b×b個小的子矩陣bij，則每一個子矩陣bij的大小為b×b。2 個相乘的矩陣B同新矩陣B一樣，劃分為b×b個小的子矩陣bij，且每一個子矩陣bij的大小同為b×b，則傳遞閉包并行算法中采用分塊矩陣乘法的定義為計算原理如圖4 所示。

圖4 傳遞閉包算法中的分塊矩陣相乘Fig.4 Multiplication of block matrix in transitive closure algorithm

在傳遞閉包并行算法的分塊矩陣乘法中，采用靜態(tài)方式定義大小為BLOCK_SIZE×BLOCK_SIZE的本地存儲器數(shù)組，用于存儲矩陣B子塊數(shù)據。

__local float as[BLOCK_SIZE][BLOCK_SIZE]

__local float bs[BLOCK_SIZE][BLOCK_SIZE]

為從全局存儲器預取計算子矩陣到本地存儲器，根據工作組的ID 和工作項的ID 確定B的計算子矩陣的位置，并將B中用于計算的2 個計算子矩陣分別預取至本地數(shù)組as和bs中。每個工作項負責計算一對元素的乘積和PPvalue+=as[ty][k]×bs[k][tx]。原來矩陣的一行或一列數(shù)據需要從全局存儲器讀取m次，現(xiàn)在只需要讀取m/BBLOCK_SIZE次，這樣在新矩陣B的計算過程中矩陣數(shù)據需要從全局存儲器讀取m×m次，優(yōu)化后只需要讀取m2/BBLOCK_SIZE次。因此，通過對GPU的存儲帶寬進行充分的利用，減少從全局存儲器中重復讀取數(shù)據。使用本地存儲器不僅可以降低訪問延遲以此提高訪問速率，同時節(jié)約了對全局存儲器的訪問帶寬。

4 實驗測試與結果分析

本節(jié)將給出所描述的傳遞閉包方法的測試結果。由于單精度浮點運算針對現(xiàn)代計算機，特別是在GPU 上進行了高度優(yōu)化，因此本文選擇單精度數(shù)據類型實現(xiàn)算法。

4.1 測試環(huán)境和實驗結果

實驗軟硬件平臺如下：

1）硬件平臺

平臺1：CPU 為AMD Ryzen5 1600X 3.6 GHz（六核心），24.0 GB 的系統(tǒng)內存。GPU 型號是NVIDIA GeForce GTX 1070，CUDA 核心1 920 顆，1 506 MHz的核心頻率，1 683 MHz 的流處理器頻率，8 GB GDDR5 的顯存，256 bit 的顯存位寬，256 Gb/s 的顯存帶寬，顯存存取速率為8 Gb/s。

平臺2：CPU 為AMD Ryzen5 1600X 3.6 GHz（六核心），24.0 GB 的系統(tǒng)內存。GPU 型號是AMD Radeon RX 570，其中，計算單元32 組，每組計算單元具有64 個處理單元，總計2 048 顆流處理單元，1 168 MHz 的核心頻率，256 bit 的顯存位寬，8 GB GDDR5 顯存。

2）軟件平臺：操作系統(tǒng)采用微軟Windows 8.1 64位；集成開發(fā)環(huán)境為微軟Visual Studio 2017；系統(tǒng)編譯環(huán)境為CUDA Toolkit 8.0，OpenCL 1.2 標準被支持。

有向圖的頂點集合大小n分別取為20、40、50、70、200、300、500、1 024，作為rand()隨機數(shù)函數(shù)的隨機數(shù)種子分別生成一組隨機數(shù)，構成布爾矩陣A。根據本文的傳遞閉包算法的描述，基于OpenMP 平臺和基于CUDA 平臺的傳遞閉包并行算法均在文中實現(xiàn)。

傳遞閉包算法運行在基于OpenMP系統(tǒng)、基于CUDA系統(tǒng)、基于AMD GPU 的OpenCL系統(tǒng)和基于NVIDIA GPU 的OpenCL系統(tǒng)的上處理時間，如表1所示。處理時間包括傳遞閉包算法的所有處理步驟。在OpenCL 中實現(xiàn)GPU 并行算法時，必須執(zhí)行額外的步驟，如內核創(chuàng)建（讀取、創(chuàng)建和構建最終內核對象）、主機內存和GPU 全局存儲器之間的數(shù)據傳輸以及數(shù)據結構初始化。

表1 傳遞閉包算法執(zhí)行時間Table 1 Execution time of transitive closure algorithm

用加速比作為加速效果的衡量標準，可以直觀地驗證各種架構下并行算法的效率，其定義如下：

CPU 串行算法執(zhí)行時間與并行算法執(zhí)行時間的比值即為加速比：

其中：Tserial是在CPU 上單個線程的順序運算時間；Tparallel是在多核CPU 或CPU+GPU 上多線程實現(xiàn)的并行運算時間。

相對加速比1基于OpenMP 的并行算法運算時間與基于NVIDIA GPU 的OpenCL 并行算法運算時間的比值：

其中：Tparallel-OpenMP是在多核CPU 上多線程的并行運算時間；Tparallel-NOpenCL是在NVIDIA GPU 上OpenCL 的并行運算時間。

相對加速比2基于NVIDIA GPU 平臺的CUDA 并行算法運算時間與基于NVIDIA GPU 平臺的OpenCL 并行算法運算時間的比值：

其中：Tparallel-CUDA是在CUDA 上的并行執(zhí)行時間；Tparallel-NOpenCL是在NVIDIA GPU 上OpenCL 并行實現(xiàn)的并行執(zhí)行時間。Tparallel-CUDA和Tparallel-NOpenCL定義如下：

其中：Tkernel為OpenCL 內核在CPU 和GPU 上總的執(zhí)行時間；Tovehead為在CPU 和GPU 上數(shù)據傳輸時間開銷的總和；Tother為數(shù)據結構初始化等操作總的運行時間。

為了更好地對應用系統(tǒng)速度進行客觀評價，采用加速比指標來反映在一定的計算架構下的并行算法相較串行算法的效率提升幅度。使用相對加速比1 指標來反映基于NVIDIA GPU 的OpenCL 并行算法相比基于多核CPU 的OpenMP 并行算法的效率提升情況，相對加速比2 指標則反映出基于NVIDIA GPU 的OpenCL 并行算法相比基于GPU 的CUDA 并行算法的效率提升情況，如表2 所示。

表2 傳遞閉包并行算法性能對比Table 2 Performance comparison of transitive closure parallel algorithm

4.2 實驗數(shù)據分析

4.2.1 系統(tǒng)性能瓶頸分析

在存儲器讀寫操作時，需要鄰接矩陣數(shù)據的m×m×m次存儲器讀取，有向圖的傳遞閉包矩陣數(shù)據的m×m次存儲器寫入操作。設一個m=200 點的有向圖，每個像素值分配存儲空間大小是4 Byte，所以，存儲器存取數(shù)據總量約為0.032 GB，除以kernel 實際執(zhí)行的時間0.000 257 s，得到的帶寬數(shù)值是約124.51 GB/s，這已經接近GeForce Tesla C2075 顯示存儲器的150.34 GB/s 帶寬。因此，可以很明顯地看出，基于OpenCL 架構的傳遞閉包并行算法的效率受限于全局存儲器帶寬。

從表2 可以看出，基于CPU+GPU 的算法加速效果明顯，但GPU 并行算法的加速比隨著有向圖頂點數(shù)的增加呈現(xiàn)緩慢下降的趨勢。主要原因是在OpenCL 并行算法操作中，CPU 負責讀取和輸出圖的鄰接矩陣數(shù)據，而這一過程并沒有加速。隨著被處理鄰接矩陣規(guī)模的增加，讀取和輸出鄰接矩陣數(shù)據所花費的時間也在增加。因此，OpenCL 架構下的傳遞閉包并行算法的性能瓶頸是顯存帶寬和主存與顯存之間數(shù)據傳輸?shù)膸挕?/p>

4.2.2 傳遞閉包并行算法性能分析

不同并行計算平臺下的傳遞閉包并行算法加速比對比曲線如圖5 所示。在多核CPU 平臺上，傳遞閉包算法的運算速度得到加速。然而，限于核心數(shù)，系統(tǒng)的加速比相對較小且變化不大，但由于GPU 具有較豐富的計算資源，在CUDA 架構和OpenCL 架構下的傳遞閉包算法就可以擁有足夠的工作項來進行大量數(shù)據的并行處理。1 920 個處理單元通過時間分割機制分配到一定數(shù)量的工作項，加速比得到較大提高且增幅明顯。通過表2 分析，在對計算密集型特征明顯的大規(guī)模數(shù)據集計算時，GPU系統(tǒng)運算時間有小量增幅，體現(xiàn)了GPU 用于計算密集型的任務運算不如CPU 敏感，顯現(xiàn)出GPU 強大的運算能力。

圖5 傳遞閉包并行算法的加速比對比Fig.5 Comparison of acceleration ratios of transitive closure parallel algorithm

由圖5 可知，隨著布爾矩陣規(guī)模的增加，GPU 加速下的加速比曲線斜率急劇變大，曲線變得十分陡峭。加速比呈現(xiàn)出快速增加的趨勢，比較明顯地體現(xiàn)出并行處理的性能提升效果。然而當布爾矩陣大小超過70×70 繼續(xù)增大時，曲線呈現(xiàn)出一種下降趨勢。雖然隨著布爾矩陣規(guī)模的增大，工作空間中包含的工作組數(shù)也隨之增多，系統(tǒng)中可同時執(zhí)行更多的子矩陣，對于提高訪問全局存儲器和本地存儲器的效率有益，也越容易隱藏存儲器延時，但是布爾矩陣規(guī)模的增大，主機端和設備端存儲器之間交互數(shù)據的時間成本變大，較大程度地抵消了GPU 并行計算的優(yōu)勢，導致GPU系統(tǒng)加速性能下降，整體系統(tǒng)性能受到制約。

4.2.3 傳遞閉包并行算法跨平臺性分析

可移植性不但要求源碼能夠在不同的平臺上成功地編譯、運行，而且還需要算法應當有相當?shù)男阅?。運算結果表明，在CUDA 架構下的傳遞閉包并行算法受到單一硬件平臺的限制，而基于OpenCL 的傳遞閉包并行算法則在多種硬件平臺上獲得了較好的可移植性和兼容性，其最大加速比為593.14 倍，如圖6 所示。

圖6 OpenCL 加速比趨勢Fig.6 OpenCL acceleration ratio trend

由于采用離線編譯內核讀寫數(shù)據文件的OpenCL加速的傳遞閉包并行算法，相比在線編譯內核讀寫數(shù)據文件的CUDA 加速的傳遞閉包并行算法減少了應用初始化時間。在同等數(shù)據集規(guī)模下，基于OpenCL 的傳遞閉包并行算法的運算耗時更少，與CUDA 計算平臺上的算法性能相比略有提升，最大獲得了1.05 倍加速比。而OpenCL 加速的傳遞閉包并行算法性能較之OpenMP 計算平臺下的算法性能則有很大的提高，加速比最大獲得了208.62 倍，如圖7 所示。

圖7 相對加速比趨勢Fig.7 Relative acceleration ratio trend

5 結束語

在許多應用系統(tǒng)中傳遞閉包是必要的基本部件，且為系統(tǒng)中較為耗時的部分，而矩陣乘對整個系統(tǒng)實時性能則有較大影響。本文針對傳遞閉包算法串行性能低下的不足，提出適合于OpenCL 架構的計算模式，并設計實現(xiàn)了傳遞閉包GPU 并行算法。實驗結果表明，基于OpenCL 架構的傳遞閉包并行算法的性能相比CPU 串行算法、基于CPU 的OpenMP 并行算法和基于GPU 的CUDA 并行算法，分別取得了593.14 倍、208.62 倍和1.05 倍的加速比。在算法的GPU 實現(xiàn)過程中配置適當?shù)膬群藚?shù)和合理的分塊參數(shù)，能有效提高處理效率，且實現(xiàn)同等計算量的GPU 相比CPU，性價比更高。因此，采用本文GPU異構計算模式對大規(guī)模數(shù)據運算且系統(tǒng)實時性要求較高的應用，將是一條新的思路。