fermat.c

/* Eisagogi ston programatismo - 2014-15 */
/* Askisi 1 */

#include <stdio.h>
#include <time.h>

#define MINNUM 1990000001
#define MAXNUM 2000000000
#define MAXTRIES 10

/* H sinartisi isprime elegxei deterministika an o arithmos n pou dinete */
/* san orisma ine protos. Epistrefei 1 gia alithia 0 gia psema */
char isprime(int n) {
	int d;
	/* Arxika elegxei an o arithmos ine polaplasio tou 2 i tou 3 */
	if (n%2==0 || n%3==0)
		return 0;
	/* Epita elegxoume olous tous arithmous tis morfis 6k+1 ke 6k-1 */
	for (d=6;(d-1)*(d-1)<=n;d+=6)
		/* An o arithmos ine polaplasiou aftis tis morfis tote  */
		if ((n%(d-1)==0 || n%(d+1)==0))
			/* epistrefete 0 */
			return 0;
	/* An to programa ftasei edo simenei oti ine prime ara epistrefei 1 */
	return 1;
}

/* Sinartisi ipsosis se dinami opos dinete sto wikipedia */
long long modular_pow(long long base, long long exponent, int modulus) {
	long long result=1;
	while (exponent>0) {
		if (exponent%2==1)
			result=(result*base)%modulus;
		exponent=exponent / 2;
		base=(base*base)%modulus;
	}
	return result;
}

/* Sinartisi elegxou an enas arithmos p ine protos me tin methodo tou fermat */
/* to orisma tries dilonei ton arithmo prospathion tou elegxou */
char fermatcheck(int p, int tries) {
	int i,a;
	/* an o arithmos ine to 1 epistrefoume oti den ine protos */
	if (p==1)
		return 0;
	/* broxos pou dokimazei gia arithmo (tries) fores an o arithmos p */
	/* ikanopoiei to theorima tou fermat */
	for (i=0;i<tries;++i) {
		/* pernoume ena tixeo a me 1<=a<p */
		a=rand()%(p-1)+1;
		/* an a^(p-1) mod p ine diaforetiko apo 1 */
		if (modular_pow(a,p-1,p)!=1)
			/* tote o arithmos den ine protos */
			return 0;
	}
	/* an o arithmos perasei ton elegxo gia tries fores tote ine protos */
	return 1;
}

int main (int argc, char *argv[]) {

	/* EROTIMA 1 - Deterministikos algorithmos elegxou proton arithmon */
	int metritisproton; /* metritisproton einai o arithmos ton proton pou brethikan */
	int n; /* prosorini metabliti pou periexei ton ekastote arithmo gia elegxo */
	clock_t stt; /* katagrafei ti ora ksekinise i epanalipsi gia ton deterministiko algorithmo */
	printf("Checking range [%d,%d] for primes\n\n",MINNUM,MAXNUM);
	/* Gia kathe arithmo apo MINNUM mexri MAXNUM */
	for (stt=((double)clock()),metritisproton=0,n=MINNUM;n<=MAXNUM;++n) {
		/* An ine protos */
		if (isprime(n))
			/* Ayxise ton metriti kata 1 */
			++metritisproton;
	}
	printf("Deterministic algorithm: Found %d primes in %.2f secs\n\n",metritisproton,(double)(clock()-stt)/CLOCKS_PER_SEC);
	
	/* EROTIMA 2 */
	int timeseed=time(NULL); /* Metabliti gia tin genitria tixeon arithmon ( 1414258936 gia epibebeosi orthotitas ) */
	int i; /* prosorini metabliti pou periexei ton ekastote arithmo ton dokimon fermat */
	srand(timeseed);
	printf("Trying Fermat test with seed:%d\n\n",timeseed);
	/* Gia kathe arithmo dokimon apo 1 mexri MAXTRIES */
	for (i=1;i<=MAXTRIES;++i) {
		/* Gia kathe arithmo apo MINNUM mexri MAXNUM */
		for (stt=((double)clock()),metritisproton=0,n=MINNUM;n<=MAXNUM;++n)
			/* An perasei ton elegxo fermat */
			if (fermatcheck(n,i))
				/* Ayxise ton metriti kata 1 */
				++metritisproton;
		printf("Probabilistic algorithm: Found %d primes in %.2f secs (tries = %d)\n",metritisproton,(double)(clock()-stt)/CLOCKS_PER_SEC,i);
	}
}