/****************************************************************************/
/* Produce el n-simo número de Fermat e intenta factorizarlo.               */
/* Los factores de Fn=2^(2^n)+1 son de la forma 2^(n+2)*k + 1. Usa LIP.     */
/*                                                                          */
/* Jaime Suarez <mcripto@bigfoot.com> 2003                                  */
/* en http://elparaiso.mat.uned.es                                            */
/****************************************************************************/

#include <stdio.h>
#include "lip.h"

main(int argc, char *argv[])
{
	verylong Fn=0,dos=0,tmp=0,tmp1=0,base=0,div=0,k=0;
	long n;

	if (argc!=2) {
		printf("%s <n> calcula el n-simo número de Fermat\n",argv[0]);
		printf("       e intenta encontrar factores.\n");
		return 1;
	}
	n=atol(argv[1]);
	zintoz(2, &dos);                                /* Calcula 2^(2^n)+1 */
	zsexp(dos,n,&tmp);
	zexp (dos,tmp,&tmp1);
	zsadd(tmp1,1,&Fn);
	printf("F%ld = ",n);                            /* Escribir Fn */
	zwriteln(Fn);
	zsexp(dos,n+2,&base);                           /* base = 2^(n+2) */
	zone(&k);                                       /* k = 1 */
	for (;;) {
		zmul(base,k,&tmp);                          /* div = base*k + 1 */
		zsadd(tmp,1,&div);
		zsq(div,&tmp);								/*  Es div^2 > Fn ? */
		zsub(tmp,Fn, &tmp1);
		if (zsign(tmp1)>0) {                        /* Divisor muy grande */
			printf("El número es primo.\n");
			break;
		}
		zdiv(Fn,div,&tmp,&tmp1);					/*  División exacta ? */
		if (ziszero(tmp1)) {
			printf("Un factor es ");
			zwriteln(div);
			break;
		}
		zsadd(k,1,&tmp);                            /* k = k + 1 */
		zcopy(tmp,&k);

	} 
	return 0;
}