double的数据，进行数值计算的误差为什么不同？

jiewenji

提示：代码在最下方！matlab2015b  cpui7-6700  win10 64位家庭版

A矩阵
1 1
2 -1
-2 4

b向量
1
2
7

现在要求Ax=b的近似解
第一种方法x=(\(A^{T }\)A)\(^{-1}\)\(A^{T}\)b
求出x1=
1
2

第二种方法对A进行QR分解得到Q=
1/3  2/3
2/3  1/3
-2/3  2/3
Qx=b
x=Q'b
x2=
-3
6

两种方法得到的近似解x虽然不同，但是得到的近似值(b在A空间的投影=b在Q空间的投影)应该一样！

但是实际在matlab中计算Ax1=
ans =


       3      
       1/4503599627370496
       6      



Qx2=
ans =


       3      
       0      
       6     



以下是各变量在内存中的数值和数据类型。其中y是对应的上文中的x2 ， x对应上文中的x1  ，不必在意A和Q本身是宽矩阵还是窄矩阵，在matlab中运算时都进行了适当的转置。要不然投影计算结果也不会仅仅是细微的误差。


计算之前使用了format rat命令，我以为这样可以提高计算精度。但是没想到是这样的结果。很想知道产生误差的原因，为什么Ax1 计算结果的第二项会得到 1/4503599627370496  而不是0，但是第二种方法就可以得到0。我看数据类型都是double
说起来Q矩阵里面存的还都是小数，直觉上计算误差应该更大一些。为什么？

说起来A‘x1的第二项进行的就是1*2+2*(-1)这样的运算。我在matlab中直接计算也可以得到0。就很奇怪！



下方代码中用的变量于截图中的变量名不同。但是可以复现出同样的误差。！
方法1代码

1. A=[1 2 -2; 1 -1 4];
   
2.   b=[1 2 7];
   
3.   b=b'
   
4. A=A';
   
5. B=inv(A'*A);
   
6. x=B*A'*b
   
7. p=A*x

复制代码

方法2代码

1. Q=[1/3 2/3 -2/3;2/3 1/3 2/3];
   
2. b=[1 2 7];
   
3. b=b';
   
4. x=Q*b
   
5. p=Q'*x

复制代码

.·.·.

> factor(4503599627370496)
cpu time = 1 ms, real time = 7 ms.
%1 =
[2 52]

明显是计算时候产生了舍入误差

jiewenji

.·.·. 发表于 2021-10-19 14:50
> factor(4503599627370496)
cpu time = 1 ms, real time = 7 ms.
%1 =

确实产生了舍入误差。但是为什么两种算法误差不一样。其实是可以说清楚的。

ShuXueZhenMiHu

可能问题出在inv上。