基于FPGA的復(fù)數(shù)長方陣SVD算法

閻 波，張威威，林水生

(電子科技大學(xué)通信與信息工程學(xué)院 成都 611731)

·通信與信息工程·

基于FPGA的復(fù)數(shù)長方陣SVD算法

閻 波，張威威，林水生

(電子科技大學(xué)通信與信息工程學(xué)院 成都 611731)

在OFDM和MIMO系統(tǒng)中普遍使用長方形矩陣復(fù)數(shù)奇異值分解運算。針對傳統(tǒng)算法運算量大，迭代次數(shù)多的問題，提出了一種基于householder和雙邊Jacobi的混合優(yōu)化算法。該算法首先通過householder變換將矩陣化解為二對角矩陣；然后提取2×2復(fù)矩陣；再進行改進型復(fù)數(shù)雙邊Jacobi變換。兼具有QR算法的高精度和Jacobi算法的低硬件實現(xiàn)成本的優(yōu)點。給出了2×8的CSVD的FPGA硬件實現(xiàn)方案并進行了板級測試。測試結(jié)果表明，該混合優(yōu)化算法較傳統(tǒng)算法在硬件資源上節(jié)省26%，延時縮短10倍，在同等位寬下計算精度至少提高了一個數(shù)量級。

復(fù)數(shù)奇異值分解; 可編程邏輯陣列; householder; Jacobi; 長方矩陣

CSVD在OFDM通信系統(tǒng)中的快衰落信道估計[1]、MIMO系統(tǒng)中的預(yù)編碼[2]、信號處理[3]、圖像處理[4]等領(lǐng)域中得到了廣泛的應(yīng)用。CSVD在工程中應(yīng)用的關(guān)鍵在于數(shù)值計算，為了保證計算精度，通常會采用迭代運算的方法。文獻[5]給出了一種直接復(fù)數(shù)雙邊Jacobi方法，該方法以一個復(fù)數(shù)22×為基本單元，并將其對角化，避免擴展矩陣維度。然而一個基本單元算法復(fù)雜，耗費的硬件資源較多，也需要多次迭代。文獻[6]提出一種復(fù)轉(zhuǎn)實的雙邊Jacobi方法，該方法將矩陣行列維度均擴大一倍，以一個實22×為基本單元，優(yōu)點是基本單元簡單，但是迭代次數(shù)增大，資源、延時和精度受限，實時性差，且這些算法對于長方矩陣具有明顯的缺點。

本文提出了基于householder和雙邊Jacobi的混合優(yōu)化算法，能有效地解決長方復(fù)矩陣奇異值分解的硬件實現(xiàn)和精度等問題。首先介紹了CSVD傳統(tǒng)的算法即QR算法和雙邊Jacobi算法，并指出其對于處理長方陣的局限性。在此基礎(chǔ)上提出一種混合優(yōu)化算法，給出其具體運算步驟，并將該方法與現(xiàn)有常見算法性能比較，具體說明混合優(yōu)化算法在硬件實現(xiàn)上的優(yōu)勢。最后以28×復(fù)矩陣為例，給出了混合優(yōu)化算法在FPGA上的實現(xiàn)及其測試結(jié)果。

1 復(fù)數(shù)矩陣奇異值分解

設(shè)M是一個mn×的復(fù)矩陣，如果存在一個分解有：

式中，U是m×m酉矩陣；S是半正定mn×對角矩陣；V是n×n酉矩陣；VH為V的共軛轉(zhuǎn)置。這樣的分解就稱為復(fù)矩陣M的奇異值分解。S對角上的值為M的奇異值。在工程實際中，只需要求得V矩陣即可。

1.1 傳統(tǒng)長方陣CSVD算法

目前傳統(tǒng)復(fù)數(shù)長方陣的奇異值分解主要有兩種較為常用：QR迭代算法和雙邊Jacobi算法。QR迭代算法首先采用householder變換將一般復(fù)矩陣化為二對角矩陣；然后采用Givens變換的方法交叉將非對角非零元素消零，如圖1所示，,ST分別代表行和列Givens矩陣，這一過程又稱為“驅(qū)逐出境”；經(jīng)過多次迭代將非對角化為零，最終達到收斂條件。

圖1 QR迭代“驅(qū)逐出境”示意圖

雙邊Jacobi算法常見有：復(fù)轉(zhuǎn)實雙邊Jacobi算法和直接復(fù)數(shù)雙邊Jacobi算法。兩種算法的本質(zhì)都是采用Brent-Luk-Van Loan(BLV)[7]脈動陣列迭代算法。所不同的是前者是將n×n復(fù)矩陣轉(zhuǎn)換為2×N2n的實矩陣，以2×2矩陣為基本單元；而后者是直接對復(fù)2×2矩陣進行對角化。

對于復(fù)轉(zhuǎn)實算法，首先將一般復(fù)矩陣M轉(zhuǎn)化為共軛對稱矩陣C，即C=MHM，令C=A+Bi ，(u+iv)為C的奇異值σ所對應(yīng)的奇異向量，則有：

將復(fù)矩陣轉(zhuǎn)成實矩陣后，以一個實2×2為基本單元。

對于直接復(fù)數(shù)雙邊Jacobi算法，以一個復(fù)數(shù)為基本運算單元進行對角化。其運算過程主要為以下兩步：

這兩種算法各有優(yōu)缺點：QR迭代算法計算復(fù)雜度低，消耗更少的乘除法等硬件資源，但不穩(wěn)定，當矩陣維度較大時，使得某次Givens矩陣為單位陣，“驅(qū)逐出境”會出現(xiàn)不收斂情況[8]，導(dǎo)致QR算法失敗；雙邊Jacobi[9]算法是基于脈動陣列結(jié)構(gòu)的迭代算法，優(yōu)點是結(jié)構(gòu)簡單對稱，易于硬件實現(xiàn)，缺點是計算量比QR大，且精度與迭代次數(shù)相關(guān)。

針對上述兩種算法的優(yōu)缺點，本文提出一種混合優(yōu)化算法，即混合householder和雙邊Jacobi的計算方法，優(yōu)化了傳統(tǒng)復(fù)數(shù)2×2復(fù)數(shù)雙邊Jacobi的計算方法。

1.2 CSVD混合優(yōu)化算法

設(shè)vec_in=(x1,x2,,xn)，則該向量的householder變換[10-12]步驟如下：

1) 計算輸入向量平方和得：

2) 構(gòu)造行向量u，并將其單位化得：

3) 構(gòu)造householder變換矩陣為：

4) 向量變換為：

對于矩陣的householder變換，則需要逐行逐列調(diào)用向量householder變換，并在式(11)用變換矩陣H右乘，更新其他行向量。對于2×4n的復(fù)矩陣只需進行行變換，將一般矩陣化簡為形如式(13)的二對角矩陣Mk，并提取左上角非零復(fù)數(shù)2×2矩陣M，

則有：

對M進行改進型復(fù)數(shù)雙邊Jacobi運算，其計算步驟如下：

1) 將復(fù)數(shù)22×矩陣轉(zhuǎn)成共軛對稱復(fù)矩陣：

2) 將Msym(1,2)轉(zhuǎn)成幅角形式：

3) 將復(fù)矩陣化為實矩陣：

4) 計算雙邊Jacobi變換旋轉(zhuǎn)角為：

5) 雙邊Jacobi變換為：

6) V矩陣計算為：

2×8CSVD的Vsvd矩陣計算：將2×2復(fù)數(shù)V矩陣補充到8×8單位陣的1,2行列得到復(fù)酉矩陣V1，然后將householder產(chǎn)生的H1，H2與V1矩陣相乘即得到最終輸出酉矩陣Vsvd，即Vsvd=H1×H2×V。

2 3種算法的比較與分析

以2×4復(fù)矩陣為例，比較直接復(fù)數(shù)雙邊Jacobi(算法1)，復(fù)轉(zhuǎn)實雙邊Jacobi(算法2)和混合優(yōu)化算法(算法3)在硬件資源消耗、迭代次數(shù)、誤差精度之間的差異。

表1給出了3種算法的性能比較，資源消耗為單次迭代次數(shù)下的運算量。圖2的精度比較結(jié)果是基于MATLAB中CORDIC函數(shù)浮點仿真結(jié)果。

圖2 復(fù)矩陣2×4在3種算法下平均絕對誤差精度比較

表1 復(fù)矩陣2×4奇異值分解3種算法性能比較

更為一般地，對于2×N或2×N的復(fù)長方矩陣，算法2和算法3的復(fù)雜度比較如圖3所示。

在圖3中將CORDIC單元、4個實數(shù)乘法統(tǒng)一等效為復(fù)數(shù)乘法。算法3中雙邊Jacobi運算在算法一基礎(chǔ)上進行了優(yōu)化，極大地減少硬件資源消耗。在單次迭代的情況下，算法2的資源消耗略低于算法3，通過硬件結(jié)構(gòu)上的優(yōu)化，可以使兩種算法在資源消耗上相當。但是增加迭代次數(shù)即增加了處理延時，導(dǎo)致矩陣吞吐率(單位時間處理矩陣個數(shù))急劇下降。且隨著矩陣列數(shù)的增長，算法2和算法3之間的資源消耗差異將立方增長。

表1中迭代次數(shù)并不固定，理論上可由BLV方法總結(jié)公式得到，即(log2N+1)(N?1)，而復(fù)轉(zhuǎn)實的方法為(log2(2N)+1)(2N?1)，N為復(fù)方陣階數(shù)。

圖3 算法2、算法3復(fù)雜度比較

計算誤差的方法如下：設(shè)M=USVH，消去U，有error=MHM?VSHSVH，分別取誤差矩陣error各元素的實部和虛部絕對值并取其平均。圖2中算法1和算法3精度高于算法2。誤差的產(chǎn)生主要取決于CORDIC核，當增加迭代次數(shù)后，使得在對復(fù)數(shù)2×2對角化時非對角元素極小，此時θb應(yīng)作0或π處理。而CORDIC核仍求兩極小數(shù)的比值來計算角度，進而導(dǎo)致計算的角度bθ誤差增大，最終導(dǎo)致輸出V矩陣誤差增大。特別是針對24n×的復(fù)矩陣，算法2需要先取其共軛轉(zhuǎn)置并與其本身相乘，再擴展為維度為88n×n的實矩陣，該矩陣至少包含8n個零元素。在迭代過程中容易產(chǎn)生極小的非對角元素。當然可以采取設(shè)置門限的方法，若非對角元素絕對值小于某一常數(shù)，直接置零。由此，對于奇異值分解的各種算法中迭代算法具有一定的不穩(wěn)定性。

對于相同維度的長方復(fù)矩陣，算法3在資源和精度上要好于算法1和算法2，由于householder變換涉及乘、除、開方，且步驟并不像其他兩種算法規(guī)則。算法3在實現(xiàn)的復(fù)雜度上要略高于其他兩種算法。

3 2×8 CSVD的硬件實現(xiàn)結(jié)構(gòu)

CSVD的硬件實現(xiàn)主要包括兩個模塊：向量householder變換和改進型22×復(fù)數(shù)雙邊Jacobi模塊。這兩個模塊的計算是順序的，因此可以設(shè)計為流水線結(jié)構(gòu)。根據(jù)具體FPGA芯片資源的多少，可以靈活選擇全流水或部分復(fù)用實現(xiàn)結(jié)構(gòu)。

3.1 向量householder變換的實現(xiàn)結(jié)構(gòu)

一方面充分考慮資源的可復(fù)用程度，合理調(diào)配運算順序；另一方面設(shè)計的結(jié)構(gòu)應(yīng)盡可能簡單、模塊化。圖4設(shè)計了一種全流水的行householder變換結(jié)構(gòu)，對于大型矩陣可以復(fù)用該模塊。

圖4 行householder變換流水線硬件框圖

3.2 改進型2×2復(fù)數(shù)雙邊Jacobi的實現(xiàn)結(jié)構(gòu)

本文提出的復(fù)數(shù)22×雙邊Jacobi方法是在文獻[5,7]的基礎(chǔ)上針對硬件結(jié)構(gòu)提出的改進。文獻[5,7]方法對于一般復(fù)矩陣計算量偏大，特別是對于FPGA、DSP等硬件平臺。但如果在對角化之前，先將復(fù)矩陣轉(zhuǎn)換為共軛對陣矩陣，則計算量會大大減小。額外的資源消耗僅僅是2個復(fù)數(shù)22×矩陣相乘。圖5的硬件框圖相比于文獻[5]的硬件框圖更省資源。由于不用迭代運算，在精度上可以考慮減少乘法、CORDIC運算的位寬。如Xilinx平臺下，一個實數(shù)乘法用1825×的位寬，可以保證最大化利用乘法器。其他運算以此為基準進行定點。圖5中CORDIC的3種模式均通過Xilinx的CORDIC IP核實現(xiàn)，本文將其設(shè)置為全并行模式，以保證整體流水線設(shè)計。

圖5 復(fù)數(shù)2×2對角化流水線硬件框圖

3.3 基于FPGA的實現(xiàn)驗證

對于輸出V矩陣維度較大的復(fù)數(shù)奇異值分解，會耗費更多的FPGA資源如乘法器、除法器和CORDIC核。當DSP資源超過50%以后，布局布線后的最高時鐘頻率較綜合后的最高時鐘頻率會大幅下降，主要原因在于乘法器核在布局布線過程中會產(chǎn)生較大的線延時，可通過減小乘法器輸入輸出的扇出解決，也可以通過更改綜合工具設(shè)置。由布局布線后的結(jié)果才能比較準確反映設(shè)計的準確性和可靠性。本文設(shè)計采用Xilinx的FPGA硬件實現(xiàn)方案，其型號為xc6vlx240t-3ff1156，基本滿足設(shè)計要求。表2給出了算法2和算法3的資源占用比較，算法2實現(xiàn)平臺為Altera Stratix IV。算法2為迭代結(jié)構(gòu)，時鐘頻率為105 MHz，延時為3 800個時鐘；而算法3為全并行，流水結(jié)構(gòu)，布局布線后的最大時鐘頻率200.240 MHz，延時為330個時鐘。算法3在資源和吞吐率上均優(yōu)于算法2。

表2 復(fù)數(shù)2×8奇異值分解FPGA資源占用

4 結(jié) 束 語

本文提出了一種主要針對24n×或42×N的CSVD混合優(yōu)化算法，若矩陣維度為奇數(shù)，需要將矩陣維度向上擴充至偶數(shù)。該算法通過對多種傳統(tǒng)算法的部分運算整合、改進，極大地減小了CORDIC核的使用且不需要迭代。通過與傳統(tǒng)算法對比，該算法比傳統(tǒng)算法至少在資源上節(jié)省26%，延時縮短10倍，精度提高一個數(shù)量級。最后對28×CSVD進行FPGA實現(xiàn)，算法上的優(yōu)勢在硬件上得以體現(xiàn)。下一步將對更高維度的長方陣CSVD作探討。

[1] AU E K S, JIN S, MCKAY M R, et al. Analytical performance of MIMO-SVD systems in ricean fading channels with channel estimation error and feedback delay[J]. IEEE Transactions on Wireless Communications,2008, 7(4): 1315-1325.

[2] SRINIVASAN J, RAJARAM S. FPGA implementation of precoding using low complexity SVD for MIMO-OFDM systems[C]//Information Communication and Embedded Systems (ICICES). [S.l.]: IEEE, 2013.

[3] CHAKROBORTY S, SAHA G. Feature selection using singular value decomposition and QR factorization with column pivoting for text-independent speaker identification [J]. Speech Communication, 2010, 52(9): 693-709.

[4] 胡謀法, 董文娟, 王書宏, 等. 奇異值分解帶通濾波背景抑制和去噪[J]. 電子學(xué)報, 2008, 36(1): 111-116. HU Mou-fa, DONG Wen-juan, WANG Shu-hong, et al. Singular value decomposition band-pass-filter for image background suppression and denoising[J]. Acta Electronica Sinica, 2008, 36(1): 111-116.

[5] WANG Y, CUNNINGHAM K, NAGVAJARA P, et al. Singular value decomposition hardware for MIMO: State of the art and custom design[C]//Reconfigurable Computing and FPGAs (ReConFig). [S.l.]: IEEE, 2010.

[6] HAN Q, ZENG L. FPGA Implementation for low-rank channel estimation of OFDM[J]. Journal of Networks, 2012, 7(10): 1631-1638.

[7] HEMKUMAR N D, CAVALLARO J R. A systolic VLSI architecture for complex SVD[C]//Circuits and Systems, ISCAS'92. [S.l.]: IEEE, 1992.

[8] 趙學(xué)智, 葉邦彥. 單向收縮QR算法在奇異值分解中的收斂特性[J]. 電子科技大學(xué)學(xué)報, 2010, 39(5): 762-767. ZHAO Xue-zhi, YE Bang-yan. Convergence characteristic of single direction shrink QR algorithm in the singular value decomposition[J]. Journal of University of Electronic Science and Technology of China, 2010, 39(5): 762-767.

[9] MA W, KAYE M E, LUKE D M, et al. An FPGA-based singular value decomposition processor[C]//Electrical and Computer Engineering. [S.l.]: IEEE, 2006.

[10] LIU J, ZHANG J. A new maximum simplex volume method based on householder transformation for endmember extraction[J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(1): 104-118.

[11] PEDRAM A, GERSTLAUER A, GEIJN R A V D. Floating point architecture extensions for optimized matrix factorization[C]//Proceedings of the 2013 IEEE 21st Symposium on Computer Arithmetic. [S.l.]: IEEE, 2013: 49-58.

[12] 張賢達. 矩陣分析與應(yīng)用[M]. 北京: 清華大學(xué)出版社有限公司, 2004. ZHANG Xian-da. Matrix analysis and applications[M]. Beijing: Tsinghua and Springer Publishing House, 2004.

編輯稅 紅

Singular Value Decomposition Algorithm of Rectangular Complex Matrix Based on FPGA

YAN Bo, ZHANG Wei-wei, and LIN Shui-sheng
(School of Communication and Information Engineering, University of Electronic Science and Technology of China Chengdu 611731)

Rectangular matrix complex singular value decomposition (CSVD) is widely used in orthogonal frequency division multiplexing (OFDM) and multiple input and multiple output (MIMO) systems. In view of large iteration computation of traditional algorithms, a householder and Jacobi based mixed optimized algorithm is proposed which diagonalizes a general complex matrix and carry out an improved complex two-sided Jacobi transform. This method combines the advantages of high precision of QR and the simple hardware structure of Jacobi. A 2×8 CSVD design is implemented on field programmable gate array (FPGA) by using MATLAB simulation and Xilinx platform. Compared with traditional algorithms, the mixed optimized algorithm saves 26% hardware resources, shortens delay time by 10 and improve the accuracy of calculation at least one order of magnitude under the same bit width.

CSVD; FPGA; householder; Jacobi; rectangular complex matrix

TP33; TN4

A doi:10.3969/j.issn.1001-0548.2015.04.001

2014 ? 03 ? 11 ；

2015 ? 03 ? 11

國家自然科學(xué)基金(61301155,61176025)；中央高校基本科研業(yè)務(wù)費專項資金(ZYGX2012J003)

閻波(1973 ? )，女，教授，主要從事通信信號處理，無線通信系統(tǒng)、通信集成電路設(shè)計等方面的研究.