Algoritmi za izračun števila Pi

Algoritmi za izračun števila Pi

Avtor: Eva Janša

Število Pi

Število je matematična konstanta, ki je prestavljena s številom 3,14159... Predstavimo ga z obsegom kroga deljenim z njegovim premerom: .

(http://www2.nauk.si/files/447/krog_pi.jpg)
pi je enak obseg/premer

Vir slike : Wikipedija

Je iracionalno število, predstavimo ga lahko tudi z ulomkom, kot na primer 22/7.

V tej predstavitvi si bomo ogledali nekaj algoritmov za izračun števila . Te algoritme bomo nato primerjali glede na hitrost in zahtevnost. Na koncu pa si bomo ogledali še nekaj kratih formul za izračun števila .

Pogledali si bomo še, katerega od teh algoritmov ali formul za izračun , najpogosteje uporabljamo v praksi.

(http://www2.nauk.si/files/447/pi_d.jpg)
pi z nekaj decimalkami

Vir slike : Google

Zgodovina

  • Število naj bi bilo prvič omenjeno v Svetem pismu. Tam je zaokroženo na 3. Kasneje naj bi s računali že Babilonci in Egipčani, predvsem naj bi število uporabljali za izgradnjo piramid.
  • Prvi algoritem za izračun števila , je podal Arhimed v stari Grčiji. Ta je izračunal s pomočjo krogu, z radijem 1, vrisanega in orisanega pravilnega n kotnika. Za n=96 je izračunal vrednost . Zaradi tega izračuna številu velikokrat pravimo tudi Arhimedovo število. Kasneje je Ptolomej izračunal na 3,1416.
  • Po tem obdobju so pi poskušali bolj natančno izračunati šele v renesansi. Kasneje so si matematiki namesto z geometrijo pomagali z neskončnimi formulami. S formulami so se začeli ukvarjati že v Indiji, prvi Evropski matematik, ki je objavil neskončni produkt, ki konvergira k je bil François Viète. Kasneje so neskončne vsote, ki predstavljajo približek za objavili še Wallis, Newton, ter Leibniz.
  • Pomembnejša formula, ki jo bomo predstavili tudi kasneje pa je formula Leibniz-Gregory. Svoj prispevek sta bazirala na starejših izračunih indijskega matematika Madhave. Angleški matematik Sharp je s pomočjo tega algoritma izračunal na 71 decimalk natančno.
  • Leta 1706 je matematik Machin s pomočjo algoritma Leibniz-Gregory sestavil formulo, ki izračuna na 100 decimalk natančno. To formulo si bomo pogledali kasneje. Variacij te formule je veliko, uporablja pa se še vedno za izračune v računalništvu.
  • Euler je svoj algoritem za dobil po tem, ko je dokazal problem Basel. S tem sicer ni pridobil na bolj natančnem izračunu števila, je pa podal nov algoritem, ki si ga bomo tudi ogledali kasneje.
  • Leta 1761 je švicarski matematik Lambert dokazal, da je iracionalno število. Leta 1794 je matematik Legendre pokazal še, da je tudi iracionalno število.
  • V dvajsetem stoletju so se pojavili računalniki, zato so izračuni decimalk števila hitro naraščali. Temu pa so pripomogla tudi odkritja rekurzivnih formul, ki so v nekaterih primerih hitrejša, kot neskončne formule, ter množenje večjih števil hitrejše, ker je pripomoglo h hitrosti izračuna. Ena od rekurzivnih formul, ki jih bom predstavila kasneje je tudi algoritem Gauss-Legendre.
  • Leta 1914 je indijski matematik Ramanujan dokazal, da njegova neskončna vsota konvergira k številu . Njegov algoritem je konvergiral do tega števila zelo hitro. Tudi ta algoritem si bomo bolj podrobno ogledali kasneje.
  • Do današnjih dni, se novi algoritmi in formule za izračun še kar pojavljajo. Njihov namen je predvsem čim hitrejši in čim bolj natančen izračun decimalk, kar bi pohitrilo tudi računanje na računalnikih.

O algoritmih

Kot rečeno si bomo ogledali nekaj pomembnejših algoritmov in formul za izračun števila . Prvi trije algoritmi so izpeljani iz neskončnih vrst. Razporejeni so po času, kdaj so bili odkriti. Četrti algoritem pa je podan rekurzivno. Na koncu si pogledamo še nekaj formul za izračun števila pi. Tudi formule so razporejene po času, kdaj so bile odkrite.

Vse algoritme smo napisali tudi v obliki metode v programskem jeziku Python. Ideja posamezne metode je, da s pomočjo vsakega od algoritmov izračunamo na dano število decimalk natančno. Poleg tega pa v vsaki metodi merimo čas izvajanja in pa število iteracij. Na koncu vsake metode torej vrnemo naračunan , trajanje izvajanja metode in pa število iteracij. Ker to naredimo za vse metode na koncu lažje primerjamo njihovo učinkovitost in hitrost.

Leibniz-Gregory

Prvi algoritem je algoritem Leibniz-Greory, nekje poimenovan tudi Madhava-Leibniz-Gregory. Tako ime ima zato, ker bazira na neskončni vsoti podani spodaj.


Prvi je vsoto zastavil že Madhava v 14. stoletju. Ker pa ni bila znana v zahodnem svetu, sta jo kasneje formurirala, najprej Leibniz, potem pa še Gregory. Zato se tudi algoritem imenuje po teh treh matematikih.
V algoritmu na vsakem sodem koraku prištejemo, na vsakem lihem pa odštejemo ulomek z lihim števcem. Členi vrste absolutno konvergirajo proti 0, celotna vsota pa konvergira proti . Konvergira počasi, saj za pravilnih decimalk potrebujemo iteracij.
Da se vsota res sešteje v dokažemo s pomočjo izračuna .

Razlaga algoritma(Python)

Metoda s pomočjo algoritma Leibniz-Gregory izračuna na dano število decimalk natančno. Da smo naračunali pravilno število decimalk vemo, ko se prejšnji in trenutno naračunani približek razlikujeta za 10^(decimalke+1) ali manj. Napako izračunamo vsakič znova na tekočem koraku. Najprej začnemo meriti čas in določimo začetne približke ter pomožne spremenljivke v katere bomo shranjevali , ter število iteracij. Potem izvajamo zanko toliko časa, dokler ne naračunamo danega števila decimalk. V zanki najprej naračunamo i-ti člen vsote. Členu najprej spremenimo predznak, potem izračunamo vrednost ulomka in ga prištejemo vsoti. Ker računamo in ne rezultat pomnožimo s 4. Potem ponastavimo še približke. Liho število povečamo za 2, števec za iteracije pa za 1. Preden pa zanko nadaljujamo z naslednjim korakom, izračunamo napako med zdajšnjim in prejšnjim približkom. V pomožno spremenljivko shranimo trenutno naračunan približek. Ko se zanka izteče izračunamo čas izvajanja ter vrnemo naračunan , čas izvajanja in število iteracij.

Kodo metoda najdete na naslednji prosojnici. Koda

Leibniz-Gregory - Python

from math import*
import time
def LeibnizGregory(decimalke):
    '''izračunamo pi do danega števila decimalk po algoritmu Leibniz-Gregory'''
    #začnemo meriti čas izvajanja
    zac = time.clock()
    #začetne spremenljivke
    trPi = 1
    mn = 1
    nas = 3
    rez = 1
    pom = 1
    napaka = 1
    #števec iteracij
    i = 1
    #dokler pi ne izračunamo do danega števila decimalk
    while napaka > 10**((decimalke + 1) * (-1)):
        #ponavljamo korake vsote
        mn *= -1
        rez += mn * (1 / nas)
        trPi = 4 * rez
        nas += 2
        #povečamo števec iteracij
        i+= 1
        #izračunamo napako in shranimo trenutno vrednost
        napaka = abs(trPi - pom)/abs(pom)
        pom = trPi
    kon = (time.clock() - zac)
    #vrnemo naračunan pi, čas in število iteracij
    return trPi, kon, i

Euler

Drugi algoritem je iznašel Euler. Tudi tukaj algoritem bazira na neskončni vrsti, ki je prikazana spodaj.


Euler je na izračun števila naletel leta 1735, ko se je ukvarjal s problemom Basel. Ta problem je zahteval natančen seštevek zgornje vsote in pa dokaz, da je ta vsota pravilna. Pri dokazu si je pomagal s Taylorjevo vrsto za sinus.

Algoritem na vsakem koraku vsoti prišteje ulomek, ki ima v imenovalcu 1, v števcu pa kvadrat števila. Členi konvergirajo proti 0, vsota pa se v neskončnosti približuje kvadratu deljenemu s 6.

Razlaga algoritma(Python)

Metoda s pomočjo Eulerjevega algoritma izračuna na dano število decimalk natančno. V metodi najprej začnemo meriti čas, določimo začetne približke in pomožne spremenljivke v katere bomo shranjevali , ter število iteracij. Potem izvajamo zanko toliko časa, dokler ne naračunamo danega števila decimalk. Da smo naračunali pravilno število decimalk vemo, ko se prejšnji in trenutno naračunani približek razlikujeta za 10^(decimalke+1) ali manj. Napako izračunamo vsakič znova na tekočem koraku. V zanki najprej naračunamo kvadrat trenutne iteracije. K vsoti nato prištejemo ulomek s tem števcem. Ker računamo in ne . Nato rezultat najprej pomnožimo s 6, potem pa ga še korenimo. Na koncu povečamo števec iteracij in izračunamo napako, kot pri prejšnjem algoritmu, ter shranimo vrednost v pomožno spremenljivko. Ko se zanka izteče izračunamo čas izvajanja ter vrnemo naračunan , čas izvajanja in število iteracij.

Kodo metoda najdete na naslednji prosojnici. Koda

Euler - Python

from math import*
import time
def Euler(decimalke):
    #začnemo meriti čas izvajanja
    zac = time.clock()
    #začetne spremenljivke
    trPi = 1
    pom = 1
    napaka = 1
    rez = 0
    #števec iteracij
    i = 1
    #dokler pi ne izračunamo do danega števila decimalk
    while napaka > 10**((decimalke + 5) * (-1)):
        #ponavljamo korake vsote
        k = i*i
        rez += 1/k
        trPi = sqrt(rez* 6)
        #števec povečamo
        i += 1
        #izračunamo napako in shranimo trenutno vrednost
        napaka = abs(trPi - pom)/abs(pom)
        pom = trPi
    kon = (time.clock() - zac)
    #vrnemo naračunan pi, čas in število iteracij
    return trPi, kon, i

Ramanujan

Najhitrejši in najučinkovitejši med predstavljenimi algoritmi je baziran na spodnji neskončni vsoti. Pravilnost vsote je okoli 1910 dokazal indijski matematik Srinivasa Ramanujan(1887 - 1920).


Vsota konvergira zelo hitro proti . Na vsakem koraku iteracije dobimo 8 dodatnih pravilnih decimalk.

Leta 1985 je to neskončno vsoto uporabil William Gosper za izračun števila na 17 miljonov decimalk natančno.

Razlaga algoritma(Python)

Metoda s pomočjo Ramanujangovega algoritma izračuna na dano število decimalk natančno. V metodi najprej začnemo meriti čas, nato določimo začetne približke in pomožne spremenljivke v katere bomo shranjevali , ter na koncu še število iteracij. Izračunamo tudi faktor s katerim pomnožimo vsoto. Potem izvajamo zanko toliko časa, dokler ne naračunamo danega števila decimalk. Da smo naračunali pravilno število decimalk vemo, ko se prejšnji in trenutno naračunani približek razlikujeta za 10^(decimalke+1) ali manj. Napako izračunamo vsakič znova na tekočem koraku. Ker je člen te vsote kar zapleten za izračun v koraku zanke izračunamo najprej posamezne operacije, ki jih na koncu združimo v vrednost člena. Najprej izračunamo obe fakulteti, nato množenje, potem še poteciranje, na koncu pa še seštevanje. Ko imamo vse delne rezultate, jih zložimo v ulomek s katerim nato izačunamo vrednost člena, ki ga prištejemo k vsoti. Vsoto nato pomnožimo še s faktorjem, povečamo števec za iteracije in izračunamo . Trenutni rezultat v vsoti je namreč enak , zato ta izračun shranimo v pomožno spremenljivko. Tudi tukaj kot v prejšnjih dveh algoritmih za naračunan izračunamo napako in ga shranimo v pomožno spremenljivko. Ko se zanka izteče izračunamo čas izvajanja ter vrnemo naračunan , čas izvajanja in število iteracij.

Kodo metoda najdete na naslednji prosojnici. Koda

Ramanujan - Python

from math import*
import time
def Ramanujan(decimalke):
    #začnemo meriti čas izvajanja
    zac = time.clock()
    #začetne spremenljivke
    trPi = 0
    F1 = 1
    F2 = 1
    E2 = 1
    pom = 1
    napaka = 1
    fk = (2*sqrt(2))/9801
    #števec iteracij
    i = 0
    #dokler pi ne izračunamo do danega števila decimalk
    while napaka > 10**((decimalke+1) * (-1)):
        #ponavljamo korake vsote
        F1 = factorial(4 * i)
        F2 = factorial(i)
        P1 = 23690 * i
        E1 = F2**4
        E2 = 396**(4 * i)
        S1 = 1103 + P1
        trPi += (F1 * S1)/(E1 * E2)
        #povečamo števec iteracij
        i += 1
        #izračunamo pi
        trPi1 = 1 / (fk*trPi)
        #izračunamo napako in shranimo trenutno vrednost
        napaka = abs(trPi1 - pom)/abs(pom)
        pom = trPi1

    kon = (time.clock() - zac)
    #vrnemo naračunan pi, čas in število iteracij
    return trPi1, kon, i

Gauss-Legendre

Algoritem Gauss-Legendre je bil prvi pomembnješi algoritem, ki je bil podan z rekurzivno formulo. Ker je formula podana rekurzivno, algoritem proti konvergira veliko hitreje, kot pred tem dokazane neskončne vsote ali formule. Tako lahko z 25 iteracijami izračunamo 45 milijonov decimalk.

Spodnji algoritem temelji na delu matematikov Carla Friedricha Gaussa in Adriena-Marie Legendre-ja, ter na nadgradnji, ki sta jo izvedla Richard Brent in Eugene Salamin.

Algoritem

Začetne vrednosti:





Rekurzija:





Izračun pi:



Razlaga algoritma(Python)

Metoda s pomočjo algoritma Gauss-Legendre izračuna na dano število decimalk natančno. V metodi najprej začnemo meriti čas, določimo začetne približke in spremenljivko za štetje iteracij. Potem izvajamo zanko toliko časa, dokler ne naračunamo danega števila decimalk. Da smo naračunali pravilno število decimalk vemo, ko se prejšnji in trenutno naračunani približek razlikujeta za 10^(decimalke+1) ali manj. Napako izračunamo vsakič znova na tekočem koraku. V zanki najprej naračunamo naslednje vrednosti začetnih spremenljivk, kot so podane zgoraj. Tako naračunamo nove a, b, t in p. Le te sedaj shranimo v stare spremenljivke. Iz teh vrednosti nato naračunamo trenutno vrednost . Nato povečamo števec iteracij, izračunamo napako in naračunani shranimo v pomožno spremenljivko. Ko se zanka izteče izračunamo čas izvajanja ter vrnemo naračunan , čas izvajanja in število iteracij.

Kodo metoda najdete na naslednji prosojnici. Koda

Gauss-Legendre - Python

from math import*
import time
def GaussLegendre(decimalke):
    #začnemo meriti čas izvajanja
    zac = time.clock()
    #začetne spremenljivke
    trPi = 1
    pom = 1
    napaka = 1
    a = 1
    b = 1/sqrt(2)
    t = 1/4
    p = 1
    i = 0
    #dokler pi ne izračunamo do danega števila decimalk
    while napaka > 10**((decimalke+1) * (-1)):
        a1 = (a+b)/2
        b1 = sqrt(a*b)
        t1 = t - p*(a-a1)**2
        p1 = 2*p
        a,b,t,p = a1,b1,t1,p1
        trPi = ((a+b)**2)/(4*t)
        #števec iteracij
        i += 1
        #izračunamo napako in shranimo trenutno vrednost
        napaka = abs(trPi - pom)/abs(pom)
        pom = trPi
    kon = (time.clock() - zac)
    #vrnemo naračunan pi, čas in število iteracij
    return trPi, kon, i

Izračuni z enačbami-Liu Hui

Liu Hui

Kitajski matematik je s to enačbo izračunal na 5 decimalk natančno. Pomagal si je s pravilnim včrtanim 96-kotnikom, v krog z radijem 1. To je dokazal s premislekom, da če krogu vrišemo pravilni n-kotnik, pri čemer n pošljemo proti neskončno, se njegova ploščina približa ploščini kroga. Za 5 decimalk je torej dovolj ploščina 96-kotnika iz katere nato izrazimo .



Koda

from math import*
import time
def LiuHui():
    '''izračun pi s formulo Liu Hui-ja'''
    z = time.clock()
    r = 768*sqrt(2-sqrt(2+sqrt(2+sqrt(2+sqrt(2+sqrt(2+sqrt(2+sqrt(2+sqrt(2+1)))))))))
    k = time.clock() - z
    return r,k

Izračuni z enačbami - Euler

Euler

Euler je poleg algoritma poiskal tudi formulo za izračun . Enake oblike kot ta formula, je tudi formula na naslednji prosojnici, torej Machinova formula.




Koda

from math import*
import time
def Euler():
    '''izračun pi s formulo Euler-ja'''
    z = time.clock()
    r = 20 * atan(1/7) + 8 * atan(3/79)
    k = time.clock() - z
    return r,k

Izračuni z enačbami - Machin

Machin

S to enačbo je John Machin število izračunal na 100 decimalk natančno. Formula izhaja iz algoritma Leibniz-Gregory. Vse formule so oblike
,
in zhajajo iz Taylorjeve vrste za arctangens.
Osnovna oblika je:

Večkrat pa je uporabljena tudi:


Koda

from math import*
import time
def Machin1():
    '''izračun pi z originalno formulo Machin-a'''
    z = time.clock()
    r = 4*(4*atan(1/5)-atan(1/239))
    k = time.clock() - z
    return r,k

def Machin2():
    '''izračun pi z eno izmed formul, podobnih Machin-ovi'''
    z = time.clock()
    r = 4*(atan(1/2)+atan(1/3))
    k = time.clock() - z
    return r,k

Primerjava hitrosti

Kviz-Prvo vprašanje

Kateri je bil prvi algoritem za izračun podan z neskončno vsoto?



Ramanujan
Gauss - Legendre
Euler
Leibniz - Gregory

Right

Pravilno!

Next

Wrong

Napačno!

Repeat

Kviz-Drugo vprašanje

S katerim algoritmom smo pi izračunali najhitreje?

Gauss-Legendre
Leibniz-Gregory
Euler
Ramanujan

Right

Pravilno!

Next

Wrong

Napačno!

Repeat

Kviz-Tretje vprašanje

Katero formulo za izračun števila pi uporabljajo v računalništvu?

Eulerjevo
Lui Hui-jevo
originalno Machinovo

Right

Pravilno!

Next

Wrong

Napačno!

Repeat

Viri

0%
0%