快速稀疏矩阵乘法

Question

快速稀疏矩阵乘法

c++matrixlinear-algebrasparse-matrixmatrix-multiplication

22

为了课程，我必须自己编写用于稀疏矩阵的线性方程求解器。我可以使用任何类型的数据结构来存储稀疏矩阵，并且我必须实现几个解算方法，包括共轭梯度法。

我想知道是否有一种著名的方式可以存储稀疏矩阵，使得矩阵与向量的乘积相对较快。

目前，我的稀疏矩阵基本上是使用std :: map< std :: pair <int，int>， double>来实现的，它存储数据（如果有）。这将把矩阵与向量的乘积从O(n²)复杂度转换为O(n²log(n))，因为我必须为每个矩阵元素执行查找。我研究了Yale Sparse矩阵格式，似乎检索元素也是在O(log(n))时间内进行的，所以我不确定它是否会更快。

作为参考，我有一个800x800的矩阵，其中填充了5000个条目。使用共轭梯度法解决这样的系统大约需要450秒。

你认为使用另一种数据结构可以更快地完成吗？

谢谢！

- lezebulon

先阅读维基百科。http://en.wikipedia.org/wiki/Sparse_matrix 它有一个很好的常见存储方法列表，可以为您提供高效的操作。 - Daniel Williams

@Song Wang：这个类的目的基本上是为了自己编写有限元方法求解器。 - lezebulon

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- David Heffernan · Accepted Answer

最常见的选择是CSC或CSR存储，对于矩阵向量乘法都非常高效。如果必须自己编写乘法程序，则这些算法也很容易编码。

尽管如此，Yale存储也可以产生非常高效的矩阵向量乘积。如果您正在执行矩阵元素查找，则可能误解了如何使用格式。建议您学习一些标准稀疏库以了解如何实现矩阵向量乘法。

即使使用当前的存储方式，您也可以在O(n)复杂度中执行矩阵乘法。我曾经见过的所有稀疏矩阵向量乘法算法都归结为相同的步骤。例如，请考虑y = Ax。

将结果向量y清零。
初始化矩阵A的非零元素的迭代器。
获取矩阵A的下一个非零元素，比如A[i,j]。请注意，i,j的模式不遵循规则模式。它只反映了A的非零元素存储顺序。
y[i] += A[i,j]*x[j]
如果还有更多的A元素，则转到第3步。

我怀疑您正在编写经典的双重for循环密集乘法代码：

for (i=0; i<N; i++)
    for (j=0; j<N; j++)
        y[i] += A[i,j]*x[j]

这就是导致您执行查找操作的原因。

但我不建议您坚持使用std::map存储。那样不会非常高效。我推荐使用CSC，因为它是最广泛使用的。