Matlab中有一种快速求逆矩阵的方法吗？

Question

Matlab中有一种快速求逆矩阵的方法吗？

algorithmmatlabmatrixperformancematrix-inverse

13

我有很多大型矩阵(约为5000 x 5000)需要在Matlab中求逆。实际上，我需要的是矩阵的逆，所以不能使用mldivide，因为它只适用于解决Ax=b这个方程，但对于单个b的求解速度要快得多。

我的矩阵来自一个问题，这意味着它们具有一些良好的性质。首先，它们的行列式为1，因此它们可逆。但是它们不可对角化，否则我会尝试将其对角化、求逆，然后再恢复成原始形式。它们的元素都是实数(实际上是有理数)。

我正在使用Matlab获取这些矩阵，并进行与它们的逆相关的操作，因此我希望能够加快Matlab的运行速度。但如果有其他更快的语言可用，请告诉我。我不懂很多其他编程语言(只会一点C和一点Java)，所以如果在其他语言中实现比较复杂，我可能无法使用它。但请仍然提供建议。

- Daniel

3

Matlab的线性代数例程通常已经被优化得很好了，所以您可能无法通过在另一种语言中手动重新实现来获得任何速度提升。不过，您可以考虑使用Eigen（C++）库。对于Matlab，您可以尝试使用并行处理工具箱。 - Jonas

2

此外，如果您告诉我们矩阵具有哪些属性或如何获得这些属性，可能会有所帮助。在将其插入MATLAB之前进行一些线性代数操作可能比任何库都更有帮助。 - user616736

2

你可能需要更详细地阐述你的特定情况。你确定你需要求逆吗？谢谢。 - eat

矩阵有多稀疏？矩阵是否具有任何特殊形状或条目之间的关系？ - nibot

既然您说您需要多次使用每个矩阵来解决方程，那么LU分解怎么样？http://www.mathworks.com/help/techdoc/ref/lu.html - Luka Rahne

显示剩余2条评论

3个回答

7

反转一个任意的5000 x 5000矩阵无论使用什么语言都不是计算上容易的。我建议研究近似方法。如果您的矩阵是低秩的，可以尝试低秩逼近M = USV'。

以下是来自math-overflow的更多想法：

https://mathoverflow.net/search?q=matrix+inversion+approximation

- BlessedKey

2

仅因矩阵的det(A)==1并不意味着低秩逼近对于工程目的不够好。例如，看看A = diag([1e5,1,1e-5])。A的行列式仍然为1，但如果您将其用于某些目的，则A也可以是diag([1e5 0 0])。 - BlessedKey

事实上，det A = 1并不意味着A在数值上是可逆的。对于大矩阵而言，条件数非常大的概率足以使大多数矩阵不可逆（除非它们来自可逆算子的离散化）。SVD方法可能是唯一明智的选择，除非OP提供了问题的背景信息。 - Alexandre C.

@AlexandreC.：SVD方法是什么？SVD代表什么？ - HelloGoodbye

@hellogoodbye：这是来自帖子（奇异值分解）的M=USV'方法。 - Alexandre C.

1

首先假设所有特征值都为1。设A是矩阵的Jordan标准型。然后您可以仅使用矩阵乘法和加法计算A^{-1}

A^{-1} = I + (I-A) + (I-A)^2 + ... + (I-A)^k

其中k < dim(A)。为什么这样做有效呢？因为生成函数很棒。回想一下展开式。

(1-x)^{-1} = 1/(1-x) = 1 + x + x^2 + ...

这意味着我们可以使用无限和来反转 (1-x)。如果您想要反转矩阵 A，那么您需要执行以下操作：

A = I - X

解出X得到X=I-A。因此，通过代入，我们有

A^{-1} = (I - (I-A))^{-1} = 1 + (I-A) + (I-A)^2 + ...

在这里，我刚刚使用了单位矩阵I代替数字1。现在我们需要处理收敛问题，但实际上这并不是一个问题。根据假设，A处于Jordan形式，并且所有特征值都等于1，我们知道A是上三角矩阵，对角线上都是1。因此，I-A是上三角矩阵，对角线上都是0。因此，I-A的所有特征值都是0，因此其特征多项式是x^dim(A)，其最小多项式是x^{k+1}，其中k < dim(A)。由于矩阵满足其最小（和特征）多项式，这意味着(I-A)^{k+1} = 0。因此，上述级数是有限的，最大的非零项是(I-A)^k。所以它收敛。

现在，对于一般情况，将您的矩阵转换为Jordan形式，这样您就会得到一个块三角矩阵，例如：

A 0 0
0 B 0
0 0 C

每个块在对角线上有一个单一的值。如果A的值为a，则使用上述技巧将1/a * A反转，然后再通过乘以a恢复。由于完整的矩阵是块三角形，因此其逆矩阵将是

A^{-1} 0      0
0      B^{-1} 0
0      0      C^{-1}

拥有三个块并没有什么特别之处，因此无论你有多少个块都可以使用这种方法。

请注意，只要您有一个约旦形式的矩阵，这个技巧就可以起作用。在Matlab中，这种情况下求逆运算非常快，因为它只涉及矩阵乘法，而且您甚至可以使用技巧来加速计算，因为您只需要单个矩阵的幂次。但是，如果将矩阵转换为约旦形式真的很昂贵，那么这可能对您没有帮助。

- PengOne

1

如果行列式为1且元素为实数，为什么特征值应该是±1？尝试 A=[1,3;2,7]。其特征多项式为 1 - 8*x + x^2。它的根是 4±sqrt(15) 而不是 ±1。 - user616736

@Yoda：某种程度上，我在考虑整数...说得好。我会编辑我的回答。 - PengOne

@Alexandre C: 什么部分？可能我在一般情况下有点仓促，但当所有特征值都为“1”时，这个解决方案确实可行。 - PengOne

@PengOne：你在此期间编辑了答案。你说：“由于det（A）= 1，特征值为+-1”。除此之外，计算Jordan形式比计算逆矩阵要慢得多。我不知道Jordan形式的稳定性如何，但我怀疑它非常依赖于A的条件数，而A的条件数可能非常高，因为矩阵很大。但是，如果所有特征值都为1，则A = I + N，其中N为幂零矩阵，并且级数技巧有效（尽管在最坏情况下需要5000个项，这与普通求逆矩阵具有相同的O（N ^ 3）复杂度）。 - Alexandre C.

@Alexandre C：您的评论是在我的编辑之后发表的，而且它指出整个答案（“这个”）都是错误的。是的，计算Jordan形式可能是代价高昂的，但基于问题中“它们不可对角化，否则我会尝试对角化它们、求逆然后再把它们放回去。”这一句话，我认为值得提及这种方法。 - PengOne

@PengOne：好的，我同意。我怀疑简单的求逆操作会比矩阵缩放更优，而且可能存在稳定性问题（对于大矩阵，det A = 1并不保证A在数值上可逆，因为这取决于条件数）。在我看来，正确的方法是使用SVD，但这取决于OP想要做什么。 - Alexandre C.

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

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

- user616736 · Accepted Answer

实际上我需要求逆矩阵，所以无法使用mldivide，...

这并不是真的，因为你仍然可以使用mldivide来获取逆矩阵。请注意A^-1 = A^-1 * I。在MATLAB中，这相当于：

invA = A\speye(size(A));

在我的机器上，对于一个 5000x5000 的矩阵，这需要大约 10.5 秒钟的时间。请注意，MATLAB 有一个 inv 函数用于计算矩阵的逆。虽然这将需要大约相同的时间，但从数值精度方面来说效率较低（链接中提供更多信息）。

首先，它们的行列式是1，所以它们肯定可逆。

与其说是 det(A)=1，不如说是你的矩阵的条件数决定了逆矩阵的精度或稳定性。请注意，det(A)=∏_i=1:n λ_i。因此，只需设置 λ₁=M，λ_n=1/M 和 λ_i≠1,n=1 即可得到 det(A)=1。但是，随着 M → ∞，cond(A) = M² → ∞，λ_n → 0，意味着你的矩阵逼近奇异性，在计算逆矩阵时会出现较大的数值误差。

我的矩阵来自一个问题，这意味着它们具有某些好的性质。

当然，如果你的矩阵是稀疏的或具有其他优越的属性，就可以采用其他更有效的算法。但是，在没有关于你特定问题的其他信息的情况下，就无法做出更多的说明。

我希望找到一种加速 MATLAB 的方法。

MATLAB 使用高斯消元法通过 mldivide 计算一般矩阵的逆（满秩、非稀疏、没有任何特殊属性），时间复杂度为 Θ(n³)，其中 n 是矩阵的大小。因此，在你的情况下，n=5000，需要进行 1.25 x 10¹¹ 次浮点运算。因此，在拥有约 10 Gflops 计算能力的合理机器上，计算逆矩阵至少需要 12.5 秒钟，并且除非利用“特殊性质”（如果可利用），否则无法解决这个问题。