我需要解决数千次 Ax=b 类型的小线性系统。这里 A 是一个不小于 3x3 且最大为 8x8 的矩阵。我知道这个http://www.johndcook.com/blog/2010/01/19/dont-invert-that-matrix/所以我认为即使矩阵很小，反转矩阵也不明智吗？那么最有效的方法是什么？我正在用 Fortran 编程，所以我应该使用 lapack 库对吗？我的矩阵是完整的并且通常是非对称的。谢谢。
问问题
2558 次
2 回答
2
警告：我没有对此进行广泛研究，但我很乐意分享一些经验。
以我的经验，解决 3x3 系统的最快方法基本上是使用Cramer 规则。如果您需要求解具有相同矩阵 A 的多个系统，则需要预先计算 A 的逆矩阵。这仅适用于 2x2 和 3x3。
如果您必须用相同的矩阵求解多个 4x4 系统，那么再次使用逆运算明显快于 LU 的正向和反向替换。我似乎记得它使用的操作更少，实际上差异更大（再次，根据我的经验）。随着矩阵大小的增长，差异缩小，并且渐近差异消失。如果您正在求解具有差异矩阵的系统，那么我认为计算逆矩阵没有优势。
在所有情况下，与使用 LU 分解相比，使用逆解系统的准确度可能会低得多，因为 A 是相当病态的。因此，如果准确性是一个问题，那么 LU 分解绝对是要走的路。
于 2013-04-02T09:05:22.327 回答
1
The LU factorization sounds like just the ticket for you, and the lapack routine dgetrf will compute this for you, after which you can use dgetrs to solve that linear system. Lapack has been optimized to the gills over the years, so in all likelihood you are better using that than writing any of this code yourself.
The computational cost of computing the matrix inverse and then multiplying that by the right-hand side vector is the same if not more than computing the LU-factorization of the matrix and then forward- and back-solving to find your answer. Moreover, computing the inverse exhibits even more bizarre pathological behavior than computing the LU-factorization, the stability of which is still a fairly subtle issue. It can be useful to know the inverse for small matrices, but it sounds like you don't need that for your purpose, so why do it?
Moreover, provided there are no loop-carried dependencies, you can parallelize this using OpenMP without too much trouble.
于 2013-04-02T06:12:00.007 回答