/*
The algorithm is very simple:
For a given n, the numbers <= n is not a prime if it can be factor by a 
number, say p, and p must be less than or equal to sqrt of n.
*/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <mpi.h>

#define n 1000000

int main(int argc, char *argv[])
{
	long i, j, xs, ps, ms;
	double clock0, clock1;
	float time;
	long m, s, count;
	long *p, *x;

	MPI_Init(&argc, &argv);
	clock0=MPI_Wtime();

	m = floor(sqrt(n));
	s = floor((m-3)/2) + 1;

	p = calloc(s, sizeof(long));

	ps = s;

	for (i = 0; i < s; i ++)
		p[i] = 0;

	j = 0;
	for (i = 3; i <= m; i += 2)
	{
		p[j] = i;
		j ++;
	}

	ms = p[ps - 1];

	s = floor((n-ms-2)/2)+1;
	x = calloc(s, sizeof(long));

	for (i = 0; i < s; i ++)
		x[i] = 0;

	j = 0;
	for (i = ms + 2; i <= n; i += 2)
	{
		x[j] = i;
		j ++;
	}

	xs = s;

	for (i = 0; i < ps - 1; i++)
	{
		for (j = i + 1; j < ps; j++)
		{
			if (p[i] != 0 && p[j] % p[i] == 0)
			{
				p[j] = 0;
			}
		}
	}

	for (i = 0; i < ps; i++)
	{
		for (j = 0; j < xs; j++)
		{
			if (p[i] != 0 && x[j] % p[i] == 0)
			{
				x[j] = 0;
			}
		}
	}
	clock1=MPI_Wtime();

	count = 0;
	for (i = 0; i < ps; i ++)
	{
		if (p[i] != 0)
		{
/*			printf("hi p %d\n", p[i]);
*/			count++;
		}
	}
	for (i = 0; i < s; i ++)
	{
		if (x[i] != 0)
		{
/*
			printf("hi x %d\n", x[i]);
*/			count++;
		}
	}

	printf("n :%d\n", n);
	printf("no. of prime :%d\n", count + 1);
	printf("time used: %e\n", clock1 - clock0);
	MPI_Finalize();

	free(p);
	free(x);
	return(0);
}