1.本发明涉及向量处理器技术领域，特别是涉及一种面向向量处理器的稀疏矩阵向量乘方法、装置及介质。
背景技术：
2.在工业界和学术界包含大量基于结构网格的数值软件，一些常用的稀疏线性系统的求解方法在数值软件中占用大量的计算时间和资源。
3.向量处理器具有多层次的存储空间和多级计算部件，能够同步进行综合数据的运算操作，若将基于稀疏线性系统的求解方法和对应的核心函数直接应用至向量处理器内，简单的移植单核版本的核心函数程序不能充分发挥向量处理器的性能，导致向量处理器的处理效率较低。
4.因此，如何充分利用向量处理器以发挥其高效性能是本领域技术人员亟需要解决的。
技术实现要素：
5.本发明的目的是提供一种面向向量处理器的稀疏矩阵向量乘方法、装置及介质，基于向量处理器的优势，并行处理各个核内的向量参数和矩阵参数，使得向量处理器的处理效率提高，充分发挥向量处理器的性能，提高计算力。
6.为解决上述技术问题，本发明提供一种面向向量处理器的稀疏矩阵向量乘方法，包括：
7.获取所述向量处理器的各核内的当前块内对应的向量参数和矩阵参数，其中，所述矩阵参数基于结构网格的稀疏线性系统且经过逻辑划分后得到的矩阵；
8.将所述向量参数和所述矩阵参数从双倍速率同步动态随机存储器分别加载至全局共享存储器空间和各所述核内的阵列存储器空间；
9.通过所述矩阵参数中的局部索引矩阵和所述向量参数从所述全局共享存储器空间索引取址至各所述核的所述阵列存储器空间；
10.在所述阵列存储器空间内通过向量部件对各所述核内的所述矩阵参数进行向量化计算处理得到向量化计算结果；
11.将所述向量化计算结果传回至所述双倍速率同步动态随机存储器，并返回所述获取所述向量处理器的各核内的当前块内对应的向量参数和矩阵参数的步骤以获取下一个当前块内的向量参数和矩阵参数直至各所述核内的所有块计算完成。
12.优选地，所述矩阵参数包括稀疏矩阵和索引矩阵，获取所述向量处理器的各核内的当前块内对应的向量参数和矩阵参数，包括：
13.获取原始向量参数和原始矩阵参数；
14.通过ell格式将所述原始向量参数和所述原始矩阵参数存储至所述双倍速率同步动态随机存储器以便于获取所述向量参数和所述矩阵参数。
15.优选地，获取所述向量参数和逻辑划分前的矩阵参数，包括：
16.根据所述全局共享存储器的空间数据确定块的维度，其中所述块的维度至少为有限元、有限差分、差分格式、三维网格或二维网格的一种；
17.根据所述块的维度对所述原始向量参数和所述原始矩阵参数进行分块排序以得到所述向量参数、所述逻辑划分前的矩阵参数、所述当前块的块向量buff层索引和所述局部索引矩阵。
18.优选地，所述向量化计算处理方式为高斯-塞戴尔迭代方式时，在获取到所述向量参数和所述逻辑划分之前的矩阵参数之后，在所述将所述向量参数和所述矩阵参数从双倍速率同步动态随机存储器分别加载至全局共享存储器空间和各所述核内的阵列存储器空间之前，还包括：
19.获取所述当前块内的各节点参数；
20.根据各所述节点参数之间的离散格式和网格维度确定所述当前块的颜色种类；
21.根据各颜色种类对各节点进行多色重排序以得到各所述节点对应的颜色矩阵参数。
22.优选地，所述逻辑划分后的所述矩阵参数的获取过程，包括：
23.根据所述向量处理器的核数对所述当前块内的各所述颜色种类下的节点对应的所述逻辑划分前的矩阵参数、所述局部索引矩阵和所述向量参数进行均匀划分以得到各所述核内的逻辑划分后的第一矩阵参数；
24.根据所述阵列存储器的空间数据对划分后的所述第一矩阵参数进行划分以得到所述矩阵参数。
25.优选地，所述向量化计算处理方式为稀疏矩阵向量乘方式时，所述逻辑划分后的所述矩阵参数的获取过程，包括：
26.根据所述向量处理器的核数对所述当前块内的节点对应的所述逻辑划分前的矩阵参数、所述局部索引矩阵和所述向量参数进行均匀划分以得到各所述核内的逻辑划分后的第一矩阵参数；
27.根据所述阵列存储器的空间数据对划分后的所述第一矩阵参数进行划分以得到所述矩阵参数。
28.优选地，所述将向量化计算后的结果传回至所述双倍速率同步动态随机存储器，包括：
29.将所述当前块内的当前颜色对应的结果传回至所述全局共享存储器空间；
30.获取所述当前颜色的下一个所述当前颜色的结果并传回至所述全局共享存储器空间直至所述当前块内的所有颜色种类对应的结果全部计算完成；
31.将各所述结果从所述全局共享存储器空间传输至所述双倍速率同步动态随机存储器。
32.为解决上述技术问题，本发明还提供一种面向向量处理器的稀疏矩阵向量乘装置，包括：
33.获取模块，用于获取所述向量处理器的各核内的当前块内对应的向量参数和矩阵参数，其中，所述矩阵参数基于结构网格的稀疏线性系统且经过逻辑划分后得到的矩阵；
34.加载模块，用于将所述向量参数和所述矩阵参数从双倍速率同步动态随机存储器
分别加载至全局共享存储器空间和各所述核内的阵列存储器空间；
35.索引取址模块，用于通过所述矩阵参数中的局部索引矩阵和所述向量参数从所述全局共享存储器空间索引取址至各所述核的所述阵列存储器空间；
36.处理模块，用于在所述阵列存储器空间内通过向量部件对各所述核内的所述矩阵参数进行向量化计算处理得到向量化计算结果；
37.传回模块，用于将所述向量化计算结果传回至所述双倍速率同步动态随机存储器，并返回所述获取所述向量处理器的各核内的当前块内对应的向量参数和矩阵参数步骤以获取下一个当前块内的向量参数和矩阵参数直至各所述核内的所有块计算完成。
38.为解决上述技术问题，本发明还提供一种面向向量处理器的稀疏矩阵向量乘装置，包括：
39.存储器，用于存储计算机程序；
40.处理器，用于执行所述计算机程序时实现如上述所述的面向向量处理器的稀疏矩阵向量乘方法的步骤。
41.为解决上述技术问题，本发明还提供一种计算机可读存储介质，所述计算机可读存储介质上存储有计算机程序，所述计算机程序被处理器执行时实现如上述所述的面向向量处理器的稀疏矩阵向量乘方法的步骤。
42.本发明提供的一种面向向量处理器的稀疏矩阵向量乘方法，该方法通过与向量处理器结合，获取到各个核内的当前块内的向量参数和矩阵参数，对于矩阵参数是基于结构网格的稀疏线性系统且经过逻辑划分后得到的矩阵，根据当前块的向量参数和矩阵参数加载至对应的全局共享存储器空间和各核内的阵列存储器空间，在各个核内索引取址至阵列存储器空间，可以节省取址时间，再通过向量处理器的向量部件对各个核的矩阵参数并行进行向量化计算处理，以将结果传回。基于向量处理器的优势，并行处理各个核内的向量参数和矩阵参数，使得向量处理器的处理效率提高，充分发挥向量处理器的性能，提高计算力。
43.另外，本发明还提供了一种面向向量处理器的稀疏矩阵向量乘装置及介质，具有如上述面向向量处理器的稀疏矩阵向量乘方法相同的有益效果。
附图说明
44.为了更清楚地说明本发明实施例，下面将对实施例中所需要使用的附图做简单的介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
45.图1为本发明实施例提供的一种面向向量处理器的稀疏矩阵向量乘方法，的流程图；
46.图2为本发明实施例提供的一种向量处理器的结构图；
47.图3为本发明实施例提供的一种4
×
4的结构网格点模型的示意图；
48.图4为本发明提供的一种矩阵原始格式存储示意图；
49.图5为本发明提供的一种分块示意图；
50.图6为本发明实施例提供的一种高斯-塞戴尔迭代方式下的分块排序的示意图；
51.图7为本发明实施例提供的一种稀疏矩阵向量乘方式下的分块排序的示意图；
52.图8为本发明实施例提供的一种多色重排序的结构网格点模型的示意图；
53.图9为本发明实施例提供的一种高斯-塞戴尔迭代方式下的多色重排序后的示意图；
54.图10为本发明实施例提供的一种高斯-塞戴尔迭代方式下的逻辑划分示意图；
55.图11为本发明实施例提供的一种稀疏矩阵向量乘方式下的逻辑划分示意图；
56.图12为本发明实施例提供的一种面向向量处理器的稀疏矩阵向量乘装置的结构图；
57.图13为本发明实施例提供的另一种面向向量处理器的稀疏矩阵向量乘装置的结构图；
58.图14为本发明实施例提供的一种高斯-塞戴尔迭代方式下的拷贝索引向量的示意图；
59.图15为本发明实施例提供的一种高斯-塞戴尔迭代方式下的加载gsm存储的示意图；
60.图16为本发明实施例提供的一种高斯-塞戴尔迭代方式下的加载至am空间的示意图；
61.图17为本发明实施例提供的一种高斯-塞戴尔迭代方式下的索引取址的示意图；
62.图18为本发明实施例提供的一种高斯-塞戴尔迭代方式下的加载至am空间的示意图；
63.图19为本发明实施例提供的一种高斯-塞戴尔迭代方式下的向量参数计算过程示意图；
64.图20为本发明实施例提供的一种高斯-塞戴尔迭代方式下的传回gsm的示意图；
65.图21为本发明实施例提供的一种高斯-塞戴尔迭代方式下的传回示意图；
66.图22为本发明实施例提供的一种稀疏矩阵向量乘方式下的加载至am空间的示意图；
67.图23为本发明实施例提供的一种稀疏矩阵向量乘方式下的索引取址的示意图；
68.图24为本发明实施例提供的最红稀疏矩阵向量乘方式下的向量参数计算过程示意图；
69.图25为本发明实施例提供的一种稀疏矩阵向量乘方式下的传回示意图。
具体实施方式
70.下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下，所获得的所有其他实施例，都属于本发明保护范围。
71.本发明的核心是提供一种面向向量处理器的稀疏矩阵向量乘方法、装置及介质，基于向量处理器的优势，并行处理各个核内的向量参数和矩阵参数，使得向量处理器的处理效率提高，充分发挥向量处理器的性能，提高计算力。
72.为了使本技术领域的人员更好地理解本发明方案，下面结合附图和具体实施方式对本发明作进一步的详细说明。
73.需要说明的是，向量处理器是一类能够执行向量指令的处理器，而不仅仅局限于标量。向量操作对于同样数目的数据进行操作的时候比相应的一系列标量指令要快得多。向量处理器尤其对于大规模的科学工程计算较为有效，包括汽车碰撞模拟和天气预报。这些应用程序通常需要一台超级计算机跑上几打个小时来处理gigabyte级别的数据。高速的标量处理器依赖于cache来减少访问主存的延迟，但是大规模长时间运行的科学计算程序通常有很大规模的工作集，并且通常局部性非常低，这导致memory hierarchy的性能非常糟糕。所以标量处理器会提供旁路cache的机制如果软件发现访存的局部性很差。但是使得主存饱和需要硬件跟踪数百上千条的in flight的标量访存操作，而这在标量处理器的指令集架构(instruction set architecture，isa)中已被证实开销是非常大的。相反，向量isa可以只使用一条向量指令就可以发起对于一整个向量中元素的访存操作，所以非常简单的逻辑就可以提供很高的带宽。
74.图1为本发明实施例提供的一种面向向量处理器的稀疏矩阵向量乘方法，的流程图，如图1所示，包括：
75.s11：获取向量处理器的各核内的当前块内对应的向量参数和矩阵参数；
76.其中，矩阵参数基于结构网格的稀疏线性系统且经过逻辑划分后得到的矩阵；
77.s12：将向量参数和矩阵参数从双倍速率同步动态随机存储器分别加载至全局共享存储器空间和各核内的阵列存储器空间；
78.s13：通过矩阵参数中的局部索引矩阵和向量参数从全局共享存储器空间索引取址至各核的阵列存储器空间；
79.s14：在阵列存储器空间内通过向量部件对各核内的矩阵参数进行向量化计算处理得到向量化计算结果；
80.s15：将向量化计算结果传回至双倍速率同步动态随机存储器，并返回步骤s11，以获取下一个当前块内的向量参数和矩阵参数直至各核内的所有块计算完成。
81.图2为本发明实施例提供的一种向量处理器的结构图，如图2所示，该向量处理器集成若干个dsp核，cn表示dsp核数，例如一个6mb的全局共享存储器(global share memory，gsm)和一个dma引擎。向量处理器内核由标量处理部件(scalar processing unit，spu)和向量处理部件(vector processing unit，vpu)组成，其中，spu负责标量计算和流控，vpu提供主要的向量计算能力，其集成了16个向量处理单元(vector processing unit，vpe)，vpe的所有相同编号的局部寄存器在逻辑上构成一个1024位的向量寄存器。每个dsp核提供64kb的标量存储器(scalar memory，sm)和768kb的阵列存储器(array memory，am)用于标量访问和向量访问。直接存储访问(direct memory access，dma)模式可以实现各级存储空间(sm、am、gsm、ddr)的快速访问。
82.dma模式提供点对点传输指令dma_p2p和间接索引取址传输指令dma_sg。其中dma_p2p表示点对点传输，可以对具有固定偏移的连续内存进行数据传输；dma_sg表示索引取址传输，可以通过索引对存储进行取址并进行数据传输。
83.步骤s11中的获取向量处理器的各核内的当前块内对应的向量参数和矩阵参数，需要说明的是，本步骤获取到的向量参数和矩阵参数是已经在ddr中按照分块的方式重新排序后的向量参数和矩阵参数。对于矩阵参数，还是在重新排序后再经过逻辑划分得到的。也就是说均是经过预先处理后的矩阵参数和向量参数。向量处理器包括多个核，一个核内
包括多个块。由于并行计算的特点，多个核可以并行处理，每个核内的块进行循环迭代。
84.对于在何种存储空间获取步骤s11中的向量参数和矩阵参数，需要将原始向量参数和原始矩阵参数经过整理存储至向量处理器的ddr中，对应地存储方式可以是ell存储，也可以是其他存储，在此不做限定。只要是从ddr中获取到向量参数和矩阵参数即可。对于分块过程，其为了将存储在ddr中的矩阵参数和向量参数与向量处理器更好的匹配，便于后续的计算过程。对应匹配过程需要设置块的属性，对于块的属性不做限定，可以根据维度等因素设置，将矩阵参数和向量参数按照分块的方式进行重新排序。
85.需要说明的是，由于向量化计算方式不同，其对应的分块方式需要根据向量化计算方式本身的特点进行排序分块。例如，向量化计算方式的高斯-塞戴尔迭代方式，使得矩阵的节点具有强烈的依赖性，一个节点处理完后才可以进行下一个节点，在矩阵中由于各个节点的相关性不同，有的节点与其他节点不相关，则需要将不相关的节点同色处理，即划分为同一种颜色，同一种颜色的节点可以并行处理以适应向量处理器的并行处理特点，需要在根据块的方式排序基础上再进行多色重排序，提高计算效率。对于向量化计算方式的稀疏矩阵向量乘方式则无需考虑依赖性特点，故不会考虑多色重排序的问题，仅需要按照块的方式排序即可。
86.本实施例中的逻辑划分是在分块排序的基础上，根据向量处理器自身的特点(dsp核数和am空间的大小)进行逻辑划分，以便于后续的am块的存储。
87.在每一块矩阵进行计算时，将分块后的块向量和buff向量从ddr中加载到gsm空间中，将逻辑划分后的矩阵从ddr中加载至对应核上的am空间中。通过局部索引矩阵并行地从全局共享存储器(gsm)中的向量取址到每个dsp核上的阵列存储器(am)，每个dsp核并行地在阵列存储器(am)中通过向量部件进行对应的向量化计算处理，并把计算结果传回ddr，直到所有的块被计算完成。
88.本发明实施例提供的一种面向向量处理器的稀疏矩阵向量乘方法，该方法通过与向量处理器结合，获取到各个核内的当前块内的向量参数和矩阵参数，对于矩阵参数是基于结构网格的稀疏线性系统且经过逻辑划分后得到的矩阵，根据当前块的向量参数和矩阵参数加载至对应的全局共享存储器空间和各核内的阵列存储器空间，在各个核内索引取址至阵列存储器空间，可以节省取址时间，再通过向量处理器的向量部件对各个核的矩阵参数并行进行向量化计算处理，以将结果传回。基于向量处理器的优势，并行处理各个核内的向量参数和矩阵参数，使得向量处理器的处理效率提高，充分发挥向量处理器的性能，提高计算力。
89.作为一种实施例，在获取过程中对应的在何种存储空间中进行获取作为限定，矩阵参数包括稀疏矩阵和索引矩阵，获取向量处理器的各核内的当前块内对应的向量参数和矩阵参数，包括：
90.获取原始向量参数和原始矩阵参数；
91.通过ell格式将原始向量参数和原始矩阵参数存储至双倍速率同步动态随机存储器以便于获取向量参数和矩阵参数。
92.具体地，预先将基于结构网格生成的稀疏矩阵、索引矩阵、待求解向量、右端向量、稀疏矩阵对角向量和中间向量存储在ddr中，以便于后续进行向量化计算处理。这里的向量化计算处理方式可以是高斯-塞戴尔迭代方式或者是稀疏矩阵向量乘方式，还可以是其他
处理方式，在此不做限定，在选取对应的向量化计算处理方式时需要根据对应的处理方式特征进行存储、加载等。在进行稀疏矩阵向量乘方式时，其向量参数为待乘向量和结果向量，在进行高斯-塞戴尔迭代方式时，向量参数为待求解向量、右端向量、稀疏矩阵对角向量和中间向量。
93.在本实施例中的存储过程无论采用何种向量化计算处理方式其存储过程相同，需要将基于结构网格的稀疏矩阵和索引矩阵(原始矩阵参数)、待求解向量、右端向量、稀疏矩阵对角向量和中间向量(原始向量参数)通过ell格式进行存储。本实施例以稀疏矩阵为例进行说明，索引矩阵相同。
94.基于结构网格生成的稀疏矩阵的特点是：对应网格内部的每行的非零项是固定的，在网格边界上的非零项较少。图3为本发明实施例提供的一种4
×
4的结构网格点模型的示意图，如图3所示，以5点差分格式，通常情况下一个节点存在多个自由度，一个自由度寓意数据代表内容不同，例如一个节点代表速度、压力、温度等参数的自由度，同时一个自由度中还可以根据不同的坐标系对应不同的维度参数。为了说明，以单节点上单自由度分量为例，5点差分格式中的中心点只与自身和上下左右的点相关联。图4为本发明提供的一种矩阵原始格式存储示意图，如图4所示，深色表示非零元素，浅色表示零元素。
95.图4中的ell格式，采用两个n_row
×
nnz_row的数组a和idx分别存储矩阵元素和列索引，n_row表示原始矩阵的维度，等于总节点数nn乘以单结点上的自由度分量dn，即n_row＝nn
×
dn。nnz_row表示每行最大非零元数，与离散格式相关联，在没有值的部分即灰色的部分用零进行填充。nnz_row＝dn
×
n_neib表示，其中n_neib表示与自身相关的点数，在5点差分格式中，n_neib＝5，dn＝1。对于最大非零元数，主要保留非零元素，需要根据各个节点的相关联的节点确定整个矩阵或者索引数组的最大行和列，例如，图3中的节点6关联4个节点再加自己的一个节点，一共5个节点，是16个节点中最多的关联个数，则定为五列，由于为一个自由度，其一个节点对应一行，一共16个节点，对应16行。
96.结构网格中的每个节点中的每个分量对应矩阵中的一行。在同步动态随机存储器(ddr)中的采用列连续的方式对矩阵数组进行存储，如图4中箭头所示，由于单个节点只有一个分量，先对node_1—node_16的第一个非零元连续存储，再对node_1—node_16的第二个非零元连续存储，直到结束。如果单个节点包含dn个分量，存储顺序为node_1的1-dn分量的第一个非零元、node_2的1-dn分量的第一个非零元、
……
、node_nn的1-dn分量的第一个非零元、node_1的1-dn分量的第二个非零元，直到结束。
97.索引数组采用与矩阵数组完全相同的方式，按列连续存储，索引数组存储的值是与矩阵元素相对应的向量x的索引。待求解向量x、右端向量b和稀疏矩阵对角向量diaga的均采用n_row长的数组进行存储。
98.上述实施例说明了在何种存储器中获取向量参数和矩阵参数，对于怎么获取向量参数和矩阵参数，需要对其进行分块处理过程，作为一种实施例，获取向量参数和逻辑划分前的矩阵参数，包括：
99.根据全局共享存储器的空间数据确定块的维度，其中块的维度至少为有限元、有限差分、差分格式、三维网格或二维网格的一种；
100.根据块的维度对原始向量参数和原始矩阵参数进行分块排序以得到向量参数、逻辑划分前的矩阵参数、当前块的块向量buff层索引和局部索引矩阵。
101.需要说明的是，由于分块排序后并未得到步骤s11中的矩阵参数，还需要再进行逻辑划分才可以得到，故本实施例仅是针对分块处理，即分块处理后得到的是向量参数和逻辑划分前的矩阵参数。
102.在分块之前需要确定块的维度问题，由于块向量需要存放在gsm中，需要满足以下条件：
103.在二维双精度的情况下，(nx 2)
×
(ny 2)
×
sizeof(double)《＝6mb；
104.在三维双精度的情况下，(nx 2)
×
(ny 2)
×
(nz 2)
×
sizeof(double)《＝6mb。
105.其中，nx、ny和nz分别代表块在x、y和z方向上的节点数。sizeof表示以字节为单位时数据类型的大小。对应的加2仅是一种实施例，加2是因为在每个方向上具有2层buff层，图5为本发明提供的一种分块示意图，如图5所示，因为每个块的边界点进行向量化处理计算时，都需要用到周围节点对应的x值。
106.可以理解的是，块的维度可以是有限元、有限差分、差分格式、三维网格或二维网格的一种，但是不限于提到的维度方式，还可以包括其他维度方式，在此不做限定。
107.根据块的维度对原始向量参数和原始矩阵参数进行分块排序以得到向量参数、逻辑划分前的矩阵参数、当前块的块向量buff层索引和局部索引矩阵。
108.具体地，根据设置的块的维度，对矩阵matrxi_a、向量vector_b、向量vector_x和向量vector_diaga进行分块和重新排序，bn表示总块数。图5中的原始节点数为nn，原始网格被分为了bn个block(图中bn＝9)，每个block的维度完全相等，每个block的节点数为bnn，用block_row进行表示每个block的矩阵行，n_row＝nn
×
dn，block_row＝n_row/bn＝bnn
×
dn。图6为本发明实施例提供的一种高斯-塞戴尔迭代方式下的分块排序的示意图，如图6所示，在ddr空间内，重排序的实质是根据节点的编号对向量x和矩阵a进行重排，如果原始节点号为1，重排序后节点号位2，原始节点号1对应的x1-xdn，b1-bdn，diaga1-diagadn和a1，1—and，nnz_row在新的向量和数组中将对应节点号2的位置。在新的排序中，每个block内的节点连续排序，block之间依次排序。在图6中，dof为自由度，一个节点对应多个自由度，一个自由度对应一个向量vector_b、一个向量vector_x和一个向量vector_diaga，一个矩阵matrxi_a内的数组。在一个块内包括多个节点。
109.同样在进行稀疏矩阵向量乘方式时，图7为本发明实施例提供的一种稀疏矩阵向量乘方式下的分块排序的示意图，如图7所示，根据设置的块的维度，对矩阵a、向量x进行分块和重新排序，bn表示总块数。图5中的原始节点数为nn，原始网格被分为了bn个block(图中bn＝9)，每个block的维度完全相等，每个block的节点数为bnn，用block_row进行表示每个block的矩阵行，n_row＝nn
×
dn，block_row＝n_row/bn＝bnn
×
dn。如图7所示，在ddr空间内，重排序的实质是根据节点的编号对向量x和矩阵a进行重排，如果原始节点号为1，重排序后节点号位2，原始节点号1对应的x1-xdn和a1，1—and，nnz_row在新的向量和数组中将对应节点号2的位置。在新的排序中，每个block内的节点连续排序，block之间依次排序。向量ax是不需要进行重排序的，因为初始所有值都为零，只需要存储最终结果。
110.如图5所示，对分块后的线性系统进行向量化计算处理是，在块的边界点对应计算gauss-seidel迭代需要与块相邻点的x值，也就是buff层。用数组buff_idx存储buff层的索引，单个块的buff层的数组的长度，在二维时等于(nx 2)
×
(ny 2)-nx
×
ny，在三维时等于(nx 2)
×
(ny 2)
×
(nz 2)-nx
×
ny
×
nz。ddr中需要用数组buff_idx存储所有块的buff层索
引，总长度为buff_size
×
block_n。经过重排序后，每个块的维度完全相同，且每个块有着完全相同的排序方式，因此每个块的局部列索引矩阵完全相同，但是由于重排序后向量x的位置发生了变化，需要重新生成新的索引，局部块索引数组的维度为block_row
×
nnz_row。
111.在上述实施例的基础上，向量化计算处理方式为高斯-塞戴尔迭代方式时，在获取到向量参数和逻辑划分之前的矩阵参数之后，在将向量参数和矩阵参数从双倍速率同步动态随机存储器分别加载至全局共享存储器空间和各核内的阵列存储器空间之前，还包括：
112.获取当前块内的各节点参数；
113.根据各节点参数之间的离散格式和网格维度确定当前块的颜色种类；
114.根据各颜色种类对各节点进行多色重排序以得到各节点对应的颜色矩阵参数。
115.具体地，现有的gauss-seidel迭代本质是对稀疏线性系统ax＝b进行迭代求解，给出ellpack-itpack(ell)格式下的一次gauss-seidel迭代计算流程。
116.步骤1、令i＝0；
117.步骤2、令ax[i]＝0；
[0118]
步骤3、令j＝0；
[0119]
步骤4、计算ax[i] ＝a[i n_row*j]*x[idx[i n_row*j]]；
[0120]
步骤5、计算j＝j 1；如果j《nnz_row跳转到步骤3，如果j＝nnz_row跳转到步骤6；
[0121]
步骤6、计算x[i] ＝(b[i]-ax[i])/diaga[i]；
[0122]
步骤7、计算i＝i 1；
[0123]
如果8、如果i《n_row跳转到步骤2，如果i＝n_row结束。
[0124]
可以看出gauss-seidel迭代具有强依赖性，点i的x的计算需要点i之前的x完成计算，不能直接并行计算。需要多色重排序去除依赖，再进行并行计算。图8为本发明实施例提供的一种多色重排序的结构网格点模型的示意图，如图8所示，对点进行重排序，先对深色的点进行编号，再对浅色的点进行编号，按照编号次序进行计算，先计算深色的点，深色的点之间计算没有依赖性，可以并行地进行计算，再计算浅色的点，浅色的点之间计算没有依赖性，可以并行地进行计算。但是两种颜色之间有依赖性。根据网格的维数(2维或者3维)和离散格式(5点格式，7点格式，9点格式，27点格式)可以进行不同的多色排序(2色，4色，8色)。在ddr空间内，图9为本发明实施例提供的一种高斯-塞戴尔迭代方式下的多色重排序后的示意图，如图9所示，用con表示颜色数量，color_row表示单块内一种颜色的矩阵行数，cnn表示单块内一种颜色的结点数，有如下关系color_row＝block_row/con＝cnn
×
dn，bnn＝cnn
×
con。
[0125]
在上述实施例的基础上，逻辑划分后的矩阵参数的获取过程，包括：
[0126]
根据向量处理器的核数对当前块内的各颜色种类下的节点对应的逻辑划分前的矩阵参数、局部索引矩阵和向量参数进行均匀划分以得到各核内的逻辑划分后的第一矩阵参数；
[0127]
根据阵列存储器的空间数据对划分后的第一矩阵参数进行划分以得到矩阵参数。
[0128]
具体地，图10为本发明实施例提供的一种高斯-塞戴尔迭代方式下的逻辑划分示意图，如图10所示，对矩阵matrxi_a的每个block中每一种颜色和局部索引矩阵local_idx中的每一种颜色，对向量vector_b、向量vector_x和向量vector_diaga的每个block中每一种颜色，根据dsp核数进行均匀划分，如图10所示。单个核分到的矩阵块的维度为core_row
×
nnz_row，core_row表示单个核分到的矩阵行，core_row＝block_row/cn，其中block_row表示每一块矩阵的行数，cn表示dsp核数，局部索引矩阵local_idx也作了完全相同的逻辑划分得到单个dsp核内的第一矩阵参数(矩阵块matrxi_a和局部索引矩阵local_idx)。
[0129]
对单个dsp核的矩阵块matrxi_a和局部索引矩阵local_idx根据阵列存储器(am)空间大小进行划分，每个am块表示am空间一次计算的矩阵块的大小，受到am存储空间(768kb)的影响，am块具体的大小为am_row
×
nnz_row，am_row表示单个am块的矩阵行数，需要满足如下条件：
[0130]
am_row
×
(nnz_row
×
2 4)
×
sizeof(double)《＝768kb；
[0131]
am_row应当取得尽可能的大才能获得较好的效果，am_row
×
nnz_row
×2×
sizeof(double)的空间用来存放am_row
×
nnz_row的矩阵行和对应的取址后的x向量，4
×
am_row
×
sizeof(double)的空间用来存储向量vector_x，右端向量vector_b、对角向量vector_diaga和中间向量vector_ax。an表示单个核上需要计算的am块数，有an＝core_row/am_row。
[0132]
对应地，作为一种实施例，在向量化计算处理方式为稀疏矩阵向量乘方式时对应的逻辑划分过程，逻辑划分后的矩阵参数的获取过程，包括：
[0133]
根据向量处理器的核数对当前块内的节点对应的逻辑划分前的矩阵参数、局部索引矩阵和向量参数进行均匀划分以得到各核内的逻辑划分后的第一矩阵参数；
[0134]
根据阵列存储器的空间数据对划分后的第一矩阵参数进行划分以得到矩阵参数。
[0135]
具体地，图11为本发明实施例提供的一种稀疏矩阵向量乘方式下的逻辑划分示意图，如图11所示，对矩阵matrxi_a、向量ax的每个block和局部索引矩阵local_idx根据dsp核数进行均匀划分，如图11所示，单个核分到的矩阵块的维度为core_row
×
nnz_row，core_row表示单个核分到的矩阵行，core_row＝block_row/cn，其中block_row表示每一块矩阵的行数，cn表示dsp核数，局部索引矩阵local_idx也作了完全相同的逻辑划分，得到了单个dsp核的第一矩阵参数(矩阵块matrxi_a和局部索引矩阵local_idx)。
[0136]
对单个dsp核的矩阵块matrxi_a和局部索引矩阵local_idx根据阵列存储器(am)空间大小进行划分，每个am块表示am空间一次计算的矩阵块的大小，受到am存储空间(768kb)的影响，am块具体的大小为am_row
×
nnz_row，am_row表示单个am块的矩阵行数，需要满足以下条件：
[0137]
am_row
×
(nnz_row
×
2 1)
×
sizeof(double)《＝768kb；
[0138]
am_row应当取得尽可能的大才能获得较好的效果，am_row
×
nnz_row
×2×
sizeof(double)的空间用来存放am_row
×
nnz_row的矩阵行和对应的取址后的x向量，am_row
×
sizeof(double)的空间用来存储计算结果ax。an表示单个核上需要计算的am块数，有an＝core_row/am_row。
[0139]
本实施例根据不同的向量化计算处理方式分别进行了分块以及逻辑划分，基于不同的向量化计算处理方式本身的特点以及向量处理器自身的优势，进行的各参数获取过程，为后续向量处理器的计算奠定了基础，提供了便利。
[0140]
在上述实施例的基础上，在向量化计算处理方式为高斯-塞戴尔迭代方式时，将向量化计算后的结果传回至双倍速率同步动态随机存储器，包括：
[0141]
将当前块内的当前颜色对应的结果传回至全局共享存储器空间；
processing unit，cpu)；协处理器是用于对在待机状态下的数据进行处理的低功耗处理器。在一些实施例中，处理器22可以集成有图像处理器(graphics processing unit，gpu)，gpu用于负责显示屏所需要显示的内容的渲染和绘制。一些实施例中，处理器22还可以包括人工智能(artificial intelligence，ai)处理器，该ai处理器用于处理有关机器学习的计算操作。
[0158]
存储器21可以包括一个或多个计算机可读存储介质，该计算机可读存储介质可以是非暂态的。存储器21还可包括高速随机存取存储器，以及非易失性存储器，比如一个或多个磁盘存储设备、闪存存储设备。本实施例中，存储器21至少用于存储以下计算机程序211，其中，该计算机程序被处理器22加载并执行之后，能够实现前述任一实施例公开的面向向量处理器的稀疏矩阵向量乘方法的相关步骤。另外，存储器21所存储的资源还可以包括操作系统212和数据213等，存储方式可以是短暂存储或者永久存储。其中，操作系统212可以包括windows、unix、linux等。数据213可以包括但不限于面向向量处理器的稀疏矩阵向量乘方法所涉及到的数据等等。
[0159]
在一些实施例中，面向向量处理器的稀疏矩阵向量乘装置还可包括有显示屏23、输入输出接口24、通信接口25、电源26以及通信总线27。
[0160]
领域技术人员可以理解，图13中示出的结构并不构成对面向向量处理器的稀疏矩阵向量乘装置的限定，可以包括比图示更多或更少的组件。
[0161]
处理器22通过调用存储于存储器21中的指令以实现上述任一实施例所提供的面向向量处理器的稀疏矩阵向量乘方法。
[0162]
对于本发明提供的一种面向向量处理器的稀疏矩阵向量乘装置的介绍请参照上述方法实施例，本发明在此不再赘述，其具有上述面向向量处理器的稀疏矩阵向量乘方法相同的有益效果。
[0163]
进一步的，本发明还提供了一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序，计算机程序被处理器22执行时实现如上述面向向量处理器的稀疏矩阵向量乘方法的步骤。
[0164]
可以理解的是，如果上述实施例中的方法以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的全部或部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：u盘、移动硬盘、只读存储器(read-only memory，rom)、随机存取存储器(random access memory，ram)、磁碟或者光盘等各种可以存储程序代码的介质。
[0165]
对于本发明提供的一种计算机可读存储介质的介绍请参照上述方法实施例，本发明在此不再赘述，其具有上述面向向量处理器的稀疏矩阵向量乘方法相同的有益效果。
[0166]
作为一种实施例，以高斯-塞戴尔迭代方式为例，整体的计算过程如下：
[0167]
步骤1、令i＝1；
[0168]
步骤2、开始计算block_i块；
[0169]
步骤3、图14为本发明实施例提供的一种高斯-塞戴尔迭代方式下的拷贝索引向量的示意图，如图14所示，根据buff_idx的block_i块在ddr中从向量vector_x中拷贝得到索
引向量buff_x的block_i块；
[0170]
步骤4、dma_p2p操作，把buff_x向量和vector_x向量对应的block_i块从ddr中加载到向量处理器的gsm存储上，在gsm中两个部分共同组成了local_xv，图15为本发明实施例提供的一种高斯-塞戴尔迭代方式下的加载gsm存储的示意图，如图15所示；
[0171]
步骤5、令m＝1；
[0172]
步骤6、开始计算color_m块；
[0173]
步骤7、令j＝1；
[0174]
步骤8、开始计算am_j块；
[0175]
步骤9、采用dma_p2p操作，将am_j块对应的矩阵matrix_a、向量vector_b和向量vector_diaga并行地从ddr加载到每个dsp核的am空间上，图16为本发明实施例提供的一种高斯-塞戴尔迭代方式下的加载至am空间的示意图，如图16所示。
[0176]
步骤10、图17为本发明实施例提供的一种高斯-塞戴尔迭代方式下的索引取址的示意图，如图17所示，使用dma_sg操作，通过am_j块局部索引矩阵local_idx从gsm中的local_xv向量，并行地索引取址至每个dsp核的am空间，索引后得到的向量用vector_xi表示。
[0177]
步骤11、采用dma_p2p操作，将am_j块对应的vector_x并行地从gsm加载到每个dsp核的am空间上，图18为本发明实施例提供的一种高斯-塞戴尔迭代方式下的加载至am空间的示意图，如图18所示。
[0178]
步骤12、所有的dsp核并行地在am空间上通过向量部件对matrix_a、vector_x、vector_b、vector_diaga和vector_xi进行gauss-seidel迭代计算。图19为本发明实施例提供的一种高斯-塞戴尔迭代方式下的向量参数计算过程示意图，如图19所示，图19描述了单个核上的情况，向量部件长度为1024位，vec_length表示向量长度，以双精度为例，一条指令可以同时操作16个双精度数的，即vec_length＝16。am_j块矩阵行的长度为am_row，vec_n表示am_row的向量个数，vec_n＝am_row/vec_length，在逻辑划分时，am_row已经被提前划分成了vec_length的倍数。用vec_ax、vec_b、vec_x、vec_a和vec_xi分别表示1024位长度的vector_ax、vector_b、vector_x、matrix_a和vector_xi向量。vector_ax只作为中间变量不需要进行传输。
[0179]
步骤13、令k＝0；
[0180]
步骤14、令vec_ax[k]＝0；
[0181]
步骤15、令n＝0；
[0182]
步骤16、vec_ax[k]＝vec_ax[k] vec_a[k vec_n*n]*vec_xi[k vec_n*n]；
[0183]
步骤17、令n＝n 1，如果n小于nnz_row则跳转到步骤16，如果n等于nnz_row则跳转到步骤18；
[0184]
步骤18、令vec_x[k]＝vec_x[k] (vec_b[k]-vec_ax[k])/vec_diaga[k]；
[0185]
步骤19、令k＝k 1，如果k小于vec_n则跳转到步骤14，如果k等于vec_n则跳转到步骤20；
[0186]
步骤20、图20为本发明实施例提供的一种高斯-塞戴尔迭代方式下的传回gsm的示意图，如图20所示，采用dma_p2p操作，把计算得到的结果vector_x传回gsm中，gsm中的vector_x完成了迭代更新。
[0187]
步骤21、令j＝j 1，如果j小于等于an则跳转到步骤8，如果j大于an则跳转到步骤22；
[0188]
步骤22、令m＝m 1，如果m小于等于con则跳转到步骤6，如果m大于con则跳转到步骤23；
[0189]
步骤23、图21为本发明实施例提供的一种高斯-塞戴尔迭代方式下的传回示意图，如图21所示，采用dma_p2p操作，把gsm中block_i的vector_x传回ddr中，ddr中block_i的vector_x完成了迭代更新。
[0190]
步骤24、令i＝i 1，如果i小于等于bn则跳转到步骤2，如果i大于bn则跳转到步骤25；
[0191]
步骤25、结束。
[0192]
作为一种实施例，以稀疏矩阵向量乘方式为例，整体的计算过程如下：
[0193]
步骤1、令i＝1；
[0194]
步骤2、开始计算block_i块；
[0195]
步骤3、该拷贝过程和高斯-塞戴尔迭代方式相同，可参照图14，根据buff_idx的block_i块在ddr中从向量vector_x中拷贝得到索引向量buff_x的block_i块。
[0196]
步骤4、dma_p2p操作，把buff_x向量和vector_x向量对应的block_i块从ddr中加载到向量处理器的gsm存储上，在gsm中两个部分共同组成了local_xv，该加载gsm存储过程和高斯-塞戴尔迭代方式相同，可参照图15；
[0197]
步骤5、令j＝1；
[0198]
步骤6、开始计算am_j块；
[0199]
步骤7、采用dma_p2p操作，将am_j块矩阵对应的matrix_a并行地从ddr加载到每个dsp核的am空间上，图22为本发明实施例提供的一种稀疏矩阵向量乘方式下的加载至am空间的示意图，如图22所示。
[0200]
步骤8、图23为本发明实施例提供的一种稀疏矩阵向量乘方式下的索引取址的示意图，如图23所示，使用dma_sg操作，通过am_j块局部索引矩阵local_idx从gsm中的local_xv向量，并行地索引取址至每个dsp核的am空间，索引后得到的向量用vector_xi表示。
[0201]
步骤9、所有的dsp核并行地在am空间上通过向量部件对matrix_a和vector_xi进行spmv操作。图24为本发明实施例提供的最红稀疏矩阵向量乘方式下的向量参数计算过程示意图，图24描述了单个核上的情况，向量部件长度为1024位，vec_length表示向量长度，以双精度为例，一条指令可以同时操作16个双精度数的，即vec_length＝16。am_j块矩阵行的长度为am_row，vec_n表示am_row的向量个数，vec_n＝am_row/vec_length，在逻辑划分时，am_row已经被提前划分成了vec_length的倍数。用vec_ax、vec_a和vec_xi分别表示1024位长度的vector_ax、matrix_a和vector_xi向量。
[0202]
步骤10、令k＝0；
[0203]
步骤11、令vec_ax[k]＝0；
[0204]
步骤12、令n＝0；
[0205]
步骤13、vec_ax[k]＝vec_ax[k] vec_a[k vec_n*n]*vec_xi[k vec_n*n]；
[0206]
步骤14、令n＝n 1，如果n小于nnz_row则跳转到步骤13，如果n等于nnz_row则跳转到步骤15；
[0207]
步骤15、令k＝k 1，如果k小于vec_n则跳转到步骤11，如果k等于vec_n则跳转到步骤16；
[0208]
步骤16、图25为本发明实施例提供的一种稀疏矩阵向量乘方式下的传回示意图，如图25所示，采用dma_p2p操作，把计算得到的结果vector_ax传回ddr中。
[0209]
步骤17、令j＝j 1，如果j小于等于an则跳转到步骤6，如果j大于an则跳转到步骤18；
[0210]
步骤18、令i＝i 1，如果i小于等于bn则跳转到步骤2，如果i大于bn则跳转到步骤19；
[0211]
步骤19、结束。
[0212]
对于本发明提供的一种分别对应高斯-塞戴尔迭代方式和稀疏矩阵向量乘方式的面向向量处理器的稀疏矩阵向量乘方法的介绍请参照上述方法实施例，本发明在此不再赘述，其具有上述面向向量处理器的稀疏矩阵向量乘方法相同的有益效果。
[0213]
以上对本发明所提供的一种面向向量处理器的稀疏矩阵向量乘方法、面向向量处理器的稀疏矩阵向量乘装置及介质进行了详细介绍。说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的装置而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以对本发明进行若干改进和修饰，这些改进和修饰也落入本发明权利要求的保护范围内。
[0214]
还需要说明的是，在本说明书中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个
……”
限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。