Segmented Sieve of Eratosthenes

seg_sieve.c

/* Title : Segmented Sieve of Eratosthenes
 *
 * Find all primes between two large numbers and say 'low' and 'high'.
 * We first find all primes in range [2, sqrt(high)+1].
 * Then, we use these for finding all primes in each 'segment sieve'
 * whose size is sqrt(high) + 1.
 * */

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

int main() {
	int low, high;
	scanf("%d %d", &low, &high);

	int prime_limit = (int)sqrt(high) + 1;
	int* prime = (int*) malloc(sizeof(int) * prime_limit);
	bool *is_prime_sqrt = (bool*) malloc(sizeof(bool) * prime_limit);
	
	/* ------------------------------------------------------- */
	/* get primes in range [2, sqrt(high)+1] using simple sieve*/
	for (int i = 2; i < prime_limit; i ++)
		is_prime_sqrt[i] = true;

	for (int p = 2; p < prime_limit; p++) {
		if (is_prime_sqrt[p]) {
			for (int i = p*p; i < prime_limit; i += p)
				is_prime_sqrt[i] = false;
		}
	}

	int prime_size = 0;
	for (int p = 2; p < prime_limit; p++) {
		if (is_prime_sqrt[p])
			prime[prime_size++] = p;
	}
	free(is_prime_sqrt);
	/* End of getting primes */
	/* -------------------------------------------------------- */



	/* -------------------------------------------------------- */
	/* Segmented Sieve of Eratostenes */
	int segment_size = prime_limit; // sqrt(high) + 1
	int temp, prime_set_size = 0;
	if (high == low)
		temp = 1;
	else
		temp = (high-low)/2;
	if (segment_size < 32768) // 2^15, L1 data cache
		segment_size = 32768;
	int temp1 = temp;

	int *prime_set = malloc(temp * sizeof(int));
	bool *is_prime = malloc((segment_size+1) * sizeof(bool));
	int segment_low = low,
		segment_high = low + segment_size;

	for (; segment_low < high;
			segment_low += segment_size,
			segment_high += segment_size) {
		// init
		temp = segment_size+1;
		for (int i = 0; i < temp; i++)
			is_prime[i] = true;

		// 1 is not a prime
		if (segment_low <= 1)
			is_prime[1-segment_low] = false;

		for (int i = 0; i < prime_size; i++) {
			// loLim : the smallest composite number that is bigger than segment_low.
			int loLim = segment_low - segment_low % prime[i] + prime[i];

			if (!(segment_low % prime[i]))
				loLim = segment_low;
			if (loLim == prime[i])
				loLim += prime[i];

			for (int j = loLim; j <= segment_high; j += prime[i])
				is_prime[j-segment_low] = false;
		}
		temp = segment_high+1;
		for (int i = segment_low; i < temp && i <= high; i++) {
			if (is_prime[i-segment_low])
				prime_set[prime_set_size++] = i;
		}
	}
	/* -------------------------------------------------------- */


	for (int i = 0; i < prime_set_size; i++) {
		printf("prime_set[%d] : %d\n", i, prime_set[i]);
	}

	/* deallocation */
	free(is_prime);
	free(prime);
	free(prime_set);

	return 0;
}