面向gpdsp的矩阵lu分解向量化计算的方法
【技术领域】
[0001]本发明主要涉及通用计算数字信号处理器(General-Purpose Digital SignalProcessor,简称GPDSP)，特指一种适用于GPDSP的矩阵LU分解向量化计算的方法。
【背景技术】
[0002]稠密线性方程组求解是高性能计算和科学计算领域内最普遍的计算应用，而矩阵LU分解(LU Factorizat1n, LU)是求解稠密线性方程组最常用的一种方法，尤其是大规模稠密线性方程组求解。高性能基准测试程序(High Performance Linpack，HPL)是T0P500最重要的测试基准，HPL通过求解稠密线性代数方程组，以测试和评价高性能计算机系统的浮点性能。HPL通过矩阵分块LU分解将大规模的矩阵分割为子块矩阵的LU分解计算，矩阵LU分解后再通过两次三角方程组求解实现最终的方程组求解，其中矩阵LU分解的计算量占据整个HPL计算量的95%以上。因此，优化矩阵LU分解的性能对提高HPL的效率具有非常重要的应用价值。
[0003]在专利申请号为201310725118.6的文献(处于实审阶段)中提供了一种通用计算数字信号处理器(General-Purpose Digital Signal Processor，简称 GPDSP)，它包含CPU核单元和DSP核单元，CPU核单元主要用于负责包括存储管理、文件控制、进程调度、中断管理任务在内的通用事务管理以及提供对通用操作系统的完整支持；DSP核单元包含若干强大计算能力的64位向量处理阵列，用于支持高密集运算任务的解算。传统的面向Cache结构的矩阵LU分解方法不适合GPDSP的非Cache的向量阵列存储访存模式和向量处理阵列并发向量处理的体系结构特征，难以发挥GPDSP向量计算优势。
【发明内容】
[0004]本发明要解决的技术问题就在于:针对现有技术存在的技术问题，本发明提供一种原理简单、操作方便、能充分利用GPDSP的DSP核向量处理阵列的强大并行计算、高带宽向量数据加载能力，显著提高DSP核计算访存比的面向GPDSP的矩阵LU分解向量化计算的方法。
[0005]为解决上述技术问题，本发明采用以下技术方案:
[0006]一种面向GPDSP的矩阵LU分解向量化计算的方法，其步骤为:
[0007]S1:根据GPDSP的体系结构特征确定最佳的LU分解的矩阵规模N值；
[0008]S2:GPDSP的DSP核通过DMA从片外DDR存储器将需要处理的矩阵数据传输到片内共享存储阵列中；
[0009]S3:GPDSP的DSP核按照列选主元方法，计算第i列的主元值以及对应的列元素序号值，初始i = O ;
[0010]S4:GPDSP的DSP核根据上述计算得到的列^的列主元P 及对应的列元素序号值Vi对矩阵A内列主元P i所在的行与列首元I ti]所在的行进行交换；
[0011]S5:GPDSP的DSP核对列^按照公式I I ^1^进行列消元计算，更新列I i;
[0012]S6:GPDSP的DSP核按照L，尸L，「I’ ^ui进行更新矩阵panel的计算；
[0013]S7:判断i是否等于N-1，若不是，令i = i+Ι，转步骤S3，若是转步骤S8 ;
[0014]S8:矩阵A的LU分解计算完毕；DSP核使用DMA将矩阵A从片内共享存储阵列传输到片外DDR存储器的原存储位置。
[0015]作为本发明的进一步改进:在步骤SI中，设所处理的LU分解的矩阵规模为NXN阶非奇异矩阵，每个矩阵元素的数据为u字节；上述矩阵规模的N值的确定方法是:设条件
(I)2*u*p*s*N彡q ； (2)u*N*N彡r ;按照满足条件(I) (2)取得的最大N值为所需要的N值。
[0016]作为本发明的进一步改进:在步骤S2中，所述片内共享存储阵列中的矩阵标记为矩阵A，其规模为NXN阶，首地址记为AAcitl;沿左上至右下的对角线方向，将矩阵A划分为上三角区域，标记为C，以及下三角区域，标记为D ;将下三角矩阵D按照从左至右的方向按列分别记为L、，将每一列的在对角线上的数据，即列首元标记为1M、
I [N-1] °
[0017]作为本发明的进一步改进:在所述步骤S5中，将片内共享存储阵列矩阵A中，列首元Im右下方的矩阵记为L’ i，矩阵A中第i行位于上三角C的部分记为Ui,列首元1?」右方包括Ui部分的矩阵记SLi，列Ii不包括列首元I ti]的部分记为列I’ -
[0018]作为本发明的进一步改进:所述步骤S3的详细流程为:
[0019]S3.1:DSP核通过DMA将列^从片内共享存储阵列中传输至片内向量阵列存储器；设初始值为全O的向量寄存器Z2和向量寄存器Z5，以及初始值为{0，1，2，一,P-1}的向量寄存器Z4;
[0020]S3.2:DSP核的向量处理阵列依次通过向量LOAD指令加载p个元素，存入向量寄存器Ztl;
[0021]S3.3:DSP核的向量处理阵列对上述向量数据Ztl进行向量绝对值操作，结果存入向量寄存器Z1;
[0022]S3.4:DSP核的向量处理阵列对上述向量数据Z# Z 2的值进行向量比较操作，将比较结果存入向量寄存器Z3;
[0023]S3.5:根据上述比较结果值Z3，DSP核的向量处理阵列使用向量MOV操作将Z1*比Z2数值更大的数据值更新向量寄存器Z 2;dsp核的向量处理阵列将数值更大的数据所对应存储在&中的列元素序号替换存入向量寄存器z 5中；
[0024]S3.6:DSP核的向量处理阵列使用向量加法操作将向量寄存器Z4的向量数值加P ;
[0025]上述步骤中向量寄存器Z4存储所处理列相对应的列元素序号；
[0026]重复步骤S3.1至步骤S3.6，直到完成该列的所有元素计算，最终得到p个列元素值以及对应的P个列元素序号值；
[0027]S3.7:DSP核的向量处理阵列对上述P个列元素值进行数值大小的比较，计算出绝对值最大的列元素值作为列^的列主元，标记为P i，对应的列元素序号值标记为
[0028]作为本发明的进一步改进:所述步骤S4的详细流程为:
[0029]S4.1:根据列主元PdA列元素序号值V i，以及矩阵A在片内共享存储阵列中的首地址AA。。，计算得出列主元PjZf在行的首地址AP i= AA oo+u* (Vi+i) *N ；
[0030]S4.2:根据列首元Im的列元素序号值，以及矩阵A在片内共享存储阵列中的首地址AAcic!，计算得出列首元1[:1]所在行的首地址Al ti]= AA 00+u*i*N ；
[0031]S4.3:将片内向量阵列存储器划分为上下两片半区，半区的大小为q/2字节，DSP核使用DMA将列主元PiK在行共N个元素从片内共享存储阵列中传输至片内向量阵列存储器的上半区；
[0032]S4.4:DSP核使用DMA将列首元I ti]所在行共N个元素从片内共享存储阵列中传输至片内向量阵列存储器的下半区；
[0033]S4.5:DSP核使用DMA将上述片内向量阵列存储器的下半区的N个元素数据传输至片内共享存储阵列中以APiS首地址的位置；
[0034]S4.6:DSP核使用DMA将上述片内向量阵列存储器的上半区的N个元素数据传输至片内共享存储阵列中以Alm为首地址的位置。
[0035]作为本发明的进一步改进:所述步骤S5的详细流程为:
[0036]S5.1:DSP核使用DMA将列Ii除列首元素I⑴以外的共N_i_l个元素数据从片内共享存储阵列传输到片内向量阵列存储器；
[0037]S5.2:DSP核使用标量LOAD指令读取片内共享存储阵列内的列首元素Im，存入标量寄存器S。；
[0038]S5.3:DSP核使用标量浮点除法指令计算l/lm值，结果存入标量寄存器S 0,并使用标量广播指令将Stl的数据广播至向量寄存器Z 0；
[0039]S5.4:DSP核使用向量LOAD指令依次加载上述已经传输至片内向量阵列存储器的该列的其他元素数据，每次加载P个元素，存入向量寄存器Z1;
[0040]S5.5:DSP核使用向量乘法指令将向量数据Z1与向量数据Z ^相乘，结果存入向量寄存器Z1;
[0041]S5.6:DSP核使用向量STORE指令将Z1中的向量数据结果存入片内向量阵列存储器中的原位置；
[0042]重复步骤S5.4到步骤S