稀疏矩阵计算：介绍SpGEMM

SpGEMM在多重网格方法和图分析领域有着广泛的应用，引起了许多研究人员的持续关注。几十年来，研究人员设计和实现了SpGEMM计算的多种优化方法，并将其应用到特定的领域和计算架构中。本章的目标是提供一个关于SpGEMM研究的结构化的全面概述，首先举例介绍了SpGEMM的应用领域，其次总结了从1977年到2020年关于SpGEMM优化方法的研究进展，根据这些研究工作重点针对的问题和计算架构进行分类，分别从结果矩阵大小预测、矩阵划分和负载均衡、局部结果汇总、面向目标体系架构的优化这四个角度介绍了现有的研究工作。除此之外，本章还提供了目前主要的基于CPU和GPU的SpGEMM实现方案的对比实验结果，并在最后提出了SpGEMM目前面临的一些挑战。

SpGEMM有较高的计算时间开销，是代数多重网格求解器、三角形计数、多源宽度优先搜索、最短路径求解、交叉着色和子图等许多科学计算应用和图分析算法的计算核心。针对SpGEMM的优化往往能够同时推动多个应用领域的发展。

输入两个稀疏矩阵A和B，SpGEMM计算出一个稀疏矩阵C满足C=A×B（A∈Rp×q ，B∈Rq×r ，C∈Rp×r ），运算符 × 代表通用矩阵乘法。令aij代表矩阵A中第i行第j列的元素，ai*代表矩阵A的第i行，a*j代表矩阵A的第j列。在通用矩阵乘法（GEMM）中，结果矩阵C的第i行可以被表示为ci*=（ai*·b*1，ai*·b*2，…，ai*·b*r），其中，运算符·代表两个向量的点积。(www.xing528.com)

因为稠密形式的矩阵需要大量的内存和计算代价，相比于GEMM，SpGEMM中的两个输入矩阵和结果矩阵中包含大量的零元，采用稠密的矩阵存储形式往往会造成不必要的内存浪费和计算开销，所以需要充分挖掘和利用三个相关矩阵的稀疏特征。大多数关于SpGEMM的研究都是基于Gustavson，运算符*代表数乘运算，且此运算只发生在稀疏矩阵A和B的非零项中。Gustavson提出的SpGEMM算法的伪代码如算法4−1所示。的逐行SpGEMM算法提出的。结果矩阵C的第i行表示为

已有关于SpGEMM的综述文章出现较早，本章在此基础上对SpGEMM的应用领域及其最新的设计和实现方法进行概述，主要阐述了现有算法、数据结构和可用库的情况，同时尝试发现不同设计实现背后的设计决策，深入描述了针对SpGEMM计算的许多算法和可用库。