CN104636316A

CN104636316A - 面向gpdsp的大规模矩阵乘法计算的方法

Info

Publication number: CN104636316A
Application number: CN201510063830.3A
Authority: CN
Inventors: 刘仲; 陈书明; 万江华; 陈磊; 田希; 彭元喜; 陈虎; 扈啸; 孙永节; 陈胜刚; 孙海燕; 阳柳; 张雪萌; 马胜
Original assignee: National University of Defense Technology
Current assignee: National University of Defense Technology
Priority date: 2015-02-06
Filing date: 2015-02-06
Publication date: 2015-05-20
Anticipated expiration: 2035-02-06
Also published as: CN104636316B

Abstract

本发明公开了一种面向GPDSP的大规模矩阵乘法计算的方法，其步骤为：S1：由CPU核在片外DDR存储器为矩阵A、B和C分配存储空间，并生成初始化数据或者从其他数据来源传输本次计算所需要的数据；S2：根据GPDSP的体系结构特征确定最佳的分块矩阵乘法所需要的分块大小MB、KB和NB；S3：由CPU核根据步骤S2确定的MB、KB、NB值将矩阵A、B和C进行逻辑分块；S4：由GPDSP的CPU核调度DSP核计算子块矩阵的乘加计算：C_ik+＝A_ij*B_jk(1≤i≤mi,1≤j≤ki,1≤k≤ni)；S5：计算完毕。本发明原理简单、操作方便，能充分利用GPDSP的CPU核通用计算和DSP核向量处理阵列强大并行计算、高带宽向量数据加载能力，显著提高DSP核计算访存比。

Description

面向GPDSP的大规模矩阵乘法计算的方法

技术领域

本发明主要涉及通用计算数字信号处理器(General-Purpose Digital Signal Processor，简称GPDSP)，特指一种适用于GPDSP的大规模矩阵乘法计算的方法。

背景技术

基本线性代数算法库(Basic Linear Algebra Subroutines，BLAS)是各类科学计算最常用的核心数学算法库之一，工业界针对各自的处理器平台都推出了高度优化的BLAS实现，如IBM的ESSL、Intel的MKL、AMD的ACML等。其中，矩阵乘法(General Matrix-MatrixMultiplication，GEMM)是BLAS库的核心算法。GEMM是典型的计算密集和访存密集型应用，对处理器的运算能力、访存带宽及延迟要求非常高，相关文献研究表明GEMM计算占据高性能基准测试程序(High Performance Linpack，HPL)运算量的90％以上。因此，针对处理器的体系结构特点研究GEMM优化方法对评测该处理器的计算效率、发挥处理器的计算优势和提高应用程序的运行速度均具有很重要的参考价值。

分块矩阵乘法方法是提升大规模GEMM性能的主要方法，最典型的研究是Gunnels针对基于Cache的多级存储结构，提出的分层计算的GEMM优化方法，能够降低存储层次间搬运数据的平均开销。面向Cache结构的分块矩阵乘法方法的主要思想是将大矩阵的乘法分割为一系列子矩阵的乘法，设m为Cache的容量，则分块参数中的子矩阵块大小blocksize通常满足约束条件blocksize<＝sqrt(m/3)，使得子矩阵计算时的数据访问能够全部在Cache中命中，该部分的计算能以接近峰值的性能进行，从而提高整个大矩阵乘法的计算性能。

在专利申请号为201310725118.6的文献(处于实审阶段)中提供了一种通用计算数字信号处理器(General-Purpose Digital Signal Processor，简称GPDSP)，它包含CPU核单元和DSP核单元，CPU核单元主要用于负责包括存储管理、文件控制、进程调度、中断管理任务在内的通用事务管理以及提供对通用操作系统的完整支持；DSP核单元包含若干强大计算能力的64位向量处理阵列，用于支持高密集运算任务的解算。传统的面向Cache结构的分块矩阵乘法方法不适合GPDSP的非Cache的向量阵列存储访存模式和向量处理阵列并发向量处理的体系结构特征，难以发挥GPDSP向量计算优势。

发明内容

本发明要解决的技术问题就在于：针对现有技术存在的技术问题，本发明提供一种原理简单、操作方便、能充分利用GPDSP的CPU核通用计算和DSP核向量处理阵列强大并行计算、高带宽向量数据加载能力，显著提高DSP核计算访存比的面向GPDSP的大规模分块矩阵乘法计算的方法。

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

一种面向GPDSP的大规模矩阵乘法计算的方法，其步骤为：

S1：由GPDSP的CPU核在片外DDR存储器为矩阵A、B和C分配存储空间，并生成初始化数据或者从其他数据来源传输本次计算所需要的数据；

S2：根据GPDSP的体系结构特征确定最佳的分块矩阵乘法所需要的分块大小MB、KB和NB，即令矩阵的M、K、N三个维数对应的分块大小分别标记为MB、KB、NB；

S3：由GPDSP的CPU核根据步骤S2确定的MB、KB、NB值将矩阵A、B和C进行逻辑分块，令mi＝INT(M/MB)，ki＝INT(K/KB)，ni＝INT(N/NB)，其中INT()表示对括号中的值向上取整；矩阵A的子块标记为A_ij(1≤i≤mi,1≤j≤ki)，矩阵B的子块标记为B_jk(1≤j≤ki,1≤k≤ni)，矩阵C的子块标记为C_ik(1≤i≤mi,1≤k≤ni)；

S4：由GPDSP的CPU核调度DSP核计算子块矩阵的乘加计算：C_ik+＝A_ij*B_jk(1≤i≤mi,1≤j≤ki,1≤k≤ni)；

S5：计算完毕。

作为本发明的进一步改进：所述步骤S2中NB的确定方法是NB＝p*q。

作为本发明的进一步改进：所述步骤S2中MB、KB的确定方法为：设条件(1)MB*KB*w≤s3，(2)(MB*KB*w+KB*NB*w)≤s1，按照满足条件(1)(2)取得的最大MB、KB值为所需要的MB、KB值。

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

S4.1：由GPDSP的CPU核通过DMA将矩阵A的子块A_ij(1≤i≤mi,1≤j≤ki)加载至片内共享存储阵列中，i，j的初始值均为1；

S4.2：由GPDSP的CPU核调度各个DSP核分别计算子块A_ij与矩阵B的子块B_jk(1≤k≤ni)的矩阵乘加计算C_ik+＝A_ij*B_jk(1≤k≤ni)；

S4.3：由GPDSP的CPU核判断DSP核是否完成子块A_ij与矩阵B的所有子块B_jk(1≤k≤ni)的矩阵乘加计算C_ik+＝A_ij*B_jk(1≤k≤ni)，若是转步骤S4.4，若不是则转步骤S4.2，直至上述DSP核计算任务全部完成；

S4.4：判断j是否等于ki，若不是令j＝j+1，转步骤S4.1；若是，则继续判断i是否等于mi，若不是令i＝i+1，j＝1，转步骤S4.1；若是，转步骤S5；

作为本发明的进一步改进：所述步骤S4.2中共有k＝ni次子块矩阵的乘加计算，由CPU核调度各个DSP核依次计算各子块矩阵的乘加计算，每个DSP核接收到计算任务后按照CPU分配的矩阵B和C的逻辑分块由各自核的DMA加载至片内向量存储数据阵列，完成子块矩阵的乘加计算后通知CPU核；CPU核根据子块计算的完成情况派发新的子块矩阵乘加计算任务，直至上述k＝ni次子块矩阵的乘加计算全部完成。

作为本发明的进一步改进：所述步骤S4.2中的GPDSP的DSP核计算子块A_ij与矩阵B的子块B_jk(1≤k≤ni)的矩阵乘加计算C_ik+＝A_ij*B_jk(1≤k≤ni)包括如下流程：

S4.2.1：由GPDSP的DSP核将片内标量存储器L1D的存储访问模式配置为SRAM访问模式；

S4.2.2：将子块A_ij根据片内标量存储器的容量大小划分为Mb个长度为KB的连续计算子块，并标记为A_x(1≤x≤INT(MB/Mb))；

上述参数Mb的确定方法是，满足条件：2*Mb*KB*w≤s2取得的最大Mb值为所需要的Mb值；

S4.2.3：DSP核通过DMA依次将子块A_x加载至片内标量存储器L1D中，初始x＝1，计算子块矩阵乘加：C_xk+＝A_x*B_jk(1≤x≤INT(MB/Mb))；

S4.2.3计算中共有MB/Mb次子块矩阵的乘加计算，DSP核通过各自的DMA分别加载A_x至片内标量存储数据阵列，完成C_xk+＝A_x*B_jk单元块的矩阵乘加计算任务，直至上述MB/Mb次子块的矩阵乘加计算全部完成。

作为本发明的进一步改进：所述步骤S4中包括DSP核在片内标量存储器进一步分块的流程：DSP核将子块A_ij根据片内标量存储器的容量大小划分为Mb个长度为KB的连续计算子块，并标记为A_x(1≤x≤INT(MB/Mb))；上述参数Mb的确定方法是，满足条件：2*Mb*KB*w≤s2取得的最大Mb值为所需要的Mb值；DSP核通过DMA依次将子块A_x加载至片内标量存储器L1D中，初始x＝1，计算子块矩阵乘加：C_xk+＝A_x*B_jk(1≤x≤INT(MB/Mb))。

与现有技术相比，本发明的优点在于：

本发明的面向GPDSP的大规模矩阵乘法计算的方法，原理简单、操作方便，基于GPDSP的体系结构特征，利用其高效的CPU核和DSP核协同计算实现大规模矩阵乘法计算；本发明由CPU核运行操作系统，负责数据的分割和任务的调度、数据同步，DSP核运行负责分块矩阵乘法计算的核心程序，不间断的从CPU核获取新的计算任务并将运行结果报告CPU核。由上可知，本发明提供的任务分割方法能够充分发挥CPU核的通用计算和DSP核的强大向量化计算能力的优点，实现CPU核与DSP核间的紧密协同配合，高效的实现大规模矩阵乘法计算。

附图说明

图1是本发明所面向的GPDSP的简化访存结构模型示意图。

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

图3是本发明在具体应用实例中GPDSP的CPU核对矩阵A、B和C进行逻辑分块的流程示意图。

图4是本发明在具体应用实例中GPDSP的CPU核对矩阵A、B和C在M、K和N维进行任务调度的流程示意图。

图5是本发明在具体应用实例中GPDSP的DSP核在片内标量存储器进一步分块的流程示意图。

图6是本发明在具体应用实例中GPDSP的DSP核在片内标量存储器进行的子块矩阵乘加计算的流程示意图。

具体实施方式

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

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

设GPDSP中DSP核数量为r个；DSP核的向量处理阵列计算单元数量为p个，每个计算单元的MAC(乘加部件)数量为q个，DSP核的片内向量阵列存储器容量为s1字节，DSP核的片内标量存储器容量为s2字节；GPDSP的片内共享存储阵列容量为s3字节，GPDSP的片外共享DDR存储容量为s4字节。设计算的GEMM矩阵乘加运算为C+＝A*B，其中A为M×K阶的矩阵，B为K×N阶的矩阵，C为M×N阶的矩阵，每个矩阵元素的数据为w字节。

如图2所示，本发明面向GPDSP的大规模矩阵乘法计算的方法在具体应用实例中的流程为：

S1：由GPDSP的CPU核在片外DDR存储器为矩阵A、B和C分配存储空间，并生成初始化数据或者从其他数据来源传输本次计算所需要的数据。

S2：根据GPDSP的体系结构特征确定最佳的分块矩阵乘法所需要的分块大小MB、KB和NB。

令矩阵的M、K、N三个维数对应的分块大小分别标记为MB、KB、NB。

其中，NB的确定方法是NB＝p*q。MB、KB的确定方法为：设条件(1)MB*KB*w≤s3，(2)(MB*KB*w+KB*NB*w)≤s1，按照满足条件(1)(2)取得的最大MB、KB值为所需要的MB、KB值。

S3：参见图3所示，由GPDSP的CPU核根据步骤S2确定的MB、KB、NB值将矩阵A、B和C进行逻辑分块，令mi＝INT(M/MB)，ki＝INT(K/KB)，ni＝INT(N/NB)，其中INT()表示对括号中的值向上取整。

矩阵A的子块标记为A_ij(1≤i≤mi,1≤j≤ki)，矩阵B的子块标记为B_jk(1≤j≤ki,1≤k≤ni)，矩阵C的子块标记为C_ik(1≤i≤mi,1≤k≤ni)。

S4：参见图4所示，由GPDSP的CPU核调度DSP核计算子块矩阵的乘加计算：C_ik+＝A_ij*B_jk(1≤i≤mi,1≤j≤ki,1≤k≤ni)。

S4.1：由GPDSP的CPU核通过DMA将矩阵A的子块A_ij(1≤i≤mi,1≤j≤ki)加载至片内共享存储阵列中(i，j的初始值均为1)。

S4.2：由GPDSP的CPU核调度各个DSP核分别计算子块A_ij与矩阵B的子块B_jk(1≤k≤ni)的矩阵乘加计算C_ik+＝A_ij*B_jk(1≤k≤ni)。

上述步骤S4.2计算中共有k＝ni次子块矩阵的乘加计算，由CPU核调度各个DSP核依次计算各子块矩阵的乘加计算，每个DSP核接收到计算任务后按照CPU分配的矩阵B和C的逻辑分块由各自核的DMA加载至片内向量存储数据阵列，完成子块矩阵的乘加计算后通知CPU核，CPU核根据子块计算的完成情况派发新的子块矩阵乘加计算任务，直至上述k＝ni次子块矩阵的乘加计算全部完成。

S4.3：由GPDSP的CPU核判断DSP核是否完成子块A_ij与矩阵B的所有子块B_jk(1≤k≤ni)的矩阵乘加计算C_ik+＝A_ij*B_jk(1≤k≤ni)，若是转步骤4.4，若不是则转步骤4.2，直至上述DSP核计算任务全部完成。

S4.4：判断j是否等于ki，若不是令j＝j+1，转步骤S4.1；若是，则继续判断i是否等于mi，若不是令i＝i+1，j＝1，转步骤S4.1；若是，转步骤S5。

S5：计算完毕。

在本实施例中，上述步骤S4.2中的GPDSP的DSP核计算子块A_ij与矩阵B的子块B_jk(1≤k≤ni)的矩阵乘加计算C_ik+＝A_ij*B_jk(1≤k≤ni)采用如下计算方法：

如图5所示，为本实施例中GPDSP的DSP核在片内标量存储器进一步分块的流程，GPDSP的DSP核将子块A_ij根据片内标量存储器的容量大小划分为Mb个长度为KB的连续计算子块，并标记为A_x(1≤x≤INT(MB/Mb))。上述参数Mb的确定方法是，满足条件：2*Mb*KB*w≤s2取得的最大Mb值为所需要的Mb值。DSP核通过DMA依次将子块A_x加载至片内标量存储器L1D中(初始x＝1)，计算子块矩阵乘加：C_xk+＝A_x*B_jk(1≤x≤INT(MB/Mb))。

如图6所示，本实施例中，GPDSP的DSP核在片内标量存储器进行的子块矩阵乘加计算流程的具体步骤为：

S4.2.2：将子块A_ij根据片内标量存储器的容量大小划分为Mb个长度为KB的连续计算子块，并标记为A_x(1≤x≤INT(MB/Mb))。

上述参数Mb的确定方法是，满足条件：2*Mb*KB*w≤s2取得的最大Mb值为所需要的Mb值。

S4.2.3：DSP核通过DMA依次将子块A_x加载至片内标量存储器L1D中(初始x＝1)，计算子块矩阵乘加：C_xk+＝A_x*B_jk(1≤x≤INT(MB/Mb))。

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

Claims

1.一种面向GPDSP的大规模矩阵乘法计算的方法，其特征在于，步骤为：

S5：计算完毕。

2.根据权利要求1所述的面向GPDSP的大规模矩阵乘法计算的方法，其特征在于，所述步骤S2中NB的确定方法是NB＝p*q。

3.根据权利要求1所述的面向GPDSP的大规模矩阵乘法计算的方法，其特征在于，所述步骤S2中MB、KB的确定方法为：设条件(1)MB*KB*w≤s3，(2)(MB*KB*w+KB*NB*w)≤s1，按照满足条件(1)(2)取得的最大MB、KB值为所需要的MB、KB值。

4.根据权利要求1或2或3所述的面向GPDSP的大规模矩阵乘法计算的方法，其特征在于，所述步骤S4的具体流程为：

5.根据权利要求4所述的面向GPDSP的大规模矩阵乘法计算的方法，其特征在于，所述步骤S4.2中共有k＝ni次子块矩阵的乘加计算，由CPU核调度各个DSP核依次计算各子块矩阵的乘加计算，每个DSP核接收到计算任务后按照CPU分配的矩阵B和C的逻辑分块由各自核的DMA加载至片内向量存储数据阵列，完成子块矩阵的乘加计算后通知CPU核；CPU核根据子块计算的完成情况派发新的子块矩阵乘加计算任务，直至上述k＝ni次子块矩阵的乘加计算全部完成。

6.根据权利要求4所述的面向GPDSP的大规模矩阵乘法计算的方法，其特征在于，所述步骤S4.2中的GPDSP的DSP核计算子块A_ij与矩阵B的子块B_jk(1≤k≤ni)的矩阵乘加计算C_ik+＝A_ij*B_jk(1≤k≤ni)包括如下流程：

7.根据权利要求4所述的面向GPDSP的大规模矩阵乘法计算的方法，其特征在于，所述步骤S4中包括DSP核在片内标量存储器进一步分块的流程：DSP核将子块A_ij根据片内标量存储器的容量大小划分为Mb个长度为KB的连续计算子块，并标记为A_x(1≤x≤INT(MB/Mb))；上述参数Mb的确定方法是，满足条件：2*Mb*KB*w≤s2取得的最大Mb值为所需要的Mb值；DSP核通过DMA依次将子块A_x加载至片内标量存储器L1D中，初始x＝1，计算子块矩阵乘加：C_xk+＝A_x*B_jk(1≤x≤INT(MB/Mb))。