16

円周率の計算に助けが必要です。Pi から X 桁までを計算する Python プログラムを作成しようとしています。Pythonメーリングリストからいくつか試してみましたが、私の使用には遅すぎます。Gauss-Legendre Algorithmについて読んだことがありますが、それを Python に移植しようとしましたが、成功しませんでした。

私はここから読んでいます。私がどこで間違っているのかについての意見をいただければ幸いです!

出力: 0.163991276262

from __future__ import division
import math
def square(x):return x*x
a = 1
b = 1/math.sqrt(2)
t = 1/4
x = 1
for i in range(1000):
    y = a
    a = (a+b)/2
    b = math.sqrt(b*y)
    t = t - x * square((y-a))
    x = 2* x

pi = (square((a+b)))/4*t
print pi
raw_input()
4

3 に答える 3

28
  1. の周りの括弧を忘れました4*t:

    pi = (a+b)**2 / (4*t)
    
  2. decimalより精度の高い計算ができます。

    #!/usr/bin/env python
    from __future__ import with_statement
    import decimal
    
    def pi_gauss_legendre():
        D = decimal.Decimal
        with decimal.localcontext() as ctx:
            ctx.prec += 2                
            a, b, t, p = 1, 1/D(2).sqrt(), 1/D(4), 1                
            pi = None
            while 1:
                an    = (a + b) / 2
                b     = (a * b).sqrt()
                t    -= p * (a - an) * (a - an)
                a, p  = an, 2*p
                piold = pi
                pi    = (a + b) * (a + b) / (4 * t)
                if pi == piold:  # equal within given precision
                    break
        return +pi
    
    decimal.getcontext().prec = 100
    print pi_gauss_legendre()
    

出力:

3.141592653589793238462643383279502884197169399375105820974944592307816406286208\
    998628034825342117068
于 2008-12-07T16:29:55.790 に答える
3
pi = (square((a+b)))/4*t

する必要があります

pi = (square((a+b)))/(4*t)
于 2008-12-07T16:37:09.413 に答える
3
  1. PI を 1000 桁まで計算する場合は、1000 桁の精度をサポートするデータ型を使用する必要があります (例: mxNumber )
  2. |ab|までa、b、t、およびxを計算する必要があります。< 10** 桁、桁回繰り返しません。
  3. @JFが示唆するように、平方と円周率を計算します。
于 2008-12-07T16:38:24.790 に答える