首页 > 解决方案 > C中浮点运算中的下溢错误

问题描述

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

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

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

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

f(x) = x-1

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

同样的逻辑不适用于较小的数字。对于小数(大约 1e-7),它们会很快被截断为零,甚至在平方之前。我怀疑这与 C 中的浮点精度有关。在我的教科书中,它说浮点类型的最小可能值为 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;
}

标签: cfloating-pointunderflow

解决方案


对于小输入,执行 1+x^2 时会出现截断错误。如果x=1e-7f,x*x将很高兴地适合 32 位浮点数(由于1e-7没有精确的浮点表示,但x*x会比 1 小得多,浮点精度不足以代表1+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 科学库具有高效计算任意整数幂的功能。


推荐阅读