CN105701068A - 基于分时复用技术的cholesky矩阵求逆系统 - Google Patents

Abstract

本发明公开了一种基于分时复用技术的cholesky矩阵求逆系统，降低了在雷达信号处理中对共轭对称矩阵求逆所使用的硬件资源，同时通过提高矩阵求逆运算的并行度来保证运算速度。本发明根据公式A = L * D * L^H，通过主控状态机将cholesky分解矩阵求逆划为三个步骤进行，分别是cholesky分解、上三角矩阵求逆和矩阵相乘。在主控状态机的控制下，通过复选器，分时复用运算资源和存储资源。该技术方案是基于单端口存储器和单精度浮点复数运算单元的并行流水处理，并支持任意阶数的共轭对称矩阵求逆，有效加快了运算速度，提高了硬件利用率，且具有良好的兼容性。

Description

基于分时复用技术的cholesky矩阵求逆系统

技术领域

本发明属于雷达信号处理技术领域，具体是一种基于FPGA的正定复数矩阵求逆的硬件结构。在雷达信号处理中，常对正定复数矩阵求逆，若使用通用矩阵求逆的算法和结构，会使用额外的资源。本发明提出的硬件结构专门针对正定复数矩阵求逆，适合雷达信号处理中的求逆操作。

背景技术

矩阵求逆的方法有很多，包括伴随矩阵法、高斯消元法、初等变换法、矩阵分解法等等。

上述各种矩阵求逆的方法中，伴随矩阵法、高斯消元法、初等变换法中需要求大量的行列式，每个行列式都几乎要计算到所有的矩阵元素，计算量大，对存储空间的需求大，在FPGA上实现比较困难。

矩阵分解法是将矩阵分解为一些较简单矩阵的乘积，如三角矩阵或酉矩阵，他们都具有某种特性，较容易得到逆矩阵。对分解得到的矩阵求逆后相乘，得到原矩阵的逆矩阵。矩阵分解法分为几种：1)LU分解法，2)QR分解法，3)奇异值分解法，4)cholesky分解法。

前三种矩阵分解法，都是对所有可逆复数矩阵进行求解与计算，算法较复杂，资源使用较多；而第四种分解法是针对共轭对称矩阵求逆，算法较简单，资源使用也较少。而在雷达信号处理中，主要是对自相关后得到的共轭对称矩阵进行求逆操作，因此本发明选择cholesky分解法。

cholesky分解法常用于线性方程的求解、矩阵求逆和行列式计算，在雷达信号处理、科学计算等诸多领域也有广泛应用。根据公式A＝L*D*L^H，cholesky分解法是将原矩阵A分解为单位下三角矩阵L(对角线元素为1)、对角矩阵D，上三角矩阵L^H(是L的共轭转置矩阵)，原矩阵A对角线上的元素必须是实数，其他元素可以是复数。

因三角矩阵几乎一半数据都为0，所以求逆简单，且L^H是L的共轭转置矩阵，因此只需求得(L^H)^-1，取共轭转置即得到了L^-1。

令G_ij＝L_ij*D_ii,则cholesky分解公式可由此推得：

d₁₁＝a₁₁，

对i＝2...n有

$\left\{\begin{array}{c}{g}_{ij}=a_{ij}-\underset{k=1}{\overset{j-1}{\Σ}}{g}_{ik}*l_{jk}^{*},(j=1,\mathrm{2...}i-1);\\ l_{ij}=\frac{g_{ij}}{d_{ij}},(j=1,\mathrm{2...}i-1);\\ d_{ii}=d_{ii}-\underset{k=1}{\overset{i-1}{\Σ}}{g}_{ik}*l_{ik}^{*};\end{array}\right.$

由于d_ii出现在分母中，因此只有当d_ii≠0(i＝1...n)时，分解才能进行到底。此外，在实际的计算过程中，即使不为零，但当它很小的时候，由于计算机的精度限制，会导致下溢，也会使分解运算不稳定。

发明内容

发明目的：在FPGA芯片上提供一种硬件复杂度低，资源利用率高的正定矩阵浮点求逆结构。在尽可能减少硬件资源消耗的情况下，通过提高矩阵求逆运算的并行度来保证运算速度。

技术方案：一种基于分时复用技术的cholesky矩阵求逆系统，包括五个模块：1)主控状态机模块2)资源复用模块3)cholesky分解模块4)上三角矩阵求逆模块5)矩阵相乘模块。

所述的主控状态机模块与系统总线相连接，在求逆操作开始前，接收使能拉高，根据待求逆矩阵维度，从外部按行顺序接收原矩阵数据，存放在资源复用模块的存储器中，接收完成后，将接收使能拉低，并通过cholesky分解模块、上三角矩阵求逆模块、矩阵相乘模块分别进行cholesky分解、上三角矩阵求逆、矩阵乘法3个步骤进行cholesky矩阵求逆操作。求逆完成后向总线发出求逆结束信号，通知系统总线从存储模块读取运算结果。

所述的资源复用模块内部有3个存储器，1个单精度浮点复数乘法器，1个单精度浮点复数累加器，存储器用于存放求逆过程中各个阶段的矩阵数据，在cholesky分解，上三角求逆，矩阵乘法的运算中，都需进行复数乘法，复数累加操作，且三个运算步骤是顺序执行，因此可复用单精度浮点复数乘法器和单精度浮点复数累加器(见表1)。资源复用模块与主控状态机模块、cholesky分解模块、上三角矩阵求逆模块、矩阵相乘模块相连接。根据主控状态机的状态，资源复用内部的存储资源和运算资源受不同模块控制。

所述的cholesky分解模块是通过对原矩阵元素进行乘法，加法，除法运算将原矩阵A分解为单位下三角矩阵L(对角线元素为1)、对角矩阵D，上三角矩阵L^H(是L的共轭转置矩阵)。cholesky分解模块内部是一个地址生成状态机和两个单精度浮点除法器，外部与主控状态机和资源复用模块相连，在分解过程中将使用资源复用模块中的2个存储器和浮点乘法器，浮点累加器，2个存储器分别存放原始矩阵A和分解得到的上三角矩阵L^H。

在接收到主控状态机模块输入的cholesky分解起始信号后，cholesky分解模块的内部地址生成状态机根据矩阵维度开始按cholesky分解公式的计算顺序生成原矩阵读地址和分解后上三角矩阵的写地址。将读出的数据依次送入浮点乘法器，浮点累加器，浮点除法器，除法器结果即为分解后上三角矩阵的元素。在分解完成后，通过复用浮点除法器对分解后得到的对角矩阵的所有对角元素求倒数(即求对角阵的逆)。

所述的上三角矩阵求逆模块，因choelsky分解后得到的上三角矩阵的对角线元素全为1，所以可以对上三角矩阵求逆公式进行简化，省去除法操作。上三角矩阵求逆模块内部是一个地址生成状态机，外部与主控状态机和资源复用模块相连，在上三角矩阵求逆过程中将使用资源复用模块中的2个存储器、浮点乘法器和浮点累加器，2个存储器分别存放的是分解得到的上三角矩阵L^H和求逆得到的上三角逆矩阵(L^H)^-1，在接收到主控状态机模块输入的上三角矩阵求逆起始信号后，上三角矩阵求逆模块的内部地址生成状态机根据矩阵维度开始按上三角矩阵求逆公式的计算顺序生成原上三角矩阵读地址和上三角逆矩阵的写地址，将读出的数据依次送入浮点乘法器，浮点累加器，累加器结果即为上三角逆矩阵的元素。

所述的矩阵乘法模块内部是一个地址生成状态机，外部与主控状态机和资源复用模块相连，在矩阵乘法过程中将使用资源复用模块中的3个存储器、浮点乘法器和浮点累加器，3个存储器分别存放的是求逆得到的上三角逆矩阵(L^H)^-1、矩阵(L^H)^-1和矩阵D^-1相乘的结果(L^H)^-1*D^-1和原矩阵的逆A^-1。在接收到主控状态机模块输入的矩阵乘法模块起始信号后，通过矩阵乘法模块内部的地址生成状态机控制的存储器读写，先进行上三角矩阵(L^H)^-1与对角阵D^-1相乘，再进行上三角矩阵(L^H)^-1*D^-1与下三角矩阵L^-1相乘，得到元矩阵的逆A^-1。

工作原理：cholesky分解矩阵求逆运算划分为三个步骤进行：

(1)cholesky分解，根据公式A＝L*D*L^H，将原矩阵A分解为单位下三角矩阵L(对角线元素为1)、对角矩阵D，上三角矩阵L^H(L^H是L的共轭转置矩阵)，并且将D对角线上的各个元素分别求倒数(即D的逆)，原矩阵A对角线上的元素必须是实数，其他元素可以是复数；

(2)上三角矩阵求逆，将L^H矩阵求逆得到其逆矩阵(L^H)^-1；

(3)矩阵乘法，根据公式A^-1＝(L^H)^-1*D^-1*L^-1，先将矩阵(L^H)^-1和矩阵D^-1相乘，再将矩阵(L^H)^-1*D^-1与矩阵L^-1相乘得到原矩阵的逆A^-1。

有益效果：在雷达信号处理中，常对正定复数矩阵求逆，若使用通用矩阵求逆的算法和结构，会使用额外的资源。本发明提出的硬件结构专门针对正定复数矩阵求逆，适合雷达信号处理中的求逆操作。同时，通过分时复用存储器和运算资源，最大限度地提高硬件资源的利用率以及运算效率。

附图说明

图1为cholesky矩阵求逆硬件架构示意图；

图2为cholesky分解流程图。

具体实施方式

下面结合具体实施例，进一步阐明本发明，应理解这些实施例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。

基于分时复用技术的cholesky矩阵求逆系统，包括五个模块：1)主控状态机模块2)资源复用模块3)cholesky分解模块4)上三角矩阵求逆模块5)矩阵相乘模块。

主控状态机模块与系统总线相连接，在求逆操作开始前，接收使能拉高，根据待求逆矩阵维度，从外部按行顺序接收原矩阵数据，存放在资源复用模块的存储器中，接收完成后，将接收使能拉低，并通过cholesky分解模块、上三角矩阵求逆模块、矩阵相乘模块分别进行cholesky分解、上三角矩阵求逆、矩阵乘法3个步骤进行cholesky矩阵求逆操作。求逆完成后向总线发出求逆结束信号，通知系统总线从存储模块读取运算结果。

资源复用模块内部有3个存储器(图1中RAMA-C)，1个单精度浮点复数乘法器，2个单精度浮点累加器。资源复用模块与主控状态机模块、cholesky分解模块、上三角矩阵求逆模块、矩阵相乘模块相连接。根据主控状态机的状态，资源复用内部的存储资源和运算资源受不同模块控制。

cholesky分解模块是通过对原矩阵元素进行乘法，加法，除法运算将原矩阵A分解为单位下三角矩阵L(对角线元素为1)、对角矩阵D，上三角矩阵L^H(是L的共轭转置矩阵)。cholesky分解模块内部是一个地址生成状态机和两个单精度浮点除法器，外部与主控状态机和资源复用模块相连，在分解过程中将使用资源复用模块中的存储资源和浮点乘法器，浮点累加器。

在接收到主控状态机模块输入的cholesky分解起始信号后，cholesky分解模块的内部地址状态机根据矩阵维度开始按顺序生成原矩阵读地址和分解后上三角矩阵的写地址。在分解完成后，通过复用浮点除法器对分解后得到的对角矩阵的所有对角元素求倒数(即求对角阵的逆)。

上三角矩阵求逆模块，因choelsky分解后得到的上三角矩阵的对角线元素全为1，所以可以对上三角矩阵求逆公式进行简化，省去除法操作。上三角矩阵求逆模块内部是一个地址生成状态机，外部与主控状态机和资源复用模块相连，在上三角矩阵求逆过程中将使用资源复用模块中的存储资源、浮点乘法器和浮点累加器。在接收到主控状态机模块输入的上三角矩阵求逆起始信号后，上三角矩阵求逆模块的内部地址生成状态机根据矩阵维度开始按顺序生成原上三角矩阵读地址和上三角逆矩阵的写地址。

矩阵乘法模块内部是一个地址生成状态机，外部与主控状态机和资源复用模块相连，在矩阵乘法模过程中将使用资源复用模块中的存储资源、浮点乘法器和浮点累加器。在接收到主控状态机模块输入的矩阵乘法模块起始信号后，通过矩阵乘法模块内部状态机控制的存储器读写，先进行上三角矩阵与对角阵相乘，再进行上三角矩阵与下三角矩阵相乘。

cholesky矩阵求逆操作分为3个步骤，分别是cholesky分解、上三角矩阵求逆、矩阵相乘，由cholesky矩阵求逆顶层模块的一个主控状态机控制三步顺序执行，并分时复用资源复用模块中的3个存储器、浮点乘法器和浮点累加器，整体架构如图1。

(1)cholesky分解

根据公式A＝L*D*L^H，令G＝L*D,d₁₁＝a₁₁，对i＝2...n有

$\left\{\begin{array}{c}{g}_{ij}=a_{ij}-\underset{k=1}{\overset{j-1}{\Σ}}{g}_{ik}*l_{jk}^{*},(j=1,\mathrm{2...}i-1);\\ l_{ij}=\frac{g_{ij}}{d_{ij}},(j=1,\mathrm{2...}i-1);\\ d_{ii}=d_{ii}-\underset{k=1}{\overset{i-1}{\Σ}}{g}_{ik}*l_{ik}^{*};\end{array}\right.$

硬件实现过程如下：

Ⅰ.计算矩阵G当前列的所有元素：由cholesky分解公式，按列依次从资源复用模块的两片RAM中分别读出g_ik和l_ik的共轭(虚部符号取反)送入浮点复数乘法器，再将乘法结果送入累加器(实虚部符号位取反实现累减)得到g_ij，并将g_ij写入对应的RAM，读写g_ik和l_jk的地址，由cholesky分解模块内部的地址生成状态机产生。

Ⅱ.计算三角矩阵L当前列所有元素：由cholesky分解公式，将累加器的结果g_ij实部和虚部分别和存放在寄存器中的d_ii送入两个浮点除法器，再将除法结果l_ij写入对应的RAM，l_ij的写地址由状态机计算产生。因为输入原矩阵的特殊性，分解得到的对角矩阵D的各元素都是实数，所以浮点除法器不需要设为通用的浮点复数除法器，节省了运算资源。

Ⅲ.更新对角矩阵D当前列元素：根据公式可得，矩阵D的对角元素等于分解后得到的矩阵G的对角元素；因此，若当前列是第一列，则直接令d₁₁＝a₁₁，存入寄存器中，否则，通过计数器判断此时矩阵G当前元素是否为对角元素，若是，则将此时矩阵G元素存放在存储D矩阵元素的寄存器中。

重复上述三个步骤，直到cholesky分解完成。

Ⅳ.在cholesky分解完成后，通过复选器改变浮点除法器的输入，从存储D矩阵元素的寄存器中将D矩阵对角线上的各个元素读出并分别求倒数(即D的逆)；

(2)上三角矩阵求逆

对上三角矩阵求逆，有以下计算公式：

$b_{ij}(i\≤j)=\left\{\begin{array}{c}\frac{1}{a_{ii}},(i=j)\\ -\frac{1}{a_{ii}}\·\underset{k=i+1}{\overset{j}{\Σ}}{a}_{ik}\·b_{kj},(i\≠j)\end{array}\right.$

而我们cholesky分解得到的上三角矩阵LH对角元素都为1，因此上三角矩阵求逆公式可化简为：

$b_{ij}(i\≤j)=\left\{\begin{array}{c}1,(i=j)\\ -\underset{k=i+1}{\overset{j}{\Σ}}{a}_{ik}\·b_{kj},(i\≠j)\end{array}\right.$

求逆过程从对角线开始，按斜行进行，硬件实现过程如下：

Ⅰ.对角线元素初始化：分解得到的上三角矩阵L^H和上三角逆矩阵(L^H)^-1分别存放在资源复用模块的两片RAM中，通过计数器判断当前斜行是否为对角线，若是，则直接将1写入(L^H)^-1对角线对应的RAM地址。

Ⅱ.非对角线元素计算：若当前斜行不是对角线，依次从资源复用模块的两片RAM中读取L^H和(L^H)^-1对应位置的元素，送入浮点复数乘法器，再将乘法器结果送入浮点累加器(实虚部符号位取反实现累减)得到(L^H)^-1当前计算的元素并写入对应RAM的地址，经过迭代直到算出(L^H)^-1右上角的值后，上三角求逆结束。对L^H和(L^H)^-1对应RAM的读写地址，由上三角求逆模块内部的地址生成状态机产生。

(3)矩阵相乘

得到原矩阵A的逆要通过计算公式A^-1＝(L^H)^-1*D^-1*L^-1，硬件实现过程如下：

Ⅰ.上三角矩阵乘对角阵：从资源复用模块的RAM中按行读取(L^H)^-1的值，将其与寄存器中的D^-1的对角元素按计算顺序送入浮点乘累加单元，并将计算结果存储在资源复用模块的RAM中。RAM的读写地址和对角阵元素的读取由模块内部的地址生成状态机产生。

Ⅱ.上三角矩阵乘下三角阵：从资源复用模块的RAM中读取矩阵(L^H)^-1*D^-1和矩阵L^-1的值，将它们顺序送入浮点乘累加单元，让后把累加器的结果存入资源复用模块的存放A^-1的RAM中。RAM的读写地址由矩阵相乘模块内部的地址生成状态机产生。

本设计具体在V7-690芯片上进行验证，主频达到200MHz。基于200MHz时钟频率下，复数浮点乘法器、浮点累加器、浮点除法器均为的流水运算单元,且延时拍数可调，所有运算单元均为单精度浮点运算单元，72阶矩阵求逆的运行时间为2.15ms，与传统cpu进行矩阵求逆操作相比，大幅提高了计算速度，且结果相对误差小于10^-4。

表1cholesky矩阵求逆资源复用表

	RAM A	RAM B	RAM C	复数乘法器	复数累加器
						cholesky	A	LH	×	√	√
上三角求逆	(LH)-1	LH	×	√	√
						矩阵相乘	(LH)-1	(LH)-1*D-1	A-1	√	√

Claims (6)

1.一种基于分时复用技术的cholesky矩阵求逆系统，其特征在于，包括五个模块：主控状态机模块、资源复用模块、cholesky分解模块、上三角矩阵求逆模块和矩阵相乘模块；

所述的主控状态机模块与系统总线相连接，在求逆操作开始前，接收使能拉高，根据待求逆矩阵维度，从外部按行顺序接收原矩阵数据，存放在资源复用模块的存储器中，接收完成后，将接收使能拉低，并通过cholesky分解模块、上三角矩阵求逆模块、矩阵相乘模块分别进行cholesky分解、上三角矩阵求逆、矩阵乘法3个步骤进行cholesky矩阵求逆操作；求逆完成后向总线发出求逆结束信号，通知系统总线从存储模块读取运算结果。

2.如权利要求1所述的基于分时复用技术的cholesky矩阵求逆系统，其特征在于，所述的资源复用模块内部有3个存储器，1个单精度浮点复数乘法器，1个单精度浮点复数累加器；资源复用模块与主控状态机模块、cholesky分解模块、上三角矩阵求逆模块、矩阵相乘模块相连接；根据主控状态机的状态，资源复用内部的存储资源和运算资源受不同模块控制。

3.如权利要求1所述的基于分时复用技术的cholesky矩阵求逆系统，其特征在于，所述的cholesky分解模块是通过对原矩阵元素进行乘法，加法，除法运算将原矩阵A分解为单位下三角矩阵L、对角矩阵D，上三角矩阵L^H；cholesky分解模块内部是一个地址生成状态机和两个单精度浮点除法器，外部与主控状态机和资源复用模块相连，在分解过程中将使用资源复用模块中的存储资源和浮点乘法器，浮点累加器。

4.如权利要求1所述的基于分时复用技术的cholesky矩阵求逆系统，其特征在于，在接收到主控状态机模块输入的cholesky分解起始信号后，cholesky分解模块的内部地址状态机根据矩阵维度开始按顺序生成原矩阵读地址和分解后上三角矩阵的写地址；在分解完成后，通过复用浮点除法器对分解后得到的对角矩阵的所有对角元素求倒数。

5.如权利要求1所述的基于分时复用技术的cholesky矩阵求逆系统，其特征在于，所述的上三角矩阵求逆模块，因choelsky分解后得到的上三角矩阵的对角线元素全为1，所以可以对上三角矩阵求逆公式进行简化，省去除法操作；上三角矩阵求逆模块内部是一个地址生成状态机，外部与主控状态机和资源复用模块相连，在上三角矩阵求逆过程中将使用资源复用模块中的存储资源、浮点乘法器和浮点累加器；在接收到主控状态机模块输入的上三角矩阵求逆起始信号后，上三角矩阵求逆模块的内部地址生成状态机根据矩阵维度开始按顺序生成原上三角矩阵读地址和上三角逆矩阵的写地址。

6.如权利要求1所述的基于分时复用技术的cholesky矩阵求逆系统，其特征在于，所述的矩阵乘法模块内部是一个地址生成状态机，外部与主控状态机和资源复用模块相连，在矩阵乘法模过程中将使用资源复用模块中的存储资源、浮点乘法器和浮点累加器；在接收到主控状态机模块输入的矩阵乘法模块起始信号后，通过矩阵乘法模块内部状态机控制的存储器读写，先进行上三角矩阵与对角阵相乘，再进行上三角矩阵与下三角矩阵相乘。