1

我有以下代码,其中我使用命令 scipy.linalg.lu() 计算给定方阵的 L 矩阵,然后我再次做同样的事情,除非然后使用 scipy 应用于给定矩阵的稀疏形式。 sparse.linalg.slu()。这是代码:

import numpy as np
from scipy.sparse.linalg import splu
from scipy.sparse import csc_matrix
import scipy.linalg
A1 = csc_matrix([[1., 0, 0.], [5., 0, 2], [0, -1., 0]])
A2 = np.array([[1., 0, 0.], [5., 0, 2], [0, -1., 0]])
B = splu(A1)
P,L,U = scipy.linalg.lu(A2)
print(L);print(csr_matrix.todense(B.L))

它返回以下内容:

[[ 1.   0.   0. ]
 [ 0.   1.   0. ]
 [ 0.2 -0.   1. ]]
[[1. 0. 0.]
 [0. 1. 0.]
 [0. 0. 1.]]

正如我们所见,这些矩阵并不相同。我是误解了这两个命令的作用还是有其他问题?任何帮助表示赞赏。谢谢!

4

1 回答 1

0

I think the answer here is that the "SuperLU" decomposition of a sparse matrix requires permutations of both the rows and columns (see the docs):

Pr * A * Pc = L * U

These are provided by the mapping of indices in the perm_r and perm_c attributes. So,

Pr = csc_matrix((3,3))
Pr[B.perm_r, np.arange(3)] = 1
Pc = csc_matrix((3,3))
Pc[np.arange(3), B.perm_c] = 1

(Pr.T @ B.U @ B.L @ Pc.T).A

gives, as required:

array([[ 1.,  0.,  0.],
       [ 5.,  0.,  2.],
       [ 0., -1.,  0.]])

The same as the non-sparse result which requires only permutation of the L-matrix, P @ L @ U.

于 2018-12-07T09:12:11.800 回答