Notebook

これは日々の作業を通して学んだことや毎日の生活で気づいたことをを記録しておく備忘録である。

HTML ファイル生成日時: 2024/11/24 14:07:03.493 (台灣標準時)

円周率を求めるための新しい公式

インド人の研究者によって、円周率を求めることに使える新しい公式が見つけ られたそうでござる。

この公式を使って、円周率の近似値を計算するプログラムを書いてみたでござ る。


#!/usr/pkg/bin/python3.12

#
# Time-stamp: <2024/08/02 13:06:52 (UT+8) daisuke>
#

# importing math module
import math

# a function to calculate approximate value of pi
def new_pi (l, n):
    # initialising value of pi
    pi = 4.0
    # calculating n terms
    for i in range (1, n+1):
        a = 1.0 / (i + l) - 4.0 / (2.0 * i + 1.0)
        b = (2.0 * i + 1.0)**2 / (4.0 * (i + l) ) - i
        c = i - 1.0
        d = math.factorial (i)
        pi += a * math.gamma (b + c) / math.gamma (b) / d
    # returning calculated approximate value of pi
    return (pi)

# calculation of approximate value of pi
pi = new_pi (10.0, 100)

# printing result
print (f'approximate value of pi')
print (f'50 digits:       3.14159265358979323846264338327950288419716939937510')
print (f'new_pi (10,100): {pi}')

実行結果は以下の通りでござる。


% ./pi_2024.py 
approximate value of pi
50 digits:       3.14159265358979323846264338327950288419716939937510
new_pi (10,100): 3.1415926535897927

λ=10 で、 100 項も計算すれば、 15 桁くらいまで合っているようでござる。

fig_202408/pi_2024_code.png
fig_202408/pi_2024_result.png

上の方法だと、倍精度の範囲内でしか円周率の近似値を求められないでござる。 更に、桁数を増やすにはどうすればよいのか考えてみたでござる。 Python に 標準で含まれる Decimal モジュールを使えばよいかと思ったのでござるが、 ガンマ関数の計算をする部分を自分で書かねばならぬので、面倒そうでござっ た。調べてみると、 mpmath というパッ ケージがあることがわかったでござる。このパッケージにはガンマ関数を計算 してくれる関数も用意されているようなので、これを使えば簡便に円周率の近 似値を求められそうでござる。以下のようなコードを準備してみたでござる。


#!/usr/pkg/bin/python3.12

#
# Time-stamp: <2024/08/02 14:50:52 (UT+8) daisuke>
#

# importing mpmath module
import mpmath

# settings for mpmath module
mpmath.mp.dps    = 100
mpmath.mp.pretty = True

# a function to calculate approximate value of pi
def new_pi (l, n):
    # initialising value of pi
    mp_1  = mpmath.mpf (1)
    mp_2  = mpmath.mpf (2)
    mp_4  = mpmath.mpf (4)
    mp_l  = mpmath.mpf (l)
    mp_pi = mp_4
    # calculating n terms
    for i in range (1, n+1):
        mp_i   = mpmath.mpf (i)
        a      = mp_1 / (mp_i + mp_l) - mp_4 / (mp_2 * mp_i + mp_1)
        b      = (mp_2 * mp_i + mp_1)**mp_2 / (mp_4 * (mp_i + mp_l) ) - mp_i
        c      = mp_i - mp_1
        d      = mpmath.factorial (i)
        mp_pi += a * mpmath.gamma (b + c) / mpmath.gamma (b) / d
    # returning calculated approximate value of pi
    return (mp_pi)

# calculation of approximate value of pi
l  = 10
n  = 10000000
pi = new_pi (l, n)

# printing result
print (f'approximate value of pi:')
print (f'100 digits:')
print (f'3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679')
print (f'new_pi ({l},{n}):')
print (f'{pi}')

実行結果は、以下の通りでござる。


% ./pi_2024_100.py
approximate value of pi:
100 digits:
3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679
new_pi (10,10000000):
3.141592653589793238462643383279502884197169399375105820974944592307816406286206269139035079106063948

1000 万項の計算をすると、 78 桁くらいまでいけるようでござる。

fig_202408/pi_2024_100_code.png
fig_202408/pi_2024_100_result.png


Frequently accessed files

  1. Computer___Python/20220518_0.html
  2. Computer___Network/20230726_00.html
  3. Misc___Taiwan/20240207_00.html
  4. Computer___Network/20230516_00.html
  5. Computer___FreeBSD/20220621_0.html
  6. Computer___Python/20220715_0.html
  7. Computer___Network/20230508_00.html
  8. Food___Taiwan/20220429_0.html
  9. Computer___NetBSD/20220817_3.html
  10. Computer___Python/20220410_0.html
  11. Computer___Network/20240416_00.html
  12. Computer___Network/20240130_00.html
  13. Computer___Debian/20210223_1.html
  14. Computer___NetBSD/20230119_00.html
  15. Computer___Python/20210124_0.html
  16. Computer___Python/20221013_0.html
  17. Computer___NetBSD/20220818_1.html
  18. Computer___NetBSD/20220428_0.html
  19. Science___Math/20220420_0.html
  20. Computer___NetBSD/20240101_02.html
  21. Computer___NetBSD/20220808_0.html
  22. Computer___TeX/20230503_00.html
  23. Computer___NetBSD/20230515_00.html
  24. Science___Astronomy/20220503_0.html
  25. Computer___NetBSD/20210127_0.html
  26. Computer___Python/20240101_00.html
  27. Computer___Network/20220413_1.html
  28. Computer___Python/20220816_1.html
  29. Computer___NetBSD/20210204_0.html
  30. Travel___Taiwan/20220809_2.html


HTML file generated by Kinoshita Daisuke.