当预期结果有限时返回“inf”

html5 • • 问答

输出不是 Python/Mathematica 之前在调用test_func(3000, 10). 相反，我的代码返回inf，我不知道为什么。

#include <iostream>
#include <iomanip>
#include <fstream>
#include <string>
#include <chrono>
#include <cmath>

using namespace std;

//I0
double I0 (double y0, double y)
{
    return(
        log ( (exp(y0+y)-1)/(exp(y0-y)-1) ) - y
        );
}

//test_func
double test_func (double k, double M)
{
    double k0 = sqrt(pow(k, 2.0)+pow(M, 2.0));
    cout << "k0: " << k0 << endl;
    double Izero = I0(k0, k);
    cout << "I0: " << Izero << endl;
    double k3 = pow(k, 3.0);
    cout << "k3: " << k3 << endl;
    return(
        Izero/k3
        );
}

int main ()
{
    cout << test_func(3000, 10) << endl;
    return 0;
}

我得到的输出是

k0: 3000.02
I0: inf
k3: 2.7e+10
inf

但I0应该是3004.1026690762033，而函数的结果应该是test_func(3000, 10)=1.1126306181763716e-07。我很困惑。你知道它有什么问题吗？我是 C++ 初学者，所以非常欢迎任何帮助。

回答

double范围有限。在大多数现代实现中，doubles 采用标准 IEEE 浮点binary64格式，范围可达约1.8e308。任何更大的（例如你的exp(6000)）四舍五入到无穷大。

切换到更大的类型long double可能不是最好的主意。尽管浮点类型的范围大致随大小呈指数增长，但您使用该exp函数的事实很容易击败额外的精度。例如，在 80 位long double（范围可达 about 1.2e4932）的实现上，修改test_func和I0使用long double仍然无法评估test_func(5700, 10)。

相反，可以重新设计I0以避免大量数字。让我们从拆分log.

double I0(double y0, double y) {
    return log(exp(y0+y)-1) - log(exp(y0-y)-1) - y;
}

当您计算 时log(exp(y0+y) - 1)，如果exp(y0+y)给出无穷大，您可以使用y0+y代替来恢复计算log。即我们忽略了-1，因为如果我们的数字那么大， 的精度double不足以实际记录最终结果中的差异。此外，您可能希望将两个exp(x) - 1s都替换为expm1。这是因为当x接近 时0，exp(x) - 1往往会失去精度。例如，exp(1e-16) - 1 == 0但expm1(1e-16) > 0假设 IEEE。我怀疑这对你来说没有那么重要。

double I0(double y0, double y) {
    double num = expm1(y0 + y), den = expm1(y0 - y);
    num = isinf(num) ? y0 + y : log(num);
    den = isinf(den) ? y0 - y : log(den); // though I suspect the domain of I0 is such that you don't actually need this
    return num - den - y;
}

这只是一个非常基本的更正。从浮点中挤出最大的正确性通常是非常困难的，并且是它自己的整个编程领域。但是，这足以使您的案例有效，甚至可以处理那些幼稚long double也失败的大型输入。（我不是专家，但我也注意到这y0-y本身是有问题的，因为y0显然被选择为接近y。它们之间的减法失去了可能（也可能不会）破坏结果的精​​度。）

如果您不想重写公式以适应浮点数的限制（完全可以理解，考虑到可能存在错误！），我建议您遵循@MM 的建议并使用像 MPFR 这样的任意精度数学库。这是类似于更换double用long double，但现在你可以保持在抛出问题位，直到他们离开，而你最终会跑出来内置浮点类型。

以上是当预期结果有限时返回“inf”的全部内容。
THE END
分享
二维码
< <上一篇
下一篇>>