c++ - dgesv 中的舍入错误？

Question

我在 C++ 中使用 dgesv 和 dgemm fortran 子例程来做一些简单的矩阵乘法和左除法。

对于随机矩阵 A 和 B，我这样做：

A\(A\(A*B));

其中 * 使用 dgemm 定义，\ 使用 dgesv 定义。显然，这个表达式应该简化为单位矩阵。我正在针对 MATLAB 测试我的答案，我在对角线上得到了或多或少的 1，但其他条目略有偏差（数字在 e-15 的数量级上，所以它们已经接近 0） .

我只是想知道这个结果是否可以预期？因为如果我做这样的事情：

C = A+B;
D = A*B;

D\(C\(C*C));

结果应该输出到 D\C。基本上，C(C*C) 非常准确（与 MATLAB 完美匹配），但第二次我做 D\CI 得到的东西与 e-1 甚至 e+00 不同。我猜这不应该发生？

Answer 1

您的问题似乎与 C/C++ 中浮点变量的有限精度有关。你可以在这里阅读更多关于它的信息。有一些技术可以最大限度地减少这种影响（其中一些在 wiki 文章中进行了描述），但在一些操作之后总会有一些准确性损失。您可能想要使用一些支持任意精度数的第三方数学库（例如 GMP）。但是仍然——只要你坚持数值方法，你的计算的准确性总是会受到污染。