我一直在研究。net反汇编和GCC源代码,但似乎在任何地方都找不到sin()和其他数学函数的实际实现…他们似乎总是在引用别的东西。

有人能帮我找到他们吗?我觉得不太可能运行C语言的所有硬件都支持三角函数,所以一定有软件算法,对吧?


我知道有几种方法可以计算函数,并且为了好玩,我写了自己的例程来使用泰勒级数计算函数。我很好奇生产语言是如何做到这一点的,因为我的所有实现总是慢几个数量级,尽管我认为我的算法非常聪明(显然不是)。


当前回答

如果您想查看这些函数在C语言中的实际GNU实现,请查看glibc的最新主干。参见GNU C库。

其他回答

在GNU libm中,sin的实现依赖于系统。因此,您可以在sysdeps的适当子目录中找到每个平台的实现。

一个目录包含一个由IBM贡献的C语言实现。自2011年10月以来,这是在典型的x86-64 Linux系统上调用sin()时实际运行的代码。它显然比汇编指令中的f_f快。源代码:sysdeps/ieee754/dbl-64/s_sin.c,查找__sin (double x)。

这段代码非常复杂。没有一种软件算法在整个x值范围内尽可能快且准确,因此库实现了几种不同的算法,它的第一项工作是查看x并决定使用哪种算法。

When x is very very close to 0, sin(x) == x is the right answer. A bit further out, sin(x) uses the familiar Taylor series. However, this is only accurate near 0, so... When the angle is more than about 7°, a different algorithm is used, computing Taylor-series approximations for both sin(x) and cos(x), then using values from a precomputed table to refine the approximation. When |x| > 2, none of the above algorithms would work, so the code starts by computing some value closer to 0 that can be fed to sin or cos instead. There's yet another branch to deal with x being a NaN or infinity.

这段代码使用了一些我以前从未见过的数值技巧,尽管据我所知,它们可能在浮点专家中很有名。有时几行代码需要几段文字来解释。例如,这两条线

double t = (x * hpinv + toint);
double xn = t - toint;

(有时)用于将x减小到接近0的值,该值与x相差π/2的倍数,特别是xn × π/2。这种没有划分或分支的方式相当聪明。但是没有任何评论!


旧的32位版本的GCC/glibc使用fsin指令,这对于某些输入是非常不准确的。有一篇精彩的博客文章用两行代码说明了这一点。

fdlibm在纯C中实现sin要比glibc简单得多,而且注释很好。源代码:fdlibm/s_sin.c和fdlibm/k_sin.c

它们通常在软件中实现,在大多数情况下不会使用相应的硬件(即汇编)调用。然而,正如Jason所指出的,这些是特定于实现的。

请注意,这些软件例程不是编译器源代码的一部分,而是可以在相应的库中找到,例如clib或GNU编译器的glibc。看到http://www.gnu.org/software/libc/manual/html_mono/libc.html三角函数

如果你想要更大的控制权,你应该仔细评估你到底需要什么。一些典型的方法是查找表的插值、程序集调用(通常很慢)或其他近似方案,如Newton-Raphson的平方根。

计算正弦/余弦/正切其实很容易通过代码使用泰勒级数来实现。自己写一个只需5秒钟。

整个过程可以用这个方程来概括:

下面是我为C语言写的一些例程:

double _pow(double a, double b) {
    double c = 1;
    for (int i=0; i<b; i++)
        c *= a;
    return c;
}

double _fact(double x) {
    double ret = 1;
    for (int i=1; i<=x; i++) 
        ret *= i;
    return ret;
}

double _sin(double x) {
    double y = x;
    double s = -1;
    for (int i=3; i<=100; i+=2) {
        y+=s*(_pow(x,i)/_fact(i));
        s *= -1;
    }  
    return y;
}
double _cos(double x) {
    double y = 1;
    double s = -1;
    for (int i=2; i<=100; i+=2) {
        y+=s*(_pow(x,i)/_fact(i));
        s *= -1;
    }  
    return y;
}
double _tan(double x) {
     return (_sin(x)/_cos(x));  
}

我将尝试在一个C程序中回答sin()的情况,该程序用GCC的C编译器在当前的x86处理器(假设是Intel Core 2 Duo)上编译。

在C语言中,标准C库包含了一些常见的数学函数,而这些函数并不包含在语言本身中(例如pow, sin和cos分别表示幂,sin和cos)。它们的头文件包含在math.h中。

现在在GNU/Linux系统上,这些库函数是由glibc (GNU libc或GNU C库)提供的。但是GCC编译器希望您使用-lm编译器标志链接到数学库(libm.so),以启用这些数学函数的使用。我不确定为什么它不是标准C库的一部分。这些将是浮点函数的软件版本,或“软浮动”。

题外话:将数学函数分开的原因由来已久,据我所知,可能是在共享库可用之前,它仅仅是为了在非常古老的Unix系统中减少可执行程序的大小。

Now the compiler may optimize the standard C library function sin() (provided by libm.so) to be replaced with an call to a native instruction to your CPU/FPU's built-in sin() function, which exists as an FPU instruction (FSIN for x86/x87) on newer processors like the Core 2 series (this is correct pretty much as far back as the i486DX). This would depend on optimization flags passed to the gcc compiler. If the compiler was told to write code that would execute on any i386 or newer processor, it would not make such an optimization. The -mcpu=486 flag would inform the compiler that it was safe to make such an optimization.

现在,如果程序执行sin()函数的软件版本,它将基于CORDIC(坐标旋转数字计算机)或BKM算法,或者更可能是现在通常用于计算此类超越函数的表格或幂级数计算。(Src: http://en.wikipedia.org/wiki/Cordic应用程序)

任何最新的gcc版本(大约2.9倍以来)也提供了内置的sin版本__builtin_sin(),作为优化,它将用于取代对C库版本的标准调用。

我相信这是非常清楚的,但希望给你更多的信息比你期望的,和许多出发点,以了解更多自己。

是的,也有计算罪恶的软件算法。基本上,用数字计算机计算这些东西通常是用数值方法来完成的,比如近似表示函数的泰勒级数。

数值方法可以将函数近似到任意精度,因为浮点数的精度是有限的,所以它们非常适合这些任务。