賈 迅 鄔貴明 謝向輝 吳 東

(數(shù)學(xué)工程與先進(jìn)計(jì)算國(guó)家重點(diǎn)實(shí)驗(yàn)室 江蘇無(wú)錫 214125)

矩陣乘運(yùn)算核心廣泛應(yīng)用于數(shù)值計(jì)算、數(shù)字信號(hào)分析、圖像處理、人工智能等領(lǐng)域，其計(jì)算形式為

[i,j]+=A[i,k]×B[k,j],

其中，A∈M×K,B∈K×N,C∈M×N.循環(huán)k的計(jì)算存在寫后讀數(shù)據(jù)相關(guān)，無(wú)法并行；循環(huán)i和循環(huán)j的計(jì)算不存在真相關(guān)，可以并行.另外，計(jì)算過(guò)程中的存儲(chǔ)訪問(wèn)均比較規(guī)整.由于計(jì)算復(fù)雜度為O(n3)，雙精度浮點(diǎn)矩陣乘往往是大規(guī)模數(shù)值計(jì)算應(yīng)用的性能瓶頸[1]，其性能直接反映了系統(tǒng)的計(jì)算能力.以性能測(cè)試程序Linpack[2]為例，程序?qū)崿F(xiàn)了基于雙精度浮點(diǎn)矩陣乘的線性方程求解，并以單位時(shí)間內(nèi)系統(tǒng)執(zhí)行浮點(diǎn)操作的次數(shù)(FLOPS)來(lái)衡量系統(tǒng)的計(jì)算能力.因此，為雙精度浮點(diǎn)矩陣乘運(yùn)算提供高性能且高效的計(jì)算平臺(tái)對(duì)大規(guī)模數(shù)值計(jì)算應(yīng)用的性能發(fā)揮至關(guān)重要.

目前，矩陣乘的計(jì)算平臺(tái)主要包括CPU，GPGPU，F(xiàn)PGA.其中，CPU和GPGPU平臺(tái)上矩陣乘的實(shí)現(xiàn)和應(yīng)用較為成熟.處理器廠商會(huì)提供針對(duì)各自產(chǎn)品架構(gòu)進(jìn)行性能優(yōu)化的計(jì)算函數(shù)接口庫(kù)，如Intel MKL[3]，IBM ESSL[4]，NVIDIA cuBLAS[5]等.由于處理器架構(gòu)面向通用計(jì)算，CPU和GPGPU上矩陣乘的計(jì)算效率無(wú)法達(dá)到最優(yōu)，如Intel Knights Corner眾核處理器矩陣乘的計(jì)算效率僅為89%[6]，NVIDIA P100的計(jì)算效率約為93%[7].矩陣乘效率對(duì)系統(tǒng)計(jì)算效率的影響幾乎是線性的，矩陣乘效率的損失必然導(dǎo)致系統(tǒng)計(jì)算效率的降低[8].

FPGA計(jì)算平臺(tái)具有靈活的可編程性、細(xì)粒度的并行能力和豐富的邏輯資源.同時(shí)，矩陣乘運(yùn)算數(shù)據(jù)并行性好和存儲(chǔ)訪問(wèn)規(guī)整的算法特點(diǎn)適合并行計(jì)算結(jié)構(gòu)的設(shè)計(jì)[9]，因此基于FPGA平臺(tái)實(shí)現(xiàn)針對(duì)矩陣乘的并行計(jì)算結(jié)構(gòu)成為學(xué)術(shù)研究的熱點(diǎn)，并且已經(jīng)取得一些研究成果[10].以浮點(diǎn)乘加運(yùn)算部件作為計(jì)算單元(processing elements， PE)的核心，多個(gè)計(jì)算單元組成廣播[11]或陣列結(jié)構(gòu)[12-14]進(jìn)行加速計(jì)算是典型的實(shí)現(xiàn)方式.其中文獻(xiàn)[12]提出的線性陣列計(jì)算結(jié)構(gòu)可以處理任意規(guī)模的矩陣乘運(yùn)算，且數(shù)據(jù)傳輸和存儲(chǔ)的效率較高.由于矩陣乘的性能受限于FPGA的可用的邏輯資源和時(shí)鐘頻率，即使以資源利用和時(shí)鐘頻率為優(yōu)化目標(biāo)的結(jié)構(gòu)設(shè)計(jì)，其計(jì)算性能也只能達(dá)到203.1 GFLOPS[15](252個(gè)PE，時(shí)鐘頻率為403 MHz).顯然，目前基于FPGA的矩陣乘無(wú)法滿足大規(guī)模數(shù)值計(jì)算應(yīng)用對(duì)計(jì)算能力的需求.

近年來(lái)，面向特定應(yīng)用的硬件定制加速成為高性能、高效計(jì)算結(jié)構(gòu)設(shè)計(jì)的一個(gè)重要趨勢(shì)[16-17].Google公開(kāi)了TPU(tensor processing unit)芯片的設(shè)計(jì)[18]，其核心是整數(shù)矩陣乘脈動(dòng)陣列計(jì)算結(jié)構(gòu).NVIDIA最新發(fā)布的V100 GPGPU[19]芯片中引入了稱為“Tensor Core”的16 b，32 b混合精度矩陣乘陣列計(jì)算結(jié)構(gòu)，單芯片集成了672個(gè)TensorCore，可提供120 TFLOPS的計(jì)算能力.由于TPU和V100芯片中矩陣乘定制加速結(jié)構(gòu)的設(shè)計(jì)面向深度學(xué)習(xí)類應(yīng)用，結(jié)構(gòu)支持的計(jì)算精度較低[20-21]，因而無(wú)法直接應(yīng)用于對(duì)計(jì)算精度要求較高的數(shù)值計(jì)算領(lǐng)域.

本文以線性陣列計(jì)算結(jié)構(gòu)為基礎(chǔ)，對(duì)其進(jìn)行結(jié)構(gòu)優(yōu)化，并采用硬件定制的方法設(shè)計(jì)實(shí)現(xiàn)了支持雙精度浮點(diǎn)運(yùn)算的矩陣乘協(xié)處理器.同時(shí)，本文建立了性能模型并深入分析了各結(jié)構(gòu)設(shè)計(jì)參數(shù)對(duì)協(xié)處理器實(shí)際計(jì)算性能和效率的影響.此外，本文還驗(yàn)證了矩陣乘協(xié)處理器的功能正確性并評(píng)估了其硬件實(shí)現(xiàn)的開(kāi)銷.本文探索了硬件定制結(jié)構(gòu)設(shè)計(jì)在雙精度浮點(diǎn)矩陣乘加速計(jì)算中的應(yīng)用，研究成果對(duì)提升現(xiàn)有計(jì)算系統(tǒng)的性能和效率有一定的借鑒意義.

1 線性陣列計(jì)算結(jié)構(gòu)及大規(guī)模矩陣乘算法

線性陣列計(jì)算結(jié)構(gòu)[12]如圖1所示，多個(gè)計(jì)算單元線性互連，每個(gè)計(jì)算單元包含局部存儲(chǔ)器c以存儲(chǔ)矩陣C的分塊數(shù)據(jù)，寄存器a和b分別存儲(chǔ)矩陣A和矩陣B的數(shù)據(jù)元素，一個(gè)雙精度浮點(diǎn)乘加運(yùn)算部件和一個(gè)計(jì)數(shù)器CNT.

Fig. 1 Linear array computing architecture圖1 線性陣列計(jì)算結(jié)構(gòu)

數(shù)據(jù)b在PE之間進(jìn)行傳輸，其有效時(shí)驅(qū)動(dòng)PE執(zhí)行一次乘加計(jì)算：c[CNT]=c[CNT]+a×b.

算法1給出了基于上述線性陣列計(jì)算結(jié)構(gòu)的大規(guī)模分塊矩陣乘算法[12].參數(shù)Sp，St分別是矩陣C行數(shù)M和列數(shù)N的分塊粒度，同時(shí)也表示線性陣列計(jì)算結(jié)構(gòu)中PE的數(shù)目和各PE集成的局部存儲(chǔ)器c的深度.算法中索引為p和t的最內(nèi)層循環(huán)對(duì)應(yīng)線性陣列計(jì)算結(jié)構(gòu)執(zhí)行的浮點(diǎn)矩乘迭代計(jì)算，即Sp個(gè)PE實(shí)現(xiàn)矩陣A中規(guī)模為Sp×1的矩陣子塊與矩陣B中規(guī)模為1×St的矩陣子塊乘計(jì)算，并將計(jì)算結(jié)果與矩陣C中規(guī)模為Sp×St的子塊相加.陣列計(jì)算結(jié)構(gòu)在單次迭代計(jì)算中共執(zhí)行2SpSt次雙精度浮點(diǎn)運(yùn)算.

Fig. 2 Data block mapping for the linear array圖2 矩陣分塊數(shù)據(jù)到陣列計(jì)算結(jié)構(gòu)的映射

算法1. 基于線性陣列結(jié)構(gòu)的大規(guī)模矩陣乘算法.

輸入：雙精度浮點(diǎn)矩陣A、矩陣B、矩陣C；

輸出：雙精度浮點(diǎn)矩陣C.

① ForTp=0 toM-1,Sp

② ForTt=0 toN-1,St

③ LoadC[Tp:Tp+Sp-1,Tt:Tt+St-1]；

④ Fork=0 toK-1

⑤ LoadA[Tp:Tp+Sp-1,k]；

⑥ LoadB[k,Tt:Tt+St-1]；

⑦ Forp=TptoTp+Sp-1

⑧ Fort=TttoTt+St-1

⑨C[p,t]+=A[p,k]×B[k,t]；

⑩ End For

算法1中索引為Tp，Tt，k的3層循環(huán)在線性陣列計(jì)算結(jié)構(gòu)進(jìn)行矩陣乘迭代計(jì)算的同時(shí)，處理矩陣分塊數(shù)據(jù)的加載與寫回，其功能在硬件上由數(shù)據(jù)傳輸控制部件DTC(data transfer controller)實(shí)現(xiàn).各矩陣分塊數(shù)據(jù)到線性陣列計(jì)算結(jié)構(gòu)的映射如圖2所示：

線性陣列結(jié)構(gòu)經(jīng)過(guò)K輪迭代計(jì)算才能得到矩陣C中一個(gè)分塊數(shù)據(jù)的計(jì)算結(jié)果.因此，在K輪迭代計(jì)算開(kāi)始前，數(shù)據(jù)傳輸控制部件需要先將矩陣C中規(guī)模為Sp×St的分塊數(shù)據(jù)C[Tp:Tp+Sp-1,Tt:Tt+St-1]分別加載到Sp個(gè)計(jì)算單元的深度為St的局部存儲(chǔ)器c中，并在K輪迭代計(jì)算完成后將對(duì)應(yīng)的分塊計(jì)算結(jié)果寫回主存.另外，每次迭代計(jì)算中各PE需要訪問(wèn)矩陣A中相同的數(shù)據(jù)元素.因此，每次迭代計(jì)算開(kāi)始前，數(shù)據(jù)傳輸控制部件還需要將矩陣A中規(guī)模為Sp×1的分塊數(shù)據(jù)A[Tp:Tp+Sp-1,k]加載到Sp個(gè)計(jì)算單元的寄存器a中.

矩陣A和矩陣C的分塊數(shù)據(jù)全部加載完成后，數(shù)據(jù)傳輸控制部件再加載矩陣B的數(shù)據(jù)元素以驅(qū)動(dòng)陣列計(jì)算結(jié)構(gòu)進(jìn)行迭代計(jì)算.矩陣乘運(yùn)算執(zhí)行完成共需K×(MSp)×(NSt)次迭代計(jì)算.對(duì)計(jì)算過(guò)程中數(shù)據(jù)傳輸控制部件傳輸矩陣A，B，C分塊數(shù)據(jù)的數(shù)據(jù)量、時(shí)間間隔以及次數(shù)的量化分析如表1所示：

Table 1 Analysis on Data Transfer for Matrix Multiplication表1 針對(duì)矩陣乘運(yùn)算過(guò)程中數(shù)據(jù)傳輸?shù)姆治?/p>

文獻(xiàn)[12]基于Xilinx公司Virtex-5系列的FPGA實(shí)現(xiàn)了上述線性陣列計(jì)算結(jié)構(gòu).PE個(gè)數(shù)Sp設(shè)置為96、局部存儲(chǔ)器深度St設(shè)置為1 024的情況下，F(xiàn)PGA的最大運(yùn)行頻率為158.08 MHz.輸入矩陣的規(guī)模M，N，K均為8 192時(shí)，陣列計(jì)算結(jié)構(gòu)完成雙精度浮點(diǎn)矩陣乘運(yùn)算所需的時(shí)間為40.48 s，即結(jié)構(gòu)的計(jì)算性能為27.16 GFLOPS，效率為89.49%.

線性陣列計(jì)算結(jié)構(gòu)雖然存儲(chǔ)效率較高，但每次迭代計(jì)算中各PE需要等待計(jì)算所需的矩陣分塊數(shù)據(jù)全部加載完成后才能開(kāi)始計(jì)算，計(jì)算與存儲(chǔ)訪問(wèn)在時(shí)間上是串行的，計(jì)算資源的閑置必然導(dǎo)致計(jì)算效率的損失，高效計(jì)算結(jié)構(gòu)的設(shè)計(jì)必須對(duì)此進(jìn)行改進(jìn).由于結(jié)構(gòu)設(shè)計(jì)面向FPGA計(jì)算平臺(tái)，文獻(xiàn)[12]僅以FPGA可用的邏輯資源作為設(shè)計(jì)約束，并未考慮局部存儲(chǔ)器的深度、輸入矩陣的規(guī)模等參數(shù)對(duì)線性陣列計(jì)算結(jié)構(gòu)實(shí)際計(jì)算性能和效率的影響，而這些參數(shù)在定制協(xié)處理器的設(shè)計(jì)中尤為重要.另外，本文協(xié)處理器的設(shè)計(jì)除了考慮核心的計(jì)算結(jié)構(gòu)外，還需要考慮其與主機(jī)端的接口，以實(shí)現(xiàn)控制信息和計(jì)算數(shù)據(jù)的傳輸，從而構(gòu)成完整的加速計(jì)算系統(tǒng).

2 協(xié)處理器的設(shè)計(jì)

2.1 線性陣列計(jì)算結(jié)構(gòu)的優(yōu)化

根據(jù)算法1給出的大規(guī)模分塊矩陣乘運(yùn)算算法，線性陣列計(jì)算結(jié)構(gòu)在進(jìn)行矩陣乘迭代計(jì)算前，需要先從主存讀入矩陣A和C的分塊數(shù)據(jù)；K輪迭代計(jì)算完成后，需要先將矩陣C的分塊計(jì)算結(jié)果寫回主存，才能開(kāi)始下一輪的迭代計(jì)算.矩陣分塊數(shù)據(jù)在主存和線性陣列計(jì)算結(jié)構(gòu)之間進(jìn)行傳輸?shù)倪^(guò)程中，結(jié)構(gòu)中各計(jì)算單元無(wú)法進(jìn)行計(jì)算，計(jì)算資源被閑置，從而導(dǎo)致線性陣列計(jì)算結(jié)構(gòu)整體的效率較低.針對(duì)這個(gè)問(wèn)題，本文采用雙緩沖技術(shù)對(duì)陣列計(jì)算結(jié)構(gòu)中各計(jì)算單元的設(shè)計(jì)進(jìn)行了優(yōu)化.優(yōu)化后，計(jì)算單元的結(jié)構(gòu)如圖3所示：

Fig. 3 Structure of PE optimized with double buffering圖3 采用雙緩沖優(yōu)化技術(shù)的PE結(jié)構(gòu)

圖3所示的PE結(jié)構(gòu)中除了已有的寄存器a，b和局部存儲(chǔ)器c之外，還包括雙緩沖優(yōu)化引入的存儲(chǔ)矩陣A，B數(shù)據(jù)元素的寄存器a_buf，b_buf以及存儲(chǔ)矩陣C分塊數(shù)據(jù)的局部存儲(chǔ)器c_buf.

應(yīng)用雙緩沖優(yōu)化技術(shù)后，陣列計(jì)算結(jié)構(gòu)可以在使用當(dāng)前緩沖中的數(shù)據(jù)進(jìn)行迭代計(jì)算的同時(shí)，提前將下一輪迭代計(jì)算所需的矩陣分塊數(shù)據(jù)裝入到另一緩沖，實(shí)現(xiàn)計(jì)算與訪存的重疊，隱藏存儲(chǔ)訪問(wèn)的開(kāi)銷.此時(shí)，高效的存儲(chǔ)訪問(wèn)調(diào)度成為應(yīng)用雙緩沖優(yōu)化提升計(jì)算結(jié)構(gòu)資源利用率的關(guān)鍵.

2.2 存儲(chǔ)訪問(wèn)調(diào)度

由表1列出的矩陣乘運(yùn)算過(guò)程中矩陣分塊數(shù)據(jù)的傳輸可知，線性陣列計(jì)算結(jié)構(gòu)在進(jìn)行矩陣乘運(yùn)算時(shí)，每個(gè)時(shí)鐘周期僅需訪問(wèn)矩陣B的一個(gè)數(shù)據(jù)元素，其對(duì)存儲(chǔ)帶寬的需求較低；而矩陣A，C的訪問(wèn)時(shí)間間隔分別為St和K×St個(gè)時(shí)鐘周期，因此可以利用讀取矩陣B數(shù)據(jù)元素剩余的存儲(chǔ)訪問(wèn)帶寬提前加載矩陣A，C的分塊數(shù)據(jù)至PE的緩沖中，或?qū)⒕仃嘋的分塊計(jì)算結(jié)果從緩沖寫回主存.

算法2給出了基于雙緩沖的存儲(chǔ)訪問(wèn)調(diào)度算法.數(shù)據(jù)傳輸控制部件在硬件上實(shí)現(xiàn)了該算法，并支持雙緩沖讀寫訪問(wèn)的控制.

算法2. 基于雙緩沖的存儲(chǔ)訪問(wèn)調(diào)度算法.

輸入：雙精度浮點(diǎn)矩陣A、矩陣B、矩陣C；

輸出：雙精度浮點(diǎn)矩陣C.

① LoadC[0:Sp-1,0:St-1]；

② LoadA[0:Sp-1,0]；

③ ForTp=0 toM-1,Sp

④ ForTt=0 toN-1,St

⑤ Fork=0 toK-1

⑥ LoadB[k,Tt:Tt+St-1]；

LoadA[Tp:Tp+Sp-1,k+1] or

StoreC[Tp:Tp+Sp-1,Tt-St:Tt-1] or

LoadC[Tp:Tp+Sp-1,Tt+St-1:Tt+2St-1]；

⑦ Forp=TptoTp+Sp-1

⑧ Fort=TttoTt+St-1

⑨C[p,t]+=A[p,k]×B[k,t]；

⑩ End For

算法2中，陣列計(jì)算結(jié)構(gòu)進(jìn)行第一次迭代計(jì)算所需的矩陣A，C的分塊數(shù)據(jù)先被讀入.讀取矩陣B的數(shù)據(jù)元素進(jìn)行迭代計(jì)算時(shí)，優(yōu)先利用剩余的存儲(chǔ)訪問(wèn)帶寬提前讀取下一次迭代計(jì)算所需的矩陣A的分塊數(shù)據(jù)；待矩陣A的分塊數(shù)據(jù)讀入完成，再利用剩余帶寬寫回前K輪迭代計(jì)算得到的矩陣C的分塊計(jì)算結(jié)果；待計(jì)算結(jié)果寫回完成，最后利用剩余帶寬裝入后K輪迭代計(jì)算所需的矩陣C的分塊數(shù)據(jù).所有迭代計(jì)算完成后，矩陣C中最后一個(gè)分塊計(jì)算結(jié)果被寫回.

理想情況下，線性陣列計(jì)算結(jié)構(gòu)在開(kāi)始新一次的迭代計(jì)算時(shí)，各PE所需的矩陣數(shù)據(jù)均已在緩沖中，即存儲(chǔ)訪問(wèn)的延遲可以有效地被計(jì)算隱藏，結(jié)構(gòu)中的計(jì)算資源可以得到充分利用.而實(shí)際的存儲(chǔ)訪問(wèn)調(diào)度受限于存儲(chǔ)帶寬、計(jì)算單元數(shù)、局部存儲(chǔ)器的深度以及矩陣規(guī)模等參數(shù)，其對(duì)矩陣乘協(xié)處理器實(shí)際性能的影響將在2.3節(jié)進(jìn)行分析.

2.3 系統(tǒng)結(jié)構(gòu)

基于矩陣乘協(xié)處理器的加速計(jì)算系統(tǒng)的整體結(jié)構(gòu)如圖4所示.系統(tǒng)由主機(jī)端處理器芯片和雙精度浮點(diǎn)矩陣乘協(xié)處理器芯片2部分組成.其中，矩陣乘協(xié)處理器通過(guò)IO總線與主機(jī)端處理器相連.系統(tǒng)采用PCIe-3.0×16作為芯片間的互連接口，雙向通信帶寬最高為31.5 GBps.

Fig. 4 System architecture based on matrix multiplication coprocessor圖4 基于矩陣乘協(xié)處理器計(jì)算系統(tǒng)的架構(gòu)

主機(jī)端和協(xié)處理器端存儲(chǔ)空間相互獨(dú)立，待計(jì)算的矩陣數(shù)據(jù)存儲(chǔ)在協(xié)處理器端的內(nèi)存.本文協(xié)處理器的設(shè)計(jì)面向大規(guī)模矩陣乘，支持的矩陣規(guī)模至少為16 384×16 384.因此，主機(jī)端與協(xié)處理器端的內(nèi)存容量配置為8 GB，并采用DDR3-2133[22]規(guī)格的存儲(chǔ)芯片，接口的數(shù)據(jù)寬度為128 b，帶寬為34.13 GBps.

系統(tǒng)對(duì)主機(jī)端處理器的類型沒(méi)有限制，只要支持PCIe傳輸協(xié)議即可，其指令集系統(tǒng)架構(gòu)可以是x86，Power，ARM，RISC-V[23]等.矩陣乘協(xié)處理器由應(yīng)用雙緩沖優(yōu)化的線性陣列計(jì)算結(jié)構(gòu)、存儲(chǔ)訪問(wèn)控制器、協(xié)處理器與主機(jī)端接口(coprocessor-host interface， CHI)3部分構(gòu)成.乘加運(yùn)算部件是計(jì)算單元的核心，本文采用全流水且兼容IEEE-754浮點(diǎn)數(shù)標(biāo)準(zhǔn)[24]的乘加運(yùn)算部件設(shè)計(jì)，以保證協(xié)處理器的計(jì)算性能.CHI接口實(shí)現(xiàn)主機(jī)端與協(xié)處理器之間控制和數(shù)據(jù)的傳輸.主機(jī)端向矩陣乘協(xié)處理器發(fā)送的主要命令如下：

1)Initialization對(duì)協(xié)處理器的控制與狀態(tài)寄存器進(jìn)行初始化；

2)DataIn將待計(jì)算矩陣數(shù)據(jù)從主機(jī)端內(nèi)存加載到協(xié)處理器端內(nèi)存；

3)MatrixMultiply啟動(dòng)協(xié)處理器上的雙精度浮點(diǎn)矩陣乘運(yùn)算；

4)DataOut將矩陣乘運(yùn)算的結(jié)果從協(xié)處理器端內(nèi)存讀出到主機(jī)端內(nèi)存.

接收到來(lái)自主機(jī)端發(fā)送的命令，CHI接口對(duì)命令進(jìn)行分析并生成線性陣列計(jì)算結(jié)構(gòu)相關(guān)的控制信號(hào).協(xié)處理器在進(jìn)行矩陣乘運(yùn)算的過(guò)程中，主機(jī)端與協(xié)處理器之間不需要進(jìn)行數(shù)據(jù)傳輸，從而簡(jiǎn)化了協(xié)處理器的設(shè)計(jì)與調(diào)試.

應(yīng)用程序在計(jì)算系統(tǒng)上加速矩陣乘運(yùn)算的過(guò)程如圖5所示.當(dāng)矩陣乘協(xié)處理器連接至主機(jī)端時(shí)，主機(jī)端發(fā)送Initialization命令對(duì)協(xié)處理器進(jìn)行初始化，如設(shè)置浮點(diǎn)計(jì)算的舍入模式、溢出控制等.應(yīng)用程序運(yùn)行在主機(jī)端，當(dāng)其調(diào)用雙精度浮點(diǎn)矩陣乘運(yùn)算核心，如基礎(chǔ)線性代數(shù)子程序庫(kù)(basic linear algebra subroutines， BLAS)中的函數(shù)blas_dgemm(M,N,K,A,B,C,…)時(shí)，主機(jī)端通過(guò)CHI接口向協(xié)處理器發(fā)送DataIn命令和函數(shù)調(diào)用傳入的矩陣規(guī)模、矩陣數(shù)據(jù)在內(nèi)存中的起始地址等參數(shù).協(xié)處理器的數(shù)據(jù)傳輸控制部件將待計(jì)算的矩陣數(shù)據(jù)從主機(jī)端內(nèi)存以直接內(nèi)存訪問(wèn)(direct memory access， DMA)方式傳輸至協(xié)處理器端內(nèi)存.數(shù)據(jù)傳輸完成后，主機(jī)端發(fā)送MatrixMultiply命令以啟動(dòng)協(xié)處理器上的矩陣乘運(yùn)算.運(yùn)算過(guò)程中，數(shù)據(jù)傳輸控制部件按算法2設(shè)計(jì)的存儲(chǔ)訪問(wèn)調(diào)度實(shí)現(xiàn)矩陣分塊數(shù)據(jù)在協(xié)處理器端內(nèi)存和線性陣列計(jì)算結(jié)構(gòu)之間的傳輸.

Fig. 5 Computation based on matrix multiplication coprocessor圖5 基于矩陣乘協(xié)處理器的計(jì)算過(guò)程

協(xié)處理器完成矩陣乘運(yùn)算后，CHI接口向主機(jī)端發(fā)出中斷信號(hào).主機(jī)端收到中斷后，向CHI接口發(fā)送DataOut命令，協(xié)處理器的數(shù)據(jù)傳輸控制部件以DMA的方式將矩陣乘運(yùn)算的結(jié)果從協(xié)處理器端內(nèi)存?zhèn)鬏敾刂鳈C(jī)端內(nèi)存.運(yùn)算結(jié)果傳輸完成后，矩陣乘運(yùn)算調(diào)用結(jié)束，應(yīng)用程序繼續(xù)執(zhí)行.

3 性能模型

為了量化分析雙精度浮點(diǎn)矩陣乘協(xié)處理器的各結(jié)構(gòu)設(shè)計(jì)參數(shù)對(duì)其實(shí)際計(jì)算性能和效率的影響，本文基于算法2設(shè)計(jì)的存儲(chǔ)訪問(wèn)調(diào)度算法構(gòu)建了協(xié)處理器的性能模型.模型中的基本參數(shù)包括線性陣列計(jì)算結(jié)構(gòu)中PE的個(gè)數(shù)Sp、PE集成的局部存儲(chǔ)器的深度St以及輸入矩陣的規(guī)模M，N，K.此外，本文將模型中協(xié)處理器的時(shí)鐘頻率記為Freq(單位為GHz)，協(xié)處理器端內(nèi)存的訪問(wèn)帶寬記為Band(單位為GBps)，矩陣乘運(yùn)算核心在協(xié)處理器上執(zhí)行完成所需的時(shí)間記為T(單位為s)，協(xié)處理器的實(shí)際計(jì)算性能和效率記為Perf(單位為GFLOPS)和Eff.

當(dāng)計(jì)算過(guò)程中的存儲(chǔ)訪問(wèn)延遲完全被計(jì)算隱藏，協(xié)處理器的計(jì)算性能可以充分發(fā)揮，此時(shí)的計(jì)算效率最高.就本文提出的存儲(chǔ)訪問(wèn)調(diào)度算法，每次迭代計(jì)算提前讀入矩陣A中規(guī)模為Sp的分塊數(shù)據(jù)的延遲與讀入矩陣B中規(guī)模為St的分塊數(shù)據(jù)的延遲之和應(yīng)當(dāng)小于等于對(duì)應(yīng)的計(jì)算延遲，即St個(gè)時(shí)鐘周期.同理，K輪迭代計(jì)算的總存儲(chǔ)訪問(wèn)延遲應(yīng)當(dāng)小于等于對(duì)應(yīng)的總計(jì)算延遲，即K×St個(gè)時(shí)鐘周期.存儲(chǔ)訪問(wèn)的數(shù)據(jù)粒度為8 B的雙精度浮點(diǎn)數(shù)，上述2個(gè)約束條件可表示為

(1)

協(xié)處理器端DDR3-2133存儲(chǔ)芯片的數(shù)據(jù)傳輸率為2.133 GTs，單次傳輸?shù)臄?shù)據(jù)寬度為128 b.若協(xié)處理器的工作時(shí)鐘頻率Freq=1.5 GHz，則其每個(gè)時(shí)鐘周期可訪問(wèn)2個(gè)雙精度浮點(diǎn)數(shù).此時(shí)，式(1)的約束條件可進(jìn)一步改寫為

(2)

若參數(shù)Sp，St，K可以滿足式(2)中的約束條件，則計(jì)算過(guò)程中的存儲(chǔ)訪問(wèn)延遲可以完全被計(jì)算隱藏.此時(shí)，矩陣乘運(yùn)算核心在協(xié)處理器上的執(zhí)行完成所需的時(shí)間T最小.其由3部分時(shí)間構(gòu)成，包括讀入矩陣A和C第1個(gè)子塊數(shù)據(jù)的時(shí)間、迭代計(jì)算的時(shí)間以及將矩陣C中最后1個(gè)子塊的計(jì)算結(jié)果寫回內(nèi)存的時(shí)間.計(jì)算時(shí)間T可表示為

(3)

若參數(shù)Sp，St，K不滿足式(2)的約束條件，存在St1-SpSt這2種可能.下文將分別分析這2種情況下矩陣乘運(yùn)算核心在協(xié)處理器上的執(zhí)行時(shí)間.

對(duì)于第1種情況，提前讀入矩陣A子塊數(shù)據(jù)的存儲(chǔ)訪問(wèn)延遲無(wú)法完全被計(jì)算時(shí)間隱藏，每次迭代計(jì)算的時(shí)間增加(Sp-St)2個(gè)時(shí)鐘周期.根據(jù)算法2，當(dāng)前迭代計(jì)算執(zhí)行過(guò)程中，寫回前K輪迭代計(jì)算中矩陣C子塊計(jì)算結(jié)果和提前讀入后K輪迭代計(jì)算中矩陣C子塊計(jì)算數(shù)據(jù)的存儲(chǔ)訪問(wèn)延遲均無(wú)法被計(jì)算隱藏，每K輪迭代計(jì)算的時(shí)間增加Sp×St個(gè)時(shí)鐘周期.此時(shí)，矩陣乘運(yùn)算核心在協(xié)處理器上執(zhí)行完成所需的時(shí)間T可表示為

(4)

對(duì)于第2種情況，當(dāng)前迭代計(jì)算的執(zhí)行時(shí)間可以完全隱藏提前讀入矩陣A子塊數(shù)據(jù)的存儲(chǔ)訪問(wèn)延遲，但只能部分隱藏寫回前K輪迭代計(jì)算中矩陣C子塊計(jì)算結(jié)果和提前讀入后K輪迭代計(jì)算中矩陣C子塊計(jì)算數(shù)據(jù)的存儲(chǔ)訪問(wèn)延遲，其導(dǎo)致每K輪迭代計(jì)算的時(shí)間增加Sp×St-K×(St-Sp)2個(gè)時(shí)鐘周期.此時(shí)，矩陣乘運(yùn)算核心在協(xié)處理器上執(zhí)行完成所需的時(shí)間T可表示為

(5)

協(xié)處理器完成矩陣乘核心的加速計(jì)算共需執(zhí)行2MNK次雙精度浮點(diǎn)操作.基于對(duì)3種情況下矩陣乘運(yùn)算核心在協(xié)處理器上執(zhí)行時(shí)間T的分析結(jié)果，可以構(gòu)建式(6)所示的協(xié)處理器的性能模型.此時(shí)，協(xié)處理器的計(jì)算效率可表示為實(shí)際計(jì)算性能與峰值計(jì)算性能的比值，即Eff=Perf(2×Sp×Freq).

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

本節(jié)以第3節(jié)提出的性能模型為基礎(chǔ)，分析矩陣乘協(xié)處理器結(jié)構(gòu)設(shè)計(jì)中雙緩沖優(yōu)化技術(shù)、局部存儲(chǔ)器的深度和輸入矩陣規(guī)模等參數(shù)對(duì)協(xié)處理器實(shí)際計(jì)算性能和效率的影響.實(shí)驗(yàn)假設(shè)線性陣列計(jì)算結(jié)構(gòu)的PE數(shù)為1 024，協(xié)處理器的時(shí)鐘頻率為1.5 GHz，待計(jì)算的矩陣均為方陣，即M=N=K=n.此外，本節(jié)還驗(yàn)證了矩陣乘協(xié)處理器功能實(shí)現(xiàn)的正確性；在主流工藝下評(píng)估了其硬件實(shí)現(xiàn)的開(kāi)銷，并與同等工藝的GPGPU芯片進(jìn)行了對(duì)比，以探索其硬件可實(shí)現(xiàn)性和在實(shí)際應(yīng)用中的優(yōu)勢(shì).

4.1 雙緩沖優(yōu)化技術(shù)對(duì)矩陣乘協(xié)處理器性能的影響

實(shí)驗(yàn)在不同的矩陣規(guī)模下，分別統(tǒng)計(jì)了應(yīng)用雙緩沖優(yōu)化技術(shù)前后矩陣乘協(xié)處理器的性能，以分析應(yīng)用雙緩沖優(yōu)化技術(shù)對(duì)協(xié)處理器性能的影響.實(shí)驗(yàn)中，各PE局部存儲(chǔ)器的深度均設(shè)置為2 048，實(shí)驗(yàn)的結(jié)果如圖6所示.

從圖6可以看出，不同矩陣規(guī)模下，雙緩沖優(yōu)化技術(shù)的應(yīng)用均可以提升矩陣乘協(xié)處理器的計(jì)算性能，協(xié)處理器的性能平均提升30.06%.矩陣規(guī)模為4 096時(shí)，性能提升最為明顯，達(dá)到了45.45%.可見(jiàn)，本文設(shè)計(jì)的存儲(chǔ)訪問(wèn)調(diào)度算法可以在計(jì)算過(guò)程中利用雙緩沖優(yōu)化有效隱藏存儲(chǔ)訪問(wèn)的延遲.

同時(shí)，雙緩沖優(yōu)化對(duì)矩陣乘協(xié)處理器的性能提升受矩陣規(guī)模的影響.對(duì)此分析如下：當(dāng)矩陣規(guī)模較小時(shí)，計(jì)算無(wú)法有效隱藏存儲(chǔ)訪問(wèn)的延遲，應(yīng)用雙緩沖優(yōu)化對(duì)性能的提升作用有限；隨著矩陣規(guī)模逐漸增加，計(jì)算可以隱藏更多的訪存延遲，此時(shí)應(yīng)用雙緩沖優(yōu)化可以顯著提升協(xié)處理器的計(jì)算性能；當(dāng)矩陣規(guī)模較大時(shí)，計(jì)算可以很好隱藏存儲(chǔ)訪問(wèn)的延遲，而實(shí)驗(yàn)中局部存儲(chǔ)器的深度保持不變，雙緩沖優(yōu)化對(duì)性能的提升作用逐漸降低.

Fig. 6 Effect of double-buffering optimization technique on coprocessor performance圖6 雙緩沖優(yōu)化技術(shù)對(duì)矩陣乘協(xié)處理器性能的提升

4.2 局部存儲(chǔ)器深度對(duì)矩陣乘協(xié)處理器性能的影響

Fig. 7 Effect of local memory depth on coprocessor performance圖7 局部存儲(chǔ)器深度對(duì)矩陣乘協(xié)處理器性能的影響

局部存儲(chǔ)器的深度直接影響協(xié)處理器硬件實(shí)現(xiàn)的開(kāi)銷，因而是矩陣乘協(xié)處理器結(jié)構(gòu)設(shè)計(jì)中一個(gè)重要的參數(shù).實(shí)驗(yàn)在128，1 024，8 192，16 384這4種矩陣規(guī)模下，分別將計(jì)算單元局部存儲(chǔ)器的深度h從128增加至8 192，統(tǒng)計(jì)矩陣乘協(xié)處理器實(shí)際計(jì)算性能的變化情況，實(shí)驗(yàn)結(jié)果如圖7所示:

從圖7可以看出，4種矩陣規(guī)模下，矩陣乘協(xié)處理器的實(shí)際計(jì)算性能均隨局部存儲(chǔ)器深度的增加而提升，局部存儲(chǔ)器深度增加大于2 048時(shí)，矩陣乘協(xié)處理器的性能基本保持不變.協(xié)處理器實(shí)際計(jì)算性能隨局部存儲(chǔ)器深度的變化可以根據(jù)式(2)中計(jì)算完全隱藏存儲(chǔ)訪問(wèn)延遲的約束條件進(jìn)行分析.當(dāng)St

4.3 矩陣規(guī)模對(duì)矩陣乘協(xié)處理器計(jì)算效率的影響

盡管本文矩陣乘協(xié)處理器的設(shè)計(jì)面向數(shù)值計(jì)算領(lǐng)域的大規(guī)模應(yīng)用，矩陣規(guī)模對(duì)協(xié)處理器實(shí)際計(jì)算性能和效率的影響仍然值得關(guān)注.實(shí)驗(yàn)在應(yīng)用雙緩沖優(yōu)化、局部存儲(chǔ)器深度設(shè)置為2 048時(shí)，統(tǒng)計(jì)了不同矩陣規(guī)模下，矩陣乘協(xié)處理器實(shí)際的計(jì)算性能和效率，以及存儲(chǔ)訪問(wèn)延遲占矩陣乘運(yùn)算核心總執(zhí)行時(shí)間的比例，實(shí)驗(yàn)結(jié)果如表2所示:

Table 2 Effect of Matrix Size on Coprocessor Performance

從表2中的數(shù)據(jù)可以看出，隨著矩陣規(guī)模的增加，矩陣乘協(xié)處理器的實(shí)際計(jì)算效率不斷提高：矩陣規(guī)模為128時(shí)，協(xié)處理器的計(jì)算效率僅為6.25%；矩陣規(guī)模增加至4 096時(shí)，協(xié)處理器的計(jì)算效率達(dá)到96.97%；矩陣規(guī)模大于等于8 192時(shí)，協(xié)處理器的計(jì)算效率可以達(dá)到99%.

對(duì)協(xié)處理器計(jì)算效率與矩陣規(guī)模相關(guān)性的分析如下：當(dāng)PE數(shù)Sp=1 024、局部存儲(chǔ)器深度St=2 048時(shí)，式(3)表示的計(jì)算完全隱藏訪存延遲的約束條件是否滿足僅取決于計(jì)算中的矩陣規(guī)模K.隨著矩陣規(guī)模的增大，計(jì)算可以更好地隱藏存儲(chǔ)訪問(wèn)的延遲，協(xié)處理器的計(jì)算效率也進(jìn)一步提升.表2中存儲(chǔ)訪問(wèn)延遲占矩陣乘總執(zhí)行時(shí)間的比例驗(yàn)證了上述分析，即較大的矩陣規(guī)模對(duì)應(yīng)較低的存儲(chǔ)訪問(wèn)延遲占比.

4.4 功能正確性驗(yàn)證與硬件實(shí)現(xiàn)開(kāi)銷的評(píng)估

本文基于Verilog硬件描述語(yǔ)言實(shí)現(xiàn)了雙精度浮點(diǎn)矩陣乘協(xié)處理器.其中，線性陣列計(jì)算結(jié)構(gòu)中各計(jì)算單元采用雙緩沖優(yōu)化技術(shù)，PE數(shù)設(shè)置為1 024、局部存儲(chǔ)器的深度設(shè)置為2 048.

對(duì)于矩陣乘協(xié)處理器功能正確性的驗(yàn)證，本文以Linpack性能測(cè)試程序HPL[25]隨機(jī)生成的雙精度浮點(diǎn)矩陣數(shù)據(jù)作為測(cè)試數(shù)據(jù)集，通過(guò)調(diào)用OpenBLAS[26-27]線性代數(shù)計(jì)算函數(shù)庫(kù)實(shí)現(xiàn)了雙精度浮點(diǎn)矩陣乘運(yùn)算核心.本文將OpenBLAS的計(jì)算結(jié)果與矩陣乘協(xié)處理器的仿真計(jì)算結(jié)果進(jìn)行了比對(duì)，從而驗(yàn)證了協(xié)處理器功能實(shí)現(xiàn)的正確性.

對(duì)于矩陣乘協(xié)處理器硬件開(kāi)銷的評(píng)估，本文基于Synopsys Design Compiler設(shè)計(jì)工具在1.5 GHz、主流工藝下對(duì)協(xié)處理器的硬件設(shè)計(jì)進(jìn)行了邏輯綜合，綜合過(guò)程中的參數(shù)配置如表3所示：

Table 3 Configurations for Synopsys Design Compiler表3 邏輯綜合的參數(shù)設(shè)置

根據(jù)Design Compiler輸出的綜合結(jié)果，矩陣乘協(xié)處理器的硬件設(shè)計(jì)完全滿足1.5 GHz頻率下的時(shí)序要求.其功耗為38.57 W，總面積為116.26 mm2.其中，組合邏輯的面積為16.78 mm2，非組合邏輯的面積為99.48 mm2.

本文將矩陣乘協(xié)處理器與采用同等工藝的NVIDIA K40 GPGPU[19]芯片進(jìn)行了對(duì)比，比較的內(nèi)容包括時(shí)鐘頻率、峰值計(jì)算性能，進(jìn)行矩陣乘運(yùn)算的實(shí)際性能和效率以及硬件開(kāi)銷等，結(jié)果如表4所示.針對(duì)雙精度浮點(diǎn)矩陣乘運(yùn)算，本文協(xié)處理器的實(shí)際性能是K40 GPGPU芯片的1.95倍，而面積開(kāi)銷僅為后者的21.05%(不考慮GPGPU芯片中PCIe、圖形顯示和存儲(chǔ)控制器等接口的硬件開(kāi)銷).

Table 4 Coprocessor versus K40 GPGPU on Performance

5 結(jié)束語(yǔ)

矩陣乘運(yùn)算廣泛應(yīng)用于科學(xué)與工程領(lǐng)域.因計(jì)算復(fù)雜度較高，且在CPU，GPGPU，F(xiàn)PGA等現(xiàn)有計(jì)算平臺(tái)上的性能和效率受限，雙精度浮點(diǎn)矩陣乘運(yùn)算往往成為大規(guī)模數(shù)值計(jì)算應(yīng)用的性能瓶頸.近年來(lái)，面向應(yīng)用的硬件定制成為高性能、高效計(jì)算結(jié)構(gòu)設(shè)計(jì)的重要趨勢(shì).本文將這一結(jié)構(gòu)設(shè)計(jì)方法應(yīng)用到矩陣乘的加速計(jì)算中，設(shè)計(jì)了雙精度浮點(diǎn)矩陣乘定制協(xié)處理器和基于協(xié)處理器的加速計(jì)算系統(tǒng).

本文雙精度浮點(diǎn)矩陣乘協(xié)處理器的設(shè)計(jì)以線性陣列計(jì)算結(jié)構(gòu)為基礎(chǔ)，應(yīng)用了雙緩沖技術(shù)并對(duì)存儲(chǔ)訪問(wèn)調(diào)度進(jìn)行了優(yōu)化.同時(shí)，本文構(gòu)建了針對(duì)協(xié)處理器的性能模型，并基于性能模型深入分析了協(xié)處理器的各結(jié)構(gòu)設(shè)計(jì)參數(shù)對(duì)其計(jì)算性能和效率的影響.此外，本文還驗(yàn)證了矩陣乘協(xié)處理器功能實(shí)現(xiàn)的正確性，并在主流工藝下評(píng)估了其硬件開(kāi)銷.本文雙精度浮點(diǎn)矩陣乘協(xié)處理器的實(shí)際計(jì)算性能可達(dá)3 TFLOPS，計(jì)算效率最高可達(dá)99%，有效實(shí)現(xiàn)了大規(guī)模數(shù)值計(jì)算應(yīng)用的加速.與同等工藝的GPGPU計(jì)算平臺(tái)相比，本文設(shè)計(jì)的矩陣乘協(xié)處理器在計(jì)算性能和硬件開(kāi)銷方面均具有明顯優(yōu)勢(shì).

本文探索了硬件定制的結(jié)構(gòu)設(shè)計(jì)在高性能計(jì)算中的應(yīng)用.后續(xù)研究將基于BLAS庫(kù)和Linpack性能測(cè)試程序評(píng)估本文雙精度浮點(diǎn)矩陣乘協(xié)處理器的應(yīng)用對(duì)大規(guī)模計(jì)算系統(tǒng)性能和效率的提升效果；同時(shí)分析協(xié)處理器實(shí)際性能隨陣列計(jì)算結(jié)構(gòu)中計(jì)算單元數(shù)增加的可擴(kuò)展性，從而進(jìn)一步優(yōu)化協(xié)處理器的結(jié)構(gòu)設(shè)計(jì).另外，協(xié)處理器的應(yīng)用適應(yīng)性也是未來(lái)研究工作的重點(diǎn).