NumPy 如何计算矩阵的逆？

Question

给定一个方阵 A 作为 NumPy 数组，比如

import numpy as np

A = np.array(
     [
         [1, 2, 3],
         [3, 4, 6],
         [7, 8, 9],
    ]
)

NumPy linalg.inv 在

时内部使用哪种算法

np.linalg.inv(A)

被调用来计算 A?

的矩阵求逆

特别地，由于矩阵求逆可能在数值上不稳定（取决于矩阵的condition number），是否有根据某些矩阵性质考虑的特殊情况？

Answer 1

根据@projjal 的评论，所有这些都等同于计算方阵的逆：

import numpy as np
from scipy.linalg import lu_factor, lu_solve

A = np.array([[1, 2, 3],[3, 4, 6],[7, 8, 9]])

A_inv_1 = np.linalg.inv(A)

A_inv_2 = np.linalg.solve(A,np.eye(A.shape[0]))

A_LU = lu_factor(A) # this way, you can potentially reuse the factorization for different RHS
A_inv_3 = lu_solve(A_LU,np.eye(A.shape[0]))

# check
np.allclose(A_inv_1,A_inv_2)
>>> True
np.allclose(A_inv_1,A_inv_3)
>>> True

Answer 2

您可能应该注意到，在 numpy 源代码的深处（参见 https://github.com/numpy/numpy/blob/master/numpy/linalg/umath_linalg.c.src），inv 例程会尝试从您的系统 LAPACK 包中调用 dgetrf 函数，然后它会对您的原始代码执行 LU 分解矩阵.

您可以查看第 1036 行代码中的注释，其中指出：

/* As in the linalg package, the determinant is computed via LU factorization
 * using LAPACK.
 * slogdet computes sign + log(determinant).
 * det computes sign * exp(slogdet).
 */
/**begin repeat
   #TYPE = FLOAT, DOUBLE, CFLOAT, CDOUBLE#
   #typ = npy_float, npy_double, npy_cfloat, npy_cdouble#
   #basetyp = npy_float, npy_double, npy_float, npy_double#
   #cblas_type = s, d, c, z#
*/

因此，正如@Yacola 所证明的那样，numpy 逆矩阵和 LU 分解的结果是相同的。这再次证明了事实。

NumPy 如何计算矩阵的逆？

How does NumPy compute the inverse of a matrix?

python

numpy

linear-algebra