C 中浮点运算的下溢错误

Question

我是 C 新手，我的任务是创建一个函数

f(x) = sqrt[(x^2)+1]-1

可以处理非常大的数字和非常小的数字。我正在在线界面上提交我的脚本来检查我的答案。

对于非常大的数字，我将表达式简化为：

f(x) = x-1

只需使用最高功率即可。这是正确的答案。

同样的逻辑不适用于较小的数字。对于小数字（大约为 1e-7），即使在平方之前，它们也会很快被截断为零。我怀疑这与 C 中的浮点精度有关。在我的教科书中，它说 float 类型的最小可能值是 1.17549e-38，具有 6 位精度。因此，虽然 1e-7 比 1.17e-38 大得多，但它的精度更高，因此四舍五入为零。这是我的猜测，如有错误请指正。

作为解决方案，我认为当 x < 1e-6 时，我应该将 x 转换为 long double 。但是，当我这样做时，我仍然遇到同样的错误。有任何想法吗？如果我能澄清的话请告诉我。代码如下：

#include <math.h>
#include <stdio.h>

double feval(double x) {
    /* Insert your code here */
    if (x > 1e299) 
    {;
        return x-1;
    }
    if (x < 1e-6)
    {
        long double g;
        g = x;
        printf("x = %Lf\n", g);
        long double a;
        a = pow(x,2);
        printf("x squared = %Lf\n", a);
        return sqrt(g*g+1.)- 1.;
    }
    else
    { 
        printf("x = %f\n", x);
        printf("Used third \n");
        return sqrt(pow(x,2)+1.)-1;
    }
}

int main(void)
{
    double x;
    printf("Input: ");
    scanf("%lf", &x);
    double b;
    b = feval(x);
    printf("%f\n", b);
    return 0;
}

Answer 1

对于小输入，执行 1+x^2 时会出现截断错误。如果x=1e-7f，将很高兴地适合 32 位浮点数（由于没有精确的浮点表示，因此会x*x产生一点误差，但会比 1 小得多，以至于浮点精度不足以代表。1e-7x*x1+x*x

对 sqrt(1+x^2) 进行泰勒展开会更合适，其最低阶为

sqrt(1+x^2) = 1 + 0.5*x^2 + O(x^4)

然后，你可以将结果写为

sqrt(1+x^2)-1 = 0.5*x^2 + O(x^4),

避免将非常小的数字加到 1 的情况。

作为旁注，您不应该使用pow整数幂。对于 x^2，你应该这样做x*x。有效地执行任意整数幂有点棘手；例如，GNU 科学库有一个有效计算任意整数幂的函数。