CN114491404A - 应用于计算设备的混合精度SpMV优化系统及方法 - Google Patents

Abstract

本发明公开了应用于计算设备的混合精度SpMV优化系统及方法，本系统包括：获取模块、第一处理模块、第二处理模块；所述获取模块用于读入稀疏矩阵到计算设备；所述第一处理模块基于稀疏矩阵中非零元素的浮点数精度，将所述稀疏矩阵划分为不同精度的子矩阵；所述第二处理模块用于计算不同精度的所述子矩阵与向量的乘法，得到混合精度SpMV计算结果。本发明利用混合精度降低SpMV的计算开销，提高访存效率，具有较高的理论意义和实际应用价值。

Description

应用于计算设备的混合精度SpMV优化系统及方法

技术领域

本发明属于稀疏线性代数求解优化领域，特别是涉及应用于计算设备的混合精度SpMV优化系统及方法。

背景技术

大规模稀疏线性方程组求解在计算电磁学、计算流体力学和空气多力学等多个领域具有非常重要的应用，且求解大型稀疏线性系统所花费的时间在求解整个问题所需的时间中占有非常高的比重，因此，如何快速有效的求解大型稀疏线性系统一直是目前的研究热点之一。直接法提供了一种数值稳定的方法来求解具有可预测浮点运算次数的线性系统。然而，对于大规模线性系统，直接分解的存储和计算成本可能不切实际。相比之下，诸如广义最小残差(GMRES)之类的迭代求解器具有良好的并行性。当与良好的预处理器结合以提高收敛性时，GMRES可以比稀疏LU求解器更快，尤其是对于大型问题。许多迭代方法使用稀疏矩阵向量乘法(SpMV)运算，例如共轭梯度法和GMRES。由于矩阵大小比向量大几个数量级，SpMV操作支配了这些迭代方法的执行时间。

随着AI技术的不断发展和进步，大规模数据处理对计算能力提出了新的需求。由于机器学习应用普遍对精度要求较低，因此有力推动了已有计算设备和新的AI加速器对16位浮点数表示的支持。CPU和GPU制造商都开始在新一代产品中支持低精度浮点计算。NVIDIA、Google和AMD的最新计算设备现在都支持16位浮点计算。Google的TPU主要是为机器学习应用定制设计的，使用bfloat16格式。AMD处理器在软件堆栈上支持bfloat16和IEEE754标准浮点格式。NVIDIA的Pascal GPU是第一代支持半精度格式的GPU。这主要是因为低精度计算不仅实现了更高的FLOPS(每秒浮点运算次数)，更重要的是，它有助于提高内存访问和数据传输的效率。计算设备在浮点计算中的新能力也为SpMV的优化提供了新的思路。然而，由于矩阵的稀疏性和SpMV计算的不规则内存访问，在GPU设备上实现SpMV的优化具有挑战性。混合精度计算降低了SpMV的数据传输和计算强度，是进一步提高SpMV性能的备选方案之一。然而，很长一段时间以来，大多数计算设备只支持单精度和双精度浮点格式。因此，在GPU设备上研究混合精度SpMV计算的研究很少。

SpMV的性能受体系结构、输入矩阵的表示和非零元结构的影响很大。目前已有的SpMV实现主要有单精度SpMV和双精度SpMV。单精度SpMV具有较高的计算性能，但是计算过程中精度损失严重；双精度SpMV具有较高精度的计算结果，但是受限于双精度算术单元的算术性能，计算开销非常高，而大多数科学计算问题都需要双精度精度。许多研究表明，以混合精度实现的程序有时可以比只使用最高精度算法更快地计算出相同精度的结果。其优势一方面是因为低精度具有更快的算术运算，另一方面更重要的是低精度相对高精度具有一定的传输数据效率的提升。混合精度计算的有效性在于必须保证最终的计算精度没有损失。

Ahmad等人提出了一种数据驱动的混合精度SpMV算法——MpSpMV，其基本思路是设定一个阈值范围(-1,1)，根据该阈值范围将矩阵中非零元素划分两个部分，其中一个为在阈值范围内的，用单精度表示，而在阈值范围之外的则用双精度表示。从理论上看，可以将浮点精度表示为：

其中u表示总的位数，据此得知p_e越小越好。通常情况下u是固定大小的，并且u＞e，所以e越小浮点精度越高。这个公式也表明数据越小，其表示的精度也越高。在做内积的时候，误差由于累积放大，所以大数的误差会更加明显。因此在计算的过程中较小数的表示精度可以适当放宽，而大数仍需用比较高的精度表示。不失一般性，作者选用了(-1,1)这个区间，区间内的用单精度表示，区间外的用双精度表示。虽然当时的GPU已经支持半精度浮点格式，但是该研究仍然是采用单精度和双精度进行混合精度数值方法设计。同时，从实验结果来看，该混合精度设计带来的性能提升较小。

长期以来，大多数计算设备仅支持单精度和双精度浮点格式。因此，很少有研究开展了混合精度SpMV计算优化工作。

在数值计算中，使用半精度或单精度浮点数据格式比使用双精度浮点数据格式具有基本的性能优势。由于更紧凑的表示，在包括寄存器文件、缓存和主存在内的内存层次结构的每一层中，可以存储更多的数据元素。同样，处理单个精度值在不同内存级别之间消耗的带宽更少，并且减少了缓存和TLB失败的数量。然而，数据移动方面主要影响内存密集型、带宽受限的应用程序，历史上并没有对混合精度算法引起太多关注。混合精度SpMV利用低精度的计算性能和高精度计算的精度，在实现性能优化的同时保持较高的精度。

发明内容

本发明提出了应用于计算设备的混合精度SpMV优化系统及方法，本发明的主要目的是利用混合精度降低SpMV的传输和计算开销，提高访存效率，具有较高的理论意义和实际应用价值。

为实现上述目的，本发明提供了如下方案：

应用于计算设备的混合精度SpMV优化系统，包括：

获取模块、第一处理模块、第二处理模块；

所述获取模块用于获取计算设备输入稀疏矩阵；

所述第一处理模块用于基于稀疏矩阵中非零元素的浮点数精度，将所述稀疏矩阵划分为不同精度的子矩阵；

所述第二处理模块用于计算不同精度的所述子矩阵与向量的乘法，得到混合精度SpMV计算结果。

优选的，所述第一处理模块对所述稀疏矩阵的划分，在所述混合精度SpMV计算之前只执行一次。

优选的，所述浮点数包括：符号、指数和尾数。

优选的，不同精度的所述子矩阵包括：半精度子矩阵、单精度子矩阵和双精度子矩阵。

优选的，不同精度的所述子矩阵均为CSR格式。

优选的，所述混合精度SpMV计算过程为:

获取所述半精度子矩阵、所述单精度子矩阵和所述双精度子矩阵以及单精度和双精度的向量；

基于划分后有多个子矩阵，设置计算矩阵一行与向量相乘的线程数，同时将线程重新编号，并合并各所述子矩阵相同行的计算结果；

基于最小重新编号的线程，将合并的所述计算结果写回所述结果向量，同时各所述子矩阵的行数自增，进行下一次的线程处理，当各所述子矩阵的行数大于所述稀疏矩阵行数时，停止混合精度SpMV计算。

本发明还公开了应用于计算设备的混合精度SpMV优化方法，包括以下步骤：

获取计算设备输入稀疏矩阵；

基于稀疏矩阵中非零元素的浮点数的精度，将所述稀疏矩阵划分为不同精度的子矩阵；

计算不同精度的所述子矩阵与向量的乘法，得到混合精度SpMV计算结果。

本发明的有益效果为：

本发明提供了应用于计算设备的混合精度SpMV优化系统及方法，本发明通过获取计算设备输入稀疏矩阵；基于稀疏矩阵中非零元素的浮点数的精度，将稀疏矩阵划分为不同精度的子矩阵；计算不同精度的所述子矩阵与向量的乘法，得到混合精度SpMV计算结果。本发明利用混合精度SpMV，采用了存储与计算分离的多精度研究思路，能够有效降低算法中的访存开销、数据传输开销和计算开销，具有较高的理论意义和实际应用价值。

附图说明

为了更清楚地说明本发明的技术方案，下面对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例一应用于计算设备的混合精度SpMV优化系统模块结构示意图；

图2为本发明实施例二应用于计算设备的混合精度SpMV优化方法流程示意图；

图3为本申请实施例三IEEE 754半精度、单精度、双精度浮点格式之间的转换示意图；

图4为本申请实施例三矩阵无损划分流程示意图；

图5为本申请实施例三矩阵无损划分的具体示例。按非零元精度将原始矩阵(a)，划分为半精度子矩阵(b)、单精度子矩阵(c)、双精度子矩阵(d)，所有矩阵均以CSR格式压缩存储；

图6为本申请实施例三GPU每一线程的混合精度SpMV计算流程示意图。

具体实施方式

下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本申请保护的范围。

为使本申请的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本申请作进一步详细的说明。

实施例一

如图1所示，应用于计算设备的混合精度SpMV优化系统，包括：

获取模块、第一处理模块、第二处理模块。

获取模块用于获取计算设备输入数据；

第一处理模块用于基于稀疏矩阵中非零元素的浮点数精度，将稀疏矩阵划分为不同精度的子矩阵；

第二处理模块用于计算不同精度的所述子矩阵与向量的乘法，得到混合精度SpMV计算结果；

具体的，所述第一处理模块对所述稀疏矩阵的划分，在所述混合精度SpMV计算之前只执行一次。

具体的，所述浮点数包括：符号、指数和尾数。

具体的，不同精度的所述子矩阵包括：半精度子矩阵、单精度子矩阵和双精度子矩阵。

具体的，不同精度的所述子矩阵均为CSR格式。

具体的，所述混合精度SpMV计算过程为:

获取所述半精度子矩阵、所述单精度子矩阵和所述双精度子矩阵以及单精度和双精度的向量；

基于划分后有多个子矩阵，设置计算矩阵一行与向量相乘的线程数为32，同时将线程重新编号，并合并各所述子矩阵相同行的计算结果；

基于最小重新编号的线程，将合并的所述计算结果写回所述结果向量，同时各所述子矩阵的行数自增，进行下一次的线程处理，当各所述子矩阵的行数大于所述稀疏矩阵行数时，停止混合精度SpMV计算。

实施例二

如图2所示，应用于计算设备的混合精度SpMV优化方法，包括以下步骤：

获取计算设备输入稀疏矩阵；

基于稀疏矩阵中非零元素的浮点数的精度，将稀疏矩阵划分为不同精度的子矩阵；

计算不同精度的所述子矩阵与向量的乘法，，得到混合精度SpMV计算结果。

实施例三

(一)矩阵无损划分

本发明的关键是基于非零元的精度将稀疏矩阵以混合精度表示，以便可以减少其内存访问开销以及SpMV计算中的计算强度。混合精度表示即根据非零元素的实际浮点精度，将稀疏矩阵划分为最多3个矩阵(半精度、单精度和双精度子矩阵)。将更多的非零元以低精度存储，意味着更少的内存占用和CPU内存和GPU全局内存之间的数据传输，以及GPU端的内存访问事务。此外，低精度计算有望获得比高精度计算更高的性能。因此，不同的存储和计算精度可能有助于提高SpMV计算的性能。

浮点数在内存中的表示分为三部分：符号s、指数e和尾数m，如图3所示。

给定一个双精度浮点数f，如果e∈[-15,15]，以及尾数m中前10位之后的所有其他位为0，则f可以无损表示为半精度。同理，当e∈[-127,127]，且尾数所需的位只包含前23位时，则浮点数f可以无损地用单精度表示。否则，这个数字只能用双精度表示。如图4所示，输入的稀疏矩阵为CSR格式，由于其存储开销低，这是最常用的稀疏矩阵压缩格式。在本发明中，也使用CSR格式来存储划分后的矩阵。图4中vals数组中以双精度存储着所有非零元，nnz表示非零元的数量。图5展示了无损划分的一个示例。无损划分在SpMV计算之前只需要执行一次。由于稀疏矩阵在迭代算法的运行过程中经常保持不变，这种一次性转换在多次迭代中摊销。矩阵划分可能会影响稀疏矩阵的稀疏性，进而影响SpMV的性能。例如，对于具有长行的不规则矩阵，矩阵划分会带来一定的性能提升，因为将长行进行划分后使负载更加均衡。但对于将相邻的非零元素划分为不同的子矩阵的情况，可能会带来性能损失，因为相邻非零元划分到不同子矩阵后不再能受益于访存局部性。通过后续的实验，证明了这种无损分区对SuiteSparse矩阵集合上的大多数矩阵都有积极的影响。

(二)混合精度SpMV实现

在CSR-Scalar内核中，每个线程处理稀疏矩阵的一行。这种方式的问题是每行非零元素分布不均导致负载不均衡。CSR-Vector内核通过为矩阵行分配一个warp的线程来解决这个问题。CUSP是一个广泛使用的稀疏线性代数和图计算库。它的SpMV实现具有更好的性能，并且整个库是开源的。本发明为不同的精度组合实现了混合精度SpMV内核。HSD(半精度、单精度和双精度)混合精度SpMV计算流程在图6中给出，其他混合精度组合的SpMV计算流程与HSD版本类似。需要注意的是，TeslaV100计算内核需要使用同步原语来保证计算的正确性。这是因为默认情况下，TeslaV100上同一warp中的线程不再同步。与一般的浮点数计算不同的是，必须使用NVIDIA提供的半精度内部函数来执行向偶舍入模式的半精度加法。

在图6中，输入是单精度和双精度向量X以及包含行偏移数组、列标数组和值数组的半精度、单精度和双精度子矩阵的CSR表示。其中M表示未划分矩阵的行数，TpV(不超过32)个线程处理矩阵相同的行，这些线程将重新编号。每个线程的结果存储在寄存器中。最后，共享内存用于在Tesla P100上合并矩阵一行的部分和，而Tesla V100用寄存器和CUDAwarp级原语进行合并。具有最小重新编号的线程将合并的结果即矩阵某一行的计算结果写回结果向量。

考虑到低精度运算可能会因溢出造成一定精度损失，半精度子矩阵的SpMV计算采用单精度，单、双精度子矩阵的SpMV计算采用双精度。由于数据以低精度格式存储，即使使用更高精度的操作，本方法仍然可以提高SpMV性能。不同稀疏矩阵的TpV的设置不一定相同，要根据矩阵的稀疏性来选择。常用的设置由下式给出：

nnz表示稀疏矩阵非零元数量

numRows表示稀疏矩阵行数。

本发明提出了一种新的基于IEEE 754标准浮点格式的稀疏矩阵划分方法。双精度输入稀疏矩阵被划分并转换为几个不同精度的子矩阵，而不会损失精度。基于稀疏矩阵的这种新的混合精度表示，设计了基于GPU的混合精度SpMV，验证了该混合精度方法的有效性。对于SuiteSparse矩阵集合中2,288个稀疏矩阵，以双精度SpMV为基准，在Tesla P100上实现了1.50倍的几何平均加速比，在Tesla V100上实现了1.44倍的几何平均加速比。

以上所述的实施例仅是对本发明优选方式进行的描述，并非对本发明的范围进行限定，在不脱离本发明设计精神的前提下，本领域普通技术人员对本发明的技术方案做出的各种变形和改进，均应落入本发明权利要求书确定的保护范围内。

Claims (7)

1.应用于计算设备的混合精度SpMV优化系统，其特征在于，包括:获取模块、第一处理模块、第二处理模块；

所述获取模块用于获取计算设备输入稀疏矩阵；

所述第一处理模块用于基于稀疏矩阵中非零元素的浮点数精度，将所述稀疏矩阵划分为不同精度的子矩阵；

所述第二处理模块用于计算不同精度的所述子矩阵与向量的乘法，得到混合精度SpMV计算结果。

2.根据权利要求1所述的应用于计算设备的混合精度SpMV优化系统，其特征在于，所述第一处理模块对所述稀疏矩阵的划分，在所述混合精度SpMV计算之前只执行一次。

3.根据权利要求1所述的应用于计算设备的混合精度SpMV优化系统，其特征在于，所述浮点数包括：符号、指数和尾数。

4.根据权利要求1所述的应用于计算设备的混合精度SpMV优化系统，其特征在于，不同精度的所述子矩阵包括：半精度子矩阵、单精度子矩阵和双精度子矩阵。

5.根据权利要求4所述的应用于计算设备的混合精度SpMV优化系统，其特征在于，不同精度的所述子矩阵均为CSR格式。

6.根据权利要求5所述的应用于计算设备的混合精度SpMV优化系统，其特征在于，所述混合精度SpMV计算过程为：

获取所述半精度子矩阵、所述单精度子矩阵和所述双精度子矩阵以及单精度和双精度向量；

基于划分后有多个子矩阵，设置计算矩阵一行与向量相乘的线程数，同时将线程重新编号，并合并各所述子矩阵相同行的计算结果；

基于最小重新编号的线程，将合并的所述计算结果写回所述结果向量，同时各所述子矩阵的行数自增，进行下一次的线程处理，当各所述子矩阵的行数大于所述稀疏矩阵行数时，停止混合精度SpMV计算。

7.应用于计算设备的混合精度SpMV优化方法，其特征在于，包括以下步骤：

获取计算设备输入稀疏矩阵；

基于稀疏矩阵中非零元素的浮点数精度，将所述稀疏矩阵划分为不同精度的子矩阵；

计算不同精度的所述子矩阵与向量的乘法，得到混合精度SpMV计算结果。

Images