π is a truly magical number followed by countless people. I am not very clear about the fascinating point of an irrational number that is always repeated. In my opinion, I am happy to calculate π, which is the value of pi. Because π is an irrational number, it is infinite. This means that any calculation of π is just an approximation. If you calculate 100 bits, I can calculate 101 bits and be more precise. So far, some people have chosen supercomputers to try to calculate the most accurate π. Some extrema include the 500 million bits that calculate π. You can even find a 10 billion-bit text file containing π from the Internet (note!). It may take a while to download this file, and it won't open with the Notepad app you normally use. )。 For me, how to calculate π with a few lines of Python is my interest.
You can always use the Math.PI variable. It is included in the standard library, and you should use it before you try to calculate it yourself. In fact, we will use it to calculate the accuracy. As a starting point, let's look at a very straightforward method of calculating pi. As usual, I will use Python 2.7, and the same ideas and code might apply to different versions. Most of the algorithms we will be using are from the Pi WikiPedia page and implemented. Let's take a look at the following code:
importsys
importmath
defmain(argv):
iflen(argv) !=1:
sys.exit('Usage: calc_pi.py')
print'\nComputing Pi v.01\n'
a=1.0
b=1.0/math.sqrt(2)
t=1.0/4.0
p=1.0
foriinrange(int(sys.argv[1])):
at=(a+b)/2
bt=math.sqrt(a*b)
tt=t-p*(a-at)**2
pt=2*p
a=at;b=bt;t=tt;p=pt
my_pi=(a+b)**2/(4*t)
accuracy=100*(math.pi-my_pi)/my_pi
print"Pi is approximately: "+str(my_pi)
print"Accuracy with math.pi: "+str(accuracy)
if__name__=="__main__":
main(sys.argv[1:])
This is a very simple script that you can download, run, modify, and share freely with others. You can see output similar to the following:
-
You will find that even though N is greater than 4, we are approaching Pi accuracy without much improvement. We can guess that even if the value of n is large, the same thing (pi approximation accuracy does not improve) will still occur. Fortunately, there is more than one way to uncover this mystery. Using the Python decimal (decimal) library, we can get a higher precision value to approximate the pi. Let's look at how the library function is used. With this simplified version, you can get more than 11 digits of the number typically given by the small Python floating point precision. Here is an example from the Python Decimal library:
-
See these numbers. Wrong! We entered only 3.14, why did we get some rubbish (junk)? This is memory garbage (junk). In the nut shell, Python gives you the decimal number you want, plus a little extra value. As long as the precision is less than the start of the garbage number, it does not affect any calculations as long as the precision is less than the previous garbage number (junk numbers). You can specify the number of digits you want by setting GetContext (). Prec. Let's try.
Very good. Now let's try to use this to see if we can get a better approximation of our previous code. Now, I usually object to using "From library import *", but in this case it will make the code look more beautiful.
importsys
importmath
fromdecimalimport*
defmain(argv):
iflen(argv) !=1:
sys.exit('Usage: calc_pi.py')
print'\nComputing Pi v.01\n'
a=Decimal(1.0)
b=Decimal(1.0/math.sqrt(2))
t=Decimal(1.0)/Decimal(4.0)
p=Decimal(1.0)
foriinrange(int(sys.argv[1])):
at=Decimal((a+b)/2)
bt=Decimal(math.sqrt(a*b))
tt=Decimal(t-p*(a-at)**2)
pt=Decimal(2*p)
a=at;b=bt;t=tt;p=pt
my_pi=(a+b)**2/(4*t)
accuracy=100*(Decimal(math.pi)-my_pi)/my_pi
print"Pi is approximately: "+str(my_pi)
print"Accuracy with math.pi: "+str(accuracy)
if__name__=="__main__":
main(sys.argv[1:])
Output Result:
All right. We're more accurate, but it looks like there's some rounding. from n = 100 and n = 1000, we have the same precision. Now what? OK, now let's turn to the formula. So far, the way we calculate pi is by adding a few parts together. I found some code from Dan's article on calculating Pi. He suggested that we use the following 3 formulas:
Bailey–borwein–plouffe formula
The formula of Bellard
Chudnovsky algorithm
Let's start with the Bailey–borwein–plouffe formula. It looks like this:
In the code we can write it this way:
import sys
import math
from decimal import *
def bbp(n):
pi=Decimal(0)
k=0
while k < n:
pi+=(Decimal(1)/(16**k))*((Decimal(4)/(8*k+1))-(Decimal(2)/(8*k+4))-(Decimal(1)/(8*k+5))-(Decimal(1)/(8*k+6)))
k+=1
return pi
def main(argv):
if len(argv) !=2:
sys.exit('Usage: BaileyBorweinPlouffe.py')
getcontext().prec=(int(sys.argv[1]))
my_pi=bbp(int(sys.argv[2]))
accuracy=100*(Decimal(math.pi)-my_pi)/my_pi
print"Pi is approximately "+str(my_pi)
print"Accuracy with math.pi: "+str(accuracy)
if __name__=="__main__":
main(sys.argv[1:])
Throw away the "wrapper" code, the BBP (N) function is what you really want. You give it the greater the N and give the GetContext (). Prec the larger the value, the more accurate the calculation will be. Let's take a look at some of the code results:
This has many digit digits. As you can see, we are not more accurate than before. So we need to go forward to the next formula, Beira the formula, hoping to get better accuracy. It looks like this:
We will only change our transformation formula, and the rest of the code will remain unchanged. Click here to download the Beira formula implemented by Python. Let's take a look at Bellards (n):
def bellard(n):
pi=Decimal(0)
k=0
while k < n:
pi+=(Decimal(-1)**k/(1024**k))*( Decimal(256)/(10*k+1)+Decimal(1)/(10*k+9)-Decimal(64)/(10*k+3)-Decimal(32)/(4*k+1)-Decimal(4)/(10*k+5)-Decimal(4)/(10*k+7)-Decimal(1)/(4*k+3))
k+=1
pi=pi*1/(2**6)
return pi
Output Result:
Oh, no, we got the same precision. Well, let's try the third formula, the Chudnovsky algorithm, which looks like this:
Again, let's take a look at this formula (suppose we have a factorial formula). Click here to download the Chudnovsky formula implemented in Python.
Here are the program and output results:
def chudnovsky(n):
pi=Decimal(0)
k=0
while k < n:
pi+=(Decimal(-1)**k)*(Decimal(factorial(6*k))/((factorial(k)**3)*(factorial(3*k)))*(13591409+545140134*k)/(640320**(3*k)))
k+=1
pi=pi*Decimal(10005).sqrt()/4270934400
pi=pi**(-1)
return pi
So what do we have to conclude? Fancy algorithms do not make the machine floating-point world reach a higher standard. I really look forward to having a better precision than we can get with the summation formula. I guess that's too much to ask. If you really need pi, just use the Math.PI variable. However, as fun and test how fast your computer can really be, you can always try the first one to calculate the PI million dollar or more bits that are several.