基于神威·太湖之光的非結(jié)構(gòu)網(wǎng)格計(jì)算加速算法

許樂，安虹，陳俊仕，張鵬飛，武錚

（中國(guó)科學(xué)技術(shù)大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，合肥 230026）

0 概述

近幾十年來，計(jì)算流體力學(xué)（Computational Fluid Dynamics，CFD）呈現(xiàn)蓬勃發(fā)展的態(tài)勢(shì)，在學(xué)術(shù)界與各學(xué)科相互交叉，科研成果日新月異，在工業(yè)界與各領(lǐng)域深度結(jié)合，成為輔助工程設(shè)計(jì)的新興技術(shù)。非結(jié)構(gòu)網(wǎng)格是有限元和有限體積等數(shù)值方法中最常用的離散網(wǎng)格，在CFD 計(jì)算中具有重要意義，需要基于硬件架構(gòu)精細(xì)調(diào)優(yōu)以保證其良好的性能優(yōu)勢(shì)。非結(jié)構(gòu)網(wǎng)格計(jì)算一般表示為稀疏矩陣運(yùn)算，其更新時(shí)的隨機(jī)訪存加劇了數(shù)據(jù)存儲(chǔ)的離散性。隨著網(wǎng)格規(guī)模的增加，離散訪存的規(guī)模也在成倍增加，在帶寬受限的計(jì)算系統(tǒng)中成為主要性能瓶頸。離散訪存的計(jì)算特點(diǎn)使得并行化時(shí)也會(huì)出現(xiàn)多個(gè)計(jì)算任務(wù)同時(shí)訪問相同元素的寫后讀沖突和寫寫沖突問題。

神威·太湖之光計(jì)算機(jī)系統(tǒng)［1］是世界上第一臺(tái)峰值性能超過十億億次量級(jí)的超算系統(tǒng)，使用完全自主研制的申威26010 異構(gòu)眾核處理器［2］。但非結(jié)構(gòu)網(wǎng)格在眾核處理器上的計(jì)算普遍存在離散訪存、數(shù)據(jù)依賴等問題，并行化難度較高。此外，有依賴關(guān)系的算子和對(duì)稱矩陣的計(jì)算進(jìn)一步提高了在眾核處理器上實(shí)現(xiàn)高性能非結(jié)構(gòu)網(wǎng)格計(jì)算的難度。

為了解決非結(jié)構(gòu)網(wǎng)格計(jì)算中有依賴關(guān)系算子在眾核處理器上的優(yōu)化問題，本文對(duì)大量稀疏網(wǎng)格數(shù)據(jù)進(jìn)行分析，從網(wǎng)格本身結(jié)構(gòu)和數(shù)據(jù)之間的關(guān)系出發(fā)，提出自適應(yīng)和無依賴的任務(wù)劃分策略，使任務(wù)劃分方法與具體算子不產(chǎn)生綁定關(guān)系，從而提高對(duì)不同類型算子的普適性。根據(jù)主從核架構(gòu)的特點(diǎn)，本文提出N 階對(duì)角染色算法平衡主從核計(jì)算，并在從核計(jì)算時(shí)摒棄傳統(tǒng)的寄存器通信操作，便于擴(kuò)展到新一代神威平臺(tái)。此外，考慮到計(jì)算訪存重疊技術(shù)是申威處理器的常見優(yōu)化策略，本文利用該技術(shù)進(jìn)一步提升計(jì)算效率。

1 研究背景

1.1 申威異構(gòu)眾核處理器

神威·太湖之光計(jì)算機(jī)系統(tǒng)是我國(guó)首臺(tái)完全自主研發(fā)的世界第一超算系統(tǒng)，也是我國(guó)目前使用最廣泛的高性能計(jì)算平臺(tái)之一，為經(jīng)濟(jì)和社會(huì)發(fā)展提供了有力支撐。神威·太湖之光超算系統(tǒng)由高速計(jì)算系統(tǒng)和輔助計(jì)算系統(tǒng)及配套的互連網(wǎng)絡(luò)和存儲(chǔ)系統(tǒng)組成，配備精準(zhǔn)的資源調(diào)度系統(tǒng)和豐富的并行編程環(huán)境。系統(tǒng)由40 960 塊申威26010 處理器構(gòu)成，內(nèi)存空間為1 024 TB，訪存總帶寬為4 473 TB/s，峰值運(yùn)算速度達(dá)到125PFLOPs，比其他同量級(jí)超算系統(tǒng)節(jié)能60%以上。

申威26010 處理器是我國(guó)通過自主核心技術(shù)研制的全新異構(gòu)眾核處理器，支持64 位申威指令集。申威處理器由4 個(gè)同構(gòu)核組構(gòu)成，每個(gè)核組內(nèi)有1 個(gè)控制核心（主核）和64 個(gè)計(jì)算核心（從核），共享統(tǒng)一編址的8 GB 主存，如圖1 所示。

圖1 申威26010 處理器架構(gòu)Fig.1 The architecture of SW26010 processor

在申威26010 處理器中：主核負(fù)責(zé)任務(wù)分發(fā)和調(diào)度，工作頻率為1.45 GHz，L1 Cache 為32 KB，L2 Cache（數(shù)據(jù)和指令Cache 混合）為256 KB；從核負(fù)責(zé)稠密計(jì)算，工作頻率為1.45 GHz，指令Cache 為16 KB。從核采用64 KB 的局部存儲(chǔ)（Local Data Memory，LDM）代替硬件Cache，需要用戶手動(dòng)完成數(shù)據(jù)的換入換出，有利于充分利用片上存儲(chǔ)空間，但也給編程帶來極大挑戰(zhàn)。

申威處理器98%的計(jì)算性能來源于從核陣列，因此，挖掘從核架構(gòu)特性、充分利用從核計(jì)算資源十分重要。稀疏矩陣運(yùn)算的計(jì)算強(qiáng)度遠(yuǎn)低于稠密矩陣，為達(dá)到較高的性能，需要更高訪存帶寬的支持，而申威處理器訪存性能不佳，因此，必須充分利用從核訪存機(jī)制來盡可能降低開銷。

從核可以通過gld/gst 指令直接對(duì)主存進(jìn)行訪問，但基準(zhǔn)測(cè)試顯示gld/gst 的延遲很高，達(dá)到上百節(jié)拍數(shù)，帶寬低于1.5 GB/s，因此，在密集訪存時(shí)一般不作考慮。在通常情況下，從核使用DMA 操作連續(xù)訪問主存數(shù)據(jù)可獲得明顯的性能提升。DMA 操作帶寬如表1 所示。

表1 DMA 操作帶寬Table 1 The bandwidth of DMA operation

DMA 操作是指從核LDM 和主存之間的數(shù)據(jù)傳輸，它只能由從核線程發(fā)起，有單從核模式（PE_MODE）、行模式（ROW_MODE）、廣播模式（BCAST_MODE）、廣播行模式（BROW_MODE）和行集合模式（RANK_MODE）5 種，一般情況下常用單從核模式。此外，從核還支持跨步DMA，即按照一定跨步間隔連續(xù)訪問主存。

在每個(gè)核組中，64 個(gè)從核構(gòu)成8×8 的從核陣列，從核間的數(shù)據(jù)交換通過寄存器通信機(jī)制（Register Level Communication，RLC）進(jìn)行，該機(jī)制以生產(chǎn)者/消費(fèi)者模式運(yùn)行，各從核以1 個(gè)向量長(zhǎng)度為單位在其行或列上進(jìn)行數(shù)據(jù)廣播和接收。如圖2 所示，源從核首先將256 位對(duì)齊的數(shù)據(jù)加載到寄存器中，然后通過發(fā)送緩沖區(qū)（Send Buffer）將它們發(fā)送到從核通信網(wǎng)格；目的從核通過接收緩沖區(qū)（Receive Buffer）從通信網(wǎng)格中獲取數(shù)據(jù)，并將其加載到本地寄存器。

圖2 寄存器通信機(jī)制原理Fig.2 RLC principle

寄存器的通信延遲通常低至幾個(gè)周期，如表2所示，這使得從核間的數(shù)據(jù)可以快速交換，但每個(gè)從核只能通過行廣播向同一行中的一個(gè)或多個(gè)從核發(fā)送數(shù)據(jù)，或通過列廣播在同一列中發(fā)送數(shù)據(jù)，這給數(shù)據(jù)的自由傳輸帶來極大限制，不利于開發(fā)有復(fù)雜依賴關(guān)系的從核程序。

表2 寄存器通信延遲Table 2 Register communication latency

1.2 相關(guān)研究

近年來，CFD 數(shù)值模擬軟件作為高性能計(jì)算領(lǐng)域中的重要應(yīng)用，已經(jīng)廣泛部署于眾多超算平臺(tái)上。商業(yè)軟件因其平臺(tái)適配性強(qiáng)和穩(wěn)定性高，曾是非結(jié)構(gòu)網(wǎng)格計(jì)算軟件的主流，如Fluent［3］、UMS3D［4-5］、FUN3D［6-7］、TAU［8］、CFD++［9］、NSU3D［10］等，但其內(nèi)部源代碼不對(duì)外公開，很難精準(zhǔn)解決用戶需求。相比而言，開源CFD軟件可以滿足用戶不同的開發(fā)需求，如OpenFOAM（Open source Field Operation And Manipulation）［11］、Featflow［12］、Gerris［13］、Code_Saturne［14］等，其中，OpenFOAM是目前應(yīng)用范圍最廣、可擴(kuò)展性最強(qiáng)、解算器最全的開源軟件包。

CFD 軟件中的核心部分是高精度數(shù)值求解器，越來越多的研究人員開始使用異構(gòu)加速方法來改進(jìn)線性方程組的求解。BOLZ等［15］首次在GPU 上實(shí)現(xiàn)了高效的SpMV 算子，此后，基于ELLPACK格式，BELL等［16］設(shè)計(jì)HYB存儲(chǔ)格式實(shí)現(xiàn)了SpMV，而VáZQUEZ等［17］、MONAKOV等［18］和CHOI等［19］分別設(shè)計(jì)了ELLPACK-R、sliced-ELLPACK 和blocked ELLPACK 存儲(chǔ)格式。

基于CSR格式，KOZA等［20］提出了CSMR格式以提高數(shù)據(jù)重用，GREATHOUSE等［21］和ASHARI等［22］分別提出CSR-adaptive和ACSR格式以解決負(fù)載均衡問題。此外，MERRILL等［23］和LIU等［24］還分別提出MCSR和CSR5格式，取得了良好的性能提升。對(duì)于分塊問題，BULU?等［25］、ASHARI等［26］、LIANG等［27］和YAN等［28］分別提出了CSB、BRC、HCC和BCCOO存儲(chǔ)格式。

在國(guó)產(chǎn)申威處理器上，文獻(xiàn)［29］基于CSR 格式提出動(dòng)靜態(tài)優(yōu)化方法，其相比主核實(shí)現(xiàn)取得了6 倍的加速效果，文獻(xiàn)［30］進(jìn)一步提出雙邊多級(jí)劃分方法，相比主核實(shí)現(xiàn)取得了12 倍以上的加速效果。文獻(xiàn)［31］基于非結(jié)構(gòu)網(wǎng)格實(shí)現(xiàn)了稀疏下三角方程求解器，文獻(xiàn)［32］提出基于排序思想的通用眾核優(yōu)化算法以減少非結(jié)構(gòu)網(wǎng)格計(jì)算中的隨機(jī)訪存，隨后，文獻(xiàn)［33］提出兩階段優(yōu)化方法以克服大規(guī)模不規(guī)則訪存和帶寬利用率低的問題。

2 非結(jié)構(gòu)網(wǎng)格計(jì)算

OpenFOAM 是一款對(duì)連續(xù)介質(zhì)力學(xué)問題進(jìn)行數(shù)值計(jì)算的開源C++類庫，因其模塊化和可定制程度高，目前已成為超算平臺(tái)上主流的CFD 軟件。OpenFOAM基于C++語言開發(fā)，利用操作符重載、繼承和模板等面向?qū)ο筇匦裕С謹(jǐn)?shù)據(jù)預(yù)處理、數(shù)據(jù)后處理和自定義偏微分方程求解，框架內(nèi)提供網(wǎng)格生成、有限體積法、線性方程組求解、輸入輸出處理等功能。

如圖3 所示，OpenFOAM 中非結(jié)構(gòu)網(wǎng)格一般通過鄰接矩陣表示，而非結(jié)構(gòu)網(wǎng)格的稀疏性又使得鄰接矩陣為稀疏矩陣。稀疏運(yùn)算計(jì)算強(qiáng)度低，在眾核處理器上仍有很大的優(yōu)化空間，因此，本文基于申威處理器提出非結(jié)構(gòu)網(wǎng)格計(jì)算的通用加速框架。

圖3 非結(jié)構(gòu)網(wǎng)格及其矩陣表示Fig.3 The unstructured gird and its matrix representation

2.1 基礎(chǔ)數(shù)據(jù)結(jié)構(gòu)

OpenFOAM 中稀疏矩陣按照LDU 格式存儲(chǔ)，矩陣上除對(duì)角線元素外各非零元素用一個(gè)三元組（行號(hào)，列號(hào)，數(shù)值）表示，如圖4 所示。

圖4 LDU 存儲(chǔ)格式Fig.4 LDU storage format

整個(gè)稀疏矩陣分為對(duì)角部分（Diag）、上三角部分（Upper）和下三角部分（Lower），對(duì)角線元素即網(wǎng)格頂點(diǎn)數(shù)據(jù)，上三角和下三角元素為網(wǎng)格邊數(shù)據(jù)。上三角元素的行索引（Row）對(duì)應(yīng)下三角元素的列索引（Col），而上三角元素的列索引對(duì)應(yīng)下三角元素的行索引，數(shù)據(jù)可以使用相同的索引數(shù)組存儲(chǔ)。在從核并行計(jì)算時(shí)，雖然矩陣上三角每行元素按序排列，但是下三角每行元素并不連續(xù)，訪問下三角元素時(shí)很難滿足空間局部性。

2.2 計(jì)算特征分析

非結(jié)構(gòu)網(wǎng)格能夠適應(yīng)各種復(fù)雜結(jié)構(gòu)的面網(wǎng)格劃分，是流體力學(xué)仿真軟件中最主要的空間離散化解決方案。非結(jié)構(gòu)網(wǎng)格計(jì)算存在3 種模式，即頂點(diǎn)狀態(tài)更新（點(diǎn)更新）、邊狀態(tài)更新（邊更新）和通過鄰居頂點(diǎn)與鄰接邊更新頂點(diǎn)狀態(tài)（邊點(diǎn)更新）。

如圖5 所示，邊更新計(jì)算特征為S(e)=f(S(e)，S'(v1)，S'(v2))，典型計(jì)算函數(shù)為S(e)+=∑(S'(vi))。

圖5 邊更新模式Fig.5 The edge update mode

如圖6 所示，點(diǎn)更新計(jì)算特征為S(v)=f(S(v)，S'(e1)，S'(e2)，S'(e3)，S'(e4)，S'(e5))，典型計(jì)算函數(shù)為S(v)+=∑(S'(ei))。

圖6 點(diǎn)更新模式Fig.6 The point update mode

如圖7 所示，邊點(diǎn)更新計(jì)算特征為S(v)=f(S(v)，S'(e1)，S'(e2)，S'(e3)，S'(e4)，S'(e5)，S''(v1)，S''(v2)，S''(v3)，S''(v4)，S''(v5))，典型計(jì)算函數(shù)為S(v)+=∑(S'(ei) ·S''(vi))。

圖7 邊點(diǎn)更新模式Fig.7 The edge and point update mode

3 種模式的共同特征在于狀態(tài)信息在頂點(diǎn)與邊之間流動(dòng)。本文將邊視為基本單元，頂點(diǎn)視為連接單元，計(jì)算過程則為遍歷非結(jié)構(gòu)網(wǎng)格的所有邊，獲取每條邊上左右頂點(diǎn)狀態(tài)，并與邊自身狀態(tài)進(jìn)行運(yùn)算，最終用運(yùn)算結(jié)果更新左右頂點(diǎn)狀態(tài)或邊自身狀態(tài)，如表3 所示。

表3 非結(jié)構(gòu)網(wǎng)格計(jì)算模式Table 3 Unstructured grid computing mode

由于數(shù)據(jù)關(guān)系的依賴性和對(duì)數(shù)據(jù)訪問的離散性等固有特點(diǎn)，導(dǎo)致非結(jié)構(gòu)網(wǎng)格計(jì)算具有局部性差、數(shù)據(jù)相關(guān)、離散訪存、并發(fā)度低等問題，在眾核處理器上進(jìn)行優(yōu)化時(shí)難度較大，性能偏低。

3 基于申威處理器的眾核加速方法

3.1 算法思想

由于非結(jié)構(gòu)網(wǎng)格的稀疏性，導(dǎo)致算法在計(jì)算時(shí)對(duì)元素訪問并不連續(xù)，無法充分利用訪存帶寬。此外，非結(jié)構(gòu)網(wǎng)格的離散寬度（稀疏矩陣中非零元素與該行對(duì)角線元素之間的距離）較大，造成訪存間隔過大，難以滿足空間局部性。如圖8 所示，1 號(hào)邊和2 號(hào)邊的距離較大，遍歷時(shí)雖然行索引相同，但是列索引相距較遠(yuǎn)，不滿足空間局部性，訪存性能較差。因此，本文采用分塊劃分的思想，將一段時(shí)間內(nèi)的訪存數(shù)據(jù)盡可能集中存儲(chǔ)在較快的存儲(chǔ)器上，降低從下層存儲(chǔ)器讀取數(shù)據(jù)的時(shí)間開銷。

圖8 非結(jié)構(gòu)網(wǎng)格的鄰接矩陣Fig.8 The adjacent matrix of unstructured grid

此外，非結(jié)構(gòu)網(wǎng)格計(jì)算存在大量的對(duì)稱矩陣，為節(jié)省存儲(chǔ)空間，一般僅保留上三角矩陣，但是在計(jì)算時(shí)需要對(duì)稱更新，因此，在眾核處理器上的并行化存在數(shù)據(jù)相關(guān)和輸出相關(guān)（寫沖突）的問題。例如，在圖8 中，1 號(hào)邊和2 號(hào)邊位于同一行，表示其對(duì)相同目標(biāo)結(jié)果進(jìn)行更新，存在依賴關(guān)系，如果計(jì)算任務(wù)被分配至2 個(gè)不同的線程執(zhí)行，可能會(huì)發(fā)生寫后讀沖突。此外，1 號(hào)邊和6 號(hào)邊位于同一列，且1′號(hào)邊和6 號(hào)邊位于同一行，由于對(duì)稱更新，則1 號(hào)邊和1′號(hào)邊同時(shí)更新時(shí)如果其他計(jì)算任務(wù)對(duì)6 號(hào)邊更新，就會(huì)發(fā)生寫寫沖突。本文在分析非結(jié)構(gòu)網(wǎng)格的數(shù)據(jù)特點(diǎn)和計(jì)算特征后，提出并行度更高的無依賴任務(wù)劃分方法，將數(shù)據(jù)相關(guān)和輸出相關(guān)的計(jì)算分配到相同任務(wù)隊(duì)列。

本文提出一種N 階對(duì)角染色算法，非結(jié)構(gòu)網(wǎng)格邊線沿對(duì)角方向劃分為大小相同的方塊后，將有依賴關(guān)系的方塊染上同種顏色，分配到同一任務(wù)隊(duì)列中進(jìn)行并行計(jì)算。然后，染色器不斷向外擴(kuò)展并對(duì)其他類對(duì)角方塊染色。算法根據(jù)方塊內(nèi)元素密度決定染色階數(shù)，即向外擴(kuò)展對(duì)角線的次數(shù)。該算法支持非結(jié)構(gòu)網(wǎng)格下大多數(shù)的算子模型，特別是有依賴關(guān)系的算子。算法執(zhí)行步驟如下：

1）網(wǎng)格預(yù)處理及自適應(yīng)劃分。獲得當(dāng)前頂點(diǎn)數(shù)和邊數(shù)，記錄邊線所連頂點(diǎn)。根據(jù)LDM 存儲(chǔ)空間等，針對(duì)不同網(wǎng)格拓?fù)渥赃m應(yīng)確定分塊大?。ㄟ吘€所連頂點(diǎn)范圍）和染色階數(shù)，保證計(jì)算單元負(fù)載均衡。

2）分塊染色及重排整理。根據(jù)分塊大小，從對(duì)角塊向外對(duì)網(wǎng)格逐階染色，按照邊隨頂點(diǎn)走、一階一類色的原則，同時(shí)建立索引表記錄頂點(diǎn)的塊內(nèi)位置與全局位置關(guān)系，方便后續(xù)計(jì)算結(jié)果更新。

3）任務(wù)調(diào)度。同色塊分配至同一任務(wù)隊(duì)列，采用動(dòng)態(tài)調(diào)度方法管理任務(wù)隊(duì)列以維持從核負(fù)載平衡。

4）訪存及計(jì)算。從核通過DMA 操作完成網(wǎng)格邊線重排序，將當(dāng)前隊(duì)列內(nèi)染色塊加載至LDM 中并執(zhí)行計(jì)算。同時(shí)，為了隱藏DMA 操作時(shí)間，從核在進(jìn)行當(dāng)前計(jì)算時(shí)開始下一輪DMA 操作，使得計(jì)算與訪存重疊。在從核計(jì)算過程中，主核同時(shí)負(fù)責(zé)未染色塊計(jì)算，因?yàn)槲慈旧珘K更為稀疏，局部性更差，更適合通過主核計(jì)算，而從核通過DMA 對(duì)數(shù)據(jù)的換入換出往往會(huì)帶來更高的時(shí)間開銷，不利于發(fā)揮其性能優(yōu)勢(shì)。

3.2 算例分析

考慮到計(jì)算的常見性和代表性，本文以典型算子稀疏矩陣向量乘（Sparse Matrix Vector Multiplication，SpMV）為例分析眾核加速方法。

作為最常見的稀疏運(yùn)算，雙精度SpMV 的計(jì)算強(qiáng)度僅為0.080～0.125 FLOPs/Byte，在帶寬受限的眾核處理器上性能較差。SpMV 算子描述如算法1 所示。

算法1SpMV 算子

在算法1 中：V為輸入/輸出向量值，即網(wǎng)格頂點(diǎn)狀態(tài)；E.row/E.col 和E.val 分別為稀疏矩陣的行/列索引和值，即網(wǎng)格邊狀態(tài)。頂點(diǎn)狀態(tài)的更新由相連邊狀態(tài)及其連接頂點(diǎn)狀態(tài)的乘積累加得到。本文所提算法執(zhí)行過程如下：

1）在算法2 中，設(shè)最小并行單位為Δ，即分塊最小頂點(diǎn)范圍。根據(jù)LDM 空間大小、DMA 特性和讀取的網(wǎng)格拓?fù)湫畔?，自適應(yīng)確定分塊因子大小α。

算法2分塊因子判決

2）在算法3 中，根據(jù)分塊大小αΔ掃描邊線，先將主對(duì)角塊染色并建立邊索引表，例如第一塊頂點(diǎn)范圍在0～（αΔ-1）內(nèi)，第二塊頂點(diǎn)范圍在αΔ～（2αΔ-1）內(nèi)，同時(shí)將頂點(diǎn)全局位置轉(zhuǎn)換為塊內(nèi)位置，建立關(guān)系映射表，原因是塊內(nèi)計(jì)算時(shí)不能使用全局地址。按照同樣的方法向二階及以上階擴(kuò)展，皆為雙側(cè)次對(duì)角塊，即包括上三角和對(duì)稱的下三角兩部分。隨后，算法4 從對(duì)角塊中挑選較稠密對(duì)角塊進(jìn)行染色，挑選標(biāo)準(zhǔn)由塊內(nèi)節(jié)點(diǎn)密度和網(wǎng)格整體密度決定，根據(jù)大量測(cè)試后得出。未被挑選的非稠密塊則分配給主核任務(wù)隊(duì)列，依據(jù)主從核的不同特點(diǎn)實(shí)現(xiàn)任務(wù)分配。三階對(duì)角染色示意圖如圖9 所示。

圖9 三階對(duì)角染色Fig.9 Third-order diagonal dyeing

算法3分塊染色及重排整理

算法4挑選染色對(duì)角

3）建立從核任務(wù)隊(duì)列queue，將全部一階色塊分入隊(duì)列，從核從隊(duì)列中獲取任務(wù)并完成計(jì)算。類似地，在一階色塊完成計(jì)算后，其他同階色塊依次被分配到任務(wù)隊(duì)列，隊(duì)列內(nèi)從核運(yùn)行狀態(tài)一致，從而避免寫后讀沖突和寫寫沖突。

4）在算法5 中，從核通過DMA 獲取塊內(nèi)頂點(diǎn)狀態(tài)和邊狀態(tài)以及索引表并完成更新，在計(jì)算時(shí)可以預(yù)取下一輪數(shù)據(jù)，從而使得DMA 時(shí)間被有效隱藏，如圖10 所示。在從核計(jì)算的同時(shí)，主核負(fù)責(zé)其他未染色塊的計(jì)算，從而實(shí)現(xiàn)主從核異步并行，進(jìn)一步提升計(jì)算效率。

圖10 計(jì)算訪存異步重疊Fig.10 The asynchronous overlap of computation and memory access

算法5訪存及計(jì)算

4 實(shí)驗(yàn)結(jié)果與分析

本次實(shí)驗(yàn)基于申威26010 眾核處理器，硬件參數(shù)如表4 所示，采用swg++編譯器編譯全部C/C++程序。

表4 申威26010 處理器硬件參數(shù)Table 4 The hardware parameters of SW26010 processor

4.1 不同網(wǎng)格的性能分析

為保證算法性能的可靠性，本文選取典型稀疏算子SpMV 作為標(biāo)準(zhǔn)測(cè)試算子，隨機(jī)輸入非結(jié)構(gòu)網(wǎng)格實(shí)例進(jìn)行測(cè)試分析。圖11 和圖12 分別為SpMV算子在不同網(wǎng)格規(guī)模下加速算法與主核樸素算法的運(yùn)行時(shí)間及加速比，以驗(yàn)證加速算法的通用性和優(yōu)化效果。

圖11 不同網(wǎng)格規(guī)模下的優(yōu)化性能Fig.11 Optimization performance under different grid scales

圖12 不同網(wǎng)格規(guī)模下的加速效果Fig.12 Acceleration effect under different grid scales

從圖11 和圖12 中可以看出，隨著網(wǎng)格規(guī)模的增加，加速算法的加速效果基本保持穩(wěn)定，因?yàn)榫W(wǎng)格劃分根據(jù)輸入網(wǎng)格密度和拓?fù)渥赃m應(yīng)調(diào)整，染色階數(shù)也根據(jù)當(dāng)前對(duì)角線密度自動(dòng)判決，因此加速算法能在多數(shù)網(wǎng)格規(guī)模下保持穩(wěn)定的性能優(yōu)勢(shì)，通用性較強(qiáng)。相比于主核上運(yùn)行的SpMV 算子，組合加速算法獲得了平均10 倍左右的加速比，最高加速比可達(dá)24 倍。

4.2 不同算子的性能分析

本文設(shè)計(jì)非結(jié)構(gòu)網(wǎng)格計(jì)算在申威眾核處理器上的通用加速方法，因此，需要選取多種算子進(jìn)行綜合分析。SpMV 算子的加速比如圖12 所示，Integration算子和calcLudsFcc 算子的加速比分別如圖13 和圖14 所示。

圖13 Integration 算子在不同網(wǎng)格規(guī)模下的加速效果Fig.13 Acceleration effect of Integration operator under different grid sizes

圖14 calcLudsFcc 算子在不同網(wǎng)格規(guī)模下的加速效果Fig.14 Acceleration effect of calcLudsFcc operator under different grid sizes

經(jīng)過對(duì)上述2 種算子在不同網(wǎng)格規(guī)模下的測(cè)試發(fā)現(xiàn)，組合加速算法分別獲得了10.22 倍和6.82 倍的平均加速比，而且本文算法對(duì)不同算子模型的性能表現(xiàn)差異并不明顯，在有依賴和無依賴情況下都能取得穩(wěn)定的性能優(yōu)勢(shì)，說明算法在任務(wù)劃分和數(shù)據(jù)映射上并沒有以犧牲性能為代價(jià)，自適應(yīng)和無依賴任務(wù)劃分方法取得了良好效果。由于算子在從核的計(jì)算和訪存是基于染色后的數(shù)據(jù)塊，其加速效果與數(shù)據(jù)塊中的數(shù)據(jù)分布有關(guān)，在數(shù)據(jù)集中度較高的網(wǎng)格實(shí)例中，算子能獲得顯著的性能提升，可達(dá)20 多倍，但在數(shù)據(jù)非常離散的情況下效果一般。

4.3 不同優(yōu)化策略的性能分析

為了說明N 階對(duì)角染色算法和自適應(yīng)任務(wù)劃分方法的有效性，以SpMV 算子為例，分別采用非N 階對(duì)角染色的分塊算法和固定分塊大小的任務(wù)劃分方法進(jìn)行對(duì)比實(shí)驗(yàn)，并以主核樸素算法為基準(zhǔn)，實(shí)驗(yàn)結(jié)果如圖15 所示。

圖15 不同優(yōu)化方法的加速效果Fig.15 Acceleration effect of different optimization methods

N 階對(duì)角染色算法通過分析對(duì)角塊密度來選擇是否染色當(dāng)前對(duì)角塊，而普通分塊算法缺少對(duì)角塊信息，易將過于稀疏的對(duì)角塊交由從核陣列計(jì)算。將全部矩陣塊映射到從核陣列會(huì)造成極大的性能損失，本文通過固定前100 階對(duì)角塊由從核計(jì)算來模擬非染色的普通分塊算法性能。自適應(yīng)劃分方法根據(jù)LDM 空間大小、DMA 特性和網(wǎng)格拓?fù)湫畔⒋_定分塊大小，可以充分利用LDM 空間，而固定分塊大小則容易造成對(duì)LDM 空間的利用不足。實(shí)驗(yàn)結(jié)果表明，非N 階對(duì)角染色的分塊算法平均加速比為2.64 倍，固定分塊大小的任務(wù)劃分方法平均加速比僅為1.8 倍，難以發(fā)揮眾核架構(gòu)的計(jì)算能力，甚至有負(fù)優(yōu)化效果出現(xiàn)。采用自適應(yīng)任務(wù)劃分的N 階對(duì)角染色算法能有效利用LDM 空間并根據(jù)塊內(nèi)密度選擇從核或主核來執(zhí)行計(jì)算，平均加速比可達(dá)10 倍。

5 結(jié)束語

為提升非結(jié)構(gòu)網(wǎng)格計(jì)算中有依賴關(guān)系算子在眾核處理器上的性能，本文針對(duì)非結(jié)構(gòu)網(wǎng)格的計(jì)算特點(diǎn)，提出一種眾核處理器上的通用加速方法，并基于申威26010 處理器架構(gòu)對(duì)其進(jìn)行精細(xì)調(diào)優(yōu)。通過自適應(yīng)任務(wù)劃分方法將從核離散訪存組織為批量訪存，以降低訪存開銷。采用無依賴劃分策略避免計(jì)算時(shí)的數(shù)據(jù)沖突，通過N 階對(duì)角染色算法將計(jì)算任務(wù)分類調(diào)度執(zhí)行，從而有效利用主從核的架構(gòu)特點(diǎn)。此外，采用計(jì)算訪存重疊技術(shù)進(jìn)一步提升計(jì)算性能。實(shí)驗(yàn)結(jié)果表明，該方法在不同網(wǎng)格規(guī)模和不同算子模型下都取得了良好的加速效果，在有依賴和無依賴情況下均能保持穩(wěn)定的性能優(yōu)勢(shì)，證明了任務(wù)劃分方法的有效性。但是，本文所提方法對(duì)于數(shù)據(jù)分布極為分散的非結(jié)構(gòu)網(wǎng)格仍存在一定局限性，下一步將結(jié)合排序算法對(duì)網(wǎng)格數(shù)據(jù)進(jìn)行重排整理，提升數(shù)據(jù)的局部性，增加在從核陣列計(jì)算的數(shù)據(jù)塊，從而滿足更多稀疏網(wǎng)格數(shù)據(jù)的眾核計(jì)算需求。