np.linalg.inv() 给出意想不到的结果

问题描述

我正在使用 numpy 对矩阵求逆,但得到了一些意想不到的结果。

运行我的程序,我得到 p1 的最终结果:

>>> p1
array([[1.69133481e+11,3.74575030e+09,8.29681977e+07,1.83800903e+06],[3.74575030e+09,1.83800903e+06,4.07236156e+04],[8.29681977e+07,4.07236156e+04,9.02416997e+02],[1.83800903e+06,9.02416997e+02,2.00000000e+01]])

然后当我尝试使用 np.linalg.inv() 反转 p1 时,我得到:

>>> np.linalg.inv(p1)
array([[ 2.33378273e+00,-3.16566294e+02,1.43119558e+04,-2.15657094e+05],[-3.16566293e+02,4.29411791e+04,-1.94139249e+06,2.92538609e+07],[ 1.43119557e+04,8.77723669e+07,-1.32261292e+09],[-2.15657092e+05,2.92538606e+07,-1.32261291e+09,1.99302540e+10]])

这显然是不正确的:

>>> (p1 @ np.linalg.inv(p1))
array([[ 9.99968764e-01,-1.28189335e-03,-7.44723976e-01,3.60136516e+00],[-2.53043124e-06,9.99986814e-01,-3.83602548e-02,1.12390996e-01],[-8.71254524e-08,9.18930849e-06,9.99317258e-01,4.33950341e-03],[-6.98491931e-10,7.45058060e-08,-1.23977661e-05,1.00001526e+00]])
>>> (p1 @ np.linalg.inv(p1)).astype(int)
array([[0,3],[0,0],1]])

我的结果显然不是单位矩阵。

然而,这就是事情变得奇怪的地方。如果我在键入“p1”时使用打印到终端的值重新定义 p1 并计算与之前相同的命令:

>>> p1 = np.array([[1.69133481e+11,...        [3.74575030e+09,...        [8.29681977e+07,...        [1.83800903e+06,2.00000000e+01]])
>>> p1 @ np.linalg.inv(p1)
array([[ 9.99999999e-01,-4.85794104e-07,-2.14803652e-05,-4.95130838e-04],[ 2.61335278e-10,9.99999940e-01,9.92065715e-07,-3.11477404e-05],[ 1.87780333e-13,-5.48610778e-10,1.00000001e+00,-1.02864912e-07],[ 9.94759830e-14,-2.00088834e-11,5.82076609e-10,9.99999998e-01]])
>>> (p1 @ np.linalg.inv(p1)).astype(int)
array([[0,1,0]])

这里,结果与我得到单位矩阵时的预期非常一致。

在重新定义 p1 之前矩阵逆不正确的事实导致我的代码出现问题。有什么想法吗?

解决方法

结果是正确的,您的问题在于浮点 (in) 精度:

In [1]: arr = np.array([[1.69133481e+11,3.74575030e+09,8.29681977e+07,1.83800903e+06],...:        [3.74575030e+09,1.83800903e+06,4.07236156e+04],...:        [8.29681977e+07,4.07236156e+04,9.02416997e+02],...:        [1.83800903e+06,9.02416997e+02,2.00000000e+01]])                                     

In [2]: inv = np.linalg.inv(arr)                                                                                      

In [3]: np.isclose(arr @ inv,np.eye(*arr.shape),atol=1e-3)                                                          
Out[3]: 
array([[ True,True,True],[ True,True]])