计算科学 - 机器精度和局部误差 - 吾爱随笔录

机器精度和局部误差

计算科学 C++ 误差估计时间积分库达龙格库塔

2021-12-17 02:20:56

我正在使用 RKF45 集成器，我在 GPU 上使用 CUDA C++ 进行了编程，并且正在思考一些问题，因为我试图追踪我的代码的一些问题。

我在计算中使用双精度值。如果机器精度是双精度，我是否可以要求我的积分器在 4 阶和 5 阶估计之间的差异 < 1 x 10 ^-16？
如果我确实将局部误差设为 10E-16 并且我积分超过例如 20,000 个时间步，那么我是否应该期待这两个数字乘积的最大全局误差：2E-12？

确实，如果是这种情况，那还不错，但是我正在尝试跟踪泄漏的来源，其中我的数字变得不切实际，并且我正在尝试排除机器精度作为该错误的根源。

编辑

我的方程式的 RHS 看起来像这样。

double frag_term = 0;

    double flux = 0;
    if (r == ((maxlength)-1))
        flux = -km*(r)*conc[r]+2*(ka)*conc[r-1]*conc[0];

    else if ( (r > ((nc)-1)) && (r != ((maxlength)-1)) )
    {
        frag_term = conc2[maxlength-1] - conc2[r];  

        flux = -(km)*(r)*conc[r] + 2*(km)*frag_term - 2*(ka)*conc[r]*conc[0] + 2*(ka)*conc[r-1]*conc[0];
    }
    else if (r == ((nc)-1))
    {
        frag_term = conc2[maxlength-1] - conc2[r];

        flux = (kn)*pow(conc[0],(nc)) + 2*(km)*frag_term - 2*(ka)*conc[r]*conc[0];
    }
    else if (r < ((nc)-1))
        flux = 0;

其中 conc 是一个浓度数组，其值我随着时间的推移进行积分， conc2 是一个部分和数组，使得 conc2[r] = 小于并包括第 r 项的所有 conc 数组项的总和。

2个回答

我在计算中使用双精度值。如果机器精度是双精度，我是否可以要求我的积分器在 4 阶和 5 阶估计之间的差异 < 1 x 10 ^-16？

假设您在谈论绝对误差容限，一般来说，是的，这是有道理的。最后一个单位（ulp）将是科学记数法的第 16 位小数，不一定对应于。假设您有一个标量 ODE，并且对于给定的时间步长，五阶方法计算的值是，而四阶方法计算的值是 . 估计的局部截断误差为，不会因浮点运算而损失精度，因此对于此示例时间步的绝对容差，而无需不得不求助于减小时间步的大小。 $10^{-16}$ $1.01 \times 10^{-15}$ $1 \times 10^{-15}$ $1 \times 10^{-17}$ $5 \times 10^{-17}$

如果我确实将局部误差设为 10E-16 并且我积分超过例如 20,000 个时间步，那么我是否应该期待这两个数字乘积的最大全局误差：2E-12？

这一切都取决于如何实现时间步长的总和。对于简单的求和，我预计舍入误差会随着时间步数大致线性增长。当采取大量时间步长时，可以使用补偿求和方法很好地控制舍入误差。（该链接参考了数值算法的准确性和稳定性，Higham 的第 2 版，第 87 页；有关舍入误差如何在各种不同数值方法中传播的更多详细信息，请参见该书。）

然而，我不一定期望全局误差线性增长，纯粹基于局部截断误差。

在我看来，比较机器精度水平是有意义的。在最坏的情况下，由于舍入错误，您将无法满足您的标准，但这是一个“安全”的类型 I 错误（尽管时间步符合要求，但仍被拒绝）。
不，这取决于您的增量函数放大了多少错误。这个数量主要取决于你右手边的 Lipshitz 常数。

最后，由于您要追求高精度，请记住，每一步方法都存在固有的不稳定性，即 $\mathcal O(\frac{\epsilon}{\tau})$ ，即与机器精度除以时间步长大小成比例的误差。这会阻止您达到非常高的精度，并且对于非常小的时间步长，可能会破坏任何收敛。

其它你可能感兴趣的问题

上一篇来自不完全 Cholesky 分解的条件数下一篇Python：使用并行编程进行多启动优化