CN110222307B - 基于fpga的实对称矩阵的特征值分解的并行实现方法 - Google Patents

基于fpga的实对称矩阵的特征值分解的并行实现方法 Download PDF

Info

Publication number
CN110222307B
CN110222307B CN201910504034.7A CN201910504034A CN110222307B CN 110222307 B CN110222307 B CN 110222307B CN 201910504034 A CN201910504034 A CN 201910504034A CN 110222307 B CN110222307 B CN 110222307B
Authority
CN
China
Prior art keywords
matrix
elements
diagonal
data
sequence
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201910504034.7A
Other languages
English (en)
Other versions
CN110222307A (zh
Inventor
张春杰
郝英军
司伟建
邓志安
曲志昱
侯长波
汲清波
杨梦�
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin Engineering University
Original Assignee
Harbin Engineering University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Engineering University filed Critical Harbin Engineering University
Priority to CN201910504034.7A priority Critical patent/CN110222307B/zh
Publication of CN110222307A publication Critical patent/CN110222307A/zh
Application granted granted Critical
Publication of CN110222307B publication Critical patent/CN110222307B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Computational Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computing Systems (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Complex Calculations (AREA)

Abstract

本发明属于阵列信号处理领域，具体涉及基于FPGA的实对称矩阵的特征值分解的并行实现方法。具体实现步骤如下：根据阵元数目构建特征值分解的脉动阵列结构，设定所需的处理单元；对接收的阵元信号进行预处理；求解旋转角度并将其转换为角度值；查表得到对应的正弦值和余弦值；更新矩阵元素和特征向量；判断是否达到要求迭代次数；若未达到，在阵列结构中交换矩阵元素为下次迭代做准备；判断是否需要改变处理单元内部的输入输出顺序；若是，则改变输入输出数据的顺序。本方法通过处理单元之间数据的传递以及处理单元内部的数据顺序的转换，提高了迭代效率而且运算速度快，应用前景广阔。

Description

基于FPGA的实对称矩阵的特征值分解的并行实现方法
技术领域
本发明属于阵列信号处理领域，具体涉及基于FPGA的实对称矩阵的特征值分解的并行实现方法。
背景技术
近几年来，阵列信号处理理论在地震勘探，雷达通信，方位侦测等领域均发挥着重要的作用。在此领域中矩阵的特征值分解是一个应用广泛的矩阵运算。然而特征值分解算法包含大量的非线性运算，运算时间长，影响了算法的实时性。因此，如何快速的得到所需的特征值和特征向量是FPGA实现的难点。
目前，大多数设计者普遍采用了经典的jacobi算法，采用单次扫描的流程，不能够体现jacobi算法的高度并行性。一部分学者提出了一种用于计算n维实矩阵特征值分解的脉动阵列结构，充分体现了jacobi算法的高度并行性。但是，此结构着重于算法的宏观运用，在内部微观还存在改进之处。即，该结构使用固定化的元素交换方式。在迭代过程中绝对值较大的非对角元素的位置不停的变动，使其无法被有效的消除，迭代效率低。
因此，本发明改进了阵列结构的数据交换网络。并且通过在迭代过程中改变阵列结构内部处理单元中的输入输出数据的顺序，提高了迭代运算的效率，使算法能够更快速的收敛。
发明内容
本发明的目的在于提供基于FPGA的实对称矩阵的特征值分解的并行实现方法。
基于FPGA的实对称矩阵的特征值分解的并行实现方法，该方法包括以下步骤：
步骤1：设数据矩阵A∈Rn×n为实对称矩阵，根据矩阵的维数建立特征值分解的阵列结构；
步骤2：对接收的阵元信号进行预处理，将矩阵的复数域运算转换为实数域；
步骤3：根据并行排序规则在数据矩阵An×n中选取对应的数据初始化整个处理单元，并将特征向量矩阵设定为单位阵，即V＝En×n，将迭代次数计数器和交换数据计数器初始化为0；
步骤4：根据Cordic算法同时计算n/2个对角单元中的旋转角度θ；
步骤5：将旋转角度θ由弧度值转换为角度值，而后通过查表法得到sinθ和cosθ，对角单元纵向和横向传递sinθ和cosθ到非对角处理单元和特征向量计算模块；
步骤6：所有处理单元依据jacobi变换求解旋转后的矩阵元素，并且更新特征向量；
步骤7：判断迭代次数是否达到阈值，若达到阈值，V即为所需的特征向量，若没有达到阈值，进行步骤8；
步骤8：累加迭代次数计数器和交换次数计数器，交换相邻处理单元的数据；
步骤9：判断交换计数器是否大于(n/2-1)，若大于，则改变处理单元内部输入输出数据的顺序直至交换计数器等于(n-1)，且交换次数计数器等于(n-1)时归零，重复步骤4到步骤9。
步骤1所述的根据矩阵的维数建立特征值分解的阵列结构包括：若n为偶数，即构建(n/4)(n/2+1)个处理单元，其中包括n/2个对角单元，其余为非对角处理单元；每个对角单元包含3个矩阵元素，每个非对角单元包含4个矩阵元素；对角单元首先计算旋转角度，然后向非对角单元传递角度值；每次迭代都会消除n/2个非对角元素，而每轮的清扫需要(n-1)次的迭代；三轮清扫后，特征值分解达到收敛；若n为奇数，可在矩阵边缘添加1行与1列0元素后构建脉动阵列结构。
步骤6所述的jacobi变换过程为：设A∈Rn×n为实对称矩阵，则存在正交矩阵G∈Rnxn使下式成立：
Figure BDA0002091179740000021
其中，
Figure BDA0002091179740000022
为矩阵A的n个特征值，矩阵
Figure BDA0002091179740000023
中的列向量则为特征值所对应的特征向量，G∈Rnxn为jacobi变换的旋转矩阵；
Figure BDA0002091179740000024
矩阵G被称之为Givens矩阵，其对角线元素除了gpp＝cosθ，gqq＝sinθ其他元素都为1，非对角元素除了gpq＝-sinθ，gqp＝sinθ其他元素都为0；
对于矩阵A中的元素在经过每次旋转后会发生如下变化：
Figure BDA0002091179740000025
Figure BDA0002091179740000026
Figure BDA0002091179740000027
Figure BDA0002091179740000028
Figure BDA0002091179740000029
其中，
Figure BDA0002091179740000031
而构造Givens矩阵中的θ需满足
Figure BDA0002091179740000032
所以
Figure BDA0002091179740000033
根据旋转角度θ求得sinθ及cosθ，构造旋转矩阵；
在迭代过程中，不断得到新的Givens矩阵。将此矩阵累乘后就可以得到特征向量，即
Figure BDA0002091179740000034
所以每次迭代都需要更新特征向量，直到迭代停止。特征向量计算公式为：
Figure BDA0002091179740000035
Figure BDA0002091179740000036
步骤8所述的交换相邻处理单元的数据的方法为：根据jacobi矩阵变换的旋转关系，确定阵列结构中的第一行元素和对角线上的元素即可确定整个阵列结构中的元素，所以只需确定变换后的第一行元素和对角元素就可推出其他元素的变换顺序，而且每次迭代的变换顺序是相同的；数据调度规则为：第一行数据中，设位置序号为偶数的序号为c，且c＞2；对于序号为偶数的数据向左进行移动，放至序号为c-2的位置；序号为2的元素移动到序号为3的元素，对于序号为奇数数据向右移动，且序号为1的元素保持不变，序号为n-1的元素移动到序号为n的元素的位置；第一行数据的交换方式同样适用于对角元素，只是将左右移动变成了上下移动，其他的矩阵元素可以通过已交换顺序的对角元素和第一行元素分别确定行标和列标；而后对照前后的变换顺序可以快速的得出数据的调度规则，此外，相邻的非对角单元内部将会有一个矩阵元素不与外界交换。
步骤9所述的改变处理单元内部输入输出数据的变换方式为：第一个对角处理单元不需要交换，其余对角单元将输入的app和aqq交换；而与第一个对角单元相邻的同一行非对角单元将api和apj，aqi和aqj交换；其余非对角单元将api和aqj，apj和aqi交换，完成本次迭代后，更新的矩阵元素不能按照原来的顺序输出，也需要改变其输出顺序，其改变方式与输入数据的改变方式相对应，即最终的输出数据的顺序与(n/2-1)次迭代之前相同。
本发明的有益效果在于：
本发明通过处理单元之间数据的传递以及处理单元内部的数据顺序的转换，提高了迭代效率而且运算速度快，应用前景广阔。
附图说明
图1为FPGA特征值分解算法流程。
图2为jacobi算法的并行排序规则。
图3为8×8实对称矩阵并行分解步骤。
图4为8×8实对称矩阵并行分解阵列结构。
图5为8×8实对称矩阵数据调度规则。
图6为处理单元内部输入数据变换规则。
图7为FPGA和Matlab的估计误差。
具体实施方式
下面结合附图对本发明做进一步的描述。
本发明具体实现步骤如下：根据阵元数目构建特征值分解的脉动阵列结构，设定所需的处理单元；对接收的阵元信号进行预处理；求解旋转角度并将其转换为角度值；查表得到对应的正弦值和余弦值；更新矩阵元素和特征向量；判断是否达到要求迭代次数；若未达到，在阵列结构中交换矩阵元素为下次迭代做准备；判断是否需要改变处理单元内部的输入输出顺序；若是，则改变输入输出数据的顺序。本方法通过处理单元之间数据的传递以及处理单元内部的数据顺序的转换，提高了迭代效率而且运算速度快，应用前景广阔。
基于FPGA的实对称矩阵的特征值分解的并行实现方法，包括以下步骤：
步骤一、设数据矩阵A∈Rn×n为实对称矩阵，根据矩阵的维数建立特征值分解的阵列结构。
步骤二、对接收的阵元信号进行预处理，将矩阵的复数域运算转换为实数域。
步骤三、根据并行排序规则在数据矩阵An×n中选取对应的数据初始化整个处理单元，并将特征向量矩阵设定为单位阵，即V＝En×n。将迭代次数计数器和交换数据计数器初始化为0。
步骤四、根据Cordic算法同时计算n/2个对角单元中的旋转角度θ。
步骤五、将旋转角度θ由弧度值转换为角度值，而后通过查表法得到sinθ和cosθ。对角单元纵向和横向传递sinθ和cosθ到非对角处理单元和特征向量计算模块。
步骤六、所有处理单元依据jacobi变换求解旋转后的矩阵元素，并且更新特征向量。
步骤七、判断迭代次数是否达到阈值。若达到阈值，V即为所需的特征向量，否则进行步骤八。
步骤八、累加迭代次数计数器和交换次数计数器，交换相邻处理单元的数据。
步骤九、判断交换计数器是否大于(n/2-1)。若大于，则改变处理单元内部输入输出数据的顺序直至交换计数器等于(n-1)。且交换次数计数器等于(n-1)时归零。重复步骤四到九。
对于n维矩阵的特征值分解搭建脉动阵列结构：
n为偶数，即构建(n/4)(n/2+1)个处理单元，其中包括n/2个对角单元，其余为非对角处理单元。每个对角单元包含3个矩阵元素，每个非对角单元包含4个矩阵元素。对角单元首先计算旋转角度，然后向非对角单元传递角度值。每次迭代都会消除n/2个非对角元素，而每轮的清扫需要(n-1)次的迭代。三轮清扫后，特征值分解达到收敛。若n为奇数，可在矩阵边缘添加1行与1列0元素后构建脉动阵列结构。
每次迭代结束后交换相邻处理单元数据的方法：
根据jacobi矩阵变换的旋转关系，只要确定阵列结构中的第一行元素和对角线上的元素即可确定整个阵列结构中的元素。所以只需确定变换后的第一行元素和对角元素就可推出其他元素的变换顺序，而且每次迭代的变换顺序是相同的。数据调度规则为：第一行数据中，设位置序号为偶数的序号为c，且c＞2。对于序号为偶数的数据向左进行移动，放至序号为c-2的位置。序号为2的元素移动到序号为3的元素。对于序号为奇数数据向右移动，且序号为1的元素保持不变，序号为n-1的元素移动到序号为n的元素的位置。第一行数据的交换方式同样适用于对角元素，只是将左右移动变成了上下移动。其他的矩阵元素可以通过已交换顺序的对角元素和第一行元素分别确定行标和列标。而后对照前后的变换顺序可以快速的得出数据的调度规则。此外，由于对角单元内部只有三个元素，所以相邻的非对角单元内部将会有一个矩阵元素不与外界交换。
当交换次数大于(n/2-1)时，处理单元内部输入输出数据的变换方式：
第一个对角处理单元不需要交换，其余对角单元将输入的app和aqq交换；而与第一个对角单元相邻的同一行非对角单元将api和apj，aqi和aqj交换；其余非对角单元将api和aqj，apj和aqi交换。且完成本次迭代后，更新的矩阵元素不能按照原来的顺序输出，也需要改变其输出顺序，其改变方式与输入数据的改变方式相对应，即最终的输出数据的顺序与(n/2-1)次迭代之前相同。
基于FPGA的实对称矩阵的特征值分解的并行实现方法的步骤具体表述为：
a.确定jacobi变换前后矩阵A中元素的变换关系以及特征向量的计算：
a1.jacobi算法变换过程
设A∈Rn×n为实对称矩阵，则存在正交矩阵G∈Rnxn使得
Figure BDA0002091179740000051
其中
Figure BDA0002091179740000052
为矩阵A的n个特征值，矩阵
Figure BDA0002091179740000053
中的列向量则为特征值所对应的特征向量。G∈Rnxn为jacobi变换的旋转矩阵。
Figure BDA0002091179740000061
矩阵G被称之为Givens矩阵，其对角线元素除了gpp＝cosθ，gqq＝sinθ其他元素都为1，非对角元素除了gpq＝-sinθ，gqp＝sinθ其他元素都为0。
对于矩阵A中的元素在经过每次旋转后会发生如下变化：
Figure BDA0002091179740000062
Figure BDA0002091179740000063
Figure BDA0002091179740000064
Figure BDA0002091179740000065
Figure BDA0002091179740000066
其中，
Figure BDA0002091179740000067
而构造Givens矩阵中的θ需满足
Figure BDA0002091179740000068
所以
Figure BDA0002091179740000069
根据旋转角度θ求得sinθ及cosθ，构造旋转矩阵。
a2.特征向量的计算
在迭代过程中，不断得到新的Givens矩阵。将此矩阵累乘后就可以得到特征向量，即
Figure BDA00020911797400000610
所以每次迭代都需要更新特征向量，直到迭代停止。特征向量计算公式为：
Figure BDA00020911797400000611
Figure BDA00020911797400000612
b.构建jacobi并行变换模块
b1、设数据矩阵A∈Rn×n为实对称矩阵，n为偶数则根据矩阵维数n搭建特征值分解的脉动阵列结构。即构建(n/4)(n/2+1)个处理单元，其中包括n/2个对角单元。每次迭代都会消除n/2个非对角元素，而每轮的清扫需要(n-1)次的迭代。三轮清扫后，特征值分解达到收敛。若n为奇数，可在矩阵边缘添加1行与1列0元素后构建脉动阵列结构；
b2、设计数据调度规则：第一行数据中，设为位置序号为偶数的序号为c，且c＞2。对于序号为偶数数据向左进行移动，放至序号为c-2的位置。序号为2的元素移动到序号为3的元素。对于序号为奇数数据向右移动，且序号为1的元素保持不变，序号为n-1的元素移动到序号为n的元素的位置。第一行数据的交换方式同样适用于对角元素，只是将左右移动变成了上下移动。其他的矩阵元素可以通过已交换顺序的对角元素和第一行元素分别确定行标和列标。而后对照前后的变换顺序可以快速的得出数据的调度规则；
b3、依据并行排序规则将矩阵A中上三角数据赋值到所有处理单元中，完成初始赋值。赋值结束后，传递数据有效标志转入b4；
b4、数据有效后，对角单元利用Cordic算法计算旋转角度θ，此时θ为弧度值。传递数据有效标志后进入b5；
b5、将θ转换为角度值，而后利用查表法得到sinθ，cosθ。按照阵列结构并将其传递到相邻的非对角单元，直至遍历整个阵列结构。得到旋转角度后转入b6；
b6、对角单元和非对角单元根据公式GTAG，利用查表所得的sinθ和cosθ，按步骤a所示更新矩阵元素。同时完成特征向量的计算。所有处理单元计算完成后转入b7；
b7、判断迭代次数是否达到阈值。若达到阈值，V即为所需的特征向量，否则进行步骤b8；
b8、每一轮计算完成后累加迭代次数计数器和交换次数计数器，并且对相邻处理单元进行数据交换。根据b2所述的调度规则先确定第一行和对角线的矩阵元素，进而确定其他非对角元素。而且每一轮结束后的数据交换规则是相同。因此，只要事先确定好交换规则，就不需要改变。这样降低了程序的复杂度，便于程序的编写；
b9、判断交换计数器是否大于(n/2-1)。若大于，则改变处理单元内部输入数据的顺序。在当前迭代完成后，更新的矩阵元素不能按照原来的顺序输出，也需要改变其输出顺序，其改变方式与输入数据的改变方式相对应。当交换计数器等于(n-1)时，还原输入输出数据的顺序。改变顺序如图6。且交换次数计数器等于(n-1)时归零。重复步骤b4到b9。
实例
阵列模型为8阵元均匀圆阵。信号频率为3GHz窄带远场信号且信源数为1，圆阵半径为0.05m。阵元接收噪声为零均值高斯白噪声，快拍数为64。选用芯片xc7vx690tffg1761-2。
在估计性能包括计算精度、计算速度和资源消耗，具体可用下面指标评价。
1.计算精度
由FPGA完成特征值分解将计算的所得噪声子空间导入到Matlab中，由Matlab完成最终的谱峰搜索。接收的信号方位角和仰角分别为40°和30°。信噪比在2到12之间，在每个参数下仿真100次。通过
Figure BDA0002091179740000071
来分析Matlab和FPGA对入射角度估计的均方误差。
2.计算速度
计算消耗的时钟数Nclk越小表示计算消耗的时间越少，计算速度越快。
3.资源消耗
触发器消耗资源；LUT查找表消耗资源；RAM消耗资源，而RAM消耗资源反应在LUTRAM和BRAM上；运算单元DSP48E1消耗资源。
实施方式，具体为以下步骤：
a.仿真信号建模：
a1.MUSIC算法对于窄带信号的DOA估计的数学模型为：
X(t)＝AS(t)+N(t)
式中，X(t)为阵列的8×1维快拍数据矢量，N(t)为阵列的8×1维噪声数据矢量，
S(t)为空间信号的1×8维矢量，A为空间阵列的8×1维流型矩阵。
a2.对接收信号进行预处理：
将8阵元均匀圆阵分为两个子阵，其输出分别为X1(t)和X2(t)，令
M＝Re[X1(t)]+Re[X2(t)] Q＝Im[X1(t)]+Im[X2(t)]
L＝Im[X1(t)]-Im[X2(t)] G＝Re[X2(t)]-Re[X1(t)]
令
Figure BDA0002091179740000081
此时Z为实矩阵，所以协方差矩阵
Figure BDA0002091179740000082
b应用本发明对Rz进行特征值分解，计算其特征值和特征向量。
b1.数据矩阵Rz为8×8矩阵。如附图4根据矩阵维数建立10个PE(ProcessingElement)单元，其中包括6个非对角单元，4个对角单元。
b2.数据调度规则：如附图5，第一行数据中，设为位置序号为偶数的序号为c，且c＞2。对于序号为偶数数据向左进行移动，放至序号为c-2的位置。序号为2的元素移动到序号为3的元素。对于序号为奇数数据向右移动，且序号为1的元素保持不变，序号为7的元素移动到序号为8的元素的位置。
b3.根据jacobi并行排序规则，排列矩阵非对角元素的消除顺序。
b4.按照消除顺序初始化PE单元的输入数据，将特征向量初始化为单位阵，即V＝En×n。迭代次数计数器和交换数据计数器初始化为0。
b5.由公式
Figure BDA0002091179740000091
对角单元计算得到4个旋转角度θ1,θ2,θ3,θ4，通过查表后将其对应的sinθ和cosθ按照阵列结构传递到其他非对角单元。
b6.根据查表所得的sinθ和cosθ完成矩阵的旋转和特征向量的计算。并且令迭代计数器和交换数据计数器加1，并判断迭代次数是否超过阈值。未超过阈值则进行b6步骤，否则进行b8步骤。
b7.判断数据交换次数是否达到3次。若达到3次，除PE1以外改变其它的处理单元输入和输出数据的顺序。PE5，PE8，PE10将输入的app和aqq交换；PE2，PE3，PE4将api和apj，aqi和aqj交换；PE6，PE7，PE9将api和aqj，apj和aqi交换。完成本次迭代的矩阵旋转后，再将更新后的矩阵元素交换回原来的顺序。当交换次数达到7次时，停止输入输出数据顺序交换。且交换次数达到7次归零。
也就是说完成一轮迭代后归零。重复b5到b7步骤。
b8.提取特征值所对应的特征向量，分离噪声子空间和信号子空间。
计算速度：Nclk＝1060
资源消耗：
Figure BDA0002091179740000092
计算精度：
如附图7，FPGA与Matlab的估计误差基本一致。而且随着信噪比的增大，估计误差减小。
本发明提供一种基于FPGA的实对称矩阵特征值分解的并行实现方法，属于阵列信号处理领域。具体实现步骤如下：根据阵元数目构建特征值分解的脉动阵列结构，设定所需的处理单元；对接收的阵元信号进行预处理；求解旋转角度并将其转换为角度值；查表得到对应的正弦值和余弦值；更新矩阵元素和特征向量；判断是否达到要求迭代次数；若未达到，在阵列结构中交换矩阵元素为下次迭代做准备；判断是否需要改变处理单元内部的输入输出顺序；若是，则改变输入输出数据的顺序。本方法通过处理单元之间数据的传递以及处理单元内部的数据顺序的转换，提高了迭代效率而且运算速度快，应用前景广阔。

Claims (5)

1.基于FPGA的实对称矩阵的特征值分解的并行实现方法，其特征在于，该方法包括以下步骤：
步骤1：设数据矩阵A∈Rn×n为实对称矩阵，根据矩阵的维数建立特征值分解的阵列结构；
步骤2：对接收的阵元信号进行预处理，将矩阵的复数域运算转换为实数域；
步骤3：根据并行排序规则在数据矩阵An×n中选取对应的数据初始化整个处理单元，并将特征向量矩阵设定为单位阵，即V＝En×n，将迭代次数计数器和交换数据计数器初始化为0；
步骤4：根据Cordic方法同时计算n/2个对角单元中的旋转角度θ；
步骤5：将旋转角度θ由弧度值转换为角度值，而后通过查表法得到sinθ和cosθ，对角单元纵向和横向传递sinθ和cosθ到非对角处理单元和特征向量计算模块；
步骤6：所有处理单元依据jacobi变换求解旋转后的矩阵元素，并且更新特征向量；
步骤7：判断迭代次数是否达到阈值，若达到阈值，V即为所需的特征向量，若没有达到阈值，进行步骤8；
步骤8：累加迭代次数计数器和交换次数计数器，交换相邻处理单元的数据；
步骤9：判断交换计数器是否大于(n/2-1)，若大于，则改变处理单元内部输入输出数据的顺序直至交换计数器等于(n-1)，且交换次数计数器等于(n-1)时归零，重复步骤4到步骤9。
2.根据权利要求1所述的基于FPGA的实对称矩阵的特征值分解的并行实现方法，其特征在于，步骤1所述的根据矩阵的维数建立特征值分解的阵列结构包括：若n为偶数，即构建(n/4)(n/2+1)个处理单元，其中包括n/2个对角单元，其余为非对角处理单元；每个对角单元包含3个矩阵元素，每个非对角单元包含4个矩阵元素；对角单元首先计算旋转角度，然后向非对角单元传递角度值；每次迭代都会消除n/2个非对角元素，而每轮的清扫需要(n-1)次的迭代；三轮清扫后，特征值分解达到收敛；若n为奇数，在矩阵边缘添加1行与1列0元素后构建脉动阵列结构。
3.根据权利要求1所述的基于FPGA的实对称矩阵的特征值分解的并行实现方法，其特征在于，步骤6所述的jacobi变换过程为：设A∈Rn×n为实对称矩阵，则存在正交矩阵G∈Rnxn使下式成立：
Figure FDA0002091179730000011
其中，
Figure FDA0002091179730000012
为矩阵A的n个特征值，矩阵
Figure FDA0002091179730000013
中的列向量则为特征值所对应的特征向量，G∈Rnxn为jacobi变换的旋转矩阵；
Figure FDA0002091179730000021
矩阵G被称之为Givens矩阵，其对角线元素除了gpp＝cosθ，gqq＝sinθ其余元素都为1，非对角元素除了gpq＝-sinθ，gqp＝sinθ其余元素都为0；
对于矩阵A中的元素在经过每次旋转后会发生的变化表现为下式：
Figure FDA0002091179730000022
Figure FDA0002091179730000023
Figure FDA0002091179730000024
Figure FDA0002091179730000025
Figure FDA0002091179730000026
其中，
Figure FDA0002091179730000027
而构造Givens矩阵中的θ需满足
Figure FDA0002091179730000028
所以
Figure FDA0002091179730000029
根据旋转角度θ求得sinθ及cosθ，构造旋转矩阵；
在迭代过程中，不断得到新的Givens矩阵，将此矩阵累乘后就得到特征向量，即
Figure FDA00020911797300000210
所以每次迭代都需要更新特征向量，直到迭代停止，特征向量计算公式为：
Figure FDA00020911797300000211
4.根据权利要求1所述的基于FPGA的实对称矩阵的特征值分解的并行实现方法，其特征在于，步骤8所述的交换相邻处理单元的数据的方法为：根据jacobi矩阵变换的旋转关系，确定阵列结构中的第一行元素和对角线上的元素即确定整个阵列结构中的元素，所以只需确定变换后的第一行元素和对角元素就推出其余元素的变换顺序，而且每次迭代的变换顺序是相同的；数据调度规则为：第一行数据中，设位置序号为偶数的序号为c，且c＞2；对于序号为偶数的数据向左进行移动，放至序号为c-2的位置；序号为2的元素移动到序号为3的元素，对于序号为奇数时数据向右移动，且序号为1的元素保持不变，序号为n-1的元素移动到序号为n的元素的位置；第一行数据的交换方式同样适用于对角元素，只是将左右移动变成了上下移动，其余的矩阵元素通过已交换顺序的对角元素和第一行元素分别确定行标和列标；而后对照前后的变换顺序快速的得出数据的调度规则，此外，相邻的非对角单元内部将会有一个矩阵元素不与外界交换。
5.根据权利要求1所述的基于FPGA的实对称矩阵的特征值分解的并行实现方法，其特征在于，步骤9所述的改变处理单元内部输入输出数据的变换方式为：第一个对角处理单元不需要交换，其余对角单元将输入的app和aqq交换；而与第一个对角单元相邻的同一行非对角单元将api和apj，aqi和aqj交换；其余非对角单元将api和aqj，apj和aqi交换，完成本次迭代后，更新的矩阵元素不能按照原来的顺序输出，也需要改变其输出顺序，其改变方式与输入数据的改变方式相对应，即最终的输出数据的顺序与(n/2-1)次迭代之前相同。
CN201910504034.7A 2019-06-12 2019-06-12 基于fpga的实对称矩阵的特征值分解的并行实现方法 Active CN110222307B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910504034.7A CN110222307B (zh) 2019-06-12 2019-06-12 基于fpga的实对称矩阵的特征值分解的并行实现方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910504034.7A CN110222307B (zh) 2019-06-12 2019-06-12 基于fpga的实对称矩阵的特征值分解的并行实现方法

Publications (2)

Publication Number Publication Date
CN110222307A CN110222307A (zh) 2019-09-10
CN110222307B true CN110222307B (zh) 2022-10-28

Family

ID=67816403

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910504034.7A Active CN110222307B (zh) 2019-06-12 2019-06-12 基于fpga的实对称矩阵的特征值分解的并行实现方法

Country Status (1)

Country Link
CN (1) CN110222307B (zh)

Families Citing this family (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111427537B (zh) * 2020-03-17 2023-06-30 云南大学 一种基于fpga的脉动阵列并行排序方法及装置
CN111597782B (zh) * 2020-05-20 2023-10-27 比科奇微电子(杭州)有限公司 数据的排序处理方法及处理装置
CN111859035B (zh) * 2020-08-12 2022-02-18 华控清交信息科技(北京)有限公司 数据处理方法及装置
CN112015369B (zh) * 2020-08-25 2022-09-16 湖南艾科诺维科技有限公司 基于fpga的信号处理方法、电子设备和存储介质
CN112214729A (zh) * 2020-10-15 2021-01-12 上海无线电设备研究所 一种用于music算法特征值分解的多核并行处理方法
CN116028770B (zh) * 2023-01-18 2024-01-12 珠海微度芯创科技有限责任公司 适用于实、复协方差矩阵的特征值分解硬件实现方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106919537A (zh) * 2017-03-07 2017-07-04 电子科技大学 一种基于FPGA的Jacobi变换的高效实现方法
CN106940689A (zh) * 2017-03-07 2017-07-11 电子科技大学 基于Jacobi迭代算法的高精度矩阵特征值分解实现方法
CN109740114A (zh) * 2018-12-28 2019-05-10 中国航天科工集团八五一一研究所 基于fpga的实对称矩阵特征分解实时处理方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
TWI456516B (zh) * 2010-12-17 2014-10-11 Univ Nat Chiao Tung 獨立成分分析處理器

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106919537A (zh) * 2017-03-07 2017-07-04 电子科技大学 一种基于FPGA的Jacobi变换的高效实现方法
CN106940689A (zh) * 2017-03-07 2017-07-11 电子科技大学 基于Jacobi迭代算法的高精度矩阵特征值分解实现方法
CN109740114A (zh) * 2018-12-28 2019-05-10 中国航天科工集团八五一一研究所 基于fpga的实对称矩阵特征分解实时处理方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
实对称矩阵特征值分解高速并行算法的FPGA实现;王飞等;《空军工程大学学报(自然科学版)》;20081215(第06期);全文 *

Also Published As

Publication number Publication date
CN110222307A (zh) 2019-09-10

Similar Documents

Publication Publication Date Title
CN110222307B (zh) 基于fpga的实对称矩阵的特征值分解的并行实现方法
CN101763445B (zh) 一种高光谱图像降维芯片
CN110361691B (zh) 基于非均匀阵列的相干信源doa估计fpga实现方法
CN107390199B (zh) 一种雷达机动目标跟踪波形设计方法
CN106228240A (zh) 基于fpga的深度卷积神经网络实现方法
CN103516643B (zh) 一种mimo检测预处理装置及方法
CN110109050A (zh) 嵌套阵列下基于稀疏贝叶斯的未知互耦的doa估计方法
AU2021103515A4 (en) Data assimilation method based on multi-layer perceptron
CN106940689A (zh) 基于Jacobi迭代算法的高精度矩阵特征值分解实现方法
CN104793176A (zh) 一种基于fpga的doa估计快速实现方法
Jia et al. Tensorlib: A spatial accelerator generation framework for tensor algebra
CN108663654A (zh) 一种基于连续量子鸽群的360度全方位动态测向方法
Jovanović et al. A survey of hardware self-organizing maps
Li et al. Hardware acceleration of MUSIC algorithm for sparse arrays and uniform linear arrays
CN110110285B (zh) 一种用于FPGA的并行Jacobi计算加速实现方法
CN112800599B (zh) 一种阵元失配情况下基于admm的无网格doa估计方法
CN105608057A (zh) 一种分时复用硬件资源的信号子空间分解的fpga实现模块及其fpga实现方法
Xiao et al. HALOC: hardware-aware automatic low-rank compression for compact neural networks
Wang et al. Hardware efficient architectures of improved Jacobi method to solve the eigen problem
CN108008665B (zh) 基于单片fpga的大规模圆阵实时波束形成器及波束形成计算方法
Edman Digital hardware aspects of multiantenna algorithms
CN105022025A (zh) 基于稀疏处理的信号波达方向估计方法
CN109459768A (zh) 一种基于北斗卫星信号强度权重优化模型的快速选星方法
CN108388673A (zh) 一种用于经济管理分析数据的计算机系统
CN103761379B (zh) 一种基于包络对偶算法体系的对地观测卫星多学科优化系统

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant