计算矩阵的零空间
我正在尝试解决一组方程,形式是 Ax = 0。这里的 A 是一个已知的 6x6 矩阵,我用下面的代码通过奇异值分解(SVD)来得到向量 x,这个方法在一定程度上是有效的。得到的答案大致是对的,但精度不够,无法满足我的需求。我该如何提高计算的精度呢?把 eps 降到 1.e-4 以下会导致函数出错。
from numpy.linalg import *
from numpy import *
A = matrix([[0.624010149127497 ,0.020915658603923 ,0.838082638087629 ,62.0778180312547 ,-0.336 ,0],
[0.669649399820597 ,0.344105317421833 ,0.0543868015800246 ,49.0194290212841 ,-0.267 ,0],
[0.473153758252885 ,0.366893577716959 ,0.924972565581684 ,186.071352614705 ,-1 ,0],
[0.0759305208803158 ,0.356365401030535 ,0.126682113674883 ,175.292109352674 ,0 ,-5.201],
[0.91160934274653 ,0.32447818779582 ,0.741382053883291 ,0.11536775372698 ,0 ,-0.034],
[0.480860406786873 ,0.903499596111067 ,0.542581424762866 ,32.782593418975 ,0 ,-1]])
def null(A, eps=1e-3):
u,s,vh = svd(A,full_matrices=1,compute_uv=1)
null_space = compress(s <= eps, vh, axis=0)
return null_space.T
NS = null(A)
print "Null space equals ",NS,"\n"
print dot(A,NS)
2 个回答
6
注意:在使用Python和Matlab的SVD(奇异值分解)时,可能会有些混淆。在Python中,使用numpy.linalg.svd(A)这个函数会返回三个矩阵u、s和v,使得u*s*v = A。严格来说,应该是dot(u, dot(diag(s), v)) = A,因为在numpy中,s是一个向量,而不是一个二维矩阵。
最上面的回答是正确的,通常你会写成u*s*vh = A,并且返回的是vh,而这个回答讨论的是v,而不是vh。
简单来说:如果你有矩阵u、s和v,使得u*s*v = A,那么v的最后几行,而不是最后几列,描述的是零空间。
补充说明:[对于像我这样的人:] 每一行的最后一个向量v0都满足A*v0 = 0(前提是对应的奇异值为0)。
10
A
是满秩的——所以 x
是 0
看起来你需要一个最小二乘解,也就是想要找到 min ||A*x|| 使得 ||x|| = 1
。你可以进行奇异值分解(SVD),这样可以得到 [U S V] = svd(A)
,然后取 V
的最后一列(假设这些列是按奇异值从大到小排序的),这就是 x
。
也就是说,
U =
-0.23024 -0.23241 0.28225 -0.59968 -0.04403 -0.67213
-0.1818 -0.16426 0.18132 0.39639 0.83929 -0.21343
-0.69008 -0.59685 -0.18202 0.10908 -0.20664 0.28255
-0.65033 0.73984 -0.066702 -0.12447 0.088364 0.0442
-0.00045131 -0.043887 0.71552 -0.32745 0.1436 0.59855
-0.12164 0.11611 0.5813 0.59046 -0.47173 -0.25029
S =
269.62 0 0 0 0 0
0 4.1038 0 0 0 0
0 0 1.656 0 0 0
0 0 0 0.6416 0 0
0 0 0 0 0.49215 0
0 0 0 0 0 0.00027528
V =
-0.002597 -0.11341 0.68728 -0.12654 0.70622 0.0050325
-0.0024567 0.018021 0.4439 0.85217 -0.27644 0.0028357
-0.0036713 -0.1539 0.55281 -0.4961 -0.6516 0.00013067
-0.9999 -0.011204 -0.0068651 0.0013713 0.0014128 0.0052698
0.0030264 0.17515 0.02341 -0.020917 -0.0054032 0.98402
0.012996 -0.96557 -0.15623 0.10603 0.014754 0.17788
所以,
x =
0.0050325
0.0028357
0.00013067
0.0052698
0.98402
0.17788
而且,||A*x|| = 0.00027528
,相比于你之前的解 x
,那时 ||A*x_old|| = 0.079442