函数scipy.linear.lu中的上三角矩阵是否总是行阶梯形

Question

我有一个m × n的矩阵a，其中n> m，我试图通过它的行阶梯形来识别独立行。函数scipy. linear .lu返回我的矩阵的PLU分解，但U因子似乎不是阶梯形，即，枢轴不是阶梯模式。据我所知，U因子应该总是呈阶梯状。

考虑下面的例子:

from numpy import array
from scipy.linalg import lu
A = array([[1, 1, 1, 1, 0, 1, 1, 1, 1, 0],
           [1, 1, 0, 0, 1, 0, 1, 0, 1, 1],
           [1, 1, 0, 0, 0, 1, 1, 0, 0, 0],
           [0, 1, 0, 1, 1, 0, 1, 0, 0, 1],
           [1, 1, 0, 0, 1, 1, 1, 1, 1, 1]])
P, L, U = lu(A)

U因子不是行阶梯形。对于每一行k，枢轴应该总是在第k-1行枢轴的右边。注意，第五行枢轴不在第四行枢轴的右边:

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

Answer 1

我以前从未见过LUP分解，但我怀疑它并不完全像您想象的那样。我在Matlab中做了同样的分解，就像你在Python中做的那样，我得到了和你做的完全一样的结果。所以我不认为这是Python LU函数的问题。

更新:我也在R中实现了这一点(下面的代码)，并且再次U矩阵给出了与Matlab和Python相同的结果。不像其他的，它给出如下警告:

Warning message:
In .local(x, ...) :
  Exact singularity detected during LU decomposition: U[i,i]=0, i=4.

Matlab:

代码:

A = ([[1, 1, 1, 1, 0, 1, 1, 1, 1, 0],
           [1, 1, 0, 0, 1, 0, 1, 0, 1, 1],
           [1, 1, 0, 0, 0, 1, 1, 0, 0, 0],
           [0, 1, 0, 1, 1, 0, 1, 0, 0, 1],
           [1, 1, 0, 0, 1, 1, 1, 1, 1, 1]]);
[L,U,P] = lu(A);
U

Matlab结果:
U =

     1     1     1     1     0     1     1     1     1     0
     0     1     0     1     1     0     1     0     0     1
     0     0    -1    -1     0     0     0    -1    -1     0
     0     0     0     0     1    -1     0     0     1     1
     0     0     0     0     1     0     0     1     1     1

Python结果(来自你的代码):

U = 
array([[ 1.,  1.,  1.,  1.,  0.,  1.,  1.,  1.,  1.,  0.],
       [ 0.,  1.,  0.,  1.,  1.,  0.,  1.,  0.,  0.,  1.],
       [ 0.,  0., -1., -1.,  0.,  0.,  0., -1., -1.,  0.],
       [ 0.,  0.,  0.,  0.,  1., -1.,  0.,  0.,  1.,  1.],
       [ 0.,  0.,  0.,  0.,  1.,  0.,  0.,  1.,  1.,  1.]])

R:

代码:

library(Matrix)
A <- Matrix( c( 1, 1, 1, 1, 0, 1, 1, 1, 1, 0 , 1, 1, 0, 0, 1, 0, 1, 0, 1, 1 , 1, 1, 0, 0, 0, 1, 1, 0, 0, 0 , 0, 1, 0, 1, 1, 0, 1, 0, 0, 1 , 1, 1, 0, 0, 1, 1, 1, 1, 1, 1 ),nrow=5,ncol=10,byrow=TRUE )
expand(mylu <-lu(A))

结果:

class "dtrMatrix" (unitriangular)
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    .    .    .    .
[2,]    0    1    .    .    .
[3,]    1    0    1    .    .
[4,]    1    0    1    1    .
[5,]    1    0    1    0    1
$U
5 x 10 Matrix of class "dgeMatrix"
     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,]    1    1    1    1    0    1    1    1    1     0
[2,]    0    1    0    1    1    0    1    0    0     1
[3,]    0    0   -1   -1    0    0    0   -1   -1     0
[4,]    0    0    0    0    1   -1    0    0    1     1
[5,]    0    0    0    0    1    0    0    1    1     1
$P
5 x 5 sparse Matrix of class "pMatrix"
[1,] | . . . .
[2,] . . . | .
[3,] . . | . .
[4,] . | . . .
[5,] . . . . |

Answer 2

不能保证LU分解会得到行阶梯形的U矩阵。它可能是，它在这种情况下失败的原因是因为矩阵是奇异的;我强烈地怀疑这就是原因，因为奇异矩阵的特征之一是不存在完整的轴心集。

请参阅https://en.wikipedia.org/wiki/LU_decomposition的通用矩阵部分以及其中提供的参考资料了解详细信息。