什么是最快的方式来获得π的价值？
我正在寻找最快的方式来获得π的价值，作为个人的挑战。 更具体地说，我使用的方法不涉及使用像M_PI这样的#define常量，或者对数字进行硬编码。
下面的程序测试了我所知道的各种方式。 内联汇编版本理论上是最快的选择，尽管显然不是可移植的。 我已经将它作为基准与其他版本进行比较。 在我的测试中，使用内置插件， 4 * atan(1)版本在GCC 4.2上速度最快，因为它将atan(1)自动折叠为常量。 使用-fno-builtin指定， atan2(0, -1)版本是最快的。
这是主要的测试程序（ pitimes.c ）：
#include <math.h>
#include <stdio.h>
#include <time.h>
#define ITERS 10000000
#define TESTWITH(x) {
diff = 0.0;
time1 = clock();
for (i = 0; i < ITERS; ++i)
diff += (x) - M_PI;
time2 = clock();
printf("%st=> %e, time => %fn", #x, diff, diffclock(time2, time1));
}
static inline double
diffclock(clock_t time1, clock_t time0)
{
return (double) (time1 - time0) / CLOCKS_PER_SEC;
}
int
main()
{
int i;
clock_t time1, time2;
double diff;
/* Warmup. The atan2 case catches GCC's atan folding (which would
* optimise the ``4 * atan(1) - M_PI'' to a no-op), if -fno-builtin
* is not used. */
TESTWITH(4 * atan(1))
TESTWITH(4 * atan2(1, 1))
#if defined(__GNUC__) && (defined(__i386__) || defined(__amd64__))
extern double fldpi();
TESTWITH(fldpi())
#endif
/* Actual tests start here. */
TESTWITH(atan2(0, -1))
TESTWITH(acos(-1))
TESTWITH(2 * asin(1))
TESTWITH(4 * atan2(1, 1))
TESTWITH(4 * atan(1))
return 0;
}
以及仅适用于x86和x64系统的内联汇编代码（ fldpi.c ）：
double
fldpi()
{
double pi;
asm("fldpi" : "=t" (pi));
return pi;
}
构建脚本构建我正在测试的所有配置（ build.sh ）：
#!/bin/sh
gcc -O3 -Wall -c -m32 -o fldpi-32.o fldpi.c
gcc -O3 -Wall -c -m64 -o fldpi-64.o fldpi.c
gcc -O3 -Wall -ffast-math -m32 -o pitimes1-32 pitimes.c fldpi-32.o
gcc -O3 -Wall -m32 -o pitimes2-32 pitimes.c fldpi-32.o -lm
gcc -O3 -Wall -fno-builtin -m32 -o pitimes3-32 pitimes.c fldpi-32.o -lm
gcc -O3 -Wall -ffast-math -m64 -o pitimes1-64 pitimes.c fldpi-64.o -lm
gcc -O3 -Wall -m64 -o pitimes2-64 pitimes.c fldpi-64.o -lm
gcc -O3 -Wall -fno-builtin -m64 -o pitimes3-64 pitimes.c fldpi-64.o -lm
除了在各种编译器标志之间进行测试之外（我还将32位与64位进行了比较，因为优化是不同的），我也尝试切换测试的顺序。 但仍然， atan2(0, -1)版本每次仍然处于顶峰。
如上所述，蒙特卡洛方法应用了一些很好的概念，但它显然不是最快的，也不是最长的，也不是通过任何合理的措施。 此外，这一切都取决于你正在寻找什么样的准确性。 我知道的最快的π是硬编码的数字。 看Pi和Pi [PDF]，有很多公式。
这是一种快速收敛的方法 - 每次迭代约14位数字。 PiFast是当前最快的应用程序，它将这个公式与FFT结合使用。 我只写公式，因为代码很简单。 Ramanujan几乎找到了这个公式，并被Chudnovsky发现。 实际上他是如何计算出数字的几十亿位数的 - 所以这不是一种无视的方法。 该公式将很快溢出，并且由于我们正在分解阶乘因子，因此延迟此类计算以删除项目将是有利的。
哪里，
以下是Brent-Salamin算法。 维基百科提到，当a和b “足够接近”时， （a + b）2 / 4t将近似于π。 我不确定什么“足够接近”的意思，但从我的测试，一个迭代得到2位数，两个得到7，三个有15，当然这是双打，所以它可能有一个基于其表示和真正的计算可能会更准确。
let pi_2 iters =
let rec loop_ a b t p i =
if i = 0 then a,b,t,p
else
let a_n = (a +. b) /. 2.0
and b_n = sqrt (a*.b)
and p_n = 2.0 *. p in
let t_n = t -. (p *. (a -. a_n) *. (a -. a_n)) in
loop_ a_n b_n t_n p_n (i - 1)
in
let a,b,t,p = loop_ (1.0) (1.0 /. (sqrt 2.0)) (1.0/.4.0) (1.0) iters in
(a +. b) *. (a +. b) /. (4.0 *. t)
最后，一些pi高尔夫（800位数字）怎么样？ 160个字符！
int a=10000,b,c=2800,d,e,f[2801],g;main(){for(;b-c;)f[b++]=a/5;for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)for(b=c;d+=f[b]*a,f[b]=d%--g,d/=g--,--b;d*=b);}
我真的很喜欢这个程序，它通过查看它自己的区域近似于pi :-)
IOCCC 1988：westley.c
#define _ -F<00||--F-OO--;
int F=00,OO=00;main(){F_OO();printf("%1.3fn",4.*-F/OO/OO);}F_OO()
{
_-_-_-_
_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_
_-_-_-_
}
以下是我在高中学习计算pi的技术的一般描述。
我只是分享这一点，因为我认为这足够简单，任何人都可以无限期地记住它，再加上它会教给你“蒙特卡罗”方法的概念 - 这是统计方法得出的答案并不是立即出现通过随机过程可以推断出来。
绘制一个正方形，并在该正方形内部（半径等于正方形边的象限）内划一个象限（四分之一半圆），以便尽可能多地填充正方形）
现在在广场投掷一个飞镖，并记录它的落地位置 - 也就是说，在广场内的任何地方选择一个随机点。 当然，它落在广场内，但它在半圆内吗？ 记录下这个事实。
重复这个过程很多次 - 你会发现半圈内的点数与抛出的总数的比值，称这个比率为x。
由于平方的面积是r乘以r，所以可以推断出半圆的面积是r乘以r的倍数（即x乘以r的平方）。 因此x次4会给你pi。
这不是一个快速使用的方法。 但这是蒙特卡罗方法的一个很好的例子。 如果你环顾四周，你会发现许多问题，除了你的计算技能，可以通过这种方法来解决。
链接地址: http://www.djcxy.com/p/70979.html