摘要:在矩阵的奇异值分解(singular value decomposition,SVD)过程中,随着矩阵维数的增加,SVD的计算量呈指数型增长,从而降低了算法运行的实时性。针对这个问题,基于Hestenes-Jacobi数值计算方法,提出了一种改进的基于坐标旋转数字计算机(CORDIC)的逻辑设计,该逻辑设计采用并行的全流水线设计思想,能够提高Jacobi平面旋转变换的运行速度,进而加快任意维矩阵奇异值分解的计算速度。分析了基于Hestenes-Jacobi方法的SVD的数值计算过程,介绍了CORDIC算法的基本原理,并具体说明了基于CORDIC算法的Jacobi平面旋转模块的设计,利用Verilog语言实现设计并验证,在现场可编程门阵列(FPGA)上运行该逻辑设计单元,与Matlab软件的运行结果进行对比。实验测试结果表明,该结构能够减少计算时间,适应高速数据处理的要求。
加入收藏
引言
奇异值分解(SVD)是矩阵分析中对矩阵进行对角化的数值方法,是线性代数中最有效且最有用的工具之一。在许多科学和工程领域中,SVD的计算都是必不可少的。在通信信号处理、信号与图像处理以及应用越来越广泛的大数据应用中,SVD算法都有着相当重要的研究与应用价值[1,2,3,4]。
在过去的几十年中,人们一直致力于提高SVD的处理速度[5,6]和减少设计的资源消耗[7,8],如何快速高效地实现SVD已成为了一个重要的研究方向。文献[9]提出在图形处理器(GPU)上实现SVD,以提高处理速度。但是,利用GPU实现SVD存在线程同步和内存读取等问题,这会导致消耗大量的能量。文献[10]中提出了一种基于脉动阵列的架构用于SVD。但是,这种架构调度复杂,且仅适用于方阵,资源利用效率不高。
本文在已有研究的基础上,提出了一种基于坐标旋转数字计算机(CORDIC)算法的改进SVD硬件架构。该架构实现了基于定点数的Hestenes-Jacobi算法,可用于任意维矩阵的SVD。通过对Hestenes-Jacobi方法中旋转变换过程的分析,结合全流水的CORDIC运算单元,对该方法的实现过程进行改进,在提高处理速度的同时降低了资源消耗。
1、理论基础
1.1奇异值分解理论
一个m×n的矩阵的奇异值分解可以表示为
Am×n=Um×m×Σm×n×VTn×n (1)
(1)式中:U和V都是正交矩阵,它们分别满足UUT=I,VVT=I,I是单位矩阵;Σ是对角矩阵,其对角线上的元素σi(i=1,2…n)为矩阵A的奇异值,为非负数。
1.2Hestenes-Jacobi方法
Hestenes发现,通过利用平面旋转变换的方法使2个向量正交与将矩阵的元素置为零这二者之间存在一定的等价关系[11]。在基于此项研究发现的基础上,Hestenes提出了一种单边Jacobi方法,也被称作Hestenes-Jacobi方法。
单边Jacobi方法的核心思想是通过对维度为m×n的矩阵A采用一系列的Jacobi平面旋转变换[12],对其进行正交化
B=A(J1J2J3⋯) (2)
使得矩阵B中任意2列向量都满足关系biTbj=0。
Jacobi平面旋转变换如(4)式矩阵的结构
公式 (3)
(3)式中:(i,j)表示消去元素在矩阵中的位置;c表示cosθ;s表示sinθ,θ称为旋转角。不难看出,JTJ=I。因此,Jacobi平面旋转变换为正交变换。因为Jacobi平面旋转变换矩阵本身即正交矩阵,令V=J1J2J3…,因此,V也是正交矩阵,即VT=V-1。
正交化完毕后,对得到的矩阵B进行归一化,即
B=BΣ−1Σ=UΣ (4)
(4)式中,Σ=diag(σ1,σ2…σn-1,σn),σi2=biTbi。
所以,矩阵A的奇异值分解可表示为
A=UΣVT (5)
在对矩阵A进行正交化时,每一次的旋转变换过程中只有i和j2列元素受到影响。因此,一次正交变换可表示为
公式(6)
对2列向量进行正交化,即向量biT和bjT正交
公式 (7)
经等式变换后,即可得到
公式 (8)
根据(8)式,可求得旋转角度值,即可通过该角度值计算相应的正余弦函数值,对矩阵列向量对元素实现正交变换。
1.3CORDIC算法
CORDIC算法的思想最早在1959年被提出,其主要是通过将初始向量按照一个特定的角度序列不断旋转迭代以逼近预期的目标向量,每次旋转变换之后,对坐标值进行更新[13]。当迭代次数足够多时,就认为当前向量的横纵坐标值是所求目标角度的正余弦函数值输出,图1即为平面坐标旋转模型。
图1平面旋转变换
图1平面旋转变换下载原图
Fig.1Planerotationtransformation
如图1,直角坐标系中的点(x0,y0)逆时针转动θ角度,即可得另一个点(x1,y1),两坐标点的关系为
x1=x0cosθ−y0sinθy1=x0sinθ+y0cosθ (9)
通过提取并去除因数cosθ,方程可写成下面的形式,即
xˆ1=x0−y0tanθyˆ1=y0+x0tanθ (10)
(9)式的变换被称为伪旋转变换,其旋转的角度是正确的,但是变换前后向量的模值增加了cos-1θ倍。
CORDIC方法的核心是旋转角度θ的变换[14],通过将旋转变换的角度θ分解为一系列旋转角度ai的和,即∑i=0n−1diαi(di=±1,i=0,1,2,⋯n),其中,角度ai=arctan(2-i),完成角度θ的旋转逼近。引入第3个方程,称为角度累加器,用来在每次迭代过程中追踪累加的旋转角度。在每次旋转后,(10)式的变换关系变为
公式 (11)
(11)式中,符号di是一个判决算子,其值为±1,用于确定旋转的方向是顺时针还是逆时针。判决算子值的确定见表1。选择不同的工作模式,即判决算子的确定方式,可以得到不同的函数值的计算。
表1CORIC算法的操作模式
将CORDIC变换中的起点(x0,y0)设为(aiT,ajT)中的元素,终点(x1,y1)设为(biT,bjT)中的元素时,对比(6)式和(9)式可发现,CORDIC变换可以实现Jacobi平面旋转变换,因此,可以用基于CORDIC变换的硬件来实现Jacobi平面旋转变换,实现矩阵奇异值分解。
2、奇异值分解架构设计
本文奇异值分解的Hestenes-Jacobi算法具体流程如算法1。Hestenes-Jacobi方法是一种迭代方法,通过对每次迭代后任意2列向量的内积值,判断算法的收敛性,当达到设定的收敛阈值时,即认为算法达到收敛,TOL为设定的收敛阈值。由于矩阵A中任意2列向量进行正交化时,与其余列向量之间并无数据依赖关系。因此,算法1中Addr函数为控制多对列向量索引的生成,实现并行运算。Addr函数以后步骤即为正交化过程,根据(8)式求解正交化过程的角度值。利用求得的角度值,根据(6)式,对列向量对进行正交化运算。
算法1Hestenes-JacobiMethod。
其中,输入变量即为需要分解的矩阵A,首先计算旋转角度θ,经过多次正交旋转后,得到输出值,分别为左奇异矩阵U及右奇异矩阵V。
2.1奇异值分解
具体的奇异值分解实现电路的顶层设计如图2,图2中的COM_Top模块为奇异值分解的Hestenes-Jacobi正交化变换单元,实现列对向量的正交化运算过程,本例设计了4个COM_Top模块,进行并行运算。电路采用2个双端口RAM单元,分别用于存储左奇异向量和右奇异向量,Addr单元即对应算法1中Addr函数,并行读取矩阵A的列向量,输出作为RAM单元的地址。Converge模块为收敛判决模块,用于判断运算是否达到预期的收敛精度。Norm模块为归一化模块,用于对正交变换后的矩阵进行归一化,依次输出奇异值。状态机作为控制单元,提供同步和控制信号,在整个分解过程中每个模块执行不同的任务。
图2奇异值分解电路顶层设计
Hestenes-Jacobi方法的核心思想是实现一系列的列向量对的正交化变换。因此,本文主要是在对COM_Top模块以及其子模块的具体实现进行改进。因此,下面针对COM_Top的具体设计进行说明。COM_Top模块由3个部分组成:向量乘积运算模块,旋转角度运算模块以及旋转变换运算模块,如图3。
图3COM_Top模块硬件结构
COM_Top模块的运算过程如下:当初始矩阵数据存入双口RAM中完毕后,数据在状态机的控制下进入向量乘积运算运算单元,计算求得2个列向量的内积值,将向量的内积值输出至旋转角度运算单元。根据(8)式可知,旋转变换矩阵的三角函数角度值可通过向量的3个内积值计算得到。在计算得到角度值θ后,直接将角度值以及需要正交化的列向量送入旋转变换运算单元进行正交化。
输入分别为矩阵A的列向量以及单位向量。通过对单位向量的旋转变换对应式V=J1J2J3…,最后得到矩阵V。由于旋转变换过程相同,因此,本文仅以对矩阵A的变换过程为例进行说明。
2.2向量乘积运算单元
图4中,向量乘积模块主要是对矩阵列对向量的内积进行运算,分别为aiTai,aiTaj和ajTaj。
图4向量乘积运算模块硬件结构
当状态进入向量乘积运算时,Addr单元为双端口RAM生成相应的地址,以向量乘积运算提供操作数据,其提供的数据分别是矩阵A的第i列和第j列向量。在状态机的控制下,依次送入列向量的每个元素进行乘积和累加运算。在一组列向量完成运算后,输出内积值。
2.3CORDIC运算单元
图5为流水线方案的CORDIC硬件架构图。为达到设置流水线方案,在每一级迭代运算结束后插入寄存器单元,寄存中间的运算结果,供下一级运算使用,以提高整个模块的运行效率。CORDIC算法的用处很多,可以应用在众多运算单元中。
COM_Top模块中旋转角度运算模块和旋转变换运算模块则是以CORDIC算法为基础,如表1,根据具体的需求,可以设计成相应的运算单元。
旋转角度运算模块作用为计算正交化变换的角度值θ。本文通过利用CORDIC算法实现求解反正切函数运算单元,求解角度值。根据表1,具体的设置如下:当选择操作模式为矢量模式时,即设置角度z0的值为0,x0的值设置为aiTai-ajTaj,y0的值设置为2aiTaj,即可实现求解角度值θ。
旋转变换运算模块即为旋转变换运算单元,实现列向量对的正交化。由(6)式和(9)式可知,列向量的正交化变换过程与CORDIC算法的旋转变换运算过程保持一致。因此,针对旋转变换运算模块,本文采用CORDIC运算单元进行实现。由表1可知,具体的设置如下:选择操作模式为旋转模式,设置初始角度z0为θ,x0和y0分别为(i,j)列向量对中对应的元素。通过对列向量中所有元素进行运算,即可实现(i,j)列向量对的正交化。
通过利用全流水线方案的CORDIC运算单元实现旋转角度运算模块以及旋转变换运算模块,避免了算法1中求解正余弦函数之后的乘法运算过程,降低了资源消耗。同时,由于减少了乘法运算过程,大大降低了算法的时延,提高了算法的运行效率。
图5CORDIC模块硬件结构
3、仿真分析
本文分别在中央处理器(CPU)和FPGA2种不同的硬件平台上实现Hestenes-Jacobi算法。FPGA采用在AlteraDE2-115开发板的CycloneIVEP4CE115F29C7N上实现,工作频率为250MHz。而CPU选择的是Intel(R)CoreTMi5-8500处理器,主频为3.0GHz,利用MATLAB软件进行奇异值分解运算,与本文设计的电路结构性能进行对比。选择用于测试的矩阵维数为8×8,数据用32位定点数表示,且测试数据全部由MATLAB生成。为了验证模块的正确性,利用ModelSim软件作为仿真验证平台,对COM_Top运算单元中的主要模块进行仿真与分析。
3.1向量乘积运算模块功能仿真
本次设定的矩阵A维数为8×8,因此,针对向量乘积运算模块的仿真验证,通过编写test_benth测试文件,将矩阵A中2个8×8的列向量中的元素依次送入向量乘积运算模块进行运算,对该运算单元的功能测试结果如图6。
图6向量乘积运算单元仿真
图6中,A和B分别表示为2个列向量的元素,Vector_MulAA,Vector_MulAB和Vector_MulBB则为计算得到的列向量的内积,done表示为该组列向量对的内积计算完成,否则即为内积运算的中间值。输入数据和输出结果之间有一个时钟周期延时,且在该组列向量对计算完成后,直接开始下组列向量的内积运算。
3.2CORDIC运算模块功能仿真
针对CORDIC运算单元的验证,本文选择矢量模式,实现求解反正切函数,对其进行仿真、验证。本次设置的流水线深度为12级,设置输入为正切函数值,输出即为求得的角度值。仿真结果如图7。
在图7中,y_in是输入的函数值信号,result即为输出信号:角度值,而result_valid表示的是输出信号有效的状态信号。result信号的第1个输出值与y_in的第1个输入数据之间间隔12个周期的延时,第2个输出结果与第1个输出之间则存在1个时钟延时。分析其原因,因为设置的流水线深度为12级,所以运行完整条流水线需要12个周期。同时,由于流水线展开的CORDIC算法,大大提高了模块的整体运算效率。
3.3奇异值分析
完成基于FPGA的奇异值分解硬件设计及仿真后,对该设计进行奇异值分解得到的奇异值与MATLAB进行分解运算得到的奇异值进行对比,结果如表2。
图7CORDIC运算单元仿真
表2奇异值分解计算结果
从表2中可看出,基于FPGA硬件结构分解得到的奇异值与在MATLAB中运算得到的奇异值并不完全相同,分析其原因是由于数据的精度不同造成的。由于本文使用的是32位的定点数,数据在浮点-定点转换以及在进行乘法运算时,由于定点截断导致运算存在误差,但误差相较小,可忽略不计。
3.4系统运行时间分析
表3列出了同一算法在FPGA平台以及CPU平台下奇异值分解算法的运行时间,尽管FPGA的工作频率远低于CPU的主频,FPGA的运算耗时却远小于CPU的耗时,大大提高了算法的运行速率。
表3SVD运行时间对比
4、结论
针对奇异值分解过程中计算复杂、实时性差的问题,提出了一种改进的基于CORDIC算法的实现方案。通过对平面旋转变换运算进行分析,提出利用全流水线的CORDIC运算单元实现,去除乘法运算过程,实现列向量的旋转变换。在此基础上进行实验验证,实验结果验证,该设计方案合理、有效,具有一定的工程应用价值。
参考文献:
[1]郭强.并行Jacobi方法求解矩阵奇异值的研究[D].苏州:苏州大学,2011.
[2]张倩倩,詹永照,林涵阳.基于分块和NSCT-SVD的彩色图像数字水印算法[J].江苏大学学报(自然科学版),2017,38(5):582-589.
[3]余成波,张林,龙曦.基于FPGA数字PLL谐振频率的跟踪研究[J].重庆理工大学学报(自然科学),2019,33(4):141-146.
[4]邓翔宇,刘增力.基于改进的MCA和K-SVD的图像稀疏表示去噪算法[J].四川大学学报(自然科学版),2016,53(4):774-780.
[9]唐吉卓.基于GPU平台的SVD并行计算研究与实现[D].成都:电子科技大学,2014.
[12]马亚峰.基于FPGA的矩阵奇异值分解加速方案的设计与实现[D].北京:北京交通大学,2017.
应俊,朱云鹏.基于CORDIC矩阵奇异值分解的FPGA实现[J].重庆邮电大学学报(自然科学版),2020,32(03):434-440.
基金:国家自然科学基金(61771082,61871062,61801065)
分享:
1、基于核心素养视角的线性代数教学改革研究与实践2、三角形网格上形态开闭算子的研究3、浅谈线性代数教学中的直觉思维的培养4、工科数学课程群建设5、培养综合应用能力的离散数学教学改革研究6、疫情期间线性代数课程线上教学模式探索7、线性代数课程应用能力的培养途径和方法研究——以软件工程专业的教学实践为例8、基于独立学院的高等代数教学研究
2020-08-11“线性代数”是工科类本科生的专业基础课程,其主要处理线性关系问题。课程的内容主要包括基础概念—行列式、矩阵及其运算;基本变换—矩阵的初等变换与线性方程组,向量组的线性相关性;以及线性问题分析—相似矩阵及二次型等内容。采用MATLAB解决矩阵特征值问题的相关教学实践论文已有报道[3],此外,将微课应用到线性代数的教学中也取得了良好的效果[4]。
2020-08-10独立学院,也常称为独立二级学院,是由公办普通本科院校与社会力量联合举办的相对独立的二级学院.独立学院数量众多,招生人数众多,招生录取分数线普遍较低.为了适应社会发展对人才的需求,绝大多数独立学院的人才培养定位为“应用型”人才.在开办之初,很多独立学院是照搬母体高校的课程体系和教学模式,但是由于独立学院学生的基础较母体高校的学生之间存在着较大差异,教学上很难达到独立学院人才培养的要求,其中数学课程教学中出现的问题尤为突出.
2020-08-102020年春季学期,哈尔滨理工大学线性代数课程网络教学如期开课,哈尔滨理工大学理学院数学系线性代数课程教学团队结合学习通平台、微信群和腾讯QQ与学生保持互动学习,并建立基于网络、手机、电脑的教学新模式,基于学校自建课程进行混合式教学,立足本校学生学情和需要,使线上线下的混合式教学很好地衔接。教学模式信息化的改革,增添了教学活力,促进了教学互动,并且较高效率地完成教学任务,达到了教学大纲要求的教学效果。
2020-08-10外汇储备作为我国关键形式的储备资产,同时又是我国基础货币投放的主要渠道。自1994年我国实行以市场供求为基础、有管理的浮动汇率制度以来,我国外汇储备规模大幅扩大,从2006年开始,我国外汇储备规模位列世界第一,外汇储备的增加有提升国家声誉,吸引外商投资和拉动国民经济等积极影响。
2020-08-04在线性代数教学中,求矩阵的幂、判定矩阵对角化、求解特征值的反问题、判定矩阵合同关系以及判定实二次型的正定性等问题都可以借助于矩阵特征值与特征向量实现。本文对这些问题进行了归纳与分析,以便学生能够熟练掌握求矩阵特征值与特征向量的方法,熟悉特征值与特征向量在线性代数教学中的应用。
2020-06-28数学作为工科院校的一门基础课,对学生学好专业课起着至关重要的作用。通常情况下,工科院校开设的数学课程主要包括高等数学、概率与数理统计、线性代数、积分变换等,在这些课程中,普遍存在“教师难教,学生难学的现象,特别是线性代数,表现得更为突出。文章从教师的教学策略和学生的学习方法两个方面提出了工科线性代数教学策略。
2020-06-28目前,关于连续函数空间中基于框架的构造、性质、函数空间刻画及其在信息处理中的应用已经得到了完善与发展。然而,在实际应用中,由于输入/输出数据和滤波器都是离散的,所以基于框架的算法实现都是在数字环境中完成的。基于数字的框架又名为离散空间中的框架,并不一定能够通过连续环境中的框架离散化得到,即使可以得到,也必须附加许多条件。
2020-06-28众所周知,代数闭域K上的Sylvester矩阵方程AX-XB=C有唯一解当且仅当矩阵A和矩阵B无公共特征向量。为深入讨论,利用Sylvester算子研究了Sylvester矩阵方程在任意域K上的可解性,得到其有唯一解当且仅当A和B的特征多项式无公共素因式的结论。在Sylvester矩阵方程有解的情况下,给出了其多项式解。
2020-06-28广义柔度矩阵可利用低阶模态参数近似表示,因此受到广泛关注.本文基于广义柔度矩阵引入广义模态柔度矩阵,利用Nelson方法计算其对损伤参数的敏感性.本文提出的损伤定位方法,相比于现有方法,在振型数据不完整及含有噪声的情况下,该方法具有使用模态参数少,定位的损伤单元数比实际损伤单元数略多且不遗漏损伤单元等特点.
2020-06-28人气:4545
人气:3576
人气:3353
人气:3106
人气:2961
我要评论
期刊名称:数学进展
期刊人气:3695
主管单位:中国科学协术协会
主办单位:中国数学会
出版地方:北京
专业分类:科学
国际刊号:1000-0917
国内刊号:11-2312/O1
邮发代号:2-503
创刊时间:1955年
发行周期:双月刊
期刊开本:16开
见刊时间:一年半以上
影响因子:0.553
影响因子:0.322
影响因子:0.352
影响因子:0.000
影响因子:0.000
您的论文已提交,我们会尽快联系您,请耐心等待!
你的密码已发送到您的邮箱,请查看!