Pi greco. perchè non va?

[deian91]

// Algoritmo per il calcolo approssimato di pi greco attraverso l'approssimazione numerica
// dell'integrale: 3*int [0:1] (1/(1+x²))

#include
#include

int my_pow(int steps);

main() {
double *f,*a,*A;
double base,pi,a1=0,A1=0,x=0;
int i,steps;

printf( "\nAttenzione: Maggiore sarà il numero di passi di approssimazione, maggiore sarà il tempo\n"
"che il calcolatore impiegherà a calcolare pi!\n\n");

printf( "Quanti passi di approssimazione devo effetturare? (1-500) ");
scanf("%d",&steps);

if ((steps<1) || (steps>500)) {
printf( "Numero di passi non valido\n");
return -1;
}

f = (double*) malloc (my_pow(steps-1)*sizeof(double));
a = (double*) malloc (my_pow(steps-1)*sizeof(double));
A = (double*) malloc (my_pow(steps-1)*sizeof(double));

base = (float) 1/(my_pow(steps-1));
f[0] = 1;

for (i=1; i<=my_pow(steps-1); i++) {
x += base;
f = 1/(1+x*x);
a = base*f;
A = base*f[i-1];
}

for (i=1; i<=my_pow(steps-1); i++) {
A1 += A;
a1 += a;
}

pi = (A1+a1)*2;

printf("Approssimazione di pi: %.100lf\n",pi);
}

int my_pow(int steps) {
int i,res=1;

for (i=0; i res*=2;

return res;
}

fino a 27 passi di approssimazione funziona....poi, o:
1) mi dice "errore di segmentazione"

2) mi restituisce 0.00000...0

dove sta il problema?

[apatriarca]

La risposta è molto semplice, stai allocando $2^{s-1}*8 = 2^{s+2}$ bytes. Con $s = 28$, vuol dire allocare $2^{30}$ bytes che equivale a circa $1$ GB di memoria*. Inoltre, siccome la tua funzione restituisce un valore intero oltre un certo punto il numero restituito diventerà 0 o negativo e il numero passato a malloc perderà di avere senso. Tutto ciò è comunque del tutto inutile, non hai infatti alcun bisogno di allocare tutta quella memoria.

* E tu ne allochi 3 uguali.. per cui sono in tutto 3 GB...

[deian91]

come faccio per allocare meno memoria?

[apatriarca]

Non hai alcuna necessità di allocare memoria per questo programma. È infatti sufficiente avere accesso agli ultimi due valori calcolati e non tutti.

[apatriarca]

Non avevo guardato bene il tuo programma e pensavo usassi un metodo differente per il calcolo (metodo che non mi è ancora chiaro - che metodo usi per integrare la funzione?). Per farti un esempio ho scritto velocemente il codice per il metodo che avevi chiesto nell'altra discussione:

#include <stdio.h>
#include <math.h>

int main()
{
    unsigned n, i;
    double a = 1.0,
           b = 1.0 / sqrt(2.0),
           t = 0.25,
           p = 1.0,
           s, pi;

    puts("Programma per il calcolo di pi greco usando l'algoritmo di Gauss-Legendre.");
    fputs("Quante iterazioni deve fare il programma? ", stdout);
    scanf("%u", &n);

    for (i = 0; i < n; ++i) {
        double aN = 0.5 * (a + b);
        double bN = sqrt(a * b);
        double da = aN - a;
        double tN = t - p * da * da;
        double pN = 2.0 * p;

        a = aN;
        b = bN;
        t = tN;
        p = pN;
    }

    s = a + b;
    pi = (s * s) / (4.0 * t);

    printf("Il valore di pi greco dopo %u iterazioni e' %.20lf.\n\n", n, pi);
    
    return 0;
}

Come vedi, se si usano i double, è tutto abbastanza semplice e immediato. È inoltre necessaria relativamente poca memoria. Nota che, nonostante il metodo sia in grado di calcolare un numero incredibile di cifre decimali di $pi$ in poche iterazioni, questo non vale per il programma e $20$ cifre sono già oltre la precisione dei double. Il commento nella pagina di wikipedia era legata più che altro alla necessità di utilizzare numero di precisione altissima, che richiedono ovviamente quantità di memoria molto elevate. Anche volendo integrare la funzione $(1)/(1+x^2)$ in $[0,1]$ (moltiplicando poi per $4$ e non per $3$..) non è necessaria tutta questa memoria. Esistono infatti tantissimi metodi di integrazione e la quantità di memoria richiesta è per tutti abbastanza contenuta. Se vuoi posso mostrartene qualcuno.

[deian91]

ma non manca una graffa?

si, mi interessano anche altri metodi...

ma c'è un modo per calcolare con precisione, che so, 100k cifre decimali certe del prigreco attraverso un programma in c?

[apatriarca]

Sì, usando la già citata libreria MPFR, oppure implementando manualmente le funzionalità di tali libreria nel tuo codice (cosa che va oltre le tue capacità attuali).

[deian91]

oltre ciò, perche su windows compilo ed eseguo tranquillamente il tuo codice mentre su ubuntu, tramite gcc non va?

[apatriarca]

Non lo so, dovrebbe essere corretto. Che errore ti restituisce gcc su ubuntu?

[deian91]

/tmp/ccub6J8e.o: In function `main':
casuale.c:(.text+0xaf): undefined reference to `sqrt'
collect2: ld returned 1 exit status

[apatriarca]

Non si capisce niente. Ma credo di aver capito il problema. Aggiungi la seguente opzione quando compili il programma:

-lm

Dice al compilatore di linkare il programma alla libreria matematica.