如何理解scipy.linalg.lu_factor的枢纽矩阵?
问题描述:
我该如何手动重建一个矩阵A,该因素被分解为lu_factor? (甲 = PLU)如何理解scipy.linalg.lu_factor的枢纽矩阵?
我的当前的尝试都失败了由于矩阵P的设置。以下是我迄今为止:
A = np.random.rand(3,3)
lu, piv = lu_factor(A)
U = np.triu(lu)
L = np.tril(lu, -1)
L[np.diag_indices_from(L)] = 1.0
我要找的矩阵P这使得该行打印真:
print np.allclose(A, np.dot(P, np.dot(L, U)))
任何暗示/链接/建议表示赞赏!
答
置换矢量需要按顺序解释。如果在顺序进行piv=[1,2,2]
然后将下面的需要(与基于零索引):
- 行0的变化以行1
- 新行1度的变化与第2行和
- 新行2保持不变。
在这个代码会做的伎俩:
P = np.eye(3)
for i, p in enumerate(piv):
Q = np.eye(3,3)
q = Q[i,:].copy()
Q[i,:] = Q[p,:]
Q[p,:] = q
P = np.dot(P, Q)
对于piv=[1,2,2]
P是
[[ 0. 0. 1.]
[ 1. 0. 0.]
[ 0. 1. 0.]]
这可能不是计算P的一个非常快速的方式,但它确实诀窍并回答问题。
+0
+1很好的答案! – 2014-09-19 09:01:09
+1
除了交换Q行并设置P = P.dot(Q),您可以按照与行描述相同的顺序就地交换P的*列*。 – 2014-09-19 14:06:10
从'lu'和'piv'重构P,L和U的另一种方法是调用'scipy.linalg.lu'(http://docs.scipy.org/doc/scipy/reference/generated/scipy .linalg.lu.html)。 – 2014-09-19 14:22:46