The Amazing World of Calculus Programs ( for Math Fans only )

Paolo1956

LOREM IPSUM
Registrato
20/6/10
Messaggi
6.484
Punti reazioni
292

piof

Nuovo Utente
Registrato
8/5/09
Messaggi
13.839
Punti reazioni
545
Si, non è presentato molto bene, ma quel problema me lo ricordavo dai corsi di calcolo.
Oggi non dovrebbe riscontrarsi più...

Ps. bellina quella formula per pi greco, non la conoscevo

infatti ero un poco stupito di quel video, mi pareva insolito un errore di calcolo in un calcolo di quel tipo, è che quel cinese faceva proprio pensare vi fosse un errore, non che w.alpha desse una delle ben 36 soluzioni in C. 36 son tante...

Però il problema c'è, nel senso che se inserisco proprio un calcolo che va a creare cancellazione, (appositamente o per puro caso), la ottengo e l'errore c'è, credo anche in matlab. Certamente in excel, ho inserito ieri sera proprio quella successione e non converge a pigreco.

Ed è carina, un po' diversa dalle solite serie che convergono a pigreco
 

Paolo1956

LOREM IPSUM
Registrato
20/6/10
Messaggi
6.484
Punti reazioni
292
In R converge. Almeno, io ho provato fino alla decima iterazione. Poi avrei dovuto creare un ciclo e non avevo voglia :)

Comunque il cinese ha un sacco di video di curiosità matematiche, se sei appassionato merita.
 

piof

Nuovo Utente
Registrato
8/5/09
Messaggi
13.839
Punti reazioni
545
Fino alla decima iterazione credo vada tutto bene con qualsiasi programma, con Excel mi pare di ricordare che l errore aumenti dopo la 16esima (dà circa 2 invece di 3.14 in poche iterazioni oltre la 15 o 16esima) e con matlab forse pure oltre la sedicesima dà errore crescente. Controllo domani come era su excel..
E guarderò qualche video del cinese , uno lo avevo già notato
 

piof

Nuovo Utente
Registrato
8/5/09
Messaggi
13.839
Punti reazioni
545
un tipo creativo, direi, piacevole

...e per la nostra successione i problemi con excel iniziano dopo la 20esima iterazione
 

Paolo1956

LOREM IPSUM
Registrato
20/6/10
Messaggi
6.484
Punti reazioni
292
Veramente rognosa quella successione:

in R, usando una precisione di 2048 bit (che è una enormità) in 50 iterazioni finisce comunque per superare il valore reale di PI, dando una precisione di sole 14 cifre:

Codice:
require (Rmpfr)
Area1   <- function(n, old){
  rad1    <- sqrt(1-(4)^(1-n)*old^2)
  rad2    <- sqrt(1 - rad1)
  res     <- 2^(n-0.5)*rad2
  return(res)
}
n       <- 1
old     <- Rmpfr::mpfr(2.0, 2048)
for(i in 1:50) {
  n <- n+1
  old <- Area1(n,old)
  print(old)
}

Output:

Codice:
1 'mpfr' number of precision  2048   bits 
[1] 3.14159265358980405029618843975265921602255351755556208699626836673877826752455846143604095423679699984787137490887594299371983549190055833273516356699082307872666382578565364234403051144270305685435632938696341849555675696248819901007869493328957326129905274074323298467891386847386793379693782862221313317670193928463240816883451027147795030388156575441165639048708124228585053872633339992507795480497014134473166681904138856365026895319335081591861802378776965151616534126826373529677784539606538413582520039312590543298192933648305718076934728320863729384537514847688401410304825654122569312276207442055400934910261

che è sbagliato :(
 

piof

Nuovo Utente
Registrato
8/5/09
Messaggi
13.839
Punti reazioni
545
è che senza rielaborare la successione eliminando quella sottrazione incriminata nella radice (1 meno qualcosa che tende ad uno)...prima o poi si ha "cancellazione"
 

Paolo1956

LOREM IPSUM
Registrato
20/6/10
Messaggi
6.484
Punti reazioni
292
Credo di esserci riuscito (non facendo come nel papero):

- siccome z(n) sotto radice compare sempre al quadrato è inutile estrarre ogni volta la radice esterna.
- allora la quantità 2^(n-1/2) diventa 2^(2n-1), che calcolo usando precisione infinita
- analogamente nella radice più interna 4^(1-n) lo calcolo come 4^(n-1) in precisione infinita e poi lo trasformo in un float a 2048 bit e faccio il reciproco

Con 100 iterazioni trovo PI esatto alla 60ma cifra decimale (e soprattutto < del vero valore).

Ma che fatica!!!

Codice:
require(Rmpfr)
require(gmp)
Area2   <- function(n, old) {
  q     <- as.bigz(4)
  q     <- q^(n-1)
  q1    <- Rmpfr::mpfr(q,2048)
  q1    <- 1/q1
  rad1  <- sqrt(1 - q1*old)
  p     <- as.bigz(2)
  p     <- p^(2*n - 1)
  old   <- p*(1 - rad1)
  return(old)
}
n       <- 1
old     <- Rmpfr::mpfr(4.0, 2048)
for(i in 1:100) {
  n <- n+1
  old <- Area2(n,old)
  
}
print(sqrt(old))

1 'mpfr' number of precision  2048   bits 
[1] 3.14159265358979323846264338327950288419716939937510582097494378833892686005527597780928589562502813124705832451502523826046252968832175785865439300607273183404782709511953884763535343058989009302756785653111742322160268851567898066533329416273664962944177475194463528954749173338670885869472965487026896903609683597936676790493500791785765719271969913703091752982879382739432850806820543004089909004394575171426424151381842914373436487293506925933578900303558567822435389333670947574913435233046783945826492651967650725749757282487775161032091796352084196378096038204192226185457982406422798069987615850942542823194006
 

Paolo1956

LOREM IPSUM
Registrato
20/6/10
Messaggi
6.484
Punti reazioni
292
Per fare una cosa del genere quand'ero studente (anni 70) ci voleva il computer del Centro Nazionale di Calcolo a Pisa:

fai conto un piano di un palazzo con l'aria condizionata e l'abbattimento delle polveri riempito di cassoni con le memorie a nastro che giravano continuamente.

Che tempi!
 

piof

Nuovo Utente
Registrato
8/5/09
Messaggi
13.839
Punti reazioni
545
non ero neppure nato!
Ottima idea quella di evitare l'estrazione di radice ad ogni iterazione, per calcolarla solo alla fine. Però è strano non dia cancellazione, il passaggio incriminato (1 - rad1) resta, non è che la dà un po' oltre le 100 iterazioni? Andrebbe verificato se al crescere del numero di iterazioni accade che, oltre un certo numero di iterazioni, l'errore (la differenza tra pigreco ed il risultato dato dal programma) cresce invece di scendere:D. Ma credo ci si possa contentare di 60 decimali esatti
 

Paolo1956

LOREM IPSUM
Registrato
20/6/10
Messaggi
6.484
Punti reazioni
292
Estendendola a 200 iterazioni la precisione passa alla 200 ma cifra e rimane < pi.

Ci dovrei ragionare un po' sopra. La perdita di precisione avviene nell'estrazione di radice (z è irrazionale), mentre tutti gli altri numeri in gioco possono essere rappresentati come esatti (anche 1/4^(n-1) e quindi è critico l'algoritmo di estrazione della radice (che sarà un qualche Newton modificato). In questa situazione l'algoritmo delle tangenti tende a dare risultati approssimati per eccesso (e quindi sfuggirebbe via) ma può darsi usino qualche altra cosa.
 

piof

Nuovo Utente
Registrato
8/5/09
Messaggi
13.839
Punti reazioni
545
certo, per come è creata la successione, (aree triangoli interni al cerchio) deve uscire una successione di valori strettamente crescente, tendente a pigreco, con valori sempre minori di pigreco. Se così è allora è tutto a posto, e neppure importa molto che vi sia o no cancellazione numerica, che anzi dovrebbe restare poichè anche tolta una radice, rimane l'altra: 1- radice di qualcosa tendente ad 1. Poi quale algoritmo usa R, che in effetti potrei scaricare, per le radici proprio non lo so, in ogni caso essendo irrazionale il risultato delle varie radici, approssima per forza

Direi che rimane aperta la questione: proseguendo con le iterazioni si ottengono via via le cifre esatte di pigreco, oppure oltre un certo n la cancellazione produce i suoi effetti nefasti?

Indubbiamente è già un eccellente risultato ottenere 60 cifre decimali esatte, eliminando una delle radici, si elimina una fonte di approssimazione
 

Paolo1956

LOREM IPSUM
Registrato
20/6/10
Messaggi
6.484
Punti reazioni
292
E' bravo, anche se si vede che è un ragazzino
 

piof

Nuovo Utente
Registrato
8/5/09
Messaggi
13.839
Punti reazioni
545
quasi certamente destinato ad andar via presto dall'Italia, nel caso in cui non gli venga consentito di iscriversi ad un corso universitario prima di altri 5 anni