为什么 numpy.linalg.solve() 提供比 numpy.linalg.inv() 更精确的矩阵求逆? [英] Why does numpy.linalg.solve() offer more precise matrix inversions than numpy.linalg.inv()?

问题描述

我不太明白为什么 numpy.linalg.solve() 给出了更精确的答案，而 numpy.linalg.inv() 有点崩溃，给出 (我相信是)估计.

举一个具体的例子，我正在求解方程 C^{-1} * d 其中 C 表示一个矩阵，d是一个向量数组.为了便于讨论，C 的维度是形状 (1000,1000) 和 d 是形状 (1,1000).

numpy.linalg.solve(A, b) 求解方程A*x=b for x，即x = A^{-1}* b. 因此，我可以通过

解这个方程

(1)

inverse = numpy.linalg.inv(C)结果 = 逆 * d

或 (2)

numpy.linalg.solve(C, d)

方法 (2) 给出了更精确的结果.为什么是这样?

究竟发生了什么，以至于一个比另一个工作得更好"?

解决方案

np.linalg.solve(A, b) 不不计算的逆一个.相反，它调用 gesv LAPACK 例程之一，它首先使用 LU 分解分解 A，然后使用前向和后向替换求解 x(参见 此处).

np.linalg.inv 使用相同的方法通过求解 A^-1 来计算 A 的逆em> in A·A^-1 = I 其中 I 是身份*.分解步骤与上面完全相同，但需要更多的浮点运算来求解A^-1(一个n×n矩阵) 而不是 x (一个 n 长的向量).此外，如果您想通过恒等式 A^-1·b = x 获得 x 那么额外的矩阵乘法将导致更多的浮动点运算，因此性能较慢，数值误差较大.

不需要计算A^-1的中间步骤——直接获得x更快更准确.><小时>

* inv 的相关源代码是 此处.不幸的是，它有点难以理解，因为它是 C 模板化的.需要注意的重要一点是单位矩阵作为参数 B 传递给 LAPACK 求解器.

I do not quite understand why numpy.linalg.solve() gives the more precise answer, whereas numpy.linalg.inv() breaks down somewhat, giving (what I believe are) estimates.

For a concrete example, I am solving the equation C^{-1} * d where C denotes a matrix, and d is a vector-array. For the sake of discussion, the dimensions of C are shape (1000,1000) and d is shape (1,1000).

numpy.linalg.solve(A, b) solves the equation A*x=b for x, i.e. x = A^{-1} * b. Therefore, I could either solve this equation by

(1)

inverse = numpy.linalg.inv(C)
result = inverse * d

or (2)

numpy.linalg.solve(C, d)

Method (2) gives far more precise results. Why is this?

What exactly is happening such that one "works better" than the other?

解决方案

np.linalg.solve(A, b) does not compute the inverse of A. Instead it calls one of the gesv LAPACK routines, which first factorizes A using LU decomposition, then solves for x using forward and backward substitution (see here).

np.linalg.inv uses the same method to compute the inverse of A by solving for A^-1 in A·A^-1 = I where I is the identity*. The factorization step is exactly the same as above, but it takes more floating point operations to solve for A^-1 (an n×n matrix) than for x (an n-long vector). Additionally, if you then wanted to obtain x via the identity A^-1·b = x then the extra matrix multiplication would incur yet more floating point operations, and therefore slower performance and more numerical error.

There's no need for the intermediate step of computing A^-1 - it is faster and more accurate to obtain x directly.

---

* The relevant bit of source for inv is here. Unfortunately it's a bit tricky to understand since it's templated C. The important thing to note is that an identity matrix is being passed to the LAPACK solver as parameter B.

这篇关于为什么 numpy.linalg.solve() 提供比 numpy.linalg.inv() 更精确的矩阵求逆?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！