CN104615516A - 面向GPDSP的大规模高性能Linpack测试基准实现的方法 - Google Patents

Abstract

本发明公开了一种面向GPDSP的大规模高性能Linpack测试基准实现的方法，采用分块的方式来求解稠密线性方程组：较小的串行任务和小矩阵求逆和三角方程求解任务，由GPDSP中的CPU核负责；规整的计算任务由GPDSP中的DSP核负责。本发明具有原理简单、操作方便、能充分利用GPDSP中DSP核向量处理阵列的强大向量并行计算和高带宽向量数据加载能力、由CPU核和DSP核协同高效的实现大规模Linpack测试基准计算的优点。

Description

面向GPDSP的大规模高性能Linpack测试基准实现的方法

技术领域

本发明主要涉及通用计算数字信号处理器(General-Purpose Digital Signal Processor，简称GPDSP)，特指一种适用于GPDSP的大规模高性能Linpack测试基准实现的方法。

背景技术

Linpack测试基准是高性能计算中最广泛使用的性能测试基准，是TOP500最重要的测试基准。Linpack基准通过求解稠密线性代数方程组，以测试和评价高性能计算机系统的浮点性能。Linpack测试标准包括三个不同数据规模的测试：100×100，1000×1000和n×n，其中前两个测试所使用的测试基准程序是可以从网站(例如：http：//www.netlib.org/benchmark/)下载，编译运行程序后，该程序会给出相应机器的性能，这种测试不允许对测试程序进行任何修改。n×n数据规模的测试要求是Linpack测试标准中最为宽松的，用户可以对任意大小的问题规模，使用任意个数的CPU，使用基于高斯消去法的各种优化方法来执行该测试程序，寻求最佳的性能测试结果。HPL(High Performance Linpack)是第一个标准的公开发布的并行Linpack测试软件包，是n×n测试的MPI实现，可移植于多种体系结构的计算机系统，目前广泛用于TOP500测试，该测试主要针对分布式存储的大规模并行计算系统而设计。

在专利申请号为201310725118.6的文献(处于实审阶段)中提供了一种通用计算数字信号处理器(General-Purpose Digital Signal Processor，简称GPDSP)，它包含CPU核单元和DSP核单元，CPU核单元主要用于包括存储管理、文件控制、进程调度、中断管理任务在内的通用事务管理以及提供对通用操作系统的完整支持；DSP核单元包含若干强大计算能力的64位向量处理阵列，用于支持高密集运算任务的解算。由于GPDSP是包含CPU核和DSP核的异构多核处理器，包括寄存器文件、片内向量阵列存储器、片内共享存储阵列、片外DDR存储器等多层次存储架构，简单的移植HPL或采用传统的方法实现的Linpack测试基准，不能够取得GPDSP最佳的性能测试结果。

发明内容

本发明要解决的技术问题就在于：针对现有技术存在的技术问题，本发明提供一种原理简单、操作方便、能充分利用GPDSP中DSP核向量处理阵列的强大向量并行计算和高带宽向量数据加载能力，由CPU核和DSP核协同高效的实现大规模Linpack测试基准计算的方法。

为解决上述技术问题，本发明采用以下技术方案：

一种面向GPDSP的大规模高性能Linpack测试基准实现的方法，采用分块的方式来求解稠密线性方程组：较小的串行任务和小矩阵求逆和三角方程求解任务，由GPDSP中的CPU核负责；规整的计算任务由GPDSP中的DSP核负责。

作为本发明的进一步改进：具体步骤为：

S1：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值和分块大小NB值；Linpack测试基准求解的稠密线性方程组为Ax＝b，其中A为n×n阶的非奇异矩阵，x,b为n×1阶的向量；

S2：根据上述n值，由GPDSP的CPU核在片外DDR存储器为矩阵A和向量x,b分配空间，并生成初始化数据；记Linpack计算量为计算GPDSP的峰值性能为F₀，记CPU核的当前系统时间为t₀；

S3：GPDSP的CPU核和DSP核协同计算，按照NB×NB阶的分块规模对n×n阶的矩阵进行LU分解；

S4：根据步骤S3得到的LU分解结果，GPDSP的CPU核和DSP核协同计算，求解方程LUx＝b；

S5：记CPU核的当前系统时间为t₁；由GPDSP的CPU核计算出GPDSP的实测性能为F₁＝ops/(t₁-t₀)，Linpack效率为F₁/F₀，并输出上述计算结果，Linpack测试完毕。

作为本发明的进一步改进：所述步骤S4中求解方程LUx＝b分两步：(1)GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ly＝b；(2)GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ux＝y。

作为本发明的进一步改进：所述步骤S4的具体步骤为：

S4.1：GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ly＝b；

S4.2：根据GPDSP的体系结构特征确定求解上述三角方程所需的分块大小NC值；

上述求解三角方程所需的分块大小NC值的确定方法是：设条件(1)16*NC*p*q≤s；(2)8*(NC*NC+(r+2)*NC)≤t；(3)NC是(p*q)的整数倍，n是NC的整数倍；按照满足条件(1)(2)(3)取得的最大NC值为所需要的NC值；

S4.3：设下三角矩阵L按照上述NC值进行分块，记j＝n/NC；沿对角线的三角分块依次标记为T_i，其中1≤i≤j，以该分块为坐标，对应列方向的分块依次标记为S_i+1,i,S_i+2,i,...,S_j,i；将y,b分别划分为长度为NC的子向量，分别有j个子向量，分别标记为y₁,y₂,...,y_j；b₁,b₂,...,b_j；

S4.4：GPDSP的CPU核将三角分块T₁和子向量b₁加载到片内共享存储阵列中，计算得到计算结果y₁传出到片外DDR存储器中，同时保持在片内共享存储阵列中；

S4.5：GPDSP的r个DSP核通过DMA依次从片外DDR存储器将三角分块T₁所在列的矩阵分块S_2,1,S_3,1,...,S_r+1,1和子向量b₁后续的r个子向量b₂,b₃,...,b_r+1传入计算，其中S_2,1,S_3,1,...,S_r+1,1传入到DSP核的片内向量阵列存储器，子向量b₂,b₃,...,b_r+1传入到片内共享存储阵列中；分别执行计算：b₂＝b₂-S_2,1*y₁,b₃＝b₃-S_3,1*y₁,...,b_r+1＝b_r+1-S_r+1,1*y₁，计算后的结果b₂,b₃,...,b_r+1保持在片内共享存储阵列中；重复上述计算，直到该列方向所有的矩阵分块计算完毕；

S4.6：GPDSP的CPU核将下一个三角分块T₂加载到片内共享存储阵列中，计算得到计算结果y₂传出到DDR存储器中，同时保持在片内共享存储阵列中；

S4.7：重复步骤S4.5和步骤S4.6，直到完成所有的子向量y₁,y₂,...,y_j的计算。

作为本发明的进一步改进：所述步骤S3的具体步骤为：

S3.1：GPDSP的CPU核对分块矩阵K_i进行LU分解；由GPDSP的CPU核将分块矩阵K_i加载至片内共享存储阵列中，按照列选主元方法对分块矩阵K_i进行LU分解，LU分解结果为K_i＝LL_i*UU_i,其中LL_i为下三角矩阵，UU_i为上三角矩阵；将LU分解的计算结果传输到片外DDR存储器的原矩阵位置；

S3.2：GPDSP的CPU核计算步骤S3.1得到的LU分解的上三角矩阵UU_i的逆矩阵；由GPDSP的CPU核计算上三角矩阵UU_i的逆矩阵UU_i ^-1，计算结果保留在片内共享存储阵列中；

S3.3：GPDSP的DSP核对行panel的行矩阵U_i进行行交换；由GPDSP的DSP核根据步骤3.1得到的行交换信息按照数据并行的方法对行矩阵U_i并行的进行行交换计算；

其中，步骤S3.2和步骤S3.3在步骤S3.1之后由GPDSP的CPU核和DSP核并行执行；

S3.4：GPDSP的CPU核计算步骤S3.1得到的LU分解的下三角矩阵LL_i的逆矩阵；由GPDSP的CPU核计算下三角矩阵LL_i的逆矩阵LL_i ^-1，计算结果保留在片内共享存储阵列中；

S3.5：GPDSP的DSP核更新列panel的列矩阵L_i；

按照 $L_i=L_i*{\mathrm{UU}}_i^{-1}={\left({(L_i*{\mathrm{UU}}_i^{-1})}^T\right)}^T={({\left({\mathrm{UU}}_i^{-1}\right)}^T*L_i^T)}^T$ 其中为下三角矩阵，将列矩阵L_i的更新由矩阵与上三角矩阵的乘积转化为下三角矩阵与矩阵的乘积，由GPDSP的DSP核并行执行计算；DSP核执行下三角矩阵与矩阵的乘积计算时，由标量处理器读取的下三角矩阵数据为步骤S3.2计算得到的保留在片内共享存储阵列中的逆矩阵UU_i ^-1转置得到，由向量阵列处理器读取列矩阵L_i的矩阵数据；矩阵数据传入时按照矩阵转置的方式，计算后的结果矩阵按照矩阵转置的方式传输到片外DDR存储器的原矩阵位置；

步骤S3.4和步骤S3.5在步骤S3.3之后由GPDSP的CPU核和DSP核并行执行；

S3.6：GPDSP的DSP核更新行panel的行矩阵U_i；

按照其中为下三角矩阵，行矩阵U_i的更新为下三角矩阵与矩阵的乘积，由GPDSP的DSP核并行执行计算；DSP核执行下三角矩阵与矩阵的乘积计算时，由标量处理器读取的下三角矩阵数据为步骤S3.2计算得到的逆矩阵LL_i ^-1，由向量阵列处理器读取行矩阵U_i的矩阵数据，计算后的结果矩阵传出到片外DDR存储器的原矩阵位置；

S3.7：GPDSP的DSP核更新尾矩阵A_i；按照A_i＝A_i-L_i*U_i更新尾矩阵A_i；

执行更新尾矩阵计算时，采用对列矩阵L_i分块的方法，每次由GPDSP的一个DSP核的标量处理器先将L_i中的一个NB×NB阶的分块矩阵Q_j加载至片内共享存储阵列中；再由DSP核将后续计算所需的行矩阵U_i矩阵数据依次传输到片内向量阵列存储器，并且按照公式A_ij＝A_ij-Q_j*U_i执行矩阵与矩阵的乘加计算；上述计算依次进行，直到完成尾矩阵A_i的全部更新；

S3.8：重复步骤S3.1和步骤S3.7，直到完成全部分块K_i的LU分解。

作为本发明的进一步改进：所述步骤S1的具体步骤为：

S1.1：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵LU分解的分块大小NB值；上述矩阵LU分解的分块大小NB值的确定方法是：设条件(1)16*NB*p*q≤s；(2)8*NB*NB≤t；(3)NB是r的整数倍；按照满足条件(1)(2)(3)取得的最大NB值为所需要的NB值；

S1.2：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值；上述矩阵规模n值的确定方法是：设条件(1)n是NB的整数倍；(2)8*(n*n+2*n)<u；按照满足条件(1)(2)取得的最大n值为所需要的n值。

与现有技术相比，本发明的优点在于：本发明的面向GPDSP的大规模高性能Linpack测试基准实现的方法，采用分块的方式求解稠密线性方程组，求解过程中，将较小的串行任务以及不适合DSP核向量处理阵列向量化计算的小矩阵求逆和三角方程求解等计算量相对较小的计算任务分派到CPU核处理，将适合DSP核向量处理阵列向量化计算的规整矩阵乘法类计算任务分派到DSP核处理，后类任务的计算量占据Linpack计算的绝大部分。同时，CPU核与DSP核的计算尽可能地并行执行，减少总的执行时间。另外，本发明根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值和分块大小NB值，能显著提高GPDSP的性能测试结果。综上所述，本发明的方法充分发挥了GPDSP中DSP核向量处理阵列的强大向量化计算能力和CPU核的灵活控制能力，能够取得GPDSP的最佳性能测试结果。

附图说明

图1是本发明在具体应用实例中所面向的GPDSP计算系统的简化存储模型示意图。

图2是本发明方法的流程示意图。

图3是本发明在具体应用实例中将n×n阶的矩阵按照NB×NB阶的规模进行分块的示意图。

图4是本发明在具体应用实例中进行LU分解的示意图。

图5是本发明在具体应用实例中GPDSP的CPU核和DSP核协同求解三角方程Ly＝b的流程示意图。

具体实施方式

以下将结合说明书附图和具体实施例对本发明做进一步详细说明。

本发明的面向GPDSP的大规模高性能Linpack测试基准实现的方法，是采用分块的方式来求解稠密线性方程组。在求解过程中，较小的串行任务以及不适合DSP核向量处理阵列向量化计算的小矩阵求逆和三角方程求解任务，这类任务计算量相对较小，计算数据存放模式不适合高效向量化处理，由CPU核负责。对于规整的计算任务(如矩阵与矩阵乘法)，适合DSP核向量处理阵列向量化计算，这类任务的计算量占据Linpack计算的绝大部分，由DSP核负责。本发明提供的计算任务分割方法能够充分发挥DSP核向量处理阵列的强大向量化计算能力和CPU核的灵活控制能力。通过相关文献研究表明，分块NB值对Linpack的效率影响很大，本发明根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值和分块大小NB值，能够取得GPDSP的最佳性能测试结果。

如图1所示，是本发明所面向的GPDSP计算系统的简化存储模型示意图。GPDSP计算系统包括CPU核和DSP核，DSP核包含若干64位向量处理阵列计算单元，存储系统包括DSP核专用的片内向量阵列存储器，CPU核和DSP核共享的片内共享存储阵列、大容量的片外DDR存储器。

在具体应用实例中，设GPDSP中CPU核数量为m个，DSP核数量为r个；DSP核的向量处理阵列计算单元数量为p个，每个计算单元的MAC(乘加部件)数量为q个，DSP核的片内向量阵列存储器容量为s字节；GPDSP的片内共享存储阵列容量为t字节，GPDSP的片外共享DDR存储容量为u字节。Linpack测试基准求解的稠密线性方程组为Ax＝b，其中A为n×n阶的非奇异矩阵，x,b为n×1阶的向量，数据类型为双精度浮点。

如图2所示，本发明的具体流程为：

S1：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值和分块大小NB值；

S1.1：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵LU分解的分块大小NB值；

上述矩阵LU分解的分块大小NB值的确定方法是：设条件(1)16*NB*p*q≤s；(2)8*NB*NB≤t；(3)NB是r的整数倍；按照满足条件(1)(2)(3)取得的最大NB值为所需要的NB值。

S1.2：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值；

上述矩阵规模n值的确定方法是：设条件(1)n是NB的整数倍；(2)8*(n*n+2*n)<u；按照满足条件(1)(2)取得的最大n值为所需要的n值。

S2：根据上述已确定的n值，由GPDSP的CPU核在片外DDR存储器为矩阵A和向量x,b分配存储空间，并生成初始化数据。

记Linpack计算量为计算GPDSP的峰值性能为F₀，记CPU核的当前系统时间为t₀。

S3：GPDSP的CPU核和DSP核协同计算，按照NB×NB阶的分块规模对n×n阶的矩阵进行LU分解。

如图3所示，将n×n阶的矩阵按照NB×NB阶的规模进行分块，沿左上至右下的对角线方向，依次取NB×NB阶的矩阵分块，共有k＝n/NB块；设分块依次标记为K_i(1≤i≤k)，以该分块为坐标，行方向的行panel(NB×(k-i)NB阶的行矩阵)记为U_i，列方向的列panel((k-i)NB×NB阶的列矩阵)记为L_i，右下方向的尾矩阵((k-i)NB×(k-i)NB阶的尾矩阵)记为A_i。

如图4所示，LU分解的具体实施步骤如下：

S3.1：GPDSP的CPU核对分块矩阵K_i进行LU分解。由GPDSP的CPU核将分块矩阵K_i加载至片内共享存储阵列中，按照列选主元方法对分块矩阵K_i进行LU分解，LU分解结果为K_i＝LL_i*UU_i,其中LL_i为下三角矩阵，UU_i为上三角矩阵。将LU分解的计算结果传输到片外DDR存储器的原矩阵位置。

S3.2：GPDSP的CPU核计算步骤S3.1得到的LU分解的上三角矩阵UU_i的逆矩阵。由GPDSP的CPU核计算上三角矩阵UU_i的逆矩阵UU_i ^-1，计算结果保留在片内共享存储阵列中。

S3.3：GPDSP的DSP核对行panel的行矩阵U_i进行行交换。由GPDSP的DSP核根据步骤3.1得到的行交换信息按照数据并行的方法对行矩阵U_i并行的进行行交换计算。

其中，步骤S3.2和步骤S3.3在步骤S3.1之后由GPDSP的CPU核和DSP核并行执行。

S3.4：GPDSP的CPU核计算步骤S3.1得到的LU分解的下三角矩阵LL_i的逆矩阵。由GPDSP的CPU核计算下三角矩阵LL_i的逆矩阵LL_i ^-1，计算结果保留在片内共享存储阵列中。

S3.5：GPDSP的DSP核更新列panel的列矩阵L_i。

按照 $L_i=L_i*{\mathrm{UU}}_i^{-1}={\left({(L_i*{\mathrm{UU}}_i^{-1})}^T\right)}^T={({\left({\mathrm{UU}}_i^{-1}\right)}^T*L_i^T)}^T$ 其中为下三角矩阵，将列矩阵L_i的更新由矩阵与上三角矩阵的乘积转化为下三角矩阵与矩阵的乘积，由GPDSP的DSP核并行执行计算。DSP核执行下三角矩阵与矩阵的乘积计算时，由标量处理器读取的下三角矩阵数据为步骤S3.2计算得到的保留在片内共享存储阵列中的逆矩阵UU_i ^-1转置得到，由向量阵列处理器读取列矩阵L_i的矩阵数据。矩阵数据传入时按照矩阵转置的方式，计算后的结果矩阵按照矩阵转置的方式传输到片外DDR存储器的原矩阵位置。

上述DSP核计算下三角矩阵与矩阵的乘积中，DSP核采用DMA双缓冲策略在片内向量阵列存储器和片外DDR存储器之间传输数据。DSP核的片内向量阵列存储器设立两个缓冲区，缓冲区的存储容量为8*NB*p*q字节。计算时，下三角矩阵保持在片内共享存储阵列中，DSP核分别在两个缓冲区中计算下三角矩阵与矩阵的乘积，在一个缓冲区中计算时，将上一次缓冲区的计算结果传输到片外DDR存储器的原矩阵位置，同时将下一次计算需要的矩阵数据从片外DDR存储器传输到该缓冲区。

步骤S3.4和步骤S3.5在步骤S3.3之后由GPDSP的CPU核和DSP核并行执行。

S3.6：GPDSP的DSP核更新行panel的行矩阵U_i。

按照其中为下三角矩阵，行矩阵U_i的更新为下三角矩阵与矩阵的乘积，由GPDSP的DSP核并行执行计算。DSP核执行下三角矩阵与矩阵的乘积计算时，由标量处理器读取的下三角矩阵数据为步骤S3.2计算得到的逆矩阵LL_i ^-1，由向量阵列处理器读取行矩阵U_i的矩阵数据，计算后的结果矩阵传出到片外DDR存储器的原矩阵位置。

上述DSP核计算下三角矩阵与矩阵的乘积中，DSP核采用DMA双缓冲策略在片内向量阵列存储器和片外DDR存储器之间传输数据。DSP核的片内向量阵列存储器设立两个缓冲区，缓冲区的存储容量为8*NB*p*q。计算时，下三角矩阵保持在片内共享存储阵列中，DSP核分别在两个缓冲区中计算下三角矩阵与矩阵的乘积，在一个缓冲区中计算时，将上一次缓冲区的计算结果传输到片外DDR存储器的原矩阵位置，同时将下一次计算需要的矩阵数据从片外DDR存储器传输到该缓冲区。

S3.7：GPDSP的DSP核更新尾矩阵A_i。

按照A_i＝A_i-L_i*U_i更新尾矩阵A_i。

执行更新尾矩阵计算时，采用对列矩阵L_i分块的方法，每次由GPDSP的一个DSP核的标量处理器先将L_i中的一个NB×NB阶的分块矩阵Q_j加载至片内共享存储阵列中。再由DSP核将后续计算所需的行矩阵U_i矩阵数据依次传输到片内向量阵列存储器，并且按照公式A_ij＝A_ij-Q_j*U_i执行矩阵与矩阵的乘加计算。

上述DSP核执行矩阵与矩阵的乘加计算中，DSP核采用DMA双缓冲策略在片内向量阵列存储器和片外DDR存储器之间传输数据。DSP核在片内向量阵列存储器为A_ij和U_i分别设立两个缓冲区，缓冲区的存储容量为4*NB*p*q。计算时，分块矩阵Q_j保持在片内共享存储阵列中，DSP核分别在A_ij和U_i的两个缓冲区中执行矩阵与矩阵的乘加计算，在一个缓冲区中计算时，将上一次缓冲区的计算结果A_ij传输到片外DDR存储器的原矩阵位置，同时将下一次计算需要的矩阵数据从片外DDR存储器传输到该缓冲区。

上述计算依次进行，直到完成尾矩阵A_i的全部更新。

S3.8：重复步骤S3.1和步骤S3.7，直到完成全部分块K_i的LU分解。

S4：根据步骤S3得到的LU分解结果，GPDSP的CPU核和DSP核协同计算，求解三角方程Ly＝b和三角方程Ux＝y。设矩阵A的LU分解结果为：A＝LU，其中L为n×n阶的下三角矩阵，U为n×n阶的上三角矩阵。

在本实施例中，如图5所示，上述步骤S4中的三角方程求解采用如下方法：

S4.1：GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ly＝b。

S4.2：根据GPDSP的体系结构特征确定求解上述三角方程所需的分块大小NC值；

上述求解三角方程所需的分块大小NC值的确定方法是：设条件(1)16*NC*p*q≤s；(2)8*(NC*NC+(r+2)*NC)≤t；(3)NC是(p*q)的整数倍，n是NC的整数倍；按照满足条件(1)(2)(3)取得的最大NC值为所需要的NC值。

S4.3：设下三角矩阵L按照上述NC值进行分块，记j＝n/NC。沿对角线的三角分块依次标记为T_i(1≤i≤j)，以该分块为坐标，对应列方向的分块依次标记为S_i+1,i,S_i+2,i,...,S_j,i。将y,b分别划分为长度为NC的子向量，分别有j个子向量，分别标记为y₁,y₂,...,y_j；b₁,b₂,...,b_j；。

S4.4：GPDSP的CPU核将三角分块T₁和子向量b₁加载到片内共享存储阵列中，计算得到计算结果y₁传出到片外DDR存储器中，同时保持在片内共享存储阵列中。

S4.5：GPDSP的r个DSP核通过DMA依次从片外DDR存储器将三角分块T₁所在列的矩阵分块S_2,1,S_3,1,...,S_r+1,1和子向量b₁后续的r个子向量b₂,b₃,...,b_r+1传入计算，其中S_2,1,S_3,1,...,S_r+1,1传入到DSP核的片内向量阵列存储器，子向量b₂,b₃,...,b_r+1传入到片内共享存储阵列中。分别执行计算：b₂＝b₂-S_2,1*y₁,b₃＝b₃-S_3,1*y₁,...,b_r+1＝b_r+1-S_r+1,1*y₁，计算后的结果b₂,b₃,...,b_r+1保持在片内共享存储阵列中。重复上述计算，直到该列方向所有的矩阵分块计算完毕。

S4.6：GPDSP的CPU核将下一个三角分块T₂加载到片内共享存储阵列中，计算得到计算结果y₂传出到DDR存储器中，同时保持在片内共享存储阵列中。

S4.7：重复S4.5和S4.6，直到完成所有的子向量y₁,y₂,...,y_j的计算。

S4.8：GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ux＝y。

采用步骤S4.1至步骤S4.7求解三角方程Ly＝b相同的方法求解三角方程Ux＝y。完成向量x的计算。

S5：记CPU核的当前系统时间为t₁。由GPDSP的CPU核计算出GPDSP的实测性能为F₁＝ops/(t₁-t₀)，Linpack效率为F₁/F₀，并输出上述计算结果，Linpack测试完毕。

以上仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述实施例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，应视为本发明的保护范围。

Claims (6)

1.一种面向GPDSP的大规模高性能Linpack测试基准实现的方法，其特征在于，采用分块的方式来求解稠密线性方程组：较小的串行任务和小矩阵求逆和三角方程求解任务，由GPDSP中的CPU核负责；规整的计算任务由GPDSP中的DSP核负责。

2.根据权利要求1所述的面向GPDSP的大规模高性能Linpack测试基准实现的方法，其特征在于，具体步骤为：

S1：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值和分块大小NB值；Linpack测试基准求解的稠密线性方程组为Ax＝b，其中A为n×n阶的非奇异矩阵，x,b为n×1阶的向量；

S2：根据上述n值，由GPDSP的CPU核在片外DDR存储器为矩阵A和向量x,b分配空间，并生成初始化数据；记Linpack计算量为计算GPDSP的峰值性能为F₀，记CPU核的当前系统时间为t₀；

S3：GPDSP的CPU核和DSP核协同计算，按照NB×NB阶的分块规模对n×n阶的矩阵进行LU分解；

S4：根据步骤S3得到的LU分解结果，GPDSP的CPU核和DSP核协同计算，求解方程LUx＝b；

S5：记CPU核的当前系统时间为t₁；由GPDSP的CPU核计算出GPDSP的实测性能为F₁＝ops/(t₁-t₀)，Linpack效率为F₁/F₀，并输出上述计算结果，Linpack测试完毕。

3.根据权利要求2所述的面向GPDSP的大规模高性能Linpack测试基准实现的方法，其特征在于，所述步骤S4中求解方程LUx＝b分两步：(1)GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ly＝b；(2)GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ux＝y。

4.根据权利要求3所述的面向GPDSP的大规模高性能Linpack测试基准实现的方法，其特征在于，所述步骤S4的具体步骤为：

S4.1：GPDSP的CPU核和DSP核协同计算，按照分块的方法求解三角方程Ly＝b；

S4.2：根据GPDSP的体系结构特征确定求解上述三角方程所需的分块大小NC值；

上述求解三角方程所需的分块大小NC值的确定方法是：设条件(1)16*NC*p*q≤s；(2)8*(NC*NC+(r+2)*NC)≤t；(3)NC是(p*q)的整数倍，n是NC的整数倍；按照满足条件(1)(2)(3)取得的最大NC值为所需要的NC值；

S4.3：设下三角矩阵L按照上述NC值进行分块，记j＝n/NC；沿对角线的三角分块依次标记为T_i，其中1≤i≤j，以该分块为坐标，对应列方向的分块依次标记为S_i+1,i,S_i+2,i,…,S_j,i；将y,b分别划分为长度为NC的子向量，分别有j个子向量，分别标记为y₁,y₂,…,y_j；b₁,b₂,…,b_j；

S4.4：GPDSP的CPU核将三角分块T₁和子向量b₁加载到片内共享存储阵列中，计算得到计算结果y₁传出到片外DDR存储器中，同时保持在片内共享存储阵列中；

S4.5：GPDSP的r个DSP核通过DMA依次从片外DDR存储器将三角分块T₁所在列的矩阵分块S_2,1,S_3,1,...,S_r+1,1和子向量b₁后续的r个子向量b₂,b₃,…,b_r+1传入计算，其中S_2,1,S_3,1,...,S_r+1,1传入到DSP核的片内向量阵列存储器，子向量b₂,b₃,…,b_r+1传入到片内共享存储阵列中；分别执行计算：b₂＝b₂-S_2,1*y₁,b₃＝b₃-S_3,1*y₁,...,b_r+1＝b_r+1-S_r+1,1*y₁，计算后的结果b₂,b₃,…,b_r+1保持在片内共享存储阵列中；重复上述计算，直到该列方向所有的矩阵分块计算完毕；

S4.6：GPDSP的CPU核将下一个三角分块T₂加载到片内共享存储阵列中，计算得到计算结果y₂传出到DDR存储器中，同时保持在片内共享存储阵列中；

S4.7：重复步骤S4.5和步骤S4.6，直到完成所有的子向量y₁,y₂,…,y_j的计算。

5.根据权利要求2～4所述的面向GPDSP的大规模高性能Linpack测试基准实现的方法，其特征在于，所述步骤S3的具体步骤为：

S3.1：GPDSP的CPU核对分块矩阵K_i进行LU分解；由GPDSP的CPU核将分块矩阵K_i加载至片内共享存储阵列中，按照列选主元方法对分块矩阵K_i进行LU分解，LU分解结果为K_i＝LL_i*UU_i,其中LL_i为下三角矩阵，UU_i为上三角矩阵；将LU分解的计算结果传输到片外DDR存储器的原矩阵位置；

S3.2：GPDSP的CPU核计算步骤S3.1得到的LU分解的上三角矩阵UU_i的逆矩阵；由GPDSP的CPU核计算上三角矩阵UU_i的逆矩阵计算结果保留在片内共享存储阵列中；

S3.3：GPDSP的DSP核对行panel的行矩阵U_i进行行交换；由GPDSP的DSP核根据步骤3.1得到的行交换信息按照数据并行的方法对行矩阵U_i并行的进行行交换计算；

其中，步骤S3.2和步骤S3.3在步骤S3.1之后由GPDSP的CPU核和DSP核并行执行；

S3.4：GPDSP的CPU核计算步骤S3.1得到的LU分解的下三角矩阵LL_i的逆矩阵；由GPDSP的CPU核计算下三角矩阵LL_i的逆矩阵计算结果保留在片内共享存储阵列中；

S3.5：GPDSP的DSP核更新列panel的列矩阵L_i；

按照其中为下三角矩阵，将列矩阵L_i的更新由矩阵与上三角矩阵的乘积转化为下三角矩阵与矩阵的乘积，由GPDSP的DSP核并行执行计算；DSP核执行下三角矩阵与矩阵的乘积计算时，由标量处理器读取的下三角矩阵数据为步骤S3.2计算得到的保留在片内共享存储阵列中的逆矩阵UU_i ^-1转置得到，由向量阵列处理器读取列矩阵L_i的矩阵数据；矩阵数据传入时按照矩阵转置的方式，计算后的结果矩阵按照矩阵转置的方式传输到片外DDR存储器的原矩阵位置；

步骤S3.4和步骤S3.5在步骤S3.3之后由GPDSP的CPU核和DSP核并行执行；

S3.6：GPDSP的DSP核更新行panel的行矩阵U_i；

按照其中为下三角矩阵，行矩阵U_i的更新为下三角矩阵与矩阵的乘积，由GPDSP的DSP核并行执行计算；DSP核执行下三角矩阵与矩阵的乘积计算时，由标量处理器读取的下三角矩阵数据为步骤S3.2计算得到的逆矩阵LL_i ^-1，由向量阵列处理器读取行矩阵U_i的矩阵数据，计算后的结果矩阵传出到片外DDR存储器的原矩阵位置；

S3.7：GPDSP的DSP核更新尾矩阵A_i；按照A_i＝A_i-L_i*U_i更新尾矩阵A_i；

执行更新尾矩阵计算时，采用对列矩阵L_i分块的方法，每次由GPDSP的一个DSP核的标量处理器先将L_i中的一个NB×NB阶的分块矩阵Q_j加载至片内共享存储阵列中；再由DSP核将后续计算所需的行矩阵U_i矩阵数据依次传输到片内向量阵列存储器，并且按照公式A_ij＝A_ij-Q_j*U_i执行矩阵与矩阵的乘加计算；上述计算依次进行，直到完成尾矩阵A_i的全部更新；

S3.8：重复步骤S3.1和步骤S3.7，直到完成全部分块K_i的LU分解。

6.根据权利要求2～4所述的面向GPDSP的大规模高性能Linpack测试基准实现的方法，其特征在于，所述步骤S1的具体步骤为：

S1.1：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵LU分解的分块大小NB值；上述矩阵LU分解的分块大小NB值的确定方法是：设条件(1)16*NB*p*q≤s；(2)8*NB*NB≤t；(3)NB是r的整数倍；按照满足条件(1)(2)(3)取得的最大NB值为所需要的NB值；

S1.2：根据GPDSP的体系结构特征确定最佳的Linpack测试基准的矩阵规模n值；上述矩阵规模n值的确定方法是：设条件(1)n是NB的整数倍；(2)8*(n*n+2*n)<u；按照满足条件(1)(2)取得的最大n值为所需要的n值。