python小练习（047）：精确计算圆周率Pi第n位

jerryxjr1220

python小练习（046）：消失的硬盘容量，传送门

精确计算圆周率Pi几百年来是许多数学家一直想攻克的一道难题，现代我们有了计算机的帮助，计算精确计算圆周率Pi小数点后几万位已经不是什么难事了。

今天的小练习就来和大家分享，如何利用python实现精确计算圆周率Pi的第n位。

提示：
利用的公式：
1. 马青公式：（来源百度百科）
马青公式由英国天文学教授约翰·马青于1706年发现，他利用这个公式计算到了100位的圆周率。

                               
登录/注册后可看大图

2. arctan（x）展开式：
arctanx=x-(x^3)/3+(x^5)/5-(x^7)/7+...+[(-1)^(n+1)]*[x^(2n-1)]/(2n-1)

大家先动动手，试试看吧！

jerryxjr1220

把上述2个公式代入，得到：
Pi = (16/5-4/239) - [(16/3) * (1/5)**3 - (4/3) * (1/239)**3] + ... + (-1)**(n+1)*[(16/(2*n-1)) * (1/5)**(2*n-1) - (4/(2*n-1)) * (1/239)**(2*n-1)]
然后把通项公式写成迭代函数，就可以迭代计算了。

jerryxjr1220

python 2的代码：

1. def calcPi(n):
   
2.     w = 10**(n+10)
   
3.     x = int(16 * w / 5 - 4 * w / 239)
   
4.     for i in range(4,2*n,2):
   
5.         x1 = 16 * w / (i-1)
   
6.         x1 /= 5**(i-1)
   
7.         x2 = 4 * w / (i-1)
   
8.         x2 /= 239**(i-1)
   
9.         x += ((-1)**(i/2+1))*(x1-x2)
   
10.     return x // 10**10
    
11. 
    
12. print calcPi(100)

复制代码

输出Pi的前100位：
31415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679

jerryxjr1220

另外一个方法，运用python的精度设置的方法，可以直接取得高精度的计算结果，但是相对运算速度比较慢。
但是对于1000位以下的精度来说，还是没问题的。

1. # -*- coding: utf-8 -*-
   
2. """
   
3. Spyder Editor
   
4. 
   
5. This is a temporary script file.
   
6. """
   
7. import decimal
   
8. 
   
9. def calcPi(n):
   
10.     decimal.getcontext().prec = n
    
11. 
    
12.     x = 16 / decimal.Decimal(5) - 4 / decimal.Decimal(239)
    
13.     for i in range(4,2*n,2):
    
14.         x1 = 16 / decimal.Decimal(i-1)
    
15.         x1 /= 5**(i-1)
    
16.         x2 = 4 / decimal.Decimal(i-1)
    
17.         x2 /= 239**(i-1)
    
18.         x += ((-1)**(i/2+1))*(x1-x2)
    
19.     return x
    
20. 
    
21. print calcPi(1000)

复制代码

直接输出精确值：
3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513282306647093844609550582231725359408128481117450284102701938521105559644622948954930381964428810975665933446128475648233786783165271201909145648566923460348610454326648213393607260249141273724587006606315588174881520920962829254091715364367892590360011330530548820466521384146951941511609433057270365759591953092186117381932611793105118548074462379962749567351885752724891227938183011949129833673362440656643086021394946395224737190702179860943702770539217176293176752384674818467669405132000568127145263560827785771342757789609173637178721468440901224953430146549585371050792279689258923542019956112129021960864034418159813629774771309960518707211349999998372978049951059731732816096318595024459455346908302642522308253344685035261931188171010003137838752886587533208381420617177669147303598253490428755468731159562863882353787593751957781857780532171226806613001927876611195909216420205

jerryxjr1220

网上的另一种写法，和本题解法思路一致

1. #!/usr/bin/env python
   
2. # -*- coding: cp936 -*- #包含中文注释标记
   
3. n = int(raw_input('请键入想要计算到小数点后的位数n:')) #先键入字符串,再转化为整数
   
4. w = n+10 #多计算10位，防止尾数取舍的影响
   
5. b = 10**w #算到小数点后w位
   
6. x1 = b*4//5 #求含4/5的首项
   
7. x2 = b// -239 #求含1/239的首项
   
8. he = x1+x2 #求第一大项
   
9. n *= 2 #设置下面循环的终点，即共计算n项
   
10. for i in xrange(3,n,2): #循环初值=3，末值2n,步长=2
    
11.         x1 //= -25 #求每个含1/5的项及符号
    
12.         x2 //= -57121 #求每个含1/239的项及符号
    
13.         x = (x1+x2) // i #求两项之和
    
14.         he += x #求总和
    
15. pai = he*4 #求出π
    
16. pai //= 10**10 #舍掉后十位
    
17. print pai #输出圆周率π的值，

复制代码

liumaoshen

请问，当我在3.4.1的python运行时报错了

OverflowError: int too large to convert to float
是怎么回事，谢谢！

jerryxjr1220

liumaoshen 发表于 2017-3-2 21:29
请问，当我在3.4.1的python运行时报错了

OverflowError: int too large to convert to float

整形数太大导致转换浮点数时溢出了

liumaoshen

jerryxjr1220 发表于 2017-3-2 21:31
整形数太大导致转换浮点数时溢出了

谢谢你的回答，我在Python2.7里运行了没问题。请问在Python3.4里咋解决这个问题？

jerryxjr1220

liumaoshen 发表于 2017-3-7 15:42
谢谢你的回答，我在Python2.7里运行了没问题。请问在Python3.4里咋解决这个问题？

已经在这篇帖子里回复了，2种解决方案。
http://bbs.fishc.com/thread-83736-1-2.html