21

对于课程,我必须为稀疏矩阵编写自己的线性方程求解器。对于稀疏矩阵,我可以自由使用任何类型的数据结构,并且我必须实现几个解决方案,包括共轭梯度。

我想知道是否有一种著名的方法来存储稀疏矩阵,使得与向量的乘法相对较快。

现在我的稀疏矩阵基本上实现了一个包装std::map< std::pair<int, int>, double>来存储数据,如果有的话。这将矩阵的乘法从向量转换为 O(n²) 复杂度到 O(n²log(n)),因为我必须对每个矩阵元素执行查找。我研究了耶鲁稀疏矩阵格式,似乎元素的检索也在 O(log(n)) 中,所以我不确定它是否会更快。

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

您认为使用另一种数据结构可以更快地做到这一点吗?

谢谢!

4

1 回答 1

26

最常见的选择是CSC 或 CSR 存储。这些对于矩阵向量乘法都很有效。如果您必须自己编写这些乘法例程,也很容易编写代码。

也就是说,耶鲁存储也产生了非常有效的矩阵向量乘法。如果您正在执行矩阵元素查找,那么您误解了如何使用该格式。我建议你研究一些标准的稀疏库来了解矩阵向量乘法是如何实现的。

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

  1. 将结果向量 y 归零。
  2. 为矩阵 A 的非零元素初始化一个迭代器。
  3. 获取矩阵的下一个非零元素 A[i,j] 说。请注意, i,j 的模式不遵循常规模式。它只是反映了 A 的非零元素的存储顺序。
  4. y[i] += A[i,j]*x[j]
  5. 如果 A 的元素更多,则转到 3。

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

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

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

但我并不是建议您坚持使用std::map存储空间。这不会是超级有效的。我推荐 CSC 主要是因为它使用最广泛。

于 2013-03-28T22:57:17.033 回答