可以说我定义了一个大的二次矩阵(例如150x150)。 一次它是一个numpy数组(矩阵A),一次是scipy稀疏数组(矩阵B)。


import numpy as np
import scipy as sp
from scipy.sparse.linalg import spsolve
size = 150
A = np.zeros((size, size))
length = 1000
# Set some random numbers at random places
A[np.random.randint(0, size, (2, length)).tolist()] = \
np.random.randint(0, size, (length, ))
B = sp.sparse.csc_matrix(A)

现在,我计算两个矩阵的逆。 对于矩阵B,我使用两种方法来计算逆(sp.sparse.linalg.inv和spsolve)。


5epsilon = 10.**-8 # Is needed to prevent singularity of each matrix
inv_A = np.linalg.pinv(A+np.eye(size)*epsilon)
inv_B = sp.sparse.linalg.inv(B+sp.sparse.identity(size)*epsilon)
inv_B2 = spsolve(B+sp.sparse.identity(size)*epsilon, sp.sparse.identity(size))

为了检查A和B的两个逆数是否相等,我将求和的平方差相加。


# Is not equal zero, question: Why?
# Sometimes very small (~+-10**-27), sometimes very big (~+-10**5)
print("np.sum((inv_A - inv_B)**2): {}".format(np.sum((inv_A - inv_B)**2)))
# Should be zero
print("np.sum((inv_B - inv_B2)**2): {}".format(np.sum((inv_B - inv_B2)**2)))

问题是这样的:如果我使用小的矩阵,例如 10x10,numpy作为scipy逆函数之间的误差很小(大约?+ -10 **-32)。 但是我需要大型矩阵(例如500x500)的稀疏版本。

我在这里做错了吗,还是有可能在python中计算稀疏矩阵的正确逆?

相对误差有多大?

您是说我的问题可能存在的相对误差,还是两个矩阵之间的相对误差?

我的意思是,您计算出的误差除以您要比较的两个反函数之一的平方欧几里德长度。

np.allclose是检查浮点相等性的便捷工具。

对于10x10矩阵,我得到例如100次,对于np.allclose始终为True。对于150x150,有时为True,有时为False。对于较大的像素(如200x200),它可能会返回更多的False而不是True。

实际上,一个可能的问题是您的矩阵在构造上几乎是奇异的。我想,反转几乎是奇异的矩阵会使该算法非常困难。通常,在这种情况下,与"正常"操作数相比,机器运算中不可避免的微小误差会更容易累积和爆发。您为什么不尝试一个表现更好的示例?

我上传了一张矩阵图片pl.vc/1crz9p。在此图片中,蓝色代表值-2,红色阴影为+ 4,+ 6,+ 8。可以看到,它是一个非常对角的矩阵。另外,当我用numpy计算逆时,也没有问题。矩阵是不可逆的。

您对这些解释满意吗?

您的回答对我来说很有意义(也让我很高兴)。但是我真的无法理解,为什么numpy反演比from scipy更好或更精确(对我而言这没有任何意义)。

您怎么知道哪个更精确?你把它们相乘了吗?像A A ^ -1或A ^ -1 A?

标题问题的答案是:由于不幸选择了示例矩阵。 让我详细说明。

机器精度是有限的,因此浮点运算很少会100%准确。 你试一试


>>> np.linspace(0, 0.9, 10)[1:] == np.linspace(0.1, 1, 10)[:-1]
array([ True, True, True, True, True, False, True, True, True], dtype=bool)

通常,这没有问题,因为错误太小而无法察觉。

但是,对于许多计算而言,有些输入难以处理,并且可能会过度拉伸数值算法。 这当然适用于矩阵求逆,并且您很不幸无法选择如此困难的输入。

实际上,您可以通过查看矩阵的奇异值来检查矩阵是否可能处于"病态"状态。 以下是您的脚本生成的几个矩阵的矩阵条件编号(size=200;行为良好的矩阵的值更接近于1)


9971899214237.0
5.0134186641e+12
36848.0807109
958492416768.0
1.66615247737e+16
1.42435766189e+12
1954.62614384
2.35259324603e+12
5.58292606978e+12

切换到行为良好的矩阵,您的结果应该会大大改善。