英文:
Why is `sympy.Matrix.inv` slow?
问题
以下是您要求的内容的翻译：
这是一个小例子。令我惊讶的是它有多慢。
import numpy as np
import sympy as sp
X = np.empty(shape=(10,3), dtype=object)
for i in range(10):
for j in range(3):
X[i,j] = sp.Symbol(f'X{i}{j}')
X = sp.Matrix(X.T @ X)
print(X.inv()) # 这真的很慢...
我本以为对一个相对较小的矩阵求逆会很快。我是否漏掉了有关预期用例的某些内容？
英文:
Here is a small example. It surprised me at how slow it is.
import numpy as np
import sympy as sp
X = np.empty(shape=(10,3), dtype=object)
for i in range(10):
for j in range(3):
X[i,j] = sp.Symbol(f'X{i}{j}')
X = sp.Matrix(X.T @ X)
print(X.inv()) # This is really slow...
I would have thought that taking an inverse of a relatively small matrix would be fast. Am I missing something about the intended use case?
答案1
得分: 4
这是你的矩阵：
In [43]: X
Out[43]:
⎡ 2 2 2 2 2 2 2 2 2 2 ⎤
⎢ X₀₀ + X₁₀ + X₂₀ + X₃₀ + X₄₀ + X₅₀ + X₆₀ + X₇₀ + X₈₀ + X₉₀ X₀₀⋅X₀₁ + X₁₀⋅X₁₁ + X₂₀⋅X₂₁ + X₃₀⋅X₃₁ + X₄₀⋅X₄₁ + X₅₀⋅X₅₁ + X₆₀⋅X₆₁ + X₇₀⋅X₇₁ + X₈₀⋅X₈₁ + X₉₀⋅X₉₁ X₀₀⋅X₀₂ + X₁₀⋅X₁₂ + X₂₀⋅X₂₂ + X₃₀⋅X₃₂ + X₄₀⋅X₄₂ + X₅₀⋅X₅₂ + X₆₀⋅X₆₂ + X₇₀⋅X₇₂ + X₈₀⋅X₈₂ + X₉₀⋅X₉₂⎥
⎢ ⎥
⎢ 2 2 2 2 2 2 2 2 2 2 ⎥
⎢X₀₀⋅X₀₁ + X₁₀⋅X₁₁ + X₂₀⋅X₂₁ + X₃₀⋅X₃₁ + X₄₀⋅X₄₁ + X₅₀⋅X₅₁ + X₆₀⋅X₆₁ + X₇₀⋅X₇₁ + X₈₀⋅X₈₁ + X₉₀⋅X₉₁ X₀₁ + X₁₁ + X₂₁ + X₃₁ + X₄₁ + X₅₁ + X₆₁ + X₇₁ + X₈₁ + X₉₁ X₀₁⋅X₀₂ + X₁₁⋅X₁₂ + X₂₁⋅X₂₂ + X₃₁⋅X₃₂ + X₄₁⋅X₄₂ + X₅₁⋅X₅₂ + X₆₁⋅X₆₂ + X₇₁⋅X₇₂ + X₈₁⋅X₈₂ + X₉₁⋅X₉₂⎥
⎢ ⎥
⎢ 2 2 2 2 2 2 2 2 2 2 ⎥
⎣X₀₀⋅X₀₂ + X₁₀⋅X₁₂ + X₂₀⋅X₂₂ + X₃₀⋅X₃₂ + X₄₀⋅X₄₂ + X₅₀⋅X₅₂ + X₆₀⋅X₆₂ + X₇₀⋅X₇₂ + X₈₀⋅X₈₂ + X₉₀⋅X₉₂ X₀₁⋅X₀₂ + X₁₁⋅X₁₂ + X₂₁⋅X₂₂ + X₃₁⋅X₃₂ + X₄₁⋅X₄₂ + X₅₁⋅X₅₂ + X₆₁⋅X₆₂ + X₇₁⋅X₇₂ + X₈₁⋅X₈₂ + X₉₁⋅X₉₂ X₀₂ + X₁₂ + X₂₂ + X₃₂ + X₄₂ + X₅₂ + X₆₂ + X₇₂ + X₈₂ + X₉₂ ⎦
计算一个矩阵的符号逆矩阵通常不是一个好主意，除非矩阵非常小，稀疏或具有某种特殊的对称性等。如果矩阵是完全符号的，并且没有理由让逆矩阵中的表达式简化，那么你应该预计到逆矩阵的符号表达式会非常复杂。SymPy 的 inv 方法会尝试为你简化这些表达式，因为如果不能简化，那么表达式可能没有多大用处。简化也是为了确定矩阵是否可逆的一部分。在这里，慢的部分是简化过程。
如果你不想要简化，并且你愿意假设矩阵是可逆的，那么你可以使用标准的逆矩阵公式快速计算逆矩阵：
In [34]: %time Xinv = X.adjugate() / X.det()
CPU times: user 543 ms, sys: 957 µs, total: 544 ms
Wall time: 543 ms
即使在这里，也只花了0.5秒，因为 det 尝试进行一些简化。
我不会显示返回的 Xinv 结果，因为它非常庞大。相反，我将展示它的字符串表示有多长：
In [38]: len(str(Xinv))
Out[38]: 618867
In [39]: len(str(Xinv[0,0]))
Out[39]: 68219
即使是这个矩阵的第一个元素的字符串表示也太长了，无法粘贴到 SO 的回答中。
英文:
Here is your matrix:
In [43]: X
Out[43]:
⎡ 2 2 2 2 2 2 2 2 2 2 ⎤
⎢ X₀₀ + X₁₀ + X₂₀ + X₃₀ + X₄₀ + X₅₀ + X₆₀ + X₇₀ + X₈₀ + X₉₀ X₀₀⋅X₀₁ + X₁₀⋅X₁₁ + X₂₀⋅X₂₁ + X₃₀⋅X₃₁ + X₄₀⋅X₄₁ + X₅₀⋅X₅₁ + X₆₀⋅X₆₁ + X₇₀⋅X₇₁ + X₈₀⋅X₈₁ + X₉₀⋅X₉₁ X₀₀⋅X₀₂ + X₁₀⋅X₁₂ + X₂₀⋅X₂₂ + X₃₀⋅X₃₂ + X₄₀⋅X₄₂ + X₅₀⋅X₅₂ + X₆₀⋅X₆₂ + X₇₀⋅X₇₂ + X₈₀⋅X₈₂ + X₉₀⋅X₉₂⎥
⎢ ⎥
⎢ 2 2 2 2 2 2 2 2 2 2 ⎥
⎢X₀₀⋅X₀₁ + X₁₀⋅X₁₁ + X₂₀⋅X₂₁ + X₃₀⋅X₃₁ + X₄₀⋅X₄₁ + X₅₀⋅X₅₁ + X₆₀⋅X₆₁ + X₇₀⋅X₇₁ + X₈₀⋅X₈₁ + X₉₀⋅X₉₁ X₀₁ + X₁₁ + X₂₁ + X₃₁ + X₄₁ + X₅₁ + X₆₁ + X₇₁ + X₈₁ + X₉₁ X₀₁⋅X₀₂ + X₁₁⋅X₁₂ + X₂₁⋅X₂₂ + X₃₁⋅X₃₂ + X₄₁⋅X₄₂ + X₅₁⋅X₅₂ + X₆₁⋅X₆₂ + X₇₁⋅X₇₂ + X₈₁⋅X₈₂ + X₉₁⋅X₉₂⎥
⎢ ⎥
⎢ 2 2 2 2 2 2 2 2 2 2 ⎥
⎣X₀₀⋅X₀₂ + X₁₀⋅X₁₂ + X₂₀⋅X₂₂ + X₃₀⋅X₃₂ + X₄₀⋅X₄₂ + X₅₀⋅X₅₂ + X₆₀⋅X₆₂ + X₇₀⋅X₇₂ + X₈₀⋅X₈₂ + X₉₀⋅X₉₂ X₀₁⋅X₀₂ + X₁₁⋅X₁₂ + X₂₁⋅X₂₂ + X₃₁⋅X₃₂ + X₄₁⋅X₄₂ + X₅₁⋅X₅₂ + X₆₁⋅X₆₂ + X₇₁⋅X₇₂ + X₈₁⋅X₈₂ + X₉₁⋅X₉₂ X₀₂ + X₁₂ + X₂₂ + X₃₂ + X₄₂ + X₅₂ + X₆₂ + X₇₂ + X₈₂ + X₉₂ ⎦
Computing the symbolic inverse of a matrix is usually not a good idea unless the matrix is very small, sparse or has some special symmetry or something. If the matrix is fully symbolic and there is no reason for the expressions in the inverse to simplify then you should expect that symbolic expressions for the inverse will be very complicated. SymPy's inv method will try to simplify these for you because the expression is unlikely to be much use if it cannot be simplified. That simplification is also needed in order to determine whether the matrix actually is invertible. It is the simplification part that is slow here.
If you don't want the simplification and you are happy to assume that the matrix is invertible then you can calculate an expression for the inverse quickly using the standard inverse formula:
In [34]: %time Xinv = X.adjugate() / X.det()
CPU times: user 543 ms, sys: 957 µs, total: 544 ms
Wall time: 543 ms
Even here it is only taking 0.5 seconds because det attempts some simplification.
I won't show the result that is returned for Xinv because it is enormous. Instead I'll just show how long its string representation is:
In [38]: len(str(Xinv))
Out[38]: 618867
In [39]: len(str(Xinv[0,0]))
Out[39]: 68219
Even just the first element of this matrix has a string representation that is too long to be pasted into an answer on SO.
通过集体智慧和协作来改善编程学习和解决问题的方式。致力于成为全球开发者共同参与的知识库，让每个人都能够通过互相帮助和分享经验来进步。
评论