/*********************************************************
Languasco & Zaccagnini, Introduzione alla Crittografia, Hoepli.
Capitolo 6.4. pag. 166 e
Capitolo 9.6.3 pag. 276
VERSIONE AGGIORNATA (FILE AKS-CORR IN ERRATA)
*********************************************************/

{ /******** Funazione ausiliaria che verifica che l'ordine di ******* 
   *******       n modulo r sia piu' grande di ord       *************/ 
order(n, r, ord) = local (k, flag);
    flag=1;
    k=1;
    until ( k > ord, 
    			flag = flag && (lift(Mod(n,r)^k) !=1); 
    			if (flag == 0, return (flag));
    			k=k+1;
    );   
    /*******  flag = 0 significa che non e` stato raggiunto  *******
    *******  un ordine di n mod r maggiore di ord; flag =1  *******
    *******  significa che e` stato raggiunto  			*******/
    return(flag)                             
}

{aks(n)= local(b, r, good, ord, logn,log2n, l, m, k);

/***************** Controlli iniziali *******************/ 
if((n <= 3)&(n>1), print(n, " e` primo"); return(1));  
if( n<=1, error("Input non valido"));
if(!bittest(n,0), print(n, " e` pari"); return(0));


/********** Test per le potenze di primi dispari ********/ 

if(issquare(n), print(n," e` il quadrato di ",floor(sqrt(n))); 
			return(0));
logn=log(n);
l= ceil(exp(logn/ 3)); 
forprime(p=3,l,  
     m = ceil(logn/log(p)); 	
     forstep(k=3, m,2,  
           	if(n==p^k, print(n," e` uguale a ",p, " elevato a ",k);
           		return(0))));
print(n," non e` una potenza prima"); 

/**********   Ricerca del minimo r per cui n ha ordine  **********
 **********	  maggiore di log_2^2 n in Zr           ************/ 

log2n = logn/log(2);

ord= floor((log2n)^2);

/**********   se n deve avere almeno ordine ord allora   **********
 **********    phi(r) >= ord e quindi r>=phi(r)+1 >= ord+1 ************/ 

r = ord+1;

/* print("ord ", ord); */

/**********   n deve essere un elemento invertibile di Zr  ************/ 
until (gcd(n,r) == 1, r=r+1);

good = order(n,r, ord);

/**********   good=1 significa che r e` stato determinato  ************/ 

until( good ,  until ((gcd(n,r) == 1), r=r+1); 
             good = order(n,r,ord)
           ); 
     	   
print("r ", r); 
/* print("good ", good); */
      	    

/**********   condizioni per n <= r o con divisori <= r  ************/       	    
for( b=1, r,
     if( gcd(b,n) > 1 &&  gcd(b,n) < n, 
        	print(n," e` composto ed e` diviso da ",gcd(b,n)); return(0) )
   );

if (n <= r,    print(n," e` primo"); return(1) );
   
/************* Test polinomiale di AKS  ***************/ 

s=floor(sqrt(r)*log2n); 
print("inizio verifiche condizioni polinomiali di aks");
for( b=1, s,
     if( Mod(x + Mod(b,n),x^r-1)^n != Mod(x^n + Mod(b,n),x^r-1),
            print(n," e` composto"); return(0))
   );
print(n," e` primo");
return(1);
}