Sieve Of Erathosthenes

A simple implementation of the Sieve of Eratosthenes, a simple, ancient algorithm for finding all prime numbers up to a specified integer. More information can be found on wikipedia, here.

We use a bit array to efficiently store the state of each number. If a number is 0, its a prime. If a number is 1, it isn't. Simple, no? As an example, after the sieve runs, prime[0] will look like this:

Bit Index 7 6 5 4 3 2 1 0
Value 0 1 0 1 0 0 1 1
Decimal 7 6 5 4 3 2 1 0
Meaning Prime Non-prime Prime Non-prime Prime Prime Non-prime Non-prime

The next byte, prime[1], would look like this:

Bit Index 7 6 5 4 3 2 1 0
Value 0 0 0 1 0 1 1 1
Decimal 15 14 13 12 11 10 9 8
Meaning Non-prime Non-prime Prime Non-prime Prime Non-prime Non-prime Non-prime

Benchmark

Benchmarks are done with -O3, average is over 400 iterations, and no output. Printing takes another ~8s at least.

Language CPU Speed Average for 10 000 000
C 1 Ghz 0.175s

This sieve will do up to 1 000 000 000 in under 34 seconds on a decent machine.

Implementation

#include <stdlib.h>
#include <stdio.h>
#include <math.h>
 
#define NumberToFind    10000000  // The maximum number to find
#define LocalAlloc      1         // Allocate storage on the stack
 
#define SetBit(d,c)     ((d) |= (1L << (c))) // Set a bit to 1
#define TestBit(d,c)    ((d) &  (1L << (c))) // Check if a bit is set
 
#if !LocalAlloc
unsigned char *prime = NULL;
#endif
 
int main(void)
{
 
#if LocalAlloc
    unsigned char prime[(NumberToFind / 8)+1] = {0};
#else
    prime = (unsigned char*)malloc((NumberToFind / 8)+1);
#endif
 
    int x,y,t = NumberToFind;
 
    prime[0] |= 3; // Set the first two bits to 1 (00000011) - 0 and 1 are both non-prime.
 
    for(x = 2; x < (int)(sqrt(NumberToFind)+1); x++)
    {
        if(TestBit(prime[x/8],((x % 8))) == 0)
        {
            for(y = x*x; y < NumberToFind; y+=x)
                SetBit(prime[y/8],(y % 8));
        }
    }
 
    // Print out the results
    for(x = 0; x < NumberToFind+1; x++)
    {
        if(TestBit(prime[x/8],(x % 8)) == 0)
            printf("%d\n",x);
    }   
 
    // Wait for the enter key to exit
    fgetc(stdin);
 
#if !LocalAlloc
    free(prime);
#endif
 
    return 0;
}
 
sieve_eratosthenes.txt · Last modified: 2010/02/11 13:54 by TkTech
 
Except where otherwise noted, content on this wiki is licensed under the following license:Public Domain
Recent changes RSS feed Donate Powered by PHP Valid XHTML 1.0 Valid CSS Driven by DokuWiki