Python 中的 BBP 算法 - 使用任意精度算术
BBP-Algorithm in Python - Working with Arbitrary-precision arithmetic
我正在写一篇关于圆周率计算的学期论文。虽然我已经完成了理论站点,但我现在正在努力在 Python.
中实施 BBP 算法
您可以在此处找到 BBP 算法:
http://en.wikipedia.org/wiki/Bailey%E2%80%93Borwein%E2%80%93Plouffe_formula
这是我在 Python 中的实现:
from sympy.mpmath import *
pi = mpf(0)
mp.dps = 30;
for n in range(0,500):
pi = pi + mpf((1/16**n)*(4/(8*n+1)- 2/(8*n+4)- 1/(8*n+5)- 1/(8*n+6)) )
print(pi)
我的问题是,无论我将 k 设置多高或将 pi 的小数位设置多高,我都无法获得比 16 位更高的精度。
因为之前遇到了一些问题,所以我使用了mpmath来获得更高的精度。
如何改进我的代码,以便获得更多数字?
默认情况下,python 将使用 IEEE-754 定义的标准浮点数。这具有大约 12 位的精度,可以表示低至 2-1022 的数字,现在您可以通过调用 mpf
运算符在这个过程的早期,因此一个更具可读性和更精确的版本是:
from sympy.mpmath import *
pi = mpf(0)
mp.dps = 30;
for n in range(0,500):
u = 4.0/(8*n+1)-2.0/(8*n+4)-1.0/(8*n+5)-1.0/(8*n+6)
u = mpf(u)
d = mpf(16.0/1.0)**n
pi += u/d
print(pi)
然而,当您计算 u
部分时,这仍然存在问题。要做到这一点,你可以使用:
from sympy.mpmath import *
pi = mpf(0)
mp.dps = 50;
for n in range(0,50000):
u = mpf(4)/mpf(8*n+1)-mpf(2)/mpf(8*n+4)-mpf(1)/mpf(8*n+5)-mpf(1)/mpf(8*n+6)
d = mpf(16.0)**n
pi += u/d
print(pi)
前 50 位数字正确的是:
3.1415926535 8979323846 2643383279 502884197 1693993751
(已添加空格)
我正在写一篇关于圆周率计算的学期论文。虽然我已经完成了理论站点,但我现在正在努力在 Python.
中实施 BBP 算法您可以在此处找到 BBP 算法: http://en.wikipedia.org/wiki/Bailey%E2%80%93Borwein%E2%80%93Plouffe_formula
这是我在 Python 中的实现:
from sympy.mpmath import *
pi = mpf(0)
mp.dps = 30;
for n in range(0,500):
pi = pi + mpf((1/16**n)*(4/(8*n+1)- 2/(8*n+4)- 1/(8*n+5)- 1/(8*n+6)) )
print(pi)
我的问题是,无论我将 k 设置多高或将 pi 的小数位设置多高,我都无法获得比 16 位更高的精度。
因为之前遇到了一些问题,所以我使用了mpmath来获得更高的精度。
如何改进我的代码,以便获得更多数字?
默认情况下,python 将使用 IEEE-754 定义的标准浮点数。这具有大约 12 位的精度,可以表示低至 2-1022 的数字,现在您可以通过调用 mpf
运算符在这个过程的早期,因此一个更具可读性和更精确的版本是:
from sympy.mpmath import *
pi = mpf(0)
mp.dps = 30;
for n in range(0,500):
u = 4.0/(8*n+1)-2.0/(8*n+4)-1.0/(8*n+5)-1.0/(8*n+6)
u = mpf(u)
d = mpf(16.0/1.0)**n
pi += u/d
print(pi)
然而,当您计算 u
部分时,这仍然存在问题。要做到这一点,你可以使用:
from sympy.mpmath import *
pi = mpf(0)
mp.dps = 50;
for n in range(0,50000):
u = mpf(4)/mpf(8*n+1)-mpf(2)/mpf(8*n+4)-mpf(1)/mpf(8*n+5)-mpf(1)/mpf(8*n+6)
d = mpf(16.0)**n
pi += u/d
print(pi)
前 50 位数字正确的是:
3.1415926535 8979323846 2643383279 502884197 1693993751
(已添加空格)