黎雷生 ,楊文浩 ,馬文靜 ,張 婭 ,趙 慧 ,趙海濤 ,李會(huì)元,孫家昶

1(中國(guó)科學(xué)院 軟件研究所 并行軟件與計(jì)算科學(xué)實(shí)驗(yàn)室,北京 100190)

2(計(jì)算機(jī)科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室(中國(guó)科學(xué)院 軟件研究所),北京 100190)

HPL(high performance Linpack)是評(píng)測(cè)計(jì)算系統(tǒng)性能的程序,是早期Linpack 評(píng)測(cè)程序的并行版本,支持大規(guī)模并行超級(jí)計(jì)算系統(tǒng)[1],其報(bào)告的每秒浮點(diǎn)運(yùn)算次數(shù)(floating-point operations per second,簡(jiǎn)稱FLOPS)是世界超級(jí)計(jì)算機(jī)Top500 列表排名的依據(jù).

基于Linpack 的Top500 排名開(kāi)始于1993 年,每年發(fā)表兩次.2007 年11 月及之前的列表,排名前10 位的超級(jí)計(jì)算機(jī)的計(jì)算能力全部由同構(gòu)CPU提供.2008年6月Top500首臺(tái)性能超過(guò)1 PFLOPS 的超計(jì)算機(jī)Roadrunner使用異構(gòu)CPU 結(jié)構(gòu),通用CPU 執(zhí)行調(diào)度、管理等任務(wù),加速CPU 主要負(fù)責(zé)計(jì)算密集任務(wù).2009 年11 月排名第5 位的Tianhe-1 使用了CPU+GPU 的異構(gòu)架構(gòu).此后榜單上排名前10 的系統(tǒng)CPU+加速器的架構(gòu)成為趨勢(shì).2019年6 月最新的排名顯示,前10 位中有7 臺(tái)系統(tǒng)使用CPU+加速器架構(gòu),其中使用GPU 加速器的5 臺(tái),使用XEON Phi 的1 臺(tái),使用Matrix-2000 的1 臺(tái)[2].

HPL 浮點(diǎn)計(jì)算集中在BLAS 函數(shù),特別是DGEMM.對(duì)于同構(gòu)CPU 架構(gòu),優(yōu)化BLAS 函數(shù)特別是DGEMM性能是提高HPL 的浮點(diǎn)效率的關(guān)鍵.Dongarra 等人[1]總結(jié)了2002 年之前Linpack 發(fā)展歷史和BLAS 函數(shù)的優(yōu)化方法.

對(duì)于CPU+加速器架構(gòu),優(yōu)化方向集中于CPU端BLAS 函數(shù)、加速器端BLAS 函數(shù)、CPU 與加速器之間負(fù)載分配和數(shù)據(jù)傳輸?shù)?Kurzak 等人[3]優(yōu)化多核CPU 多GPU 的HPL,其4 個(gè)GPU 的DGEMM 浮點(diǎn)性能達(dá)到1 200 GFLOPS.Bach 等人[4]面向AMD CPU 與Cypress GPU 架構(gòu)優(yōu)化了DGEMM 和HPL,HPL 效率達(dá)到70%.Wang 等人[5]采用自適應(yīng)負(fù)載分配方法優(yōu)化CPU 與GPU 混合系統(tǒng)HPL 性能,并調(diào)優(yōu)HPL 參數(shù),目標(biāo)系統(tǒng)Tianhe1A 浮點(diǎn)性能達(dá)到0.563 PFLOPS.Heinecke 等人[6]基于Intel Xeon Phi Knights Corner 優(yōu)化了DGEMM,使用動(dòng)態(tài)調(diào)度和改進(jìn)的look-ahead 優(yōu)化HPL,100 節(jié)點(diǎn)集群HPL 效率達(dá)到76%.Gan 等人[7]優(yōu)化基于加速器的DGEMM,并應(yīng)用于Tianhe-2 HPL 評(píng)測(cè).

隨著以GPU 為代表的加速器技術(shù)的發(fā)展,加速器浮點(diǎn)性能越來(lái)越高,CPU 與加速器的浮點(diǎn)性能差距越來(lái)越大.2019 年6 月,Top500 排名第一的Summit 系統(tǒng)[8]的1 個(gè)節(jié)點(diǎn)裝備2 個(gè)CPU,理論浮點(diǎn)性能1 TFLOPS,裝備6個(gè)GPU,理論浮點(diǎn)性能42 TFLOPS.本文研究的目標(biāo)系統(tǒng)使用CPU+GPU 異構(gòu)架構(gòu),每個(gè)節(jié)點(diǎn)裝備1 個(gè)32 核CPU,4 個(gè)GPU,CPU 浮點(diǎn)計(jì)算性能約是GPU 的1/61.同時(shí)加速器本身的結(jié)構(gòu)也變得越來(lái)越復(fù)雜,通過(guò)增加特定的硬件滿足特定領(lǐng)域的需求,如Nvidia GPU 的Tensor Core 等.已有研究使用Tensor Core 的強(qiáng)大的半精度運(yùn)算能力混合雙精度計(jì)算開(kāi)發(fā)了HPL-AI[9],報(bào)告Summit 的HPL-AI 性能是全雙精度HPL 的2.9 倍.并且已有應(yīng)用采用混合精度算法加速計(jì)算,從HPL 和應(yīng)用角度來(lái)看,混合精度都是值得研究的方向.面對(duì)這種新的計(jì)算架構(gòu),內(nèi)存、總線、網(wǎng)絡(luò)以及系統(tǒng)設(shè)計(jì)都要與之相適應(yīng),形成復(fù)雜的異構(gòu)計(jì)算系統(tǒng),這為HPL 評(píng)測(cè)帶來(lái)很多機(jī)遇與挑戰(zhàn).

本文在基礎(chǔ)HPL 代碼之上,針對(duì)目標(biāo)系統(tǒng)實(shí)現(xiàn)HPL 并研究其優(yōu)化方法.第1 節(jié)簡(jiǎn)要介紹基礎(chǔ)算法.第2 節(jié)介紹目標(biāo)系統(tǒng)基本情況.第3 節(jié)描述復(fù)雜異構(gòu)系統(tǒng)HPL 平衡點(diǎn)理論.第4 節(jié)介紹復(fù)雜異構(gòu)系統(tǒng)HPL 高效并行算法.第5 節(jié)介紹基礎(chǔ)模塊的性能優(yōu)化,包括panel 分解優(yōu)化和行交換long算法的優(yōu)化.第6 節(jié)介紹目標(biāo)系統(tǒng)HPL實(shí)驗(yàn)和結(jié)果分析.第7 節(jié)總結(jié)并展望未來(lái)的工作.

1 基礎(chǔ)算法

HPL算法使用64 位浮點(diǎn)精度矩陣行偏主元LU 分解加回代求解線性系統(tǒng).矩陣是稠密實(shí)矩陣,矩陣單元由偽隨機(jī)數(shù)生成器產(chǎn)生,符合正態(tài)分布.

線性系統(tǒng)定義為:

行偏主元LU 分解N×(N+1)系數(shù)矩陣[A,b]:

其中,Pr表示行主元交換矩陣,分解過(guò)程中下三角矩陣因子L已經(jīng)作用于b,解x通過(guò)求解上三角矩陣系統(tǒng)得到:

HPL 采用分塊LU算法,每個(gè)分塊是一個(gè)NB列的細(xì)長(zhǎng)矩陣,稱為panel.LU 分解主循環(huán)采用right-looking算法,單步循環(huán)計(jì)算panel 的LU 分解和更新剩余矩陣.基本算法如圖1 所示,其中A1,1和A2,1表示panel 數(shù)據(jù).需要特別說(shuō)明的是,圖示矩陣是行主順序,HPL 代碼中矩陣是列主存儲(chǔ)的.

Fig.1 Block LU algorithm圖1 分塊LU算法

計(jì)算公式如下:

第1 個(gè)公式表示panel 的LU 分解,第2 個(gè)公式表示求解U,一般使用DTRSM 函數(shù),第3 個(gè)公式表示矩陣更新,一般使用DGEMM 函數(shù).

對(duì)于分布式內(nèi)存計(jì)算系統(tǒng),HPL 并行計(jì)算模式基于MPI,每個(gè)進(jìn)程是基本計(jì)算單元.進(jìn)程組織成二維網(wǎng)格.矩陣A被劃分為NB×NB的邏輯塊,以Block-Cycle 方式分配到二維進(jìn)程網(wǎng)格,數(shù)據(jù)布局示例如圖2 所示.

Fig.2 Block-Cycle algorithm圖2 Block-Cycle算法

對(duì)于具有多列的進(jìn)程網(wǎng)格,單步循環(huán)只有一列進(jìn)程執(zhí)行panel 分解計(jì)算,panel 分解過(guò)程中每一列執(zhí)行一次panel 的行交換算法選擇并通信最大主元行.Panel 分解計(jì)算完成后,把已分解數(shù)據(jù)廣播到其他進(jìn)程列.HPL 基礎(chǔ)代碼包含6 類廣播算法,可以通過(guò)測(cè)試選擇較優(yōu)的算法.

HPL 采用行主元算法,單步矩陣更新之前,要把panel 分解時(shí)選出的最大主元行交換到U矩陣中,需要執(zhí)行未更新矩陣的主元行交換和廣播.主元行交換和廣播后,每個(gè)進(jìn)程獲得完整的主元行數(shù)據(jù).

矩陣更新包括兩部分計(jì)算,一是使用DTRSM 求解U,二是使用DGEMM 更新矩陣數(shù)據(jù).

LU 分解完成后,HPL 使用回代求解x,并驗(yàn)證解的正確性.

2 系統(tǒng)介紹

2.1 系統(tǒng)硬件

優(yōu)化HPL 的目標(biāo)系統(tǒng)是典型的復(fù)雜異構(gòu)系統(tǒng),計(jì)算部件由通用CPU 和計(jì)算加速器GPU 組成,網(wǎng)絡(luò)接口和互聯(lián)系統(tǒng)是高速Infiniband 網(wǎng)絡(luò).

2.2 基礎(chǔ)軟件

目標(biāo)CPU 兼容x86 指令集.編譯器使用gcc,MPI 使用OpenMPI.

HPL 基礎(chǔ)代碼使用2.3 版本.HPL 會(huì)調(diào)用多個(gè)BLAS 函數(shù).CPU端BLAS 使用針對(duì)目標(biāo)CPU 優(yōu)化的OpenBLIS,支持OpenMP 實(shí)現(xiàn)多核并行.GPU端BLAS 使用針對(duì)目標(biāo)GPU 優(yōu)化的BLAS.

3 復(fù)雜異構(gòu)系統(tǒng)HPL 平衡點(diǎn)理論

從硬件的角度來(lái)看,CPU 的每個(gè)Die 與對(duì)應(yīng)的內(nèi)存系統(tǒng)、GPU 和網(wǎng)絡(luò)接口組成相對(duì)獨(dú)立的系統(tǒng).因此,在選擇并行方案時(shí),采用1 個(gè)進(jìn)程使用1 個(gè)Die,并且與相應(yīng)的內(nèi)存、GPU 和網(wǎng)絡(luò)接口綁定的架構(gòu).HPL 基礎(chǔ)代碼中只有1 列或1 行進(jìn)程時(shí)執(zhí)行流程不同,較大規(guī)模并行一般是多行多列進(jìn)程,因此使用1 個(gè)節(jié)點(diǎn)4 進(jìn)程作為研究對(duì)象,覆蓋多行多列進(jìn)程執(zhí)行流程.

對(duì)于復(fù)雜異構(gòu)系統(tǒng)HPL 計(jì)算,需要協(xié)同調(diào)度CPU、GPU、PCIe 和通信網(wǎng)絡(luò)等資源.首先研究CPU 與GPU計(jì)算任務(wù)分配方式,進(jìn)程內(nèi)CPU 與GPU 以及進(jìn)程間的數(shù)據(jù)傳輸,并且建立線程模型控制CPU、GPU 和網(wǎng)絡(luò)接口協(xié)同計(jì)算.然后分析HPL 各部分的性能,提出復(fù)雜異構(gòu)系統(tǒng)的平衡點(diǎn)理論指導(dǎo)性能優(yōu)化.

3.1 計(jì)算任務(wù)分配

文獻(xiàn)[4,6]等研究的異構(gòu)系統(tǒng)中,CPU 執(zhí)行panel 分解計(jì)算,加速器執(zhí)行矩陣更新的DTRSM 和DGEMM 計(jì)算.CPU 計(jì)算能力可以達(dá)到加速器計(jì)算能力的10%或更高,除了panel 分解,CPU 還可以計(jì)算部分矩陣更新.并且由于加速器內(nèi)存限制,矩陣不能完全存儲(chǔ)在加速器,矩陣更新一般采用流水的方式.

經(jīng)過(guò)分析和實(shí)驗(yàn)對(duì)比,由于panel 分解控制流程復(fù)雜,并且大量調(diào)用GPU 計(jì)算效率較低的DGEMV、DSCAL等函數(shù),不適合在GPU端實(shí)現(xiàn),因此,仍然采用CPU 計(jì)算panel 分解的方案.

對(duì)于矩陣更新來(lái)說(shuō),目標(biāo)系統(tǒng)CPU 和GPU 之間的數(shù)據(jù)傳輸成為瓶頸,導(dǎo)致矩陣更新流水線方法不再適用.針對(duì)單個(gè)GPU 分析,假設(shè)矩陣規(guī)模N=58368(占用約80%的系統(tǒng)內(nèi)存),NB=384,矩陣更新的DGEMM 效率為85%.可以計(jì)算出矩陣更新計(jì)算時(shí)間0.52s.PCIe 單向傳輸矩陣數(shù)據(jù)的理論時(shí)間1.70s.數(shù)據(jù)傳輸時(shí)間超過(guò)GPU 的計(jì)算時(shí)間,導(dǎo)致GPU 超過(guò)2/3 處于空閑狀態(tài),GPU 計(jì)算能力不能充分利用.

因此,采用矩陣數(shù)據(jù)常駐GPU 內(nèi)存的方式,矩陣更新時(shí)避免傳輸整個(gè)剩余矩陣.這種模式下,矩陣的規(guī)模受限于GPU 設(shè)備內(nèi)存的容量.為了運(yùn)算方便,CPU端內(nèi)存保留GPU端數(shù)據(jù)存儲(chǔ)的鏡像.

3.2 數(shù)據(jù)傳輸

數(shù)據(jù)傳輸包括進(jìn)程內(nèi)主機(jī)內(nèi)存與GPU 設(shè)備內(nèi)存之間數(shù)據(jù)傳輸和進(jìn)程間數(shù)據(jù)傳輸.主機(jī)內(nèi)存與GPU 設(shè)備內(nèi)存數(shù)據(jù)傳輸使用PCIe 總線,一般采用DMA 異步模式,數(shù)據(jù)傳輸?shù)耐瑫r(shí)CPU 和GPU 可以執(zhí)行計(jì)算任務(wù).進(jìn)程間數(shù)據(jù)傳輸使用MPI.如果數(shù)據(jù)傳輸?shù)倪M(jìn)程位于同一個(gè)節(jié)點(diǎn),MPI 使用共享內(nèi)存等技術(shù)提高數(shù)據(jù)傳輸性能.節(jié)點(diǎn)之間的數(shù)據(jù)傳輸通過(guò)網(wǎng)絡(luò)接口.

由于矩陣數(shù)據(jù)常駐GPU 內(nèi)存,panel 分解前把當(dāng)前的panel 數(shù)據(jù)從GPU 傳輸?shù)紺PU.Panel 分解在CPU端執(zhí)行panel 的主元選取與廣播,使用基礎(chǔ)代碼的算法.已完成panel 分解的數(shù)據(jù),進(jìn)程內(nèi)部通過(guò)PCIe 傳輸?shù)紾PU.進(jìn)程之間通過(guò)MPI 傳輸,非當(dāng)前panel 分解進(jìn)程接收到panel 數(shù)據(jù)后,再通過(guò)PCIe 傳輸?shù)紾PU.

相比于基礎(chǔ)代碼,矩陣更新行交換和廣播增加了CPU 與GPU 數(shù)據(jù)傳輸過(guò)程.矩陣的行在GPU 內(nèi)存不連續(xù),數(shù)據(jù)傳輸前后需要執(zhí)行封裝和重排.數(shù)據(jù)傳輸關(guān)系如圖3 所示.

Fig.3 Data transmisson圖3 數(shù)據(jù)傳輸

3.3 多線程模型

在1 個(gè)節(jié)點(diǎn)4 個(gè)進(jìn)程的架構(gòu)下,每個(gè)進(jìn)程使用8 個(gè)CPU 核心,使用多線程管理和調(diào)度這些核心.Panel 分解和panel 廣播共享數(shù)據(jù),執(zhí)行順序存在依賴關(guān)系,放在同一個(gè)線程.矩陣行交換和矩陣更新共享數(shù)據(jù),執(zhí)行順序存在依賴關(guān)系,放在同一個(gè)線程.兩組操作之間相對(duì)獨(dú)立,分別綁定到特定核心.兩個(gè)線程之間通過(guò)信號(hào)量同步.

Panel 分解線程調(diào)用和管理OpenBLIS 函數(shù),OpenBLIS 使用OpenMP 線程模型.Panel 分解調(diào)用OpenBLIS函數(shù)時(shí),與panel 分解的其他過(guò)程沒(méi)有資源共享和沖突,OpenBLIS 的OpenMP 線程綁定的核心可以包括panel分解主線程的核心,也就是最多可綁定7 個(gè)核心.

3.4 復(fù)雜異構(gòu)系統(tǒng)協(xié)同計(jì)算平衡點(diǎn)理論

目標(biāo)系統(tǒng)是復(fù)雜異構(gòu)系統(tǒng),涉及CPU、GPU、PCIe 和網(wǎng)絡(luò)接口等部件的協(xié)同計(jì)算.為了更好分析和評(píng)測(cè)各部分對(duì)HPL 性能的影響,指導(dǎo)性能優(yōu)化,提出平衡點(diǎn)理論.

HPL算法中,矩陣更新的浮點(diǎn)計(jì)算占據(jù)絕大部分,這部分計(jì)算由GPU 完成.性能分析過(guò)程中,矩陣更新被認(rèn)為是有效計(jì)算,其他功能模塊包括panel 分解和矩陣行交換的目標(biāo)是為矩陣更新提供數(shù)據(jù)準(zhǔn)備.為了獲得較高的效率,要盡可能保證GPU 處于工作狀態(tài),發(fā)揮GPU 浮點(diǎn)計(jì)算能力.對(duì)于panel 分解和矩陣行交換,盡可能讓它們與GPU 計(jì)算并行,使用GPU 計(jì)算隱藏panel 分解和矩陣行交換.

研究沒(méi)有使用特殊優(yōu)化、單步循環(huán)各部分時(shí)間以及串行執(zhí)行的總時(shí)間,如圖4 所示.采用2×2 進(jìn)程布局,圖中是進(jìn)程0 的時(shí)間曲線.進(jìn)程0 在計(jì)算過(guò)程中處于不同的位置,包括當(dāng)前列-當(dāng)前行和非當(dāng)前行-非當(dāng)前列兩種情況,圖中”A”代表當(dāng)前列-當(dāng)前行,”B”代表非當(dāng)前行-非當(dāng)前列.兩種情況下執(zhí)行的計(jì)算有差異,每種運(yùn)算有兩條時(shí)間曲線.panel 分解時(shí)間在不同規(guī)模下稍有波動(dòng).可以看出,單步循環(huán)串行執(zhí)行時(shí)間遠(yuǎn)大于GPU 執(zhí)行有效計(jì)算的時(shí)間.那么首先必須做好矩陣更新與panel 分解和行交換的并行,GPU 計(jì)算時(shí)間隱藏panel 分解與和行交換時(shí)間.

Fig.4 Time of panel factorization,row swap,matrix update and each loop圖4 Panel 分解、行交換、矩陣更新和三部分串行執(zhí)行累加時(shí)間

另一方面,HPL 的LU 分解過(guò)程分為兩個(gè)階段,矩陣更新時(shí)間大于panel 分解和行交換的時(shí)間的階段和矩陣更新時(shí)間小于panel 分解直至小于行交換時(shí)間的階段.在panel 分解與矩陣行交換并行的情況下,兩個(gè)階段的分界點(diǎn)即矩陣更新與panel 分解時(shí)間曲線的交點(diǎn)就是CPU 與GPU 協(xié)同計(jì)算的平衡點(diǎn).

針對(duì)1 節(jié)點(diǎn)4 進(jìn)程進(jìn)行分析,N=88320 情況下,從計(jì)算規(guī)模來(lái)看,平衡點(diǎn)在37.82%,LU 計(jì)算量占76.00%,之后GPU 不能充分發(fā)揮性能.有必要繼續(xù)優(yōu)化panel 分解和行交換,減少執(zhí)行時(shí)間,平衡點(diǎn)后移,延長(zhǎng)GPU 高效率計(jì)算時(shí)間,進(jìn)而提高整體的計(jì)算效率.平衡點(diǎn)理論可以很好地指導(dǎo)HPL 的優(yōu)化工作.

4 復(fù)雜異構(gòu)系統(tǒng)HPL 高效并行算法

根據(jù)第3.4 節(jié)的分析,首先做好GPU 計(jì)算與panel 分解和行交換的并行算法,實(shí)現(xiàn)GPU 計(jì)算重疊其他部分的時(shí)間.

4.1 GPU參與計(jì)算的look-ahead算法

基礎(chǔ)代碼用look-ahead算法實(shí)現(xiàn)panel 分解與panel 廣播重疊執(zhí)行.計(jì)算過(guò)程中保存兩個(gè)panel 的數(shù)據(jù)結(jié)構(gòu).執(zhí)行當(dāng)前panel 分解的進(jìn)程列,優(yōu)先執(zhí)行當(dāng)前panel 列數(shù)據(jù)的行交換和矩陣更新,然后執(zhí)行panel 分解計(jì)算并發(fā)起panel 廣播,廣播數(shù)據(jù)傳輸?shù)耐瑫r(shí)更新剩余矩陣.對(duì)于當(dāng)前不計(jì)算panel 分解的進(jìn)程列,執(zhí)行上一次循環(huán)的行交換和矩陣更新,行交換和矩陣更新計(jì)算過(guò)程中監(jiān)測(cè)本次panel 廣播數(shù)據(jù)通信,數(shù)據(jù)到達(dá)后,執(zhí)行廣播通信流程.

由于加入了GPU 計(jì)算資源,look-ahead算法相比于基礎(chǔ)版本有所改變.按照第3.3 節(jié)的線程模型,在執(zhí)行當(dāng)前panel 分解計(jì)算的進(jìn)程中,panel 分解線程CPU 優(yōu)先執(zhí)行panel 列的行交換,然后GPU 優(yōu)先更新panel 分解的列.CPU 等待panel 數(shù)據(jù)更新后,發(fā)起panel 數(shù)據(jù)從GPU 內(nèi)存到CPU 內(nèi)存的傳輸,傳輸結(jié)束后執(zhí)行panel 分解計(jì)算,最后發(fā)起panel 廣播.矩陣更新線程在panel 分解列行交換完成后發(fā)起剩余矩陣的行交換,然后調(diào)用GPU 矩陣更新計(jì)算.對(duì)于非當(dāng)前panel 分解計(jì)算進(jìn)程,矩陣更新線程執(zhí)行行交換和調(diào)用GPU 矩陣更新.Panel 分解線程只負(fù)責(zé)處理panel 廣播的操作.

4.2 行交換連續(xù)流水線算法

矩陣行交換的前提是已接收panel 分解計(jì)算的行交換索引數(shù)據(jù)(包含在panel 廣播數(shù)據(jù)中),并且上一次迭代矩陣更新已完成.HPL 基礎(chǔ)代碼中,行交換完成后啟動(dòng)矩陣更新計(jì)算,如果采用這種算法,GPU 在行交換執(zhí)行的過(guò)程是空閑的.文獻(xiàn)[6]提出行交換流水線算法優(yōu)化行交換流程.在使用GPU 的環(huán)境,也可以采用這一算法.使用行交換流水線算法,對(duì)矩陣更新的行分段,首先執(zhí)行第1 段的行交換,然后GPU 執(zhí)行第1 段的矩陣更新,同時(shí)執(zhí)行第2 段的行交換,依次類推.采用這種算法,GPU 計(jì)算與行交換數(shù)據(jù)傳輸重疊執(zhí)行,減少了GPU 的空閑時(shí)間.但是,單步循環(huán)執(zhí)行第1 段行交換的這段時(shí)間內(nèi)GPU 是空閑的.

為了避免GPU 等待,提出了連續(xù)流水線算法.算法中,第1 個(gè)分段完成上一次循環(huán)矩陣更新和當(dāng)前進(jìn)程接收到下一次循環(huán)的行交換信息之后,不需要等待上一次循環(huán)矩陣更新全部完成,就執(zhí)行第1 分段的下一次循環(huán)的行交換.采用這種算法,下一次循環(huán)第1 分段的行交換被上一次循環(huán)矩陣更新隱藏,下一次循環(huán)第1 分段的矩陣更新隱藏第2 分段的行交換,后續(xù)分段繼續(xù)流水.

使用連續(xù)流水算法避免單步循環(huán)的流水線啟動(dòng)過(guò)程,如果接收到行交換信息早于矩陣更新,則GPU 可以一直處于工作狀態(tài).從第3.4 節(jié)可以看出,行交換信息在平衡點(diǎn)之前是早于矩陣更新的,而這一階段正是浮點(diǎn)運(yùn)算量較大的情況,在這一階段充分發(fā)揮GPU 計(jì)算能力,可以提高整體的計(jì)算效率.

行交換分段方法選擇倍數(shù)遞增,也就是第2 段列數(shù)是第1 段列數(shù)的倍數(shù),依次類推,增長(zhǎng)因子可調(diào).連續(xù)行交換流水算法下實(shí)現(xiàn)look-ahead算法,對(duì)于當(dāng)前panel 分解,增加panel 分解數(shù)據(jù)提前行交換和矩陣更新的分段.

5 基礎(chǔ)模塊性能優(yōu)化

從平衡點(diǎn)理論來(lái)看,平衡點(diǎn)之前的部分GPU 計(jì)算時(shí)間可以隱藏panel 分解和行交換時(shí)間.隨著計(jì)算的推進(jìn),panel 分解和行交換成為系統(tǒng)性能的瓶頸.為了進(jìn)一步提高效率,需要針對(duì)panel 分解和行交換進(jìn)行優(yōu)化.

5.1 Panel分解優(yōu)化

5.1.1 基本參數(shù)調(diào)優(yōu)

Panel 分解浮點(diǎn)計(jì)算集中在BLAS 函數(shù),使用針對(duì)目標(biāo)系統(tǒng)優(yōu)化的Hygon OpenBLIS庫(kù).Panel 分解使用遞歸算法,中間遞歸層次的浮點(diǎn)運(yùn)算集中在BLAS 的DTRSM 和DGEMM 兩個(gè)函數(shù).當(dāng)遞歸層次包含的列數(shù)小于等于閾值時(shí),使用非遞歸算法,浮點(diǎn)運(yùn)算集中在BLAS 的DGEMV、DTRSV、DSCAL、IDAMAX 等函數(shù).對(duì)于每一列的LU 分解,需要通過(guò)交換和廣播通信主元所在的panel 行,并記錄主元行交換信息.交換和廣播通信使用binary-exchange算法.

Panel 分解計(jì)算有關(guān)參數(shù)有NB、NBMIN、PFACT、RFACT、DIV等.NB取決于GPU 執(zhí)行矩陣更新的效率,同時(shí)考慮CPU 與GPU 計(jì)算的平衡.當(dāng)節(jié)點(diǎn)規(guī)模較小時(shí),NB＝384;當(dāng)節(jié)點(diǎn)規(guī)模較大時(shí),NB＝256.通過(guò)參數(shù)調(diào)優(yōu),選擇優(yōu)化的參數(shù)組合,見(jiàn)表1.

Table 1 Parameters of panel factorization表1 Panel 分解參數(shù)

進(jìn)一步分析panel 分解各部分時(shí)間,0 號(hào)進(jìn)程主要計(jì)算函數(shù)時(shí)間見(jiàn)表2.

Table 2 Time of BLAS functions in panel factorization表2 Panel 分解BLAS 函數(shù)時(shí)間

5.1.2 GPU 加速panel 分解DGEMM

從第5.1.1 節(jié)可以看出,panel 分解中,DGEMM 時(shí)間占有最大比例,需要進(jìn)一步優(yōu)化.DGEMM 在panel 分解遞歸層次調(diào)用.根據(jù)panel 分解的left-looking算法,從左到右執(zhí)行LU 分解,左側(cè)的subpanel 完成分解之后,執(zhí)行DGEMM,更新相同層次的右側(cè)的subpanel,GPU 可以加速這部分DGEMM.

GPU 加速DGEMM 的算法過(guò)程,首先把左側(cè)subpanel 和相應(yīng)的U數(shù)據(jù)傳輸?shù)紾PU 設(shè)備內(nèi)存,然后執(zhí)行DGEMM 更新GPU 內(nèi)存的右側(cè)subpanel,更新后數(shù)據(jù)傳輸?shù)紺PU端內(nèi)存,繼續(xù)執(zhí)行后續(xù)的panel 分解.算法增加了CPU 和GPU 之間的數(shù)據(jù)傳輸,但由于GPU 執(zhí)行DGEMM 的速度遠(yuǎn)高于CPU,在一定規(guī)模DGEMM 的情況下,這一過(guò)程總的執(zhí)行時(shí)間少于CPU端執(zhí)行DGEMM 時(shí)間,整體性能提高.

并不是所有的panel 分解DGEMM 都可以使用GPU 加速,如果傳輸?shù)臅r(shí)間開(kāi)銷大于DGEMM 加速的效果,就不能采用這種方式.實(shí)際計(jì)算中,使用參數(shù)控制采用GPU 加速DGEMM 的閾值,通過(guò)調(diào)優(yōu)獲得最佳的性能.

5.1.3 Panel 廣播優(yōu)化

5.1.3.1 避免數(shù)據(jù)封裝

基礎(chǔ)代碼中,panel 分解數(shù)據(jù)使用CPU端矩陣A的內(nèi)存,A是列存儲(chǔ),panel 分解的列是分段連續(xù)的.使用MPI傳輸panel 數(shù)據(jù)之前,需要執(zhí)行數(shù)據(jù)封裝,把數(shù)據(jù)復(fù)制到發(fā)送緩沖區(qū),而復(fù)制操作會(huì)帶來(lái)一定的時(shí)間開(kāi)銷.

為此提出一種避免數(shù)據(jù)封裝的方法.待分解的panel 數(shù)據(jù)是從GPU 內(nèi)存復(fù)制的,通過(guò)使用二維復(fù)制接口把panel 數(shù)據(jù)復(fù)制到連續(xù)存儲(chǔ)區(qū)域.Panel 分解計(jì)算完成后,MPI 接口直接使用緩沖區(qū)數(shù)據(jù),避免了數(shù)據(jù)封裝.

5.1.3.2 Panel 廣播流水

基礎(chǔ)代碼中,panel 分解完成后才會(huì)執(zhí)行panel 廣播,計(jì)算與數(shù)據(jù)傳輸串行,并且一次性數(shù)據(jù)傳輸較大,傳輸時(shí)間較長(zhǎng).通過(guò)觀察可以看出,left-looking算法中,左側(cè)subpanel計(jì)算完成后,主體部分不再變化,后續(xù)只有發(fā)生行交換的數(shù)據(jù)變化.

因此,提出一種廣播流水算法,對(duì)已分解的subpanel 數(shù)據(jù)提前發(fā)起廣播,這種情況下,panel 數(shù)據(jù)廣播與后續(xù)panel 分解計(jì)算并行.采用廣播流水算法,panel 計(jì)算完成后,只需要傳輸最后的subpanel 和行交換發(fā)生變化的數(shù)據(jù),相比于傳輸整個(gè)panel 數(shù)據(jù),縮短了傳輸時(shí)間.

廣播流水算法與GPU 加速panel 分解DGEMM 協(xié)同使用,subpanel 的數(shù)據(jù)廣播與panel 加速DGEMM 過(guò)程中的CPU 與GPU 之間數(shù)據(jù)傳輸、GPU 執(zhí)行DGEMM 計(jì)算并行,充分利用系統(tǒng)的CPU、GPU、PCIe 和網(wǎng)絡(luò)接口資源.

HPL 包含的6 種廣播算法都可以使用廣播流水優(yōu)化.

第5.1.2 節(jié)描述的GPU 加速panel 分解DGEMM 和panel 廣播流程緊密耦合,表3 給出了優(yōu)化前后偽代碼對(duì)比.

Table 3 Pseudocode of accelerated panel DGEMM and broadcast pipeline表3 Panel DGEMM 加速和廣播流水偽代碼

5.2 行交換long算法優(yōu)化

矩陣更新行交換有兩種算法可選,bianary-exchange算法和long算法,經(jīng)過(guò)分析和評(píng)測(cè)選擇long算法[10].第4.2 節(jié)從行交換與GPU 并行計(jì)算的角度針對(duì)行交換流程做了優(yōu)化,本節(jié)優(yōu)化具體的行交換算法.

Long算法包括兩個(gè)步驟,spread 和roll.基礎(chǔ)代碼中,首先是spread 過(guò)程,當(dāng)前行進(jìn)程將需要換出的行分發(fā)給其他所有進(jìn)程.然后,非當(dāng)前行進(jìn)程把接收的數(shù)據(jù)與本地需要換出的數(shù)據(jù)交換.接著是roll 過(guò)程,總共需要進(jìn)行P–1 次傳輸.每次傳輸所有的進(jìn)程都要參與,每個(gè)進(jìn)程只與自己的鄰居交換.第1 步,先把自己所擁有的一部分與左(右)鄰居交換,然后,把自己剛剛得到的數(shù)據(jù)與右(左)鄰居交換.這樣左右交替,直到所有的數(shù)據(jù)都交換完畢,每個(gè)進(jìn)程都擁有了全部的U.最后所有進(jìn)程還需要把所有的數(shù)據(jù)進(jìn)行重排,放置到矩陣A相應(yīng)的位置.

基礎(chǔ)代碼中,spread 和roll 使用同一塊數(shù)據(jù)緩沖區(qū),為了避免數(shù)據(jù)沖突,spread 和roll 是嚴(yán)格串行的.通過(guò)觀察發(fā)現(xiàn),spread 的接收數(shù)據(jù)與roll 的發(fā)送數(shù)據(jù)并沒(méi)有依賴關(guān)系,spread 和roll算法上是可并行性.基于以上觀察與分析,將spread 的接收緩沖與roll 的發(fā)送緩沖分離.這樣,非當(dāng)前進(jìn)程行在等待接收spread 數(shù)據(jù)時(shí),就可以執(zhí)行roll數(shù)據(jù)封裝,并把數(shù)據(jù)傳輸?shù)桨l(fā)送緩沖區(qū),這部分時(shí)間就可以被隱藏.

在使用上述的緩沖區(qū)分離算法之后,spread 和roll 已經(jīng)沒(méi)有依賴關(guān)系,因此提出了一種新的行交換方法,即交換spread 和roll 的執(zhí)行順序,如圖5(b)所示.對(duì)當(dāng)前行進(jìn)程來(lái)說(shuō),在roll 操作之前,只需將本地需要交換的數(shù)據(jù)拷貝到發(fā)送緩沖區(qū)即可,也就是減少了行交換的啟動(dòng)時(shí)間.在roll 執(zhí)行網(wǎng)絡(luò)傳輸?shù)耐瑫r(shí),當(dāng)前行進(jìn)程可以將spread 所需的數(shù)據(jù)封裝,并異步傳輸?shù)桨l(fā)送緩沖區(qū).roll 執(zhí)行完成后,開(kāi)啟執(zhí)行spread.同時(shí),把roll 過(guò)程接收的數(shù)據(jù)異步傳輸?shù)紾PU,并執(zhí)行數(shù)據(jù)交換.在spread 過(guò)程結(jié)束后,非當(dāng)前行進(jìn)程再將換入的數(shù)據(jù)傳輸?shù)紾PU,并交換到相應(yīng)位置,當(dāng)前行進(jìn)程對(duì)需要進(jìn)行內(nèi)部交換的數(shù)據(jù)執(zhí)行GPU 上的本地交換.這種算法可以實(shí)現(xiàn)網(wǎng)絡(luò)通信傳輸、GPU 數(shù)據(jù)封裝和數(shù)據(jù)交換計(jì)算以及CPU 與GPU 數(shù)據(jù)傳輸相互重疊,減少行交換執(zhí)行時(shí)間.

Fig.5 Original row swap algorithm (a),and modified row swap procedure (b)圖5 原始行交換算法(a)及新的行交換算法(b)

6 實(shí)驗(yàn)與分析

6.1 單節(jié)點(diǎn)性能分析

首先針對(duì)單節(jié)點(diǎn)實(shí)驗(yàn),分析評(píng)估各種優(yōu)化方法的效果.根據(jù)平衡點(diǎn)理論,優(yōu)化的目標(biāo)是GPU 盡可能地處于高效工作狀態(tài),執(zhí)行有效的計(jì)算,也即GPU 盡可能地連續(xù)執(zhí)行矩陣更新.圖6 對(duì)比了優(yōu)化后單步循環(huán)的時(shí)間和GPU 矩陣更新時(shí)間.圖6 給出了1 個(gè)節(jié)點(diǎn)內(nèi)0 號(hào)進(jìn)程的時(shí)間曲線.”A”表示當(dāng)前行-當(dāng)前列,”B”表示非當(dāng)前行-非當(dāng)前列.可以看出,N=88320 情況下,LU 分解矩陣規(guī)模56.09%之前單步循環(huán)時(shí)間接近GPU 矩陣更新時(shí)間.56.09%之后,當(dāng)前進(jìn)程執(zhí)行panel 分解,單步循環(huán)時(shí)間仍然接近GPU 矩陣更新時(shí)間,GPU 大部分時(shí)間處于工作狀態(tài),說(shuō)明算法優(yōu)化已經(jīng)達(dá)到較好的效果.對(duì)于非當(dāng)前panel 分解進(jìn)程,增加了等待panel 廣播數(shù)據(jù)的時(shí)間開(kāi)銷,GPU 有一段時(shí)間處于空閑狀態(tài).

Fig.6 GPU time vs.total time of each loop圖6 GPU 時(shí)間與單步循環(huán)總時(shí)間

可以看出,圖6 優(yōu)化版本的平衡點(diǎn)是56.09%,圖4 未優(yōu)化版本版本是37.82%.從LU 分解的計(jì)算量來(lái)看,前56.09%占總的LU 分解計(jì)算量的91.53%,大部分時(shí)間充分發(fā)揮了GPU 的計(jì)算能力.優(yōu)化后單節(jié)點(diǎn)HPL 效率達(dá)到了79.51%.

和優(yōu)化有關(guān)的參數(shù)主要是GPU 加速DGEMM 和panel 廣播流水的層次,這里都是3.行交換非當(dāng)前列第1段列數(shù)是待更新列數(shù)的1/8,增長(zhǎng)因子是2.Panel 分解參數(shù)和表1 一致.

6.2 可擴(kuò)展分析

多節(jié)點(diǎn)HPL 較大規(guī)?？蓴U(kuò)展,目前最大規(guī)模4 096 個(gè)節(jié)點(diǎn).HPL 較大規(guī)模并行有關(guān)的參數(shù)和效率見(jiàn)表4,表中列出了針對(duì)相應(yīng)節(jié)點(diǎn)數(shù)執(zhí)行的實(shí)驗(yàn)中效率最高的配置.

Table 4 Parameters and efficiency of multi-nodes HPL表4 多節(jié)點(diǎn)HPL 參數(shù)和效率

參數(shù)N是較大規(guī)模節(jié)點(diǎn)HPL 評(píng)測(cè)中比較重要的參數(shù),每個(gè)節(jié)點(diǎn)計(jì)算的矩陣列數(shù)一般是NB的倍數(shù),多節(jié)點(diǎn)的N正比于節(jié)點(diǎn)數(shù)平方根.N值越大,計(jì)算效率較高的矩陣更新占比越大,一般來(lái)說(shuō)HPL 效率越高,但N受限于GPU 的內(nèi)存容量,同時(shí)還要考慮系統(tǒng)的負(fù)載.NB使用384 或256.P和Q確定二維進(jìn)程網(wǎng)格的組織方式,經(jīng)過(guò)驗(yàn)證,P=Q時(shí)目標(biāo)系統(tǒng)性能最優(yōu).Panel 廣播和矩陣行交換分別使用優(yōu)化的long算法和廣播流水算法.

從HPL 效率來(lái)看,4 節(jié)點(diǎn)效率高于單節(jié)點(diǎn)效率,原因在于4 節(jié)點(diǎn)矩陣規(guī)模大于單節(jié)點(diǎn),效率較高的計(jì)算占比較大,并且可以抵消4 節(jié)點(diǎn)規(guī)模下網(wǎng)絡(luò)傳輸路徑增加帶來(lái)的性能損失.隨著規(guī)模的擴(kuò)大,網(wǎng)絡(luò)傳輸對(duì)性能的影響增大,HPL 效率逐步降低,但降低的趨勢(shì)比較平緩.

7 結(jié)論與未來(lái)的工作

針對(duì)復(fù)雜異構(gòu)系統(tǒng),分析HPL算法的特點(diǎn),提出CPU 與GPU 協(xié)同計(jì)算方法,提出平衡點(diǎn)理論指導(dǎo)各方面優(yōu)化.實(shí)現(xiàn)CPU 與GPU 協(xié)同計(jì)算的look-ahead、行交換連續(xù)流水線算法隱藏panel 分解和行交換時(shí)間,并優(yōu)化panel 分解和行交換算法,延長(zhǎng)GPU 計(jì)算占優(yōu)的時(shí)間,最終提高整個(gè)系統(tǒng)的HPL 效率,單節(jié)點(diǎn)效率79.51%.

當(dāng)前的超級(jí)計(jì)算機(jī)正向百億億級(jí)邁進(jìn),目標(biāo)系統(tǒng)采用的節(jié)點(diǎn)內(nèi)通用CPU 裝備加速器,配合高速的節(jié)點(diǎn)內(nèi)總線被認(rèn)為是較有潛力的一種架構(gòu).未來(lái)工作主要是針對(duì)百億億級(jí)計(jì)算架構(gòu),開(kāi)展加速器DGEMM 優(yōu)化、CPU 與加速器協(xié)同計(jì)算和混合精度等相關(guān)研究.同時(shí),百億億級(jí)超級(jí)計(jì)算機(jī)節(jié)點(diǎn)規(guī)模繼續(xù)增加,節(jié)點(diǎn)數(shù)量有可能達(dá)到數(shù)萬(wàn)甚至更多,節(jié)點(diǎn)間的互連也更加重要.需要進(jìn)一步研究HPL 大規(guī)模可擴(kuò)展能力,評(píng)測(cè)行交換算法和廣播算法的性能,并進(jìn)行相應(yīng)的優(yōu)化.