圆周率的计算

圆周率的计算

整理三种计算圆周率的方法，体验下数学的美妙和编程的神奇

1.传说中计算π的超牛的c程序3行代码计算π后面800位)1:传说中计算π的超牛的C程序，3行代码，计算π后面2800位…

原c代码如下：

#include "stdio.h"
long a=10000, b, c=2800, d, e, f[2801], g;
void 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) ;
}

表示根本不理解背后的数学推理，看了几个大神的帖子，还是很难理解，后续再继续研读吧。
根据c代码，换成python代码，效率和准确度都很高，如下：

c = 2800
f = np.zeros(2801)
b = 0
for i in range(c):  # 除最后一位，其余全部是2
	f[i] = 2
e = 0
while c > 0:
	d = 0
	b = c
	while b > 0:
		d *= b
		d += (f[b] * 10)
		f[b] = d % (b * 2-1)
		d = int(d/(b * 2 - 1))  # 如果不取整，源代码结果不对
		b -= 1
	c -= 1
	print("%d" % ((e+int(d)/10) % 10), end="")
	e = d % 10

2.割圆术

我们伟大祖先勤劳和智慧的结晶，徒手开根号计算这么复杂的式子，内心深感骄傲与自豪！
上原理图：

计算到后面，精度不是很好。

i = 0
n = 15
an = 1
while i < n:
    an1Pow = 2 - math.sqrt(4-math.pow(an, 2))
    an1 = math.sqrt(an1Pow)
    an = an1
    i += 1
print("由正%d边形求出的PI为%.50f" % (6*math.pow(2, n), (3*math.pow(2, n)*an)))

3.莱布尼兹公式

世界人民伟大的智慧结晶，具体的推导超出了本文的范围，感兴趣可以查找相关资料进行证明。
计算比较占内存，也比较慢，假如计算到1w位，估计要跑了几小时…

# 为了节省空间，定义xrange函数，py3中没有该函数
def xrange(x):
    n = 1
    while n <= x:
        yield n
        n += 2

quarterPi = 0
flag = 1
for i in xrange(100000000):
    quarterPi += (1/i)*flag
    flag = -flag
print(4*quarterPi)

4.参考及整理

mygithub
五种方法来计算圆周率PI
解读传说中计算π的超牛的C程序