Introduction

Archimedes calculated the approximate values of the circumference ratio by calculating the circumferences of a regular hexagon, a regular dodecagon, a regular 24 polygon, a regular 48 square, and a regular 96 polygon circumscribing a circle in that order. tsujimotter's notebook "Archimedean and Pi" At that time, the decimal point had not been invented yet, and it was expressed as a ratio of integers. Even with these restrictions, Archimedes wanted to be less than 3 and 1/7 and greater than 3 and 10/71. (In decimal terms, between 3.1428571428571429 and 3.1408450704225352) With the help of programming, I tried to find out how accurately the pi could be calculated by proceeding with Archimedes' method as it is.

Archimedes method

The Archimedes method itself is accurately described in tujimotter's notebook, so if you are interested, please have a look there. Here, only the proof necessary for creating the program is described. If the circumference of the regular n-sided polygon is L and the radius of the circumscribed circle is 1, the circumference is 2π, so the approximate value of the circumference ratio based on the regular n-sided polygon is L / 2. 証明.JPG

[For regular hexagons]

L＝F(0)B×2×6 F (0) B = 1 / √3 * ∠F (0) OB is 30 ° (△ OF (0) C is an equilateral triangle) The approximate value of pi based on a regular hexagon is L ÷ 2 ＝ 1 / √3 × 2 × 6 ÷ 2 ＝ 2 × √3 ≒ 3.464.

[For regular dodecagons]

L＝F(1)B×2×12 F(1)B ＝ OB× F(0)B ÷（OB＋F(0)O）= F(0)B ÷（1＋√(1＋F(0)B× F(0)B)) The approximate value of pi based on a regular dodecagon is L ÷ 2 ＝ F (1) B × 2 × 12 ÷ 2 ≒ 3.215. ** Proof ** Draw a straight line passing through F (0) parallel to F (1) O, and let A be the intersection with the extension line of OB. Since ΔOBF (1) and ΔABF (0) are similar, the following holds. F(1)B：OB = F(0)B：（OB＋OA） Furthermore, since ΔOF (0) A is an isosceles triangle, F (0) O = OA. F (0) O is obtained from √ (OB × OB ＋ F (0) B × F (0) B) using Pythagorean theorem.

[For regular n-sided polygons]

Similarly, the approximate value of the pi based on the regular n-sided polygon is as follows. L÷2＝F(n)B×2×n÷2 ＝ F(n-1)B÷(1＋√(1＋F(n-1)B× F(n-1)B))×ｎ

Implementation by Python

The above result is implemented in Python as follows.

import math
import itertools
def f(n):
  if n == 0: return math.sqrt(3)/3
  return f(n-1)/(1+math.sqrt(1+f(n-1)**2))

N = int(input())
for n in range(N):
  print(f(n)*6*2**n)

However, since it is the Mendou that accurately inputs a regular n-sided polygon, the input value N is as follows. 0… Calculate a regular hexagon 1… Calculate up to a regular dodecagon 2… Calculate up to a regular 24 square The same applies below. Since the above algorithm requires the amount of calculation of O (2 ** n), I rewrote it with DP as follows. In addition, since it is not possible to accurately represent the decimal point with a binary number (*), we made the number about 15 digits larger and made it a program to format it when outputting.

Isn't the decimal number dangerous in programming? The story.

import math
import itertools

def pi2s(pi):
  return str(pi)[0] + "." + str(pi)[1:]

def regN(n):
  return format(n, '10d') + ": "

R = 1732050807568877
B = 1000000000000000
pi = "3.14159265358979323846264338327950288"

N = int(input())
F = [0.0 for n in range(N)]
F[0] = B
for n in range(1, N):
  F[n] = F[n-1]*R/(R+math.sqrt(F[n-1]**2+R**2))

for n, v in enumerate(F):
  print(regN(6*2**n) + pi2s(int(v*6*2**n*B/R)))
print(format("pi: ", '>12') + pi[:17])

It is the result of calculation up to ** regular 3145728 square **. 実行結果.JPG It outputs π on the last line. You can see that the values are fairly accurate.

at the end

I tried to trace the method of Archimedes by programming. Contrary to my own expectations, I was surprised that it was already close to 3.14 at a fairly early stage (48-sided). Archimedes, who calculated this by hand and by the ratio of integers (fractions), just took off his hat.

[PYTHON] About the accuracy of Archimedean circle calculation method