5

sinl当参数接近 pi 的非零倍数时，为什么会给出不正确的结果？为什么sinl当参数很大时会给出不正确的结果？下面的代码说明了这一点。

请注意，用于初始化变量 pi 的数字与任何 64 位长双精度值都不完全匹配。编译器选择最接近的值，即3.14159265358979323851280895940618620443274267017841339111328125. 可以使用 libquadmath、gnu MPFR lib 或在线计算器（例如http://www.ttmath.org/online_calculator ）找到预期的正弦值。

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

int main (int argc, char *argv [])
    {
    volatile long double pi = 3.14159265358979323846L;
    volatile long double big = 9223372035086174241L;
    volatile long double expected1 = -5.0165576126683320235E-20L;
    volatile long double expected2 = -4.2053336735954077951E-10L;
    double result;
    double ex1 = expected1, ex2 = expected2;

    result = sinl (pi);
    printf("expected: %g, \nreturned: %g\n\n", ex1, result);
    result = sinl (big);
    printf("expected: %g, \nreturned: %g\n\n", ex2, result);
    return 0;
    }

我正在使用 gcc 4.7.3。使用 volatile 可以防止编译器用sinl()硬编码结果替换调用。我的电脑有一个 Intel Core i7 处理器并运行 Windows。我将结果打印为 double 而不是 long double，因为我使用的 gcc 的 mingw 端口不支持打印 long double。这是程序输出：

expected: -5.01656e-020,
returned: -5.42101e-020

expected: -4.20533e-010,
returned: -0.011874
4

3 回答 3

3

错误可以追溯到 sinl 库代码使用的 fsin 处理器指令。正如英特尔声称的那样，指令 fsin、fcos 和 fptan 不准确到 1.0 ulp：http: //notabs.org/fpuaccuracy/

于 2013-06-02T06:29:04.453 回答
2

为了对 pi 的倍数实现 1 ULP 精度，内部常量 M_PI 应该具有大约 106 位精度（或长双精度为 128）。

在归约阶段，一个完美的实现必须在减去 后以某种方式生成缺失的 53 或 64 位精度(x - M_PI)，因为一个简单的实现会将这个中间值计算为零。当参数将是非零的大整数乘法时，问题当然会变得越来越大。

M_PI 的 66 位内部精度不足以满足 1 ULP 精度。再一次，人们可以重新阅读声明并检查是否声明了 1 ULP 的准确性相对于结果或论点。

于 2013-06-02T08:53:54.427 回答
1

GNU libc 的文档（可通过运行访问）列出了x86 和“x86_64/fpu”上info libc math errors的 1 ulp“已知错误” 。cosl它没有为sinl. cosl我可以在我的 x86_64 机器上重现大约 pi/2 的类似巨大错误。

也许您应该将此作为文档错误报告给 glibc 和 Linux 手册页的人们；我无法想象实施“正确”修复是值得的。

如果你真的想要一个快速准确的sinl，我不太确定去哪里找。 CRlibm确实sin（doubles 的变体）。 MPFR会处理它，但它会比fsin.

于 2013-06-02T08:07:48.137 回答