domenica 27 dicembre 2009

Algoritmo sul gaussiano di una frazione, primalità e compositezza


In questo blog ci poniamo il problema di scrivere un algoritmo in PARI/GP per determinare il periodo di una frazione.

Nel blog precedente abbiamo detto che: 
1 - Si definisce gaussiano del numero N rispetto ad una base B, indicato come k = g(N) = T(1/N), “il più piccolo valore k tale che B^k = 1 mod N”.
2 - Si dimostra che il gaussiano k di un rapporto a/b dipende solo dal denominatore b, quindi nel calcolo del gaussiano si può porre a=1 e calcolare k = T(1/N), dove t è l'algoritmo da realizzare.
3 - Un teorema che deriva dal Piccolo Teorema di Fermat porta a dire anche che il gaussiano di n è un divisore della funzione di Eulero φ(n). In particolare se n è primo è φ(n) = n – 1.
Se dobbiamo scrivere degli algoritmi di primalità o di compositezza basati sulle frazioni continue occorre scrivere prima un algoritmo per il gaussiano.
/*
*  Period of a fraction with
*  gaussian function
*/
GaussN(n)= local(g=0, r=1); {
 while( Mod(n,2) == 0, n = n/2;);
 while( Mod(n,5) == 0, n = n/5;);
 while( g == 0 | r != 1,
        r = Mod(( 10 * r), n);
        g = g + 1;
 );
 return(g); 
}
 
Considerazioni sull'algoritmo
Scrivere un algoritmo che calcoli il gaussiano comunque è semplice, ma sicuramente il metodo non è veloce per numeri grandi (basta applicarlo ad esempio al numero di Mersenne 2^31-1)
In PARI/GP esiste anche un ulteriore metodo semplice da collocare comunque nell’ambito di un algoritmo; ad esempio sappiamo che sqrt(21)= 4,1,1,2,1,8 dove la parte sottolineata è il periodo T(sqrt(21)).Se usiamo contfrac(sqrt(21)) otteniamo il vettore [4,1,1,2,1,8,1,1,2,1,8,1,…].  Ora se a1 è il primo elemento del vettore e aj è il j-esimo elemento del vettore, allora il periodo T(sqrt(21))=j-1 se aj=2*a1. Ovviamente il metodo è poco applicabile se a1=0, ovvero la frazione non è periodica ma è finita o infinita (numero irrazionale).L’algoritmo presentato per il calcolo del gaussiano non usa però contfrac. La scelta dipende dal fatto che contfrac non sempre fornisce il risultato in modo facilmente trattabile nei casi di rapporti 1/N. Contfrac è più efficace nei casi di calcoli delle radici quadrate di un numero N.
Algoritmo di compositezza  (basato sul Teorema del blog precedente)
/*
*  IsComposite ?
*  -1 Boh
*   1 Yes
*   0 No (it's a prime number)
*/
isCompN(n)= local(g=0, status=-1); {
 if( Mod(n,2) == 0, return(status););
 g = GaussN(n);
 print("gaussiano(",n,") = ", g);
 if( g == n-1,
     status=0;
     print("No, it's a prime number!");
     return(status);
 );
 if( g != 1,
     if( frac((n-1)/g) != 0.0, status = 1;);
 );     
 if( status == 1, print("Yes, it's a composite number!"););
 return(status); 
}
/*
*  IsMyPrime ?
*   1 Yes
*   0 No
*/
isMyPrime(n)= local(status=0); {
 if( eulerphi(n) == n-1, status=1;);
 return(status); 
}
Se si guarda l’algoritmo isComposite ci si rende conto che, per numeri grandi, già ad esempio per un numero di Mersenne come 2^31-1, poiché il gaussiano è calcolato per successivi incrementi unitari ne consegue che sono necessari molti cicli. In tal caso è più veloce un algoritmo come isMyPrime o le funzionalità built-in di PARI/GP come eulerphi, isprime, ispseudoprime. Queste ultime funzioni built-in soffrono ovviamente del fatto che dipendono dal numero di primi precaricati e dalla RAM disponibile e allocata al programma.
Alla prox

0 commenti:

Posta un commento