scipy.linalg.inv函数是否检查矩阵是否是对角矩阵？

Question

获取对角矩阵的逆矩阵非常简单，不需要复杂方法。在应用更复杂的方法之前，scipy.linalg.inv 是否会检查矩阵是否为对角矩阵？还是我需要自己检查？

- Make42

Answer 1

您可以看到scipy.linalg.inv的Github代码，函数inv首先调用

getrf, getri, getri_lwork = get_lapack_funcs(('getrf', 'getri','getri_lwork'),

然后函数getrf完成了其任务并给出了LU分解等结果。现在我们需要调查getrf函数如何给出LU分解。因为如果它在处理输入矩阵之前检查它是否是对角线矩阵，那么就不需要再自己检查一遍。

getrf函数是通过调用_get_funcs获得的，但我无法从那里继续向下走（_get_funcs使用以下参数进行调用：_get_funcs(names, arrays, dtype, "LAPACK", _flapack, _clapack, "flapack", "clapack", _lapack_alias)）。

我建议您使用大型对角线矩阵运行实验，以比较使用linalg拆分输出所需的时间与手动求逆所需的时间。

---

更新（由问题作者提供）：

import numpy as np
from scipy.linalg import inv
a = np.diag(np.random.random(19999))
b = a.copy()
np.fill_diagonal(a, 1/a.diagonal())
c = inv(b)

甚至不需要时间测量工具：很明显inv要慢得多…（这令人惊讶地令人失望）。

- dallonsi

---

没错。应该在Github项目上开一个问题。也许可以考虑其他特殊情况来提高反转。 - dallonsi

Answer 2

请查看：scipy.linalg.inv。如果您在try except中使用scipy.linalg.inv，并且当矩阵a是奇异矩阵时它引发了LinAlgError异常。奇异矩阵的行列式为零。最初的回答。

try:
    # your code that will (maybe) throw  scipy.linalg.inv(your matrix)

except np.linalg.LinAlgError as err:
    # It shows your matrix is singular
    # Its determinant of a matrix is equal to zero
    # The matrix does not have an inverse.
    # You can conclude if the matrix is diagonal or not

如果一个矩阵的行列式等于零：

这个矩阵的秩小于满秩。这个矩阵是奇异的。这个矩阵没有逆矩阵。

手动翻译如下：

def is_diagonal(matrix):
    #create a dummy matrix
    dummy_matrix = np.ones(matrix.shape, dtype=np.uint8)
    # Fill the diagonal of dummy matrix with 0.
    np.fill_diagonal(dummy_matrix, 0)

    return np.count_nonzero(np.multiply(dummy_matrix, matrix)) == 0

diagonal_matrix = np.array([[3, 0, 0],
                            [0, 7, 0],
                            [0, 0, 4]])
print is_diagonal(diagonal_matrix)
>>> True

random_matrix = np.array([[3, 8, 0],
                          [1, 7, 8],
                          [5, 0, 4]])
print is_diagonal(random_matrix)
>>> False

scipy.sparse.dia_matrix.diagonal 返回矩阵的第k条对角线。

最初的回答

from scipy.sparse import csr_matrix
A = csr_matrix([[1, 2, 0], [0, 0, 3], [4, 0, 5]])
A.diagonal()
array([1, 0, 5])
A.diagonal(k=1)
array([2, 3])

此外，如果输入的数组不是方阵，则from scipy.linalg import block_diag将无法创建对角矩阵。

请在Jupyter中考虑使用%timeit yourfunctionname来查找时间复杂度。

“Original Answer”翻译成中文为“最初的回答”。

- Mahsa Hassankashi

---

我真的不明白这个评论有什么帮助。让我们再读一遍问题 :) - dallonsi

@dallonsi scipy.linalg.inv 无法帮助检查矩阵是否为对角矩阵。有一些技巧可以通过 scipy.sparse.dia_matrix.diagonal 或从 scipy.linalg 中导入 block_diag 来找到它。有任何不清楚的地方，请告诉我 :) - Mahsa Hassankashi

是的，我认为问题是：“在应用某些算法之前，scipy.linalg.inv是否检查矩阵是否为对角线矩阵。”作者询问linalg.inv是否这样做：如果是对角线矩阵，则手动通过取对角线条目d_i并将其替换为1/d_i来求逆；否则，应用一些复杂且更昂贵的算法。 - dallonsi

@dallonsi 我更新了我的答案，手动和通过我提到的方法，并让我们在Jupyter中查看它们的时间复杂度。如果您将scipy.linalg.inv放在try except中，当矩阵a是奇异的时候，它会引发LinAlgError。奇异矩阵的行列式为零。 - Mahsa Hassankashi

对角矩阵不一定是奇异的，@dallonsi说得对。 - Make42