6

我有这个问题需要解决Xin AX=B. A的数量级为 15000 x 15000，并且是稀疏且对称的。B是 15000 X 7500 并且不是稀疏的。求解 X 的最快方法是什么？

我可以想到2种方法。

  1. 最简单的方法，X = A\B

  2. 使用 for 循环，

    invA = A\speye(size(A))
for i = 1:size(B,2)
    X(:,i) = invA*B(:,i);
end

有没有比上述两种更好的方法？如果不是，我提到的两者之间哪一个最好？

4

4 回答 4

9

首先要做的事情-永远永远不会计算A的逆。除非A是对角矩阵，否则永远不会稀疏。尝试一个简单的三对角矩阵。该行本身会杀死您的代码 - 内存方面和性能方面。并且计算倒数在数值上不如其他方法准确。

一般来说，\应该适合你。MATLAB 确实认识到您的矩阵是稀疏的并执行稀疏分解。如果你给一个矩阵B作为右手边，性能比你只用一个b向量求解一个方程组要好得多。所以你这样做是正确的。您可以在这里尝试的唯一其他技术方法是显式调用lu、chol或ldl，具体取决于您拥有的矩阵，并自己执行向后/向前替换。也许你在那里节省了一些时间。

事实上，求解线性方程组的方法，尤其是稀疏系统，很大程度上取决于问题。但在我想象的几乎任何（稀疏）情况下，15k 系统的因式分解应该只需要几分之一秒。这不是当今的大型系统。如果您的代码很慢，这可能意味着您的因素不再那么稀疏了。您需要确保对矩阵进行正确重新排序，以在稀疏分解期间最小化填充（添加的非零条目）。这是关键的一步。查看此页面以了解有关如何重新排序系统的一些测试和说明。并简要查看此 SO 线程中的示例重新排序。

于 2012-10-11T08:18:04.793 回答
3

由于您可以自己回答两者中哪一个更快，因此我会尝试建议下一个选项。使用 GPU 解决它。很多细节都可以在网上找到，包括这个SO 帖子、A/b 的matlab 基准测试等。此外，还有LAMG（精益代数多重网格）的 MATLAB 插件。LAMG 是一个快速的图拉普拉斯求解器。它可以在 O(m) 的时间和存储时间内解决 Ax=b。

于 2012-10-11T04:12:23.603 回答
2

如果您的矩阵A是对称正定矩阵，那么您可以采取以下措施来高效稳定地求解系统：

  • 首先，计算 Cholesky 分解，A=L*L'。由于你有一个稀疏矩阵，并且你想利用它来加速反演，你不应该chol直接应用，这会破坏稀疏模式。相反，请使用此处描述的一种重新排序方法。
  • 然后，解决系统X = L'\(L\B)

最后，如果不处理潜在的复数值，则可以替换所有L'by L.'，这会进一步加速，因为它只是尝试转置而不是计算复共轭。

另一种选择是pcgMatlab 中的预处理共轭梯度法。这个在实践中非常流行，因为你可以用速度换取准确性，即减少迭代次数，它会给你一个（通常非常好的）近似解决方案。您也永远不需要显式存储矩阵，但如果您的矩阵不适合内存，则A只需能够计算矩阵向量乘积。A

于 2012-10-11T05:43:34.367 回答
1

如果这需要永远在您的测试中解决，您可能会进入虚拟内存来解决问题。一个 15k 方（完整）矩阵将需要 1.8 GB 的 RAM 来存储在内存中。

>> 15000^2*8
ans =
      1.8e+09

您将需要一些重要的 RAM 来解决这个问题，以及 64 位版本的 MATLAB。除非您有足够的 RAM 来解决问题，否则分解不会对您有所帮助。

如果您的矩阵是真正稀疏的，那么您是否使用 MATLAB 的稀疏形式来存储它？如果不是，则 MATLAB 不知道矩阵是稀疏的，并且不使用稀疏分解。

A 有多稀疏？很多人认为半满零的矩阵是“稀疏的”。那将是浪费时间。在这样大小的矩阵上，您需要超过 99% 的零才能真正从矩阵的稀疏分解中获益。这是因为填写。否则，得到的分解矩阵几乎总是几乎满的。

如果您无法获得更多 RAM（您知道 RAM 是 cheeeeeeeeeep，当然，一旦您考虑到尝试解决此问题所浪费的时间），那么您将需要尝试迭代求解器。由于这些工具不会分解您的矩阵，如果它真的是稀疏的，那么它们将不会进入虚拟内存。这是一个巨大的节省。

由于迭代工具通常需要预处理器才能尽可能好地工作，因此可能需要一些研究才能找到最佳的预处理器。

于 2012-10-11T08:47:59.293 回答