Skip to content

Instantly share code, notes, and snippets.

@komasaru

komasaru/calc_pi_arctan.py

Last active March 9, 2018 11:03
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save komasaru/290022bcc86f380d771e687ddc3ea5f7 to your computer and use it in GitHub Desktop.
Save komasaru/290022bcc86f380d771e687ddc3ea5f7 to your computer and use it in GitHub Desktop.
Download ZIP
Python script to compute Pi with arctan's formula.
Raw
calc_pi_arctan.py
#! /usr/local/bin/python3.6
"""
Computation of Pi with arctan's formula
Formulas: 1: Machin
2: Klingenstierna
3: Euler
4: Euler(2)
5: Gauss
6: Stormer
7: Stormer(2)
8: Takano
"""
import math
import sys
import time
class CalcPiArctan:
FORMULA = [
"Machin", "Klingenstierna", "Euler", "Euler2",
"Gauss", "Stormer", "Stormer2", "Takano"
] # Formula names
KK = [ # Coefficients
[ 16, 1, 5, -4, 1, 239 ], # 1: Machin
[ 32, 1, 10, -4, 1, 239, -16, 1, 515 ], # 2: Klingenstierna
[ 20, 1, 7, 8, 3, 79 ], # 3: Euler
[ 16, 1, 5, -4, 1, 70, 4, 1, 99 ], # 4: Euler(2)
[ 48, 1, 18, 32, 1, 57, -20, 1, 239 ], # 5: Gauss
[ 24, 1, 8, 8, 1, 57, 4, 1, 239 ], # 6: Stormer
[176, 1, 57, 28, 1, 239, -48, 1, 682, 96, 1, 12943], # 7: Stormer(2)
[ 48, 1, 49, 128, 1, 57, -20, 1, 239, 48, 1, 110443] # 8: Takano
]
# File name
FILE_PRE = "PI_" # Prefix
FILE_EXT = ".txt" # Extension
# Output strings
STR_TITLE = "** Pi Computation by Arctan method **"
STR_FORMULA = " Formula = [ {:s} ]"
STR_DIGITS = " Digits = {:d}"
STR_TIME = " Time = {:f} seconds"
# Multi-digit computation
MAX_DIGITS = 100000000 # 8-digits per 1 item
def __init__(self, x, y):
self.fname = self.FILE_PRE + self.FORMULA[x - 1] + self.FILE_EXT
print("\n[ Output file :", self.fname, "]\n")
self.formula = self.FORMULA[x - 1]
self.ks = int(len(self.KK[x-1]) / 3)
self.kk = self.KK[x - 1]
self.l = y
self.l1 = int(self.l / 8) + 1
self.n = []
for i in range(0, self.ks):
n = self.l / (math.log10(self.kk[i * 3 + 2]) \
- math.log10(self.kk[i * 3 + 1])) + 1
self.n.append(int(n / 2) + 1)
def compute(self):
""" Computation and display """
try:
t0 = time.time()
a = [[0 for _ in range(self.l1 + 2)] for _ in range(self.ks)]
q = [[0 for _ in range(self.l1 + 2)] for _ in range(self.ks)]
sk = [[0 for _ in range(self.l1 + 2)] for _ in range(self.ks)]
s = [0 for _ in range(self.l1 + 2)]
for i in range(self.ks):
a[i][0] = self.kk[i * 3] * self.kk[i * 3 + 2]
if self.kk[i * 3] < 0:
a[i][0] *= -1
if self.kk[i * 3 + 1] != 1:
a[i] = self.__long_div(a[i], self.kk[i * 3 + 1])
for i in range(self.ks):
for k in range(1, self.n[i] + 1):
a[i] = self.__long_div(a[i], self.kk[i * 3 + 2])
a[i] = self.__long_div(a[i], self.kk[i * 3 + 2])
if self.kk[i * 3 + 1] != 1:
a[i] = self.__long_mul(a[i], self.kk[i * 3 + 1])
a[i] = self.__long_mul(a[i], self.kk[i * 3 + 1])
q[i] = self.__long_div(a[i], 2 * k - 1)
if k % 2 == 0:
sk[i] = self.__long_sub(sk[i], q[i])
else:
sk[i] = self.__long_add(sk[i], q[i])
for i in range(self.ks):
if self.kk[i * 3] >= 0:
s = self.__long_add(s, sk[i])
else:
s = self.__long_sub(s, sk[i])
t1 = time.time()
tt = t1 - t0
self.__display(tt, s)
except Exception as e:
raise
def __long_add(self, a, b):
""" Computation of long + long
:param list a
:param list b
:return list z
"""
try:
z = [0 for _ in range(self.l1 + 2)]
cr = 0
for i in reversed(range(self.l1 + 2)):
z[i] = a[i] + b[i] + cr
if z[i] < self.MAX_DIGITS:
cr = 0
else:
z[i] -= self.MAX_DIGITS
cr = 1
return z
except Exception as e:
raise
def __long_sub(self, a, b):
""" Computation of long - long
:param list a
:param list b
:return list z
"""
try:
z = [0 for _ in range(self.l1 + 2)]
br = 0
for i in reversed(range(self.l1 + 2)):
z[i] = a[i] - b[i] - br
if z[i] >= 0:
br = 0
else:
z[i] += self.MAX_DIGITS
br = 1
return z
except Exception as e:
raise
def __long_mul(self, a, b):
""" Computation of long * short
:param list a
:param list b
:return list z
"""
try:
z = [0 for _ in range(self.l1 + 2)]
cr = 0
for i in reversed(range(self.l1 + 2)):
w = a[i]
z[i] = (w * b + cr) % self.MAX_DIGITS
cr = int((w * b + cr) / self.MAX_DIGITS)
return z
except Exception as e:
raise
def __long_div(self, a, b):
""" Computation of long / short
:param list a
:param list b
:return list z
"""
try:
z = [0 for _ in range(self.l1 + 2)]
r = 0
for i in range(self.l1 + 2):
w = a[i]
z[i] = int((w + r) / b)
r = ((w + r) % b) * self.MAX_DIGITS
return z
except Exception as e:
raise
def __display(self, tt, s):
""" Display
:param float tt: elapsed time
:param list s: result of calculation
"""
try:
print(self.STR_TITLE)
print(self.STR_FORMULA.format(self.formula))
print(self.STR_DIGITS.format(self.l))
print(self.STR_TIME.format(tt))
out_file = open(self.fname, "w")
out_file.write(self.STR_TITLE)
out_file.write("\n")
out_file.write(self.STR_FORMULA.format(self.formula))
out_file.write("\n")
out_file.write(self.STR_DIGITS.format(self.l))
out_file.write("\n")
out_file.write(self.STR_TIME.format(tt))
out_file.write("\n\n {:d}.\n".format(s[0]))
for i in range(1, self.l1):
if i % 10 == 1:
out_file.write("{:08d}:".format((i - 1) * 8 + 1))
out_file.write(" {:08d}".format(s[i]))
if i % 10 == 0:
out_file.write("\n")
out_file.close
except Exception as e:
raise
if __name__ == '__main__':
try:
print("1:Machin, 2:Klingenstierna, 3:Euler, 4:Euler(2)")
print("5:Gauss, 6:Stormer, 7:Stormer(2), 8:Takano : ", end="")
f = int(input())
if f < 1 or f > 8:
f = 1
print("Please input number of Pi Decimal-Digits : ", end="")
n = int(input())
obj = CalcPiArctan(f, n)
obj.compute()
except Exception as e:
print("EXCEPTION!", e.args, file=sys.stderr)
sys.exit(1)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment