多核數(shù)字信號(hào)處理器并行矩陣轉(zhuǎn)置算法優(yōu)化*

裴向東，王慶林,2，廖林玉,2，李榮春,2，梅松竹,2，劉 杰,2，龐征斌

(1. 國(guó)防科技大學(xué) 計(jì)算機(jī)學(xué)院, 湖南 長(zhǎng)沙 410073； 2. 國(guó)防科技大學(xué) 并行與分布處理國(guó)防科技重點(diǎn)實(shí)驗(yàn)室， 湖南 長(zhǎng)沙 410073)

矩陣轉(zhuǎn)置是矩陣運(yùn)算中最常見(jiàn)的一種操作，廣泛應(yīng)用于科學(xué)和工程領(lǐng)域[1-2]，如信號(hào)處理、科學(xué)計(jì)算和深度學(xué)習(xí)[3-5]。作為重要的基礎(chǔ)算子[6-7]，矩陣轉(zhuǎn)置的效率高低對(duì)應(yīng)用的性能會(huì)產(chǎn)生直接的影響。特別是矩陣轉(zhuǎn)置操作屬于訪存受限型運(yùn)算，對(duì)訪存密集型應(yīng)用的影響將更大[8-9]。

面對(duì)各領(lǐng)域?qū)Ω咝阅芫仃囖D(zhuǎn)置操作的需求，國(guó)內(nèi)外學(xué)者對(duì)矩陣轉(zhuǎn)置優(yōu)化已經(jīng)進(jìn)行了大量的研究工作。遠(yuǎn)遠(yuǎn)等提出適應(yīng)對(duì)稱多處理器系統(tǒng)的兩維均衡細(xì)粒度交織矩陣轉(zhuǎn)置算法[10]，通過(guò)提升矩陣在存儲(chǔ)器上的讀寫(xiě)效率來(lái)降低矩陣轉(zhuǎn)置開(kāi)銷。Zekri 等提出了面向Intel 處理器的矩陣轉(zhuǎn)置優(yōu)化算法[11]，通過(guò)采用向量擴(kuò)展指令A(yù)VX來(lái)加速矩陣轉(zhuǎn)置運(yùn)算。王琦等提出了面向Intel KNL融合處理器的并行矩陣轉(zhuǎn)置優(yōu)化[12]。Springer等提出了面向張量轉(zhuǎn)置的編譯器TTC[13]，以離線方式自動(dòng)產(chǎn)生高性能的C++轉(zhuǎn)置實(shí)現(xiàn)。TTC效率雖高，但不能直接應(yīng)用在運(yùn)行時(shí)才能確定張量大小與轉(zhuǎn)置需求的場(chǎng)景中，因此Springer等隨后又提出了支持在線自動(dòng)調(diào)優(yōu)的張量轉(zhuǎn)置庫(kù)HPTT[14]，集成了自動(dòng)調(diào)優(yōu)、顯式向量化和多線程并行等優(yōu)化技術(shù)，在多種中央處理器(central processing unit,CPU)硬件架構(gòu)上均展現(xiàn)了優(yōu)異的性能，是當(dāng)前面向CPU架構(gòu)的流行張量轉(zhuǎn)置開(kāi)源庫(kù)。肖漢等采用OpenCL編程模型實(shí)現(xiàn)了面向圖形處理器(graphic processing unit,GPU)的并行矩陣轉(zhuǎn)置算法優(yōu)化，相比CUDA實(shí)現(xiàn)，具有更好的可移植性[15]。Hynninen等面向NVIDIA GPU構(gòu)建了張量轉(zhuǎn)置函數(shù)庫(kù)cuTT[16]。Vedurada等為GPU提出張量轉(zhuǎn)置庫(kù)TTLG[17]，構(gòu)建了性能預(yù)測(cè)模型來(lái)選擇不同轉(zhuǎn)置內(nèi)核以及分塊大小。高捷針對(duì)“神威”國(guó)產(chǎn)超算平臺(tái)優(yōu)化了張量轉(zhuǎn)置庫(kù) SWTT[18]。

由于能源效率和功耗的限制，高性能計(jì)算的異構(gòu)計(jì)算領(lǐng)域引入了低功耗嵌入式架構(gòu)，如數(shù)字信號(hào)處理器[19](digital signal processor,DSP)。與CPU和GPU相比，DSP通常采用基于超長(zhǎng)指令字(very long instruction word,VLIW)和順序執(zhí)行的向量核，集成了直接存儲(chǔ)器存取(direct memory access,DMA)引擎用于數(shù)據(jù)訪問(wèn)[20-21]，從而已有的針對(duì)CPU和GPU的并行矩陣轉(zhuǎn)置實(shí)現(xiàn)并不能直接適用于多核DSP架構(gòu)。

FT-M7032是國(guó)防科技大學(xué)自主研發(fā)的一款異構(gòu)多核DSP[19]，主頻為1.8 GHz時(shí)雙精度浮點(diǎn)峰值性能高達(dá)5.53 Tflops/s，在信號(hào)處理、科學(xué)計(jì)算和深度學(xué)習(xí)等領(lǐng)域有巨大潛力[22]，對(duì)高性能矩陣轉(zhuǎn)置實(shí)現(xiàn)提出了強(qiáng)烈需求。本文面向FT-M7032芯片的體系結(jié)構(gòu)特征，研究高性能矩陣轉(zhuǎn)置并行實(shí)現(xiàn)特性，提出了一種適配不同數(shù)據(jù)位寬(8 B、4 B以及2 B)矩陣的并行矩陣轉(zhuǎn)置算法ftmMT，基于DSP核中向量Load/Store單元支持的取模存儲(chǔ)功能構(gòu)建了向量化矩陣轉(zhuǎn)置核心，通過(guò)矩陣分塊、隱式乒乓以及多核并行的設(shè)計(jì)實(shí)現(xiàn)片上向量化轉(zhuǎn)置與片外訪存的重疊以及對(duì)片外存儲(chǔ)帶寬的高效利用。實(shí)驗(yàn)結(jié)果顯示，ftmMT獲得了高達(dá)75.95%的DDR有效帶寬利用率；與主流矩陣轉(zhuǎn)置開(kāi)源算法庫(kù)HPTT相比，ftmMT實(shí)現(xiàn)了高達(dá)8.99倍的性能加速。因此，本文工作對(duì)于推動(dòng)FT-M7032在各領(lǐng)域的應(yīng)用具有重要的意義。

1 矩陣轉(zhuǎn)置定義與FT-M7032結(jié)構(gòu)簡(jiǎn)介

本文的研究范圍僅限于二維矩陣，故令輸入矩陣為A[M][N]，轉(zhuǎn)置輸出矩陣為B=AT，其中M和N表示矩陣的行數(shù)和列數(shù)。矩陣轉(zhuǎn)置的具體計(jì)算如式(1)所示：

B(m,n)=A(n,m)

(1)

其中，0≤m

FT-M7032是由一個(gè)多核CPU和4個(gè)GPDSP簇構(gòu)成，其體系結(jié)構(gòu)如圖1所示。多核CPU是一個(gè)簡(jiǎn)化版本的FT-2000plus處理器[4-5]，由16個(gè)兼容ARMv8架構(gòu)的CPU核構(gòu)成，主要負(fù)責(zé)進(jìn)程管理和通信，上面運(yùn)行操作系統(tǒng)。4個(gè)GPDSP簇負(fù)責(zé)提供最主要的計(jì)算性能，每個(gè)GPDSP簇由8個(gè)DSP核與全局共享內(nèi)存(global shared memory,GSM)通過(guò)片上互聯(lián)網(wǎng)絡(luò)連接而成。GSM大小為6 MB, 由單簇內(nèi)的8個(gè)DSP核共享，可以作為數(shù)據(jù)或指令的片上存儲(chǔ)空間使用。同時(shí)，CPU與GPDSP共享內(nèi)存空間。每個(gè)GPDSP簇所能訪問(wèn)的DDR理論帶寬為42.62 GB/s。每個(gè)GPDSP簇中的DSP核基于VLIW架構(gòu)實(shí)現(xiàn)，包含指令調(diào)度單元(instruction fetch unit,IFU)、標(biāo)量處理單元(scalar processing unit, SPU)、向量處理單元(vector processing unit,VPU)和DMA引擎等部分，F(xiàn)T-M7032中DSP核的微體系結(jié)構(gòu)如圖2所示。

圖1 FT-M7032體系結(jié)構(gòu)Fig.1 Architecture of FT-M7032

圖2 FT-M7032中DSP核的微體系結(jié)構(gòu)Fig.2 Micro-architecture of DSP cores in FT-M7032

SPU用于指令流控制和標(biāo)量計(jì)算，主要由標(biāo)量處理部件(scalar processing elements, SPE)和64 KB標(biāo)量存儲(chǔ)(scalar memory, SM)空間組成，每周期最多能夠處理5條標(biāo)量指令。VPU主要由768 KB的向量存儲(chǔ)(vector memory,AM)空間和16個(gè)向量處理部件(vector processing elements,VPEs)構(gòu)成，為每個(gè)DSP核提供最主要的計(jì)算能力。每個(gè)VPE包含64個(gè)寄存器、3個(gè)浮點(diǎn)乘累加(floating point multiply accumulator, FMAC)和1個(gè)位操作(bit processing, BP)單元，一次可以處理1個(gè)8 B數(shù)據(jù)(如FP64、Int64)、2個(gè)4 B數(shù)據(jù)(如FP32、Int32)，或4個(gè)2 B數(shù)據(jù)(如FP16、Int16)。16個(gè)VPE以單指令多數(shù)據(jù)(single instruction multiple data,SIMD)方式運(yùn)行，因而8 B數(shù)據(jù)類型的 SIMD 寬度為 16。AM 空間可以通過(guò)VPU中的兩個(gè)向量Load/Store單元在每個(gè)周期向寄存器提供多達(dá) 512 B的數(shù)據(jù)，從而AM與向量寄存器之間的訪問(wèn)帶寬高達(dá)921.6 GB/s。

由于功耗的限制，VPU的多個(gè)VPEs之間不支持混洗功能，同時(shí)DMA部件也不支持矩陣轉(zhuǎn)置功能。盡管如此，VPU的Load/Store單元支持取模存儲(chǔ)VSTDW0M16/VSTDW1M16指令，如圖3所示。執(zhí)行完該指令后，可以將4個(gè)向量寄存器中的數(shù)據(jù)按地址模16的方式存儲(chǔ)在AM空間中，即實(shí)現(xiàn)4×16矩陣的轉(zhuǎn)置。

圖3 VPU中模16的存儲(chǔ)功能Fig.3 Store based on Module 16 in VPU

2 面向FT-M7032的矩陣轉(zhuǎn)置實(shí)現(xiàn)分析

在FT-M7032上實(shí)現(xiàn)矩陣轉(zhuǎn)置主要有四種方法：一是基于FT-M7032 上多核CPU的實(shí)現(xiàn)，二是基于DSP 純DMA操作的實(shí)現(xiàn)，三是基于DSP SPU的實(shí)現(xiàn)，四是基于DSP VPU的實(shí)現(xiàn)。

基于FT-M7032上多核CPU的實(shí)現(xiàn)方法是直接調(diào)用已有面向CPU的開(kāi)源矩陣轉(zhuǎn)置庫(kù)，如HPTT[14]等。然而，由上一節(jié)內(nèi)容可知，F(xiàn)T-M7032的主要計(jì)算能力是由GPDSP簇提供，雖然多核CPU與GPDSP簇共享存儲(chǔ)空間，但GPDSP簇與多核CPU之間進(jìn)行交互時(shí)需要將多核CPU中Cache的內(nèi)容進(jìn)行作廢或者寫(xiě)回DDR中，當(dāng)多次交互時(shí)，交互開(kāi)銷較大。

基于DSP 純DMA操作的實(shí)現(xiàn)方法是通過(guò)將DMA操作的塊大小設(shè)置為矩陣中單個(gè)數(shù)據(jù)的大小來(lái)實(shí)現(xiàn)轉(zhuǎn)置。然而，F(xiàn)T-M7032中DMA操作最小粒度為8 B，這意味著該方法只適合于單個(gè)數(shù)據(jù)為8 B的類型，如FP64、Int64等。同時(shí)，DMA的塊大小直接影響DMA的訪存效率，這意味著使用該方法不能充分發(fā)揮FT-M7032的訪存帶寬。

基于DSP SPU的實(shí)現(xiàn)方法是先將數(shù)據(jù)加載到SPU的SM空間，然后基于標(biāo)量操作來(lái)實(shí)現(xiàn)矩陣轉(zhuǎn)置運(yùn)算。在FT-M7032的SPU中，SM僅為64 KB，使得在轉(zhuǎn)置實(shí)現(xiàn)過(guò)程中需要頻繁調(diào)用DMA，DMA操作的啟動(dòng)開(kāi)銷較大。與此同時(shí)，SPU僅有一個(gè)標(biāo)量Load/Store單元，每周期僅能向寄存器提供16 B的數(shù)據(jù)，使得基于SPU的片上轉(zhuǎn)置核心在實(shí)現(xiàn)過(guò)程中性能受限。

基于DSP VPU的實(shí)現(xiàn)方法則是將數(shù)據(jù)加載到VPU的AM空間，然后基于向量操作來(lái)實(shí)現(xiàn)矩陣轉(zhuǎn)置運(yùn)算。在FT-M7032的VPU中，AM空間大小為768 KB，VPU的兩個(gè)Load/Store單元每周期可以提供512 B的數(shù)據(jù)。相比前面兩種基于DSP 的實(shí)現(xiàn)方法，可一定程度上有效克服片上存儲(chǔ)帶寬受限、頻繁調(diào)用DMA以及訪存效率較低等問(wèn)題。

綜上所述，基于DSP VPU的方法是面向FT-M7032的高性能矩陣轉(zhuǎn)置實(shí)現(xiàn)的必然選擇。但是VPU的多個(gè)VPE之間不支持混洗功能，使得基于DSP VPU的轉(zhuǎn)置實(shí)現(xiàn)成為一個(gè)難題。本文將通過(guò)充分利用VPU中Load/Store單元支持的取模存儲(chǔ)功能來(lái)解決這一問(wèn)題?？紤]到基于純DMA操作和SPU的兩種實(shí)現(xiàn)在理論上性能就遠(yuǎn)低于本文基于VPU的實(shí)現(xiàn)，故本文在后續(xù)實(shí)驗(yàn)部分僅與基于多核CPU的實(shí)現(xiàn)進(jìn)行性能比較。

3 并行矩陣轉(zhuǎn)置算法與優(yōu)化

本節(jié)先對(duì)本文提出的面向飛騰異構(gòu)多核DSP的并行矩陣轉(zhuǎn)置算法ftmMT進(jìn)行整體介紹，然后對(duì)該算法的具體思路及實(shí)現(xiàn)過(guò)程做詳細(xì)描述。

3.1 ftmMT算法整體設(shè)計(jì)

如第2節(jié)所分析，基于DSP VPU的轉(zhuǎn)置實(shí)現(xiàn)是面向FT-M7032的高性能矩陣轉(zhuǎn)置實(shí)現(xiàn)的必然選擇?；诖耍疚奶岢隽瞬⑿芯仃囅蛄炕D(zhuǎn)置算法ftmMT，如算法1所示。ftmMT算法充分結(jié)合了FT-M7032的體系結(jié)構(gòu)特征和矩陣轉(zhuǎn)置的計(jì)算訪存特點(diǎn)，主要由三個(gè)步驟構(gòu)成。

算法1 并行矩陣轉(zhuǎn)置算法ftmMTAlg.1 Parallel matrix transpose algorithm (ftmMT)

步驟1：將CPU Cache中的內(nèi)容寫(xiě)回DDR中(Step 1)。

步驟2：調(diào)用DSP端函數(shù)__MTrn()執(zhí)行矩陣轉(zhuǎn)置運(yùn)算(Step 2)。

步驟3：作廢CPU Cache中的內(nèi)容，完成矩陣轉(zhuǎn)置運(yùn)算(Step 3)。

DSP端函數(shù)__MTrn()負(fù)責(zé)在一個(gè)GPDSP簇上完成矩陣的轉(zhuǎn)置操作，其基本思想是將A矩陣劃分為Mb×Nb大小的子矩陣，依次加載到DSP的AM空間并擴(kuò)展成Ma×Na大小(第12行)，調(diào)用trnKernel()函數(shù)完成Ma×Na子矩陣的轉(zhuǎn)置運(yùn)算(第13行)，然后從AM空間讀取Nb×Mb子矩陣存儲(chǔ)回B矩陣中(第14行)。trnKernel()函數(shù)負(fù)責(zé)將AM空間中Ma×Na的矩陣進(jìn)一步劃分成ha×wa大小的子矩陣，然后依次調(diào)用向量化矩陣轉(zhuǎn)置核心Vtkernel()基于VPU完成ha×wa子矩陣的轉(zhuǎn)置。

對(duì)于DSP端的應(yīng)用，可以直接調(diào)用DSP端函數(shù)__MTrn()，避免了CPU與DSP之間的反復(fù)交互開(kāi)銷。

總的來(lái)說(shuō)，本文提出的并行矩陣轉(zhuǎn)置算法ftmMT具有以下三個(gè)特點(diǎn)：

1)不依賴任何CPU開(kāi)源轉(zhuǎn)置函數(shù)庫(kù)，也不依賴其他DSP矩陣轉(zhuǎn)置函數(shù)庫(kù)；

2)基于AM空間實(shí)現(xiàn)，增加了DMA每次傳輸?shù)臄?shù)據(jù)量，從而減少了DMA操作的總次數(shù)，大幅降低了DMA啟動(dòng)開(kāi)銷；

3)基于DSP VPU單元實(shí)現(xiàn)，能夠有效利用VPU單元的Load/Store帶寬。

3.2 向量化矩陣轉(zhuǎn)置核心

Vtkernel()是DSP上基于VPU來(lái)實(shí)現(xiàn)矩陣轉(zhuǎn)置的關(guān)鍵核心函數(shù)。為充分利用FT-M7032的向量部件來(lái)提高矩陣轉(zhuǎn)置的性能，以單個(gè)數(shù)據(jù)大小為8 B、ha和wa均為16為例，本文提出了如圖4所示的Vtkernel()的實(shí)現(xiàn)過(guò)程，一共分為6個(gè)步驟。

圖4 基于DSP VPU的向量化矩陣轉(zhuǎn)置核心流程Fig.4 Flow chart of vectorized matrix transpose kernel based on VPU of DSP

步驟1：基于VPU的向量Load功能將16×16大小子矩陣中的16行分別讀入16個(gè)對(duì)應(yīng)的向量寄存器中。比如，第0～3行，分別讀入向量寄存器VR0/4/8/12中；第4～7行，分別讀入向量寄存器VR1/5/9/13中等。

步驟2：對(duì)向量寄存器進(jìn)行分組，將向量寄存器分為VR0～VR3、VR4～VR7、VR8～VR11和VR12～VR15四組，然后將分組后的寄存器組分別調(diào)用模16存儲(chǔ)指令VSTDW0M16/VSTDW1M16，第一次將數(shù)據(jù)存入AM臨時(shí)空間，分別完成4行的轉(zhuǎn)置。比如，基于VR0～VR3，調(diào)用模16存儲(chǔ)指令VSTDW0M16/VSTDW1M16完成了原16×16子矩陣中第0、4、8、12行組成的新子矩陣的轉(zhuǎn)置。

步驟3：基于VPU的向量Load功能將步驟2產(chǎn)生的AM臨時(shí)空間中16×16大小臨時(shí)矩陣重新讀入16個(gè)對(duì)應(yīng)的向量寄存器中，讀入規(guī)則與步驟1完全一致。

步驟4：再次將向量寄存器進(jìn)行分組，并再次基于模16存儲(chǔ)指令VSTDW0M16/VSTDW1M16將數(shù)據(jù)存入AM臨時(shí)空間，分組規(guī)則與調(diào)用模16存儲(chǔ)指令的規(guī)則與步驟2完全一致，從而完成了原16×16大小子矩陣的轉(zhuǎn)置，轉(zhuǎn)置結(jié)果存儲(chǔ)在AM臨時(shí)空間。

步驟5：從AM臨時(shí)空間中讀入轉(zhuǎn)置后的臨時(shí)矩陣到向量寄存器中。

步驟6：將向量寄存器中保存的轉(zhuǎn)置后子矩陣存入Ma×Na矩陣中對(duì)應(yīng)位置，完成ha×wa大小子矩陣的轉(zhuǎn)置。

對(duì)于單個(gè)數(shù)據(jù)大小分別為4 B和2 B的矩陣來(lái)說(shuō)，處理流程大致相似，只需在步驟1與步驟2之間加入數(shù)據(jù)打包過(guò)程，即：根據(jù)最后的轉(zhuǎn)置需求，將原始子矩陣中相鄰兩行的2個(gè) 4 B數(shù)據(jù)和相鄰四行的4個(gè)2 B數(shù)據(jù)分別合并為1個(gè)8 B數(shù)據(jù)，分別產(chǎn)生新的2行和4行8 B數(shù)據(jù)存儲(chǔ)于向量寄存器中。最后，基于所產(chǎn)生的新寄存器數(shù)據(jù)執(zhí)行步驟2～6，構(gòu)建轉(zhuǎn)置后的矩陣。以單個(gè)數(shù)據(jù)大小為4 B為例，設(shè)置ha和wa為32，每?jī)尚袛?shù)據(jù)分別進(jìn)行高低位打包操作形成新的32行數(shù)據(jù)，如圖5所示。隨后將新32行數(shù)據(jù)以行號(hào)為基準(zhǔn)進(jìn)行奇偶分組，即奇數(shù)編號(hào)16行與偶數(shù)編號(hào)16行，隨后分別執(zhí)行步驟2～6，即可完成32×32大小4 B矩陣的轉(zhuǎn)置。對(duì)于單個(gè)數(shù)據(jù)大小為2 B的矩陣轉(zhuǎn)置，設(shè)置ha和wa為64，然后分別進(jìn)行打包和轉(zhuǎn)置。

圖5 VPU中高低位數(shù)據(jù)打包功能Fig.5 Packing of data in VPU

3.3 分塊設(shè)計(jì)

本文算法在設(shè)計(jì)中主要涉及ha×wa、Mb×Nb兩組分塊參數(shù)的選擇與評(píng)估。

在矩陣轉(zhuǎn)置過(guò)程中，ha×wa表示Vtkernel一次處理的子矩陣大小，其主要涉及兩個(gè)方面因素的影響：一是對(duì)于特定元素大小，向量化轉(zhuǎn)置核心功能實(shí)現(xiàn)中所要求的最小矩陣大小。如第3.2節(jié)所示，對(duì)于8 B、4 B以及2 B的矩陣轉(zhuǎn)置來(lái)說(shuō)，ha×wa最小要求為16×16、32×32以及64×64。二是為充分利用VPU中兩個(gè)向量Load/Store單元的處理能力，需要同時(shí)執(zhí)行多個(gè)最小矩陣的轉(zhuǎn)置，盡可能填充滿流水線，實(shí)現(xiàn)VPU的最大化利用。在本文實(shí)現(xiàn)中，對(duì)于8 B、4 B以及2 B的矩陣轉(zhuǎn)置，Vtkernel分別一次完成4個(gè)、1個(gè)以及1個(gè)最小矩陣的轉(zhuǎn)置，從而ha×wa分別取為32×32、32×32以及64×64。

Mb×Nb的選擇涉及四個(gè)方面因素的影響：一是Mb和Nb應(yīng)該盡可能分別是ha和wa的整數(shù)倍。二是受到AM可用空間大小的約束，即存儲(chǔ)Aam[Ma][Na]、Bam[Na][Ma]和Vtkernel實(shí)現(xiàn)中臨時(shí)空間大小Tempam之和不能超過(guò)AM可用空間的大小。對(duì)于8 B、4 B以及2 B的Vtkernel實(shí)現(xiàn)，Tempam分別為8 KB、4 KB以及8 KB。三是DMA訪問(wèn)DDR內(nèi)存的性能受連續(xù)訪問(wèn)的塊大小與調(diào)用次數(shù)的影響較大。從DDR讀取矩陣A的子矩陣A(m,n)[Mb][Nb]過(guò)程中，連續(xù)訪問(wèn)DDR的塊大小為Nb個(gè)數(shù)據(jù)；將轉(zhuǎn)置后的子矩陣Bam[Nb][Mb]存入DDR中B的過(guò)程中，塊大小為Mb個(gè)數(shù)據(jù)。因此，在前兩個(gè)因素的約束下，Mb與Nb的取值理論上應(yīng)盡可能接近，從而才能在轉(zhuǎn)置實(shí)現(xiàn)過(guò)程中同時(shí)保證DDR讀和寫(xiě)的性能。四是根據(jù)參數(shù)的實(shí)時(shí)調(diào)整，當(dāng)輸入矩陣A的寬度N小于Nb或高度M小于Mb時(shí)，將分別動(dòng)態(tài)降低Nb或Mb的大小。同時(shí)，在AM可用空間約束下，隨之增大Mb或Nb的大小。

3.4 多核并行與乒乓設(shè)計(jì)

為充分發(fā)揮多核DSP的并行性能，本文算法設(shè)計(jì)中將實(shí)現(xiàn)面向多個(gè)DSP核的線程級(jí)并行，具體為將矩陣A劃分為多個(gè)Mb×Nb子矩陣，然后調(diào)用多個(gè)DSP核并行完成多個(gè)Mb×Nb子矩陣的轉(zhuǎn)置，直到完成矩陣A中所有Mb×Nb子矩陣的轉(zhuǎn)置。在算法實(shí)現(xiàn)上，即在M和N兩個(gè)維度進(jìn)行多核并行，如算法1中第6行和第9行所示。

乒乓算法是用來(lái)實(shí)現(xiàn)計(jì)算與數(shù)據(jù)傳輸重疊的經(jīng)典方法，分為顯式乒乓與隱式乒乓兩種。顯式乒乓是在單個(gè)DSP核內(nèi)實(shí)現(xiàn)計(jì)算與數(shù)據(jù)傳輸重疊，需要將單個(gè)DSP核內(nèi)如AM等片上空間平均劃分為兩部分，一部分用于數(shù)據(jù)傳輸(DMA操作)，一部分用于計(jì)算(trnKernel函數(shù))。隱式乒乓算法是在多個(gè)DSP核之間實(shí)現(xiàn)計(jì)算與數(shù)據(jù)傳輸?shù)闹丿B，即一個(gè)DSP核的計(jì)算與另一個(gè)DSP核的數(shù)據(jù)傳輸重疊。在FT-M7032單簇8個(gè)DSP核并行計(jì)算時(shí)，隱式乒乓算法是一種天然存在的重疊形態(tài)，不需進(jìn)行特殊設(shè)計(jì)。兩種算法相比各有優(yōu)缺點(diǎn)，如隱式乒乓算法中總存在DSP計(jì)算核在等待數(shù)據(jù)的傳輸，而顯式乒乓算法雖然可以讓單個(gè)DSP核不間斷地進(jìn)行計(jì)算，但每次DMA傳輸?shù)臄?shù)據(jù)量降為隱式算法的一半，總數(shù)據(jù)傳輸開(kāi)銷變大。鑒于矩陣轉(zhuǎn)置屬于訪存密集型操作，訪存性能或數(shù)據(jù)傳輸性能才是決定性能的關(guān)鍵。因此，ftmMT采用隱式乒乓算法來(lái)實(shí)現(xiàn)trnKernel函數(shù)與DMA數(shù)據(jù)傳輸操作的重疊，通過(guò)最大化降低數(shù)據(jù)傳輸開(kāi)銷來(lái)實(shí)現(xiàn)總時(shí)間開(kāi)銷成本的最小化。

4 性能評(píng)估

在本節(jié)中，將全面評(píng)測(cè)本文所提算法ftmMT的性能。首先，測(cè)試轉(zhuǎn)置矩陣位于片上AM空間時(shí)矩陣轉(zhuǎn)置實(shí)現(xiàn)的性能，即trnKernel實(shí)現(xiàn)的性能。其次，測(cè)試不同Mb×Nb大小對(duì)于ftmMT性能的影響。最后，測(cè)試ftmMT在不同矩陣大小下的算法性能，并與運(yùn)行在FT-M7032 多核CPU上的開(kāi)源算法HPTT進(jìn)行性能比較。

在本節(jié)中將用兩個(gè)指標(biāo)來(lái)衡量矩陣轉(zhuǎn)置的性能，一個(gè)是完成矩陣轉(zhuǎn)置的時(shí)間T，另一個(gè)是轉(zhuǎn)置帶寬F，其計(jì)算如式(2)所示。

(2)

其中，QSizeof(a)表示矩陣A或B中一個(gè)數(shù)據(jù)的字節(jié)數(shù)。本文算法主要支持8 B、4 B以及2 B大小數(shù)據(jù)的矩陣，在后續(xù)討論中對(duì)應(yīng)實(shí)現(xiàn)分別標(biāo)記為ftmMT-64/trnKernel-64、ftmMT-32/trnKernel-32以及ftmMT-16/trnKernel-16?？紤]到這些實(shí)現(xiàn)所處理的基本塊大小ha×wa不盡相同，實(shí)驗(yàn)中采用了不同的輸入矩陣規(guī)模來(lái)進(jìn)行對(duì)應(yīng)實(shí)現(xiàn)的測(cè)試。

4.1 trnKernel的性能評(píng)估

當(dāng)Ma×Na取不同大小時(shí)，處理不同數(shù)據(jù)大小的trnKernel的性能如圖6～8所示。對(duì)于所有測(cè)試的子矩陣來(lái)說(shuō)，時(shí)間開(kāi)銷均低于10 μs。隨著Ma×Na大小的增加，轉(zhuǎn)置帶寬呈逐漸增加的趨勢(shì)。對(duì)于8 B矩陣，當(dāng)Ma×Na為192×224，獲得了最高145.45 GB/s的矩陣轉(zhuǎn)置帶寬，如圖6所示；對(duì)于 4 B矩陣，當(dāng)分塊大小為256×352時(shí)，最高矩陣轉(zhuǎn)置帶寬為98.12 GB/s，如圖7所示；對(duì)于2 B矩陣，當(dāng)分塊大小為384×448時(shí)，最高矩陣轉(zhuǎn)置帶寬為102.77 GB/s，如圖8所示。

圖6 trnKernel-64實(shí)現(xiàn)的性能Fig.6 Performance of trnKernel-64 implementation

圖7 trnKernel-32實(shí)現(xiàn)的性能Fig.7 Performance of trnKernel-32 implementation

圖8 trnKernel-16實(shí)現(xiàn)的性能Fig.8 Performance of trnKernel-16 implementation

4.2 不同分塊大小下ftmMT的性能評(píng)估

如第3.3節(jié)分析所知，分塊矩陣Mb×Nb受多個(gè)因素約束，同時(shí)也對(duì)最終矩陣轉(zhuǎn)置性能有較大的影響。為評(píng)估不同分塊大小的性能，本小節(jié)測(cè)試了ftmMT實(shí)現(xiàn)在Mb×Nb不同大小下的性能，實(shí)驗(yàn)結(jié)果分別如圖9～11所示。其中，橫坐標(biāo)表示不同的分塊大小，從左到右，Mb逐漸增大，Nb逐漸減小，即分塊子矩陣從矮瘦矩陣逐漸演變?yōu)楦呤菥仃嚕玀b和Nb的乘積盡可能接近ftmMT實(shí)現(xiàn)中AM空間所允許的最大值。本節(jié)測(cè)試中，輸入矩陣A的大小M×N為固定值，對(duì)于ftmMT-64和ftmMT-32，輸入矩陣A均為16 384×16 384的方形矩陣；對(duì)于ftmMT-16，輸入矩陣A的大小為32 768×32 768。

圖9 輸入矩陣A為16 384×16 384的矩陣，不同分塊大小下ftmMT-64實(shí)現(xiàn)的性能Fig.9 Performance of ftmMT-64 under different block sizes when the input matrix A is 16 384×16 384 matrix

圖10 輸入矩陣A為16 384×16 384的矩陣，不同分塊大小下ftmMT-32實(shí)現(xiàn)的性能Fig.10 Performance of ftmMT-32 under different block sizes when the input matrix A is 16 384×16 384 matrix

圖11 輸入矩陣A為32 768×32 768的矩陣，不同分塊大小下ftmMT-16實(shí)現(xiàn)的性能Fig.11 Performance of ftmMT-16 under different block sizes when the input matrix A is 32 768×32 768 matrix

對(duì)于ftmMT-64，如圖9所示，當(dāng)輸入矩陣A為16 384×16 384的矩陣，分塊大小為192×224時(shí)，矩陣轉(zhuǎn)置性能最佳，其中耗時(shí)為152.56 ms，轉(zhuǎn)置帶寬達(dá)到了26.22 GB/s，對(duì)應(yīng)的DDR有效帶寬利用率為61.61%。

對(duì)于ftmMT-32，如圖10所示，當(dāng)輸入矩陣A為16 384×16 384的矩陣，分塊大小為256×352時(shí)，矩陣轉(zhuǎn)置性能最佳，其中耗時(shí)為72.78 ms，轉(zhuǎn)置帶寬達(dá)到了27.48 GB/s，相應(yīng)的DDR有效帶寬利用率為65.17%。

對(duì)于ftmMT-16，如圖11所示，當(dāng)輸入矩陣A為32 768×32 768的矩陣，分塊大小為448×384和576×320時(shí)，矩陣轉(zhuǎn)置性能最佳，其中耗時(shí)為155.34 ms，轉(zhuǎn)置帶寬達(dá)到了25.75 GB/s，相應(yīng)的DDR有效帶寬利用率為60.41%。

通過(guò)本小節(jié)的評(píng)估，從圖9～11呈現(xiàn)的實(shí)驗(yàn)結(jié)果可以看出，對(duì)于不同分塊大小Mb×Nb，ftmMT-64、 ftmMT-32和ftmMT-16所實(shí)現(xiàn)的性能均呈現(xiàn)中間高、兩邊低的特點(diǎn)，與第3.3小節(jié)的分析設(shè)計(jì)相符。

4.3 不同大小矩陣下ftmMT的性能評(píng)估

本小節(jié)對(duì)ftmMT處理不同大小矩陣的性能進(jìn)行評(píng)估，并與運(yùn)行在FT-M7032中多核CPU上的開(kāi)源高性能張量轉(zhuǎn)置庫(kù)HPTT進(jìn)行性能比較。本小節(jié)采用第4.2節(jié)中獲得最佳性能的Mb×Nb大小作為ftmMT實(shí)現(xiàn)中的分塊大小。同時(shí)，由于HPTT僅支持8 B和4 B大小數(shù)據(jù)矩陣的轉(zhuǎn)置，因此也僅進(jìn)行相關(guān)實(shí)現(xiàn)的比較。在HPTT測(cè)試中，線程數(shù)設(shè)置為16，即使用FT-M7032中所有的CPU核來(lái)進(jìn)行并行矩陣轉(zhuǎn)置運(yùn)算。

對(duì)于ftmMT-64，測(cè)試結(jié)果如表1所示。與HPTT相比，在輸入矩陣小于1 024×1 024規(guī)模時(shí)，ftmMT-64的性能較低。主要原因是ftmMT-64的時(shí)間開(kāi)銷里包含了刷Cache的開(kāi)銷(算法1中Step 1和Step 3的開(kāi)銷)，該部分開(kāi)銷是相對(duì)固定的，屬于ftmMT的固有開(kāi)銷，與輸入的矩陣大小無(wú)關(guān)。當(dāng)輸入的矩陣規(guī)模較小時(shí)，刷Cache的開(kāi)銷在總時(shí)間開(kāi)銷中所占比例較大，從而使得ftmMT算法的性能低于HPTT。當(dāng)輸入的矩陣規(guī)模變大時(shí)，DSP端函數(shù)__MTrn()的時(shí)間開(kāi)銷增加，刷Cache開(kāi)銷所占比例也逐漸變小。當(dāng)輸入矩陣大于等于1 024×1 024時(shí)，ftmMT-64的性能均優(yōu)于HPTT。當(dāng)輸入矩陣規(guī)模為4 096×4 096時(shí)，ftmMT-64的耗時(shí)僅為HPTT的15.14%，計(jì)算速度為HPTT的6.61倍。在ftmMT-64的測(cè)試中，轉(zhuǎn)置帶寬最高為30.98 GB/s，與DDR的理論帶寬42.62 GB/s相比，對(duì)應(yīng)DDR有效帶寬利用率達(dá)到了72.69%。除了刷Cache開(kāi)銷包含在ftmMT總開(kāi)銷內(nèi)的原因外，影響ftmMT的DDR有效帶寬利用率的另一個(gè)重要因素是DDR本身的實(shí)測(cè)帶寬與理論值存在差距。

表1 不同矩陣大小下ftmMT-64與HPTT的性能Tab.1 Performance of ftmMT-64 and HPTT implementations under different matrix sizes

對(duì)于ftmMT-32，測(cè)試結(jié)果如表2所示。當(dāng)輸入矩陣規(guī)模為8 192×8 192時(shí)，轉(zhuǎn)置帶寬達(dá)到了32.37 GB/s，相應(yīng)的DDR帶寬利用率為75.95%。與HPTT相比，當(dāng)輸入矩陣小于等于1 024×1 024時(shí)，ftmMT-32的性能較低，與ftmMT-64情況一樣，也是因?yàn)樗ache開(kāi)銷占總執(zhí)行時(shí)間的比率較大。當(dāng)輸入矩陣大于1 024×1 024時(shí)，ftmMT-32的性能也是均高于HPTT。當(dāng)輸入矩陣規(guī)模為8 192×8 192時(shí)，ftmMT-32的耗時(shí)僅為HPTT的11.13%，轉(zhuǎn)置帶寬加速比高達(dá)8.99。

表2 不同矩陣大小下ftmMT-32與HPTT的性能Tab.2 Performance of ftmMT-32 and HPTT implementations under different matrix sizes

對(duì)于ftmMT-16，測(cè)試結(jié)果如表3所示。當(dāng)輸入矩陣規(guī)模為16 384×16 384時(shí)，ftmMT-32的轉(zhuǎn)置帶寬達(dá)到了31.78 GB/s，DDR的有效帶寬利用率為74.57%。

表3 不同矩陣大小下ftmMT-16的性能Tab.3 Performance of ftmMT-16 implementation under different matrix sizes

從上述實(shí)驗(yàn)結(jié)果可知，ftmMT-64、 ftmMT-32和ftmMT-16均能有效發(fā)揮DDR存儲(chǔ)器的性能，且ftmMT-64與ftmMT-32相對(duì)運(yùn)行在多核CPU上的HPTT實(shí)現(xiàn)了較高的性能提升。

5 結(jié)論

本文針對(duì)飛騰異構(gòu)多核DSP的體系結(jié)構(gòu)特征與矩陣轉(zhuǎn)置操作的特點(diǎn)，基于DSP VPU的實(shí)現(xiàn)方法，提出了一種適配不同數(shù)據(jù)位寬(8 B、4 B以及2 B)矩陣的并行矩陣轉(zhuǎn)置算法ftmMT。ftmMT基于矩陣分塊實(shí)現(xiàn)了多個(gè)DSP核的并行處理，通過(guò)隱式乒乓設(shè)計(jì)實(shí)現(xiàn)了片上向量化轉(zhuǎn)置與片外訪存的重疊，有效克服片外存儲(chǔ)帶寬受限、頻繁調(diào)用DMA以及訪存效率較低等問(wèn)題。通過(guò)開(kāi)展相關(guān)實(shí)驗(yàn)，證明了ftmMT能夠顯著加快多核DSP上的矩陣轉(zhuǎn)置操作，DDR的帶寬利用率最高達(dá)到75.95%；與HPTT相比，可獲得高達(dá)8.99倍的性能加速。該項(xiàng)研究對(duì)于推動(dòng)國(guó)產(chǎn)化多核DSP在科學(xué)計(jì)算和深度學(xué)習(xí)領(lǐng)域的應(yīng)用具有重要的意義。