圆周率公式怎么写( 四 ) _圆周率

. (2) 很容易想到，要得到超高精度的 PI 值，实数在计算机中必须以数组的形式进行存取，数组的大小跟所需的有效位数成正比。在这个算法中， PI 的有效位数 n 随（2）的求和项数线性增加。
而为计算（2）中的每一项，需要进行超高精度实数除以小整数（52, 2392, 2k+1）的循环，循环所需次数也跟 n 成正比。所以，这个算法总的时间复杂度为 O(n2) 。
这个算法的优点是简单，而且只需要进行整数运算。下面给出我写的算 PI 程序。
在程序中，我采用了一些提高速度的措施：超高精度实数以数组的形式进行存取，数组元素的类型为 64 位整数（long long），每个元素储存 12 个十进制位；对 xk (x = 1/5, 1/239) 的头部和尾部的 0 的数量进行估计，只对非 0 的部分进行计算。pi.cpp C++ 源程序，在 Linux 下以 g++ pi.cpp -o pi -O2 编译 pi.s 在 g++ 生成的汇编程序的基础上进行修改，速度更快，在 Linux 下以 g++ pi.s -o pi 编译另外，还有许多跟（1）类似的式子，但不常用。
例如： PI/4 = arctan(1/2) + arctan(1/3) PI/4 = 8 arctan(1/10) - arctan(1/239) - 4 arctan(1/515) 第二类算法：与 1/PI 有关的级数 1/PI = (sqrt(8) / 9801) sumk=0~inf { [(4k)! (1103 + 26390k)] / [(k!)4 3964k] } (Ramanujan) 1/PI = (sqrt(10005) / 4270934400) sumk=0~inf { [(6k)! (13591409 + 545140134k)] / [(k!)3 (3k)! (-640320)3k] } (Chudnovsky) 以上两个级数（还有其它类似形式的级数，但不常用）比起 arctan 的泰勒级数要复杂得多。虽然仍然是线性收敛，总的时间复杂度也仍然是 O(n2) ，但它们的收敛速度相当快，（Ramanujan）每项可以增加 8 位有效数字，（Chudnovsky）每项可以增加 14 位。
在这个算法中，除了要进行超高精度实数（数组形式）和小整数的运算外，还有一次超高精度实数的开方和倒数的运算，这需要用到 FFT（快速傅立叶变换），在下文叙述。第三类算法：算术几何平均值和迭代法算术几何平均值（Arithmetic-Geometric Mean, AGM） M(a, b) 定义如下： a0 = a, b0 = b ak = (ak-1 + bk-1) / 2, bk = sqrt(ak-1 bk-1) M(a, b) = limk->inf ak = limk->inf bk 然后，由椭圆积分的一系列理论（抱歉，过程我不懂）可以推导出如下公式： a0 = 1, b0 = 1 / sqrt(2) 1/PI = { 1 - sumk=0~inf [2k (ak2 - bk2)] } / 2M(a0, b0)2 (AGM) 根据这条公式可以制定适当的迭代算法。
在迭代过程中，有效位数随迭代次数按 2 的指数增加，即每迭代一次有效位数乘 2 。算法中的超高精度实数的乘、除、开方等运算需要使用 FFT ，在下文叙述。
综合考虑 FFT 的时间复杂度，整个算法的时间复杂度约为 O(n log(n)2）。除了（AGM）以外，还有其它的迭代序列，它们具有同样的时间复杂度。
例如下面的这个序列将按 4 的指数收敛到 1/PI: y0 = sqrt(2) - 1, a0 = 6 - 4 sqrt(2) yk = [1 - sqrt(sqrt(1 - yk-14))] / [1 + sqrt(sqrt(1 - yk-14))], ak = (1 + yk)4 ak-1 - 22k+1 yk (1 + yk + yk2) 1/PI = limk->inf ak (Borwein) FFT 如上所述，第二和第三类算法不可避免地要涉及超高精度实数（数组形式存取的多位数）的乘、除、开方等运算。多位数乘法如果按照常规方法来计算，逐位相乘然后相加，其时间复杂度将达到 O(n2) 。
使用 FFT 可大大减少计算量。设有复数数组 a[k] 和 b[k] (k=0~n-1) ，正向和反向的离散傅立叶变换（DFT）定义如下：（i = sqrt(-1)) b = FFTforward(a) : b[k] = sumj=0~n-1 ( a[j] e-i*j*2PI*k/n ) (3) b = FFTbackward(a) : b[k] = (1/n) sumj=0~n-1 ( a[j] ei*j*2PI*k/n ) (4) (3) 和（4）中的（1/n）可以放在任何一个式子中，也可以拆成（1/sqrt(n)）同时放在两个式子中，目的是保证正向和反向傅立叶变换以后不会相差一个因子。