在cmath中sqrt、sin、cos、pow等的定义

Question

是否有如sqrt(), sin(), cos(), tan(), log()和exp()（这些函数从math.h/cmath库中获取）的函数定义?

我只是想知道它们是如何工作的。

- Pythagoras of Samos

---

fdlibm提供了所有这些东西的实现，是开源的、独立的，相当易读。它们不是最简单的实现，因为它们旨在提供良好的性能。 - Steve Jessop

可能是如何计算sin()和其他数学函数的C代码？的重复问题。 - Jason Orendorff

Answer 1

这是一个有趣的问题，但是阅读高效库的源代码并不足以让你了解使用的方法。

以下是一些提示，帮助你理解经典方法。我的信息并不十分准确。以下方法仅为经典方法，具体实现可能使用其他方法。

• 查找表经常被使用。
• 三角函数通常通过CORDIC算法（在CPU或库中）实现。注意通常同时计算正弦和余弦，我一直好奇为什么标准C库没有提供sincos函数。
• 平方根使用牛顿法以及一些聪明的实现技巧：你可以在网上找到quake源代码的摘录，其中包含了令人惊异的1 / sqrt(x)的实现。
• 指数和对数使用exp(2^n x) = exp(x)^(2^n)和log2(2^n x) = n + log2(x)使参数接近零（取对数时接近1），并使用有理函数逼近（通常Padé逼近）。请注意，完全相同的技巧可以让你得到矩阵指数和对数。根据@Stephen Canon的说法，现代实现更倾向于使用泰勒展开而不是有理函数逼近，因为除法比乘法慢得多。
• 其他函数可以从这些函数中推导出来。具体实现可能提供专门的例程。
• pow(x, y) = exp(y * log(x))，所以当y为整数时不应使用pow
• hypot(x, y) = abs(x) sqrt(1 + (y/x)^2)（如果x > y，则另外处理hypot(y, x)）以避免溢出。atan2通过调用sincos和一些小的逻辑计算。这些函数是复杂算术的构建块。

对于其他的超越函数（gamma、erf、bessel等），请参阅优秀的书籍Numerical Recipes, 3rd edition获得一些想法。而经典的Abramowitz & Stegun也是有用的。在http://dlmf.nist.gov/上有一个新版本。

像Chebyshev逼近、连分数展开（实际上与Padé逼近相关）或幂级数经济化这样的技术在更复杂的函数中使用（如果您恰好读取erf、bessel或gamma的源代码，例如）。我怀疑它们在裸机简单数学函数中没有真正的用途，但谁知道呢。请参考Numerical Recipes进行概述。

- Alexandre C.

---

+1 是因为解释了数学。当我意识到三角函数只是被截断的泰勒级数展开时，我感觉好多了。否则这些逼近就像是真正的魔法！ - Ben Jackson

Ben说：好的库通常不使用截断泰勒级数，其他多项式逼近方法（Minimax、Chebyshev、Padé）具有更理想的误差特性，可以用更少的算术运算实现相同的精度。 - Stephen Canon

@Alexandre：对于exp和log，有理逼近法大多已经过时了，因为在现代硬件上，乘法和加法比除法快得多。然而，它仍然用于更严格的函数。 - Stephen Canon

@Janus "The" Abramowitz是一个经典人物...我几乎为在SO上发现它的创始人而流泪。 - Dr. belisarius

@Janus：NIST dlmf 引用 Abramowitz & Stegun 的数字表格或数值系数，这些东西有时非常关键（或者例如您懒得设计一个好的归纳计算方案）。 - Alexandre C.

Answer 2

每个实现都可能是不同的，但你可以查看GNU C库（glibc）源代码中的一个实现。

编辑：Google Code Search已经下线，所以我之前的链接已经失效。

glibc的数学库源代码位于此处：

http://sourceware.org/git/?p=glibc.git;a=tree;f=math;h=3d5233a292f12cd9e9b9c67c3a114c64564d72ab;hb=HEAD

- wkl

Answer 3

看看glibc如何实现各种数学函数，充满了神奇、近似和汇编。

- ismail

---

+1 给 glibc sourceware 链接，但哇，该网站现在非常缓慢。 - wkl

据我所知，这些是较慢的版本，而在特定架构的子文件夹中有更快的实现。 - ismail

LOL，这里的互联网速度总是很慢，看不出任何区别 ;) - ismail

Answer 4

一定要查看fdlibm源代码。它们很好，因为fdlibm库是自包含的，每个函数都有详细的数学解释，并且代码非常易于阅读。

- Daniel Trebbien

Answer 5

经过对数学代码的大量研究，我建议不要去看glibc的代码，因为这些代码通常很难跟踪，并且在很大程度上依赖于glibc的魔法。如果有需要，可以查看FreeBSD中的数学库，它更容易阅读，但速度可能会稍慢（但不会慢太多）。

对于复杂函数，主要困难在于边界情况——正确处理nan/inf/0对实函数已经很困难了，但对于复杂函数来说就更加困难了。C99标准定义了许多边角情况，有些函数甚至有10-20个边角情况。您可以查看最新的C99标准文档的附录G以获取一个概念。此外，还有一个问题是长双精度浮点数的格式未经标准化——根据我的经验，您应该预计会出现相当多的长双精度浮点数错误。希望即将推出的IEEE754修订版扩展精度会改善情况。

- David Cournapeau

---

关于边角案例的观点很好。它们在某些情况下很容易成为瓶颈（例如，在MSVC上实现ldexp函数使该函数几乎无用）。 - Alexandre C.

Answer 6

现代大多数硬件都包括浮点单元，可以非常高效地实现这些函数。

- David Heffernan

Answer 7

用法: root(number,root,depth)

示例: root(16,2) == sqrt(16) == 4
示例: root(16,2,2) == sqrt(sqrt(16)) == 2
示例: root(64,3) == 4

C#实现:

double root(double number, double root, double depth = 1f)
{
    return Math.Pow(number, Math.Pow(root, -depth));
}

用法: Sqrt(数字,深度)

示例: Sqrt(16) == 4
示例: Sqrt(8,2) == sqrt(sqrt(8))

double Sqrt(double number, double depth = 1) return root(number,2,depth);

作者：Imk0tter

- Imk0tter

---

欢迎来到SO，请参考此链接以更好地展示您的问题/答案：https://meta.stackoverflow.com/a/251362/3519504 - Sandeep Kumar