馬鄴晨，李醒飛

(天津大學精密測試技術(shù)及儀器國家重點實驗室，天津 300072)

1 概述

作為慣導系統(tǒng)核心環(huán)節(jié)，導航解算的實時性決定著慣導系統(tǒng)能否及時跟蹤運載體的運動［1］。導航解算從計算顆粒劃分，主要包括向量加法、三維向量叉乘、四元數(shù)乘法、矩陣-矢量積和矩陣乘法5種運算，其中，矩陣乘積計算量大，耗時也最多［2］。從解算流程劃分，載體坐標系旋轉(zhuǎn)更新模塊，姿態(tài)、速度、位置更新模塊中均包含一次或多次矩陣乘法，初始對準中采用的濾波器算法也包含大量的矩陣乘法運算［3］。因此，提高浮點矩陣乘法運算速度對慣導系統(tǒng)實時性的提高具有重要意義。

之前的浮點矩陣乘積運算一般都采用PC或DSP實現(xiàn)，但這種串行處理器在應對高階數(shù)、高復雜度的算法時，更新速率并不高。伴隨超大規(guī)模集成電路技術(shù)的發(fā)展，國內(nèi)外很多學者開始研究使用具有并行處理能力的FPGA來計算浮點矩陣乘積［4］。文獻［5］提出了一種各計算單元之間不存在任何通訊的并行矩陣乘法器結(jié)構(gòu)，但其所需要的存儲空間隨矩陣維數(shù)的增加而顯著增加，且效率較低。文獻［6］在Xilinx FPGA上設(shè)計了用于雙精度浮點矩陣乘法的處理單元(Process Element，PE)陣列結(jié)構(gòu)，峰值計算性能可達到3 GFLOPS，但其并未采取任何運算優(yōu)化措施，且計算性能隨矩陣維數(shù)增加而下降。文獻［7］基于Nios II系統(tǒng)，應用浮點運算宏功能模塊設(shè)計了矩陣硬件加速器，但并未探究其加速機理;另外雖然該方法總體消耗的周期數(shù)較串行方法有所下降，但部分環(huán)節(jié)效率依舊較低。

針對上述問題，本文分析了FPGA浮點矩陣乘積運算加速原理，包括經(jīng)典脈動陣列結(jié)構(gòu)并行提速以及近兩年新提出的用于增強矩陣處理并行度的數(shù)據(jù)空間、迭代空間處理方法。基于上述原理，編寫出高性能浮點矩陣乘積模塊，用于提高計算速度;在數(shù)據(jù)交換方面，摒棄傳統(tǒng)傳輸模式，利用DMA大批量高速傳輸速度優(yōu)勢進行數(shù)據(jù)交換，從而設(shè)計出性能良好的浮點矩陣乘積加速器。最后利用嵌入式邏輯分析儀SignalTap II對運算過程進行性能實測。

2 FPGA浮點矩陣乘積加速原理

2.1 計算單元陣列結(jié)構(gòu)

矩陣乘法作為操作次數(shù)大于輸入輸出象元的運算，屬于計算受限類型。常采用具有并行-流水線特點的脈動(systolic)陣列［8］實現(xiàn)并行計算。該陣列由功能相同或相近的PE組成，運算中數(shù)據(jù)受同一個時鐘控制，在各PE間沿著預設(shè)的方向流動，在數(shù)據(jù)流流入、流出陣列的過程中完成運算。以 Cm×n=Am×p×Bp×n為例，m 個 PE 串聯(lián)而成的一維 systolic線陣運算過程如下:B陣中的數(shù)據(jù)按照行序串行流過陣列中每一個PE單元，同時A陣形成m條數(shù)據(jù)流分別輸入相應的單元，在每n個周期內(nèi)保持恒定，以列序p循環(huán)。2個數(shù)據(jù)在PE中相遇進行相應運算，并把形成的部分結(jié)果存放在加法循環(huán)回路中，這樣當B陣傳輸結(jié)束時，各PE元分別輸出結(jié)果矩陣各行計算結(jié)果。

與微處理器串行實現(xiàn)方法2×m×p×n，即ο(n3)的計算復雜性相比，采用一維systolic線性陣列結(jié)構(gòu)，計算的復雜性降低至2×p×n，即ο(n2/m)。矩陣運算時空映射而得出的二維脈動陣列能將矩陣計算的復雜性降至ο(n2/mp)。

2.2 分塊映射策略

應對不同規(guī)模的矩陣，尤其針對較大矩陣，簡單依賴增大陣列規(guī)模，以資源換速度的計算方法會受到片上資源約束，因而需將大規(guī)模問題映射到較小或(固定)規(guī)模的陣列上。現(xiàn)今常用的映射方法是循環(huán)分塊［9］，即以局部并行、全局串行的方式完成調(diào)度。其核心思想是:將計算過程分為2層循環(huán)，內(nèi)層并行完成矩陣運算，外層通過有限狀態(tài)機完成過程控制、外部存儲器與計算陣列之間的數(shù)據(jù)傳輸。

針對循環(huán)中大批量的數(shù)據(jù)訪問，對數(shù)據(jù)空間分割及迭代空間合并是增加訪問并行度的有效手段［10］。其思路是:通過建模獲取循環(huán)控制流下訪存密集數(shù)組的訪問形式及相關(guān)性，若迭代空間中矩陣數(shù)目大于1，且未遇到新矩陣，則依照迭代空間維度順序由內(nèi)向外將所有數(shù)據(jù)空間內(nèi)矩陣分為(N×所需寬度)個矩陣。為不失分塊后矩陣間相關(guān)性，再合并迭代空間，合并過程減少了循環(huán)次數(shù)，因而狀態(tài)機跳轉(zhuǎn)次數(shù)也降低，循環(huán)迭代過程所占的延時降低。以A×B，分塊數(shù)4為例，處理前后的數(shù)據(jù)訪問為:

其中，array分別為矩陣A，B的行、列向量;a1～ a4，b1～b4分別為矩陣A，B分塊后的子矩陣。

數(shù)據(jù)空間和迭代空間經(jīng)優(yōu)化后，設(shè)處理器空間的PE單元數(shù)為Sp，矩陣A和B存儲模塊深度為St1，St2。此時加載矩陣A，B，加載并存儲結(jié)果矩陣C帶來的數(shù)據(jù)傳輸量分別為:

可見St1的大小對數(shù)據(jù)傳輸量沒有影響，故不考慮其他因素時可將St1設(shè)為1，即讓A矩陣的存儲模塊退化為寄存器，給St2留更多余量，可增大其值以進一步減少數(shù)據(jù)傳輸耗時，提高運算速度。

3 浮點矩陣乘積加速器設(shè)計

基于多處理器核FPGA/SOPC的導航解算以解算模塊流水線化、模塊內(nèi)各粒度計算單元同時進行的方式［7］完成解算任務(wù)。處理器調(diào)度各運算單元完成操作數(shù)的讀取、運算、結(jié)果存儲。

對于完成浮點矩陣乘積的運算單元，其任務(wù)可分為矩陣計算和對緩沖區(qū)數(shù)據(jù)訪問兩部分。上述兩方面著手提升運算速度，據(jù)此設(shè)計的加速器總體結(jié)構(gòu)如圖1所示。作為導航解算中的協(xié)處理器，加速器采用硬件方法完成復雜耗時的矩陣乘積運算，同時還可以通過Avalon-MM總線與導航系統(tǒng)中其他模塊通信。

各模塊設(shè)計說明如下:

(1)浮點矩陣乘積模塊:基于浮點矩陣乘積加速原理制成的浮點矩陣乘積單元由矩陣A，B存儲模塊、寄存器、高速緩存器、高性能點乘、求和計算單元5個部分組成。其結(jié)構(gòu)如圖2所示。

圖1 浮點矩陣乘積加速器總體結(jié)構(gòu)

圖2 浮點矩陣乘法單元結(jié)構(gòu)

應用中，依據(jù)數(shù)據(jù)空間分割及迭代空間合并理論選擇矩陣分塊參數(shù)vectorsize，block的最優(yōu)值。矩陣乘法宏功能模塊據(jù)此參數(shù)生成相應規(guī)模的DSP block二維陣列完成點乘計算。陣列中DSP block完成18×18乘積和運算并兼具高頻流水線、反饋累加與級聯(lián)功能，使得內(nèi)層迭代僅耗(vectorsize/(2×block))個時鐘周期，再經(jīng)(A矩陣的列數(shù)/vectorsize)次迭代運算，高性能求和計算單元將中間量對應相加，得到結(jié)果矩陣C。另外優(yōu)化為部分寄存器的A矩陣存儲模塊由M144K充當，相比B矩陣帶寬較窄，使用中需要在A當?shù)丶拇嫫髦写娣哦鄠€數(shù)據(jù)以拓展帶寬，使數(shù)據(jù)加載用時不超過處理耗時，達到較高的占空比。

加速器中設(shè)定矩陣分塊參數(shù)，例化生成的硬件描述語言文件，再根據(jù)時序圖將運算過程劃分為數(shù)據(jù)加載、運算開啟、運算處理、結(jié)果輸出4個狀態(tài)，編寫狀態(tài)機控制乘積單元運算，完成浮點矩陣乘積運算。

(2)DMA［11］:FPGA/SOPC 平臺上傳統(tǒng)的Avalon-MM(Avalon Memory Mapped Interface)總線數(shù)據(jù)傳輸方式基于地址映射方式，控制信號繁多，尤其在大批量數(shù)據(jù)傳輸時，CPU中斷載荷負擔重，導致傳輸速率較低。但同樣采用Avalon-MM總線的DMA以單向點對點方式傳輸，用精簡的控制方式及高速流模式進行數(shù)據(jù)傳輸，不會影響CPU的其他工作。大批量數(shù)據(jù)傳輸時具有非常高的傳輸效率。

SOPC中DMA數(shù)據(jù)傳輸依賴于硬件系統(tǒng)搭建和DMA控制器控制。加速器硬件系統(tǒng)需配置處理器、PLL、DMA、外部存儲控制器等組件，設(shè)置DMA長度寄存器位寬、構(gòu)建方式、基址和中斷級，連接讀寫主端口與源、目的地址所在的組件。控制方式采用為Nios II系統(tǒng)提供設(shè)備驅(qū)動程序的HAL庫中功能語句，而非寄存器方式。

(3)組件接口:浮點矩陣乘積模塊只有添加MM總線讀寫從機并進行封裝后才能被Nios II系統(tǒng)調(diào)用，另外DMA傳輸單個浮點數(shù)用時不恒定，而矩陣乘積運算的數(shù)據(jù)加載要求字加載速度恒定在一個時鐘周期，這樣就出現(xiàn)數(shù)據(jù)傳送速率不匹配。故在編寫MM讀寫從機程序的基礎(chǔ)上采用FIFO實現(xiàn)數(shù)據(jù)傳輸與加載間的同步。使用存儲器處理跨域傳輸，避免了握手信號和邏輯同步處理機制在同步設(shè)計中大量的時鐘周期消耗，保證了通信雙方較高的數(shù)據(jù)吞吐率。設(shè)計中編寫 FIFO讀寫時序，使得FIFO能在DMA數(shù)據(jù)線上數(shù)據(jù)正確時將其讀入，當FIFO中暫存的數(shù)據(jù)量達到指定數(shù)目后，再依照加載時序讀出，實現(xiàn)數(shù)據(jù)傳輸與加載間的同步。

浮點矩陣乘積加速器硬件系統(tǒng)設(shè)計完成后，進行軟件部分設(shè)計。軟件部分任務(wù)是控制DMA實現(xiàn)RAM與自定義外設(shè)之間的流模式數(shù)據(jù)傳輸。為保證計算結(jié)果的及時輸出，注冊PIO中斷標志量，使其在浮點矩陣乘積模塊計算結(jié)果有效時置1，并在中斷服務(wù)子程序中完成DMA2到RAM2之間的數(shù)據(jù)輸送。系統(tǒng)主程序流程如圖3所示。

圖3 系統(tǒng)主程序流程

4 實驗分析

4.1 浮點矩陣乘積模塊功能測試

本文以非方陣乘積 A10×8×B8×12為例，對浮點矩陣乘積模塊進行功能仿真測試。選取A，B矩陣vectorsize=8，blocksize=2，testbench 腳本中時鐘信號200 MHz，在Modelsim_altera平臺上進行仿真。仿真測試結(jié)果如圖4所示。

圖4中整個運算過程消耗了1885 ns，即377個時鐘周期。浮點操作數(shù)表達式為:

其中，F(xiàn)Hz表示系統(tǒng)運行頻率;tload，tprocess，tnaive分別表示數(shù)據(jù)加載、處理、簡單算法時間。依據(jù)上述公式及仿真結(jié)果得出:時鐘信號200 MHz時浮點矩陣乘積模塊的浮點操作數(shù)達到了2.09 GFLOPS。由于浮點矩陣乘積單元的峰值頻率可以達到412 MHz［12］，因此峰值計算性能可達到4.30 GFLOPS。

使用浮點矩陣乘積模塊與在PC上使用Matlab計算矩陣乘積的部分結(jié)果對比如表1所示。從運算結(jié)果可以看出，浮點矩陣乘積模塊可以正確完成乘積計算。只是乘積模塊直接計算單精度浮點數(shù)，而Matlab中間過程以雙精度形式進行，結(jié)果轉(zhuǎn)換為單精度形式輸出，從而導致了結(jié)果的微小差異。

圖4 Modelsim仿真測試結(jié)果

表1 部分運算結(jié)果對照

4.2 浮點矩陣乘積加速器性能實測

將浮點矩陣乘積加速器的.sof文件下載至Altera DE3開發(fā)板，在Nios EDA平臺上以硬件方式運行軟件程序，使用嵌入式邏輯分析儀SignalTap II對浮點矩陣乘積加速器的計算性能進行測試，選取CPU的時鐘信號作為采集時鐘，實時采集得到的信號如圖5所示。由實測圖可以看出，浮點矩陣乘積加速器計算時序與功能測試吻合。從數(shù)據(jù)讀取、矩陣運算到結(jié)果輸出一共消耗了842個時鐘周期。在200 MHz時鐘實驗條件下，耗時4.21 μs。各種浮點矩陣乘積方法的性能分析如表2所示，本文方法較其他文獻中硬、軟件實現(xiàn)方法速度分別提升71.3%，78%。

圖5 SignalTap II實測圖

表2 浮點矩陣乘積方法的性能分析

5 結(jié)束語

為提高導航解算中矩陣乘積運算速度，本文從硬件實現(xiàn)角度設(shè)計了應用于Nios II系統(tǒng)的矩陣運算加速器。首先分析了FPGA用于浮點矩陣運算的加速原理，據(jù)此編寫性能優(yōu)越的矩陣乘積模塊并使用Modelsim進行功能仿真驗證。之后設(shè)計DMA傳輸通道及外設(shè)，并將整體封裝成自定義組件作為Nios II系統(tǒng)的協(xié)處理器。為驗證其實際特性，將所有程序下載至Altera DE3開發(fā)板進行實測，利用SignalTap II拾取的信號與其他實現(xiàn)方法對比表明，本文設(shè)計的矩陣運算硬件加法器具有明顯的速度優(yōu)勢，有利于提高導航解算的實時性。

［1］秦永元.慣性導航［M］.北京:科學出版社，2006.

［2］林燦龍，馬龍華，吳鐵軍，等.高動態(tài)條件下的捷連慣導并行算法研究［J］.彈箭與制導學報，2012，32(2):1-5.

［3］李宗濤.高動態(tài)捷連慣導系統(tǒng)的并行實現(xiàn)分析［D］.杭州:浙江大學，2011.

［4］Underwood K.FPGAsvs.CPUs:Trendsin Peak Floating-point Performance［C］//Proc.of International Symposium on Field Programmable Gate Arrays.Monterey，USA:ACM Press，2004:171-180.

［5］Campbell S J，Khatri S P.Resource and Delay Efficient Matrix Multiplication Using NewerFPGA Devices［C］//Proc.of the 16th ACM Great Lakes Symposium onVLSI.Philadelphia，USA:ACM Press，2006:308-311.

［6］田 翔，周 凡，程耀武，等.基于FPGA的實時雙精度浮點矩陣乘法器設(shè)計［J］.浙江大學學報:工學版，2008，42(9):1611-1615.

［7］許 芳，席 毅，陳 虹，等.基于 FPGA/Nios-II的矩陣運算硬件加速器設(shè)計［J］.電子測量與儀器學報，2011，25(4):377-382.

［8］Dutta H，HannigF，Ruckdeschel H，et al.Efficient Control Generation for Mapping Nested Loop Programs onto Processor Arrays［J］.Journal of Systems Architecture，2007，53(5):300-309.

［9］Wolfe M J，Shanklin C，Ortega L.High Performance Compilers for Parallel Computing［M］.Boston，USA:Addison-Wesley Longman Publishing Co.，Inc.，1995.

［10］Asher Y，Rotem N.Automatic Memory Partitioning:Increasing Memory Parallelism via Data Structure Partitioning［C］//Proc.of International Conference on Hardware/Software Codesign and System Synthesis.Scottsdale，USA:IEEE Press，2010:155-161.

［11］Altera Corporation.ALTERA Embedded Peripherals IP User Guide［EB/OL］.(2011-06-10).http://www.altera.com.cn/literature/ug/ug_embedded_ip.pdf.

［12］Altera Corporation.ALTERA Floating-point Megafunctions User Guide ［EB/OL］.(2011-06-10).http://www.altera.com.cn/literature/ug/ug_altfp_mfug.pdf.

［13］雷 晶，金心宇，王 銳.矩陣相乘的并行計算及其DSP 實現(xiàn)［J］.傳感技術(shù)學報，2006，19(3):737-740.