Generating prime numbers

This document starts by presenting a new, simple and efficient algorithm to generate prime numbers progressively. The algorithm turns out to be a non-sieve version of the sieve of Eratosthenes. It is based on generating all integers progressively starting from 2, taking O(n log n) time and O(sqrt(n)) space, where n is the largest integer generated so far. An improved version of the algorithm is then presented, which turns out to be a non-sieve version of the sieve of Euler.

The proof of the algorithm is then presented. This is where the very simple but fundamental nature of the construction of all integers using prime numbers as basic building blocks is encountered. This fundamental nature is the core idea behind everything presented in this document. For that reason it is explained to its finest details with figures for illustration. It is absolutely important to understand and agree with it. This section also mentions why prime numbers seem to occur at random whereas that is just due to an ordered but accumulated effect. With the fundamental nature haven been explained, a number of sections then follow:

First are proofs of the following (already proven) theorems:

Fundamental theorem of arithmetic. This is the most obvious as it follows directly from the fundamental nature explained.
Euclid's theorem: there is an infinite number of prime numbers.
Bertrand's postulate: there is always at least one prime between any integer n>1 and its double 2n.

Next are proofs of the following conjectures:

Twin prime conjecture: there is an infinite number of primes separated by a difference of 2.
Polignac conjecture: there is an infinite number of primes separated by a difference of 2k where k is any integer.
For any n^th prime p_n there is at least one prime between p_n² and p_np_n+1.

Lastly are detailed analyses on the following conjectures:

There are as many twin primes as there are cousin primes (primes separated by a difference of 4).
The number of twin primes (or cousin primes) between any prime p and its square p² increases with p.

... to be updated ... (the techniques are currently being used to tackle other related problems!)

Check out the math fullfloor function

Implemented algorithm

/*
	Here is the source code of the algorithm. To compile with GCC, execute:
		gcc -Wall -std=c99 -pedantic -O2 generate.c -o generate

	To generate all primes in the interval from 10^9-100 to 10^9 execute:
		time ./generate 999999900 1000000000 1 1 ", "

	Equivalently, to generate primes from index 50847533 to 50847534 execute:
		time ./generate 50847533 50847534 0 1 ", "

	This takes 1 MB of memory and 32 seconds on a 3.2 GHz processor.

	Enable print_heap() in order to visualise the algorithm in execution.
*/

#include <stdio.h>
#include <stdbool.h>

#define INT unsigned long long


void generating_prime_numbers (bool (*submit_prime)(INT i, INT p));

// The function generating_prime_numbers() will
// call the function submit_prime(), for every
// generated prime starting from i=1 with p=2.
// The i-th prime is p. i is the index of p.


INT start, stop; // the range of the printed result
bool forPnotI;   // range is for 'p' else is for 'i'
bool showIndex;  // additionally print the index 'i'
const char* delimiter;


bool submit_prime (INT i, INT p)
{
	INT n = forPnotI ? p : i;
	if(n > stop) return true;   // stop generating
	if(n < start) return false; // continue generating
	if(n > start) printf("%s",delimiter);
	if(showIndex) printf("%llu ", i);
	printf("%llu", p);
	return false; // continue generating
}


int main (int argc, char** argv)
{
	do{
		if(argc!=6){
			printf(
				"\n    Use: <program> <start> <stop>"
				" <forPnotI> <showIndex> <delimiter>\n"
				"\n    Example: <program> 1 100 0 0 \", \""
				"\n\n");
			break;
		}
		start = stop = 0;
		int i;
		sscanf(argv[1], "%llu", &start);
		sscanf(argv[2], "%llu", &stop);
		sscanf(argv[3], "%d", &i); forPnotI = i;
		sscanf(argv[4], "%d", &i); showIndex = i;
		delimiter = argv[5];

		if(start < 1){
			printf("\nInvalid input: (start = %llu) < 1.\n\n", start);
			break;
		}
		if(start > stop){
			printf("\nInvalid input: (start = %llu) > (stop = %llu).\n\n", start, stop);
			break;
		}
		generating_prime_numbers(submit_prime);
		puts("");
	}while(false);
	return 0;
}


/******************************************************************************/


#include <malloc.h>
#include <string.h>

typedef struct _Prime {
	INT pm; // product p*m
	INT m;  // prime multiplier
	int p;  // prime of this node
	int i;  // wheel-index of m
} Prime;

static inline int compare (Prime p1, Prime p2)
{
	if(p1.pm < p2.pm) return -1;
	if(p1.pm > p2.pm) return +1;
	return (p1.p < p2.p) ? -1 : +1;
}

void print_heap (const Prime* heap, int size, int l, int indent)
{
	if(l>=size) return;
	print_heap (heap, size, l*2+2, indent+3); // right subtree
	for(int i=0; i < indent; i++) printf("    ");
	printf("%d*%llu=%llu\n", heap[l].p, heap[l].m, heap[l].pm);
	print_heap (heap, size, l*2+1, indent+3); // left subtree
}


void generating_prime_numbers (bool (*submit_prime)(INT i, INT p))
{
	int i, j;
	INT c = 0;

	/* submit first few primes */
	if(submit_prime(++c, 2)) return;
	if(submit_prime(++c, 3)) return;
	if(submit_prime(++c, 5)) return;

	const int wheel[] = {2, 6, 4, 2, 4, 2, 4, 6};
	const int wheelSize = sizeof(wheel)/sizeof(*wheel);

	/* get first prime used by the algorithm */
	int k = 1;             // k is the wheel-index of n
	INT n = wheel[k]+1;    // n is the latest integer generated
	INT L = n;             // L is the largest prime in the heap
	if(submit_prime(++c, n)) return;

	/* initialise heap with the first prime */
	int heapMax = 1<<14;
	Prime* heap = (Prime*) malloc (heapMax*sizeof(*heap));
	Prime root = { .pm = L*L, .m = L, .p = L, .i = k};
	heap[0] = root;
	int heapSize = 1;

	while(true){
		root = heap[0];

		/* get the next integer 'n' using its wheel-index 'k' */
		if(++k >= wheelSize) k=0;
		n += wheel[k];

		/* if n turns out to be a prime */
		if(n != root.pm){
			if(submit_prime(++c, n)) break;
			else continue;
		}

		/* print the heap data struture */
		if(false){
			print_heap(heap, heapSize, 0, 1);
			printf("press enter...");
			getchar();
		}

		/* if a new L to be added to the heap is identified */
		if(root.p == L && root.m != L)
		{
			/* get the new largest prime L of heap */
			L = root.m;
			root.p = L;
			root.pm = L*L;

			/* add at end of heap, then update heap */
			i = heapSize++;
			while(true)
			{
				j = (i-1)/2;
				if(compare(heap[j], root)<0) break;
				heap[i] = heap[j];
				i = j;
			}
			heap[i] = root;

			/* if heap memory is full then increase it */
			if(heapSize >= heapMax)
			{
				heapMax *= 2;
				heap = realloc(heap, heapMax*sizeof(*heap));
				//printf("heapMax = 0x%x\n", heapMax);
			}
		}

		while(true){
			root = heap[0];
			if(n != root.pm) break;

			if(++root.i >= wheelSize)
				root.i = 0;
			root.m += wheel[root.i];
			root.pm = root.p * root.m;

			/* update the heap data structure */
			i=0;
			while(i*2+2 <= heapSize)
			{
				Prime left  = heap[i*2+1];
				Prime right = heap[i*2+2];

				if(i*2+2 == heapSize
				|| compare(left, right)<0)
				{
					if(compare(left, root)<0)
					{
						heap[i] = left;
						i = i*2+1;
						continue;
					}
				}
				else{
					if(compare(right, root)<0)
					{
						heap[i] = right;
						i = i*2+2;
						continue;
					}
				}
				break;
			}
			heap[i] = root;
		}
	}
	free(heap);
}