Aperçu

Le 14 mars, les dépanneurs et les grands magasins sont bondés de ventes de bonbons, et je pense que beaucoup de gens sont impressionnés que le printemps approche. Les gens dans le monde attirent également l'attention le jour du rapport de circonférence, et je pense qu'il y a un effet économique considérable. Donc, aujourd'hui, j'écrirai un programme pour calculer le rapport de circonférence.

C'est une histoire simple.

Algorithme de Gauss Legendre

[Algorithme Gauss-Legendor] de wikipedia (https://ja.wikipedia.org/wiki/%E3%82%AC%E3%82%A6%E3%82%B9%EF%BC%9D%E3%83] % AB% E3% 82% B8% E3% 83% A3% E3% 83% B3% E3% 83% 89% E3% 83% AB% E3% 81% AE% E3% 82% A2% E3% 83% AB Selon% E3% 82% B4% E3% 83% AA% E3% 82% BA% E3% 83% A0), il est relativement facile à mettre en œuvre.

Je pensais que c'était un grand génie de penser ainsi.

implémentation python

Cependant, je ne sais pas comment calculer le rapport de circonférence, alors j'ai regardé les sources de différentes personnes et c'était terminé. En tant que fonctionnalité, Decimal était utilisé pour représenter un grand nombre de chiffres.

# -*- coding:utf-8 -*-
import math
import time
from decimal import *

import numpy as np


#Nombre de caractères 1002
PI_1000 = "3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679821480865132823066470938446095505822317253594081284811174502841027019385211055596446229489549303819644288109756659334461284756482337867831652712019091456485669234603486104543266482133936072602491412737245870066063155881748815209209628292540917153643678925903600113305305488204665213841469519415116094330572703657595919530921861173819326117931051185480744623799627495673518857527248912279381830119491298336733624406566430860213949463952247371907021798609437027705392171762931767523846748184676694051320005681271452635608277857713427577896091736371787214684409012249534301465495853710507922796892589235420199561121290219608640344181598136297747713099605187072113499999983729780499510597317328160963185950244594553469083026425223082533446850352619311881710100031378387528865875332083814206171776691473035982534904287554687311595628638823537875937519577818577805321712268066130019278766111959092164201989"


def get_pi(prec=1000, verbose=False):
    '''
Renvoie le rapport de circonférence de prec avec le nombre de chiffres après la virgule décimale sous forme de chaîne de caractères
    '''
    prec = prec+1+1 #précision"3."Puisqu'il y a une partie de, augmentez la précision de un, et si vous n'avez pas d'autre chiffre, il peut se déplacer en raison de l'arrondi, alors ajoutez-en un de plus
    N = 2*math.ceil(math.log2(prec)) #Nombre de répétitions Selon wikipedia, log2(prec)Il semble que le degré soit suffisant, donc juste au cas où, double-le

    getcontext().prec = prec 

    a = Decimal(1.0)
    b = Decimal(1.0) / Decimal(2.0).sqrt()
    t = Decimal(0.25)
    p = Decimal(1.0)


    start_time = time.clock() 

    #Démarrer N essais
    for _ in range(1, N):
        a_next = (a + b)/Decimal(2)
        b_next = (a*b).sqrt()
        t_next = t - p * (a - a_next)**2
        p_next = Decimal(2)*p

        a = a_next
        b = b_next
        t = t_next
        p = p_next

    #Calculer le rapport de circonférence
    pi = (a + b)**2 / (Decimal(4)*t)
    #Calculer le temps d'exécution
    elapsed = time.clock()  - start_time

    #Afficher le résultat de l'exécution
    if verbose:
        print("précision:",prec)
        print("Pi:", pi)
        print("Temps d'exécution:%f" % (elapsed))

    return str(pi)[0:prec]

Je l'ai exécuté comme ci-dessous et j'ai confirmé qu'environ 1000 chiffres étaient corrects. (Certaines personnes calculent environ 200 millions de chiffres, mais c'est très bien.)

pi = get_pi(prec=1000, verbose=True)
print(pi)
#print(PI_1000)
#print(pi == PI_1000)

Précision: 1002
Pi: 3.14159265358979323846264338327950288419716939937510582097494459230781640628620899862803482534211706798214808651328230664709384460955058223172535940812848111745028410
270193852110555964462294895493038196442881097566593344612847564823378678316527120190914564856692346034861045432664821339360726024914127372458700660631558817488152092096282925
409171536436789259036001133053054882046652138414695194151160943305727036575959195309218611738193261179310511854807446237996274956735188575272489122793818301194912983367336244
065664308602139494639522473719070217986094370277053921717629317675238467481846766940513200056812714526356082778577134275778960917363717872146844090122495343014654958537105079
227968925892354201995611212902196086403441815981362977477130996051870721134999999837297804995105973173281609631859502445945534690830264252230825334468503526193118817101000313
783875288658753320838142061717766914730359825349042875546873115956286388235378759375195778185778053217122680661300192787661119590921642019890
Temps d'exécution:0.014581

Statistiques

[Introduction aux statistiques](https://www.amazon.co.jp/ Introduction to Statistics-Basic Statistics-Statistics Department, Faculty of Liberal Arts, University of Tokyo / dp / 4130420658) Un problème est survenu lors de la vérification du tarif, alors j'ai essayé de le résoudre.

Après avoir regardé le graphique, nous effectuerons également un test du chi carré pour voir si les nombres sont constants.

Lorsque 100 nombres sont extraits du 123e chiffre

Le 123ème chiffre est approprié.

start = 123
pi_123 = pi[2+start:2+start+100]

freq = np.zeros(10, dtype=np.int)
for a in pi_123:
    freq[int(a)] += 1
print("freq(100)=",freq)

#Demander avec chisquare
print(chisquare(freq, 10*np.ones(10,dtype=np.int8)))

plt.xlim(-0.5, 9.5)
plt.bar(np.array([0,1,2,3,4,5,6,7,8,9]), freq, align='center')

En regardant le graphique, cela semble désagréable, mais comme la valeur p est d'environ 0,37, il n'y a pas de contradiction même si le taux d'apparition est presque constant.

freq(100)= [ 9 12 11  7 15 13  7  4 12 10]
Power_divergenceResult(statistic=9.8000000000000007, pvalue=0.3669177991127523)

Au moment de 1000 chiffres

freq = np.zeros(10, dtype=np.int)
for a in pi:
    if a != ".":
        freq[int(a)] += 1
print("freq(1001)=",freq)
print(chisquare(freq, 100.1*np.ones(10,dtype=np.int8)))

plt.xlim(-0.5, 9.5)
plt.bar(np.array([0,1,2,3,4,5,6,7,8,9]), freq, align='center')

Lorsqu'environ 1000 données sont collectées, il n'y a pas de contradiction même si le taux d'apparition est considéré comme constant quel que soit le nombre. La valeur p est de 0,85, ce qui est bien.

freq(1001)= [ 93 116 103 103  93  97  94  95 101 106]
Power_divergenceResult(statistic=4.7842157842157853, pvalue=0.85269838594638103)

Être inquiet

En fait, je voulais connaître les informations statistiques de la séquence de nombres, mais c'est la suivante. Je dois étudier les nombres transcendantaux.

Le 14 mars est le jour du rapport de circonférence. L'histoire du calcul du ratio de circonférence avec python