 Hans Journal of Wireless Communications 无线通信, 2013, 3, 71-76 http://dx.doi.org/10.12677/hjwc.2013.33011 Published Online June 2013 (http://www.hanspub.org/journal/hjwc.html) One-Side Jacobi Matrix Inversion Algorithm and DSP Realization Xi Yang, Zheng Li, Shuai Fang, Tian Zhou, Bin Jiang, Qin Guo Department of Information Science and Engineering, Southeast University, Nanjing Email: seuyangxi@gmail.com, lizhengjustin@gmail.com, seufangshuai@gmail.com, seu.zhoutian@gmail.com, guoqinzju@gmail.com, binjiang@seu.edu.cn Received: May 9th, 2013; revised: May 11th, 2013; accepted: May 20th, 2013 Copyright © 2013 Xi Yang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unre- stricted use, distribution, and reproduction in any medium, provided the original work is properly cited. Abstract: Link adaptive transmission and advanced receiver are two key technologies in broadband wireless communi- cation system. The design and realization of the system both involve a large number of matrix decomposition and inver- sion. The basic way to improve transmission efficiency of broadband wireless communication system is to enhance the efficiency of matrix decompose and inverse computations. For this purpose, this paper develops a kind of one-sided Jacobi algorithm based on classic Jacobi. Since this algorithm has the characteristic of parallelism, it can increase the efficiency at least twice in terms of the instruction execution cycle numbers. This article will first focus on the im- proved one-sided Jacobi algorithm as well as internal architecture and characteristics of DSP TMS320C6474. It then elaborates on how to implement this algorithm in parallel using TI’s real-time multi-tasking operating system kernel (DSP/BIOS). Finally, this paper will compare instruction execution cycle numbers between parallel and serial algorithm under the same accuracy, proving the high efficiency of the improved one-sided Jacobi algorithm. Keywords: One-Side Jacobi Algorithm; Matrix Inversion; TMS320C6474; DSP/BIOS; Parallel Realization 单侧 Jacobi 矩阵求逆算法及其 DSP 实现 阳 析，李 峥，房 帅，周 天，江 彬，郭 骎 东南大学信息科学与工程学院，南京 Email: seuyangxi@gmail.com, lizhengjustin@gmail.com, seufangshuai@gmail.com, seu.zhoutian@gmail.com, guoqinzju@gmail.com, binjiang@seu.edu.cn 收稿日期：2013 年5月9日；修回日期：2013年5月11 日；录用日期：2013 年5月20 日 摘 要：链路自适应与先进接收机是宽带无线通信系统的核心技术，其设计与实现均涉及大量的矩阵分解以及 矩阵求逆运算，提高矩阵分解和矩阵求逆运算的效率是提高宽带无线通信系统传输效能的基本途径。针对此目 的，本文提出一种在经典 Jacobi 算法上改进的单侧 Jacobi 算法。由于该算法具有并行的特性，相比于串行（单 核）实现在指令执行周期数上可提高至少两倍的运行效率。本文首先重点介绍改进的单侧 Jacobi 算法和 TMS320C6474 DSP的内部架构与特性，然后重点阐述结合TI的实时多任务操作系统内核(DSP/BIOS)并行实现 此算法，最后在同样精度的计算结果下比较并行算法与串行算法指令执行周期数，由此验证改进的单侧 Jacobi 算法在并行实现上的高效性。 关键词：单侧 Jacobi算法；矩阵求逆；TMS320C6474；DSP/BIOS；并行 1. 引言 自21 世纪以来，移动通信飞速发展，尤其是最 Copyright © 2013 Hanspub 71  单侧 Jacobi 矩阵求逆算法及其 DSP 实现 近几年，移动通信作为一种便捷、快速的通信方式融 入人们的生活而备受青睐。然而人们对实时数据业务 需求的增长给移动通信的数据速率提出了更高的要 求，高数据速率需求的增加促使3GPP 在2004年底 启动了长期演进计划(Long Term Evolution, LTE)。 3GPP LTE的主要目标包括显著提高终端用户的吞吐 量和扇区容量、降低用户延迟、显著提高用户移动性、 降低网络功耗和用户设备复杂度等。为了达到技术指 标，LTE 引入了许多关键技术，MIMO 技术便为其中 之一。MIMO，即多输入多输出，它是利用多个天线 实现多发多收，在不需要增加频谱资源和天线发送功 率的情况下，成倍地提高信道容量，从而有效提升无 线通信系统的频谱效率，提高传输速率并改善通信质 量。但是如何从收集到的信号中滤除噪声，还原出信 号源发射出的真实信号，则是一个关键问题，MIMO 的信道可以简单地抽象成一个矩阵，因此要解决这个 关键问题就必然涉及到大量的矩阵求逆。由于 MIMO 技术能提供的信道容量和天线数是成正比的，所以矩 阵的维度也随天线数的增加而增加，此时传统的算法 逐渐表现出极低的效率、极差的数值稳定性以及极高 的实现成本，所以能否开发出更为高效的矩阵求逆算 法，成为整套系统中迫切需要解决的核心问题。 Jacobi 算法因其良好的并行性[1-3]，一直是矩阵求 逆算法研究的热点。Jacobi方法的思想就是通过不断 的正交相似变换，使矩阵的非对角线元素逐步趋于0， 从而得到对角阵，此时对角上的元素就是矩阵的特征 值，由每一次的旋转矩阵J累乘可以得到特征向量矩 阵。然而，传统的双边 Jacobi 算法，每次迭代都同时 涉及行和列的更新，行和列之间存在很大的数据相关 性，这种数据相关性在算法的并行实现时将带来很大 的核间通信代价，降低算法效率。 本文首先提出一种改进的单侧 Jacobi 矩阵求逆算 法。紧接着，在简要介绍TMS320C6474 DSP的内部 架构与特性后将重点阐述如何基于DSP/BIOS 并行实 现此算法。最后，本文将分别给出在同样计算结果精 度的情况下，并行算法与串行算法的指令执行周期 数，以证明此改进的单侧 Jacobi算法在并行实现上的 高效。 2. 改进的单侧 Jacobi 算法 单侧 Jacobi 算法由于每次迭代只涉及列的数据更 新，数据相关性大为降低，只需较少的核间通信便能 并行实现，因此各处理器之间可以更加独立的工作， 算法效率也将更高。本小节将主要介绍单侧Jacobi算 法以及对其进行的改进。 以下为改进的单侧Jacobi 算法的具体步骤： 1) 旋转矩阵的求取： 设实对称矩阵 A0 A ，逐次选择修正矩阵 ，使 的两列正交化[4,5]。为减少乘法次数，改进的 单侧 Jacobi 算法的修正矩阵将采用快速 Givens 矩 阵的形式： k R k A,ij k R 100 01 ,, 01 000 ij ij 0 0 0 1 k R (1) 式(1)为快速 Givens 矩阵的形式，主要作用为将矩阵 中元素的值集中到对角线元素上，以利于后续矩阵求 逆的计算，映射到以列为矢量元素的矢量图中相当于 对矩阵中的列矢量进行旋转。其中，参数 计算如下： 设 ,,, ,, kkk kkk iij ji j Iaa Jaaeaa ， 这里 , x y表示向量 和 的内积，若，则说明 两列已经正交，取 xy0e 0 ，否则记 2 ,sgn 2 JI ttt e 1 t (2) 这样取值可以保证相应的旋转角满足 π4 ，以提高算法的数值稳定性。注意到参数 2 sgn 1ttt ，这里涉及到开方运算，复杂度 较高，考虑进行一定的近似，从而降低复杂度。注意 到 关于 t是奇对称的，所以本文仅讨论 的情形 即 0t 21(tt t 0) (3) 利用微积分的知识以及Taylor 展开式，易得： 2 2 11 lim1limlim 2 1 tt tt t tt t (4) Copyright © 2013 Hanspub 72  单侧 Jacobi 矩阵求逆算法及其 DSP 实现 2 2 00 lim1lim 12 tt t tt t (5) 所以，在实际计算时可以作出如下近似(分段点的 选取是随本文要求的精度变化而变化的，此处本文假 设控制误差在 数量级内，得出如下的近似结果)： 4 10 2 2 10 0.25 2 10.25 8 18102 2 01024 t tt tt t t tt 4 k k j (6) 即在实际计算时，采用式(6)的近似式代替式(3)， 可以避免开方运算，有效降低运算复杂度，提高运算 速度。 2) 得到快速 Givens 矩阵 后，即可对目的矩阵 进行旋转变换： k R k A 10,1, k Ak kk AR (7) 注意到在每次旋转变换时都需要计算第 列的 模值，即 I， J的值，常规的单侧Jacobi 算法采用直接 求取列向量模值的方法，共需要 4n 个flop，n为矩阵 维数。受到文献[6,7]的启发，本文提出的改进的单侧 Jacobi 算法将采取一种效率更高的列向量模值计算的 方法，如下： ,ij 从表达式中易知 ，所以 1 1 kk iij kk ji aaa aaa 111 2 111 () 2 , 2 , 2 T kkk kkk ii ijij kkk T kkk kk k jjij ij kkk I aaaaaa IeJ Jaaaa aa IeJ k k (8) 由此，可以在每次旋转变换结束前对列的模值按 照式(1-8)进行更新，这样只需要 12 个flop，仅为原来 的3/n，在n较大的情况下，效率可以得到比较大的 提升。 3) 反复进行上述两步，选择著名的Robin Ring Ordering作为列的正交化顺序，则的各列趋向于两 两正交，即趋向于 Q，而 ，即有 k A k A012 k V=RRRLR A VQ ，即 1 A QV 。 矩阵求逆：将 Q矩阵各列化为单位向量，得 0 QQ 0 ，为单位正交阵，为对角阵，所以 0 Q 1 A QV ，则 1111 00 T A VQ VQ (9) 此即完成了矩阵求逆。综上所述，改进的单侧 Jacobi 算法主要在三个方面提高了算法效率：第一， 采用快速Givens 旋转矩阵作为修正矩阵，节省了一半 的乘法运算量；第二，采用极限和泰勒展开式近似， 减少复杂的开方运算；第三，在计算列向量的模值时 采用迭代的方法，对列向量模值进行更新。 3. TMS320C6474的特性和 DSP/BIOS 本文将主要采用多核 DSP TMS320C6474来并行 实现上一节所提出的改进的单侧Jacobi 算法，因此本 小节将根据算法具体实现时所涉及到的 DSP 具体部 件简要介绍TMS320C6474 DSP 的内部架构与特性。 由上一节算法的具体实施步骤可知，对于高维矩 阵的求逆，单侧 Jacobi算法的计算量主要集中于对矩 阵的列的排序和每两列之间的数据处理，这就要求 DSP 有较高的执行能力、足够的存储空间以及较快的 存储器访问速度。如图 1(C6474 内部结构简略图，主 要示出算法的 DSP 实现中所需模块)所示，TMS320C 6474 包括三个 TMS320C64x + TM DSP 核，总主频可 达3.6 GHz，性能高达 28800MIPS (Million Instructions 64x+ DSP Co re 20-30% Bette Code Size r 20% Better Cycle Performan 3 MB L2 32KB L1P Cache Switched Central Res ource (SCR) EDMA3 SRIO DDR2-667 TCP2 VCP2 32KB L1D Cache 128 128 128 64 64 256 128 Figure 1. Internal architecture of TMS320C6474 图1. TMS320C6474的内部架构 Copyright © 2013 Hanspub 73  单侧 Jacobi 矩阵求逆算法及其 DSP 实现 Per Second)且3个核之间相对独立，每个核可以执行 不同的程序。片上共有 8个功能单元，每个时钟周期 执行一条指令，最高情况下可同时执行八条指令；此 外C6474 DSP 还集成了大量片上存储器，称之为三级 存储系统：第一级为Cache，包括 32KB 的一级程序 存储器(L1P)和32KB 的一级数据存储器(L1D)；第二 级(L2)存储器中三个核共享 3MB 的存储空间，有两 种不同的配置方式，速度比第一级慢，但是容量大； 第三级存储器为通过DDR接口扩展外部 DDR 存储 器。这些特性为单侧 Jacobi算法的并行实现提供了必 要的物理基础。 由于算法是并行实现，这就要求同时建立多个线 程，并对多线程之间的执行顺序进行正确的调度，此 外还需对有限的资源进行合理分配。DSP/BIOS 是TI 公司为其 TMS320C6000TM 系列 DSP 平台所设计开 发的一个尺寸可裁剪的实时多任务操作系统内核，它 主要由三部分组成：多线程实时内核(抢占式多线程)， 实时分析工具和芯片支持库。它以模块化的方式提供 给用户对线程(TSK)、中断(HWI，SWI)、定时器(CLK)、 内存资源(MEM)、所有外设资源的管理能力，这些资 源都可以根据需要剪裁。 4. 单侧 Jacobi 算法的 DSP 实现 论文在前面两小节详细讲述单侧 Jacobi 算法以及 DSP TMS320C6474的内部架构与主要特性后，本小 节将重点阐述如何结合TI的实时多任务操作系统内 核(DSP/BIOS)在C6474 上并行实现此改进的单侧 Jacobi 算法。 在DSP中主要有两种并行处理模式：Master/Slave 模式和 Data Flow 模式。 Master/Slave 模式：表现为由一个处理器中央控 制多个处理器共同完成，中央控制处理器需要调度任 务的执行线程和相关资源，并将线程分配给空闲的处 理器来完成；此种模式主要适用于少量数据运算，大 量任务调度的应用程序，一般和相应的操作系统相结 合，它的缺点是控制处理器之间通信和调度的代码量 较大，导致整个程序所需存储空间大，它难点则在于 怎样保证在线程很多情况下的实时性。 Data Flow模式：当一个处理器对数据完成相应 的操作后，将数据传递给下一个处理器以便对数据进 行进一步的处理；此模式适合于需要大量的复杂的数据 运算的程序，它的难点在于怎样将程序分解成合适的 处理流程以便实现流水线操作和较高的数据吞吐率。 考虑到算法的特点，本文最终采用结合了Master/ Slave 模式和 Data Flow模式二者特点的混合模式。算 法实现中，将要被处理的数据存放在DSP 的DDR2 memory 中，主核中建立有多个线程(Task)，由 各自 的 信号量(semaphore)控制执行(running)和暂停(block- ing)；主核(核0)通过 IPC module控制其它两个副核(核 1、核 2)对不同数据同时进行操作，两个副核完成相 应操作后反馈信息给主核，主核则根据信息判断目前 所处的状态并进行下一步的操作的部署直至算法完 成后在 message_log 中输出所得结果。同时程序利用 clock()；函数对算法的执行时间进行了较精确的统计， 以用于算法之间效率的比较。 如图 2所示，首先在 DDR2 memory中开辟出一 块共属于三个核的公共数据的存储区域。然后由主核 (核0)控制系统的初始化工作，包括初始化 DDR2 memory，载入需要处理的相关数据，启动 clock()；函 数等。一切准备工作就绪后，主核通过核间通信模块 (IPC module)通知另外两个副核开始进行数据处理， 并行的单侧Jacobi 算法开始执行。 核0 TSK0、TSK1、TSK2 核1 TSK0 核2 TSK0 DSP/B IOS DDRmemory 存储公共数据 系统初始化 ① ① ② ② ③③ ④ ④ clock(); t ① 主核发送控制信号； ② 从memory 中读取数据； ③ 将数据写到memory 中； ④ 副核反馈信息给主核； Figure 2. DSP implementation schematic 图2. DSP具体实现示意图 Copyright © 2013 Hanspub 74  单侧 Jacobi 矩阵求逆算法及其 DSP 实现 Copyright © 2013 Hanspub 75 3个核的详细功能为： 行周期数， 以证明此改进的单侧Jacobi 算法在并行实 现上的高效。 核0：主核，控制整个算法的执行、计算整个算 法执行所需的 cycle 数以及输出最后的结果。在 main(); 函数使能IPC中断后执行 TSK0，TSK0 的主要功能为 提示整个算法开始执行，并抛出信号量给 TSK1进而 真正开启整个算法。TSK1 则是控制中心，在对所要 处理的矩阵的列号进行排序后，通过一定的算法每次 将特定的两列的序号分别发送给其它两个副核以进 行进一步的操作，遍历完后，对矩阵的列号再一次排 序并发送，以此类推直至下一次循环即完成一轮，经 过几轮迭代后整个算法执行完毕，激活 TSK2，TSK2 负责对最后所得结果进行一定的处理以输出最后的 结果。 串行实现是指仅使用 C6474 的一个核执行该算 法，并行实现则是指采用上一节所述的思想由C6474 的三个核相互合作共同执行该算法。用clock();函数分 别对算法的完成所需的cycle 数进行 统计，并以 cycle 数与结果的精度为主要指标对 串行与并行的效率进行比较。 在此本文以 8 × 8维的浮点数矩阵为例，所处理 的浮点数矩阵如下图 3： 根据算法的特性，当增加迭代的次数时结果的精 度会不断提高，直到迭代的次数大于某一界限时(该界 限与所处理的矩阵的维数有关，矩阵维数越大，所需 迭代次数越多)，采用控制变量法对效率进行比较，即 在保持其他条件如矩阵维数、数据精度等不变的情况 下，增加迭代次数，所得结果如下所示： 核1和核 2：两个副核的操作除所处理的数据不 同外，其它均相同。这两个核是算法执行的主体，算 法的实现由这两个核共同完成。副核的主要操作均在 TSK0 中，由信号量 SEM0进行控制。TSK0 的主要任 务为：根据接收到的主核发来的序号对相应列进行 Givens旋转变换，以使最后得到的矩阵各列两两正 交，操作完成后反馈主核。 为更加清楚地体现并行实现的优势，将表 1的数 据以直方图的形式展示出来，见图 4。从图 4中可以 看出，当迭代次数较少时，串行实现与并行实现所需 的cycle 数差距不大，且并行实现的cycle数有时会出 现比串行实现所需的 cycle数多的情况；但当迭代次 数增加时，并行实现的优势逐渐体现出来，当迭代次 数为 20 次时，串行的效率就只是并行的69.7%，随着 迭代次数的加大，并行的效率更是不断上升。 5. 串/并行效率比较 本小节将在上一节的基础上比较在同样计算结 果精度的情况下，并行算法与串行算法所需的指令执 Figure 3. Experimental floating-point 图3. 实验浮点矩阵 Table 1. Serial and parallel implementations with cycle number 表1. 串并行实现所用 cycle 数 迭代次数 串行实现所用 cycle 数 并行实现所用cycle 数 8 3278268 4247376 20 8104920 5653764 50 19944208 9191054 100 39415158 15100850 200 78025414 26880504 0 2000 000 0 4000 000 0 6000 000 0 8000 000 0 20 50 200 串行 并行 8100 Figure 4. Serial and parallel im p l em e ntations with cycle number comparison chart 图4. 串并行实现所用 cycle 数对比图  单侧 Jacobi 矩阵求逆算法及其 DSP 实现 6. 总结 综上所述，本文在提出一种改进的单侧 Jacobi 矩 阵求逆算法后，利用 TI 的TMS320C6474 DSP 和 DSP/BIOS 并行实现了此算法。最后为证明此改进的 单侧 Jacobi 算法在并行实现上的高效，本文比较了在 同样计算结果精度的情况下，并行算法与串行算法的 指令执行周期数。由本文可见，在对单侧Jacobi算法 进行改进并采用并行实现后，算法的执行效率有了大 幅度的提高，矩阵的求逆也更为高效和快速，这将为 MIMO 通信中高维数矩阵的求逆提供了一种可能的 解决方法。 但由于改进的算法采用了近似，因此当矩阵维数 增加至百维以上时，迭代次数急剧增加，近似所带来 的误差也会增加，由此误差带来的逆矩阵的非正交性 将是本文的下一步研究目标。此外，采用汇编语言来 实现并行算法，进一步提高运行效率也是本文的努力 方向。 7. 致谢 本论文是基于东南大学国家大学生创新性实验 计划项目《宽带无线通信系统中高效矩阵运算及其 DSP 实现》的实验成果撰写的，项目编号为 111028609。 感谢东南大学信息科学与工程学院在项目的研究过 程中提供了实验环境与实验平台，同时感谢项目指导 老师金石老师与江彬老师，他们在算法的改进过程中 提供了很有价值的建议和指导。 参考文献 (References) [1] J. Demmel, K. Veselic. Jacobi’s method is more accurate than QR. SIAM Journal on Matrix Analysis and Applications, 1992, 11: 1204-1246. [2] A. H. Sameh. On Jacobi and Jacobi-like algorithms for a parallel computer. Mathematics of Computation, 1971, 25(115): 579- 590. [3] R. P. Brent, F. T. Luk. The solution of singular-value and sym- metric eigenvale problems on multiprocessor arrays. SIAM Journal on Scientific and Statistical Computing, 1985, 6(1): 69- 84. [4] M. Gentleman. Least squares computations by Givens transfor- mations without square roots. Journal of the Institute of Mathe- matics and its Applications, 1973, 12(3): 329-336. [5] A. A. Anda, H. Park. Fast plane rotations with dynamic scaling. SIAM Journal on Matrix Analysis and Applications, 1994, 15(1): 162-174. [6] J. Yan, S. Osher. A local discontinuous Galerkin method for directly solving Hamilton-Jacobi equations, Journal of Compu- tational Physics, 2011, 230(1): 232-244. [7] M. Baggio, J. de Boer and K. Holsheimer. Hamilton-Jacobi renormalization for Lifshitz spacetime. Journal of High Energy Physics, 2012. http://arxiv.org/abs/1107.5562 Copyright © 2013 Hanspub 76 |