Sieve of Eratosthenes in C

The Sieve of Eratosthenes is one of the simplest and easy to understand algorithms to generate prime numbers. It is so simple, so even my Expr programming language can implement it.

let numbers = 2..1000; 
reduce(numbers, { 
    let n = #; 
    filter(#acc, # == n || # % n != 0) 
}, numbers)

But here is an implementation of the sieve in C:

#include <stdbool.h>
#include <stddef.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

static const size_t LIMIT = (size_t)1000000000u;

static inline bool bit_get(const uint8_t *bits, size_t i) {
    return (bits[i >> 3] >> (i & 7u)) & 1u;
}

static inline void bit_clear(uint8_t *bits, size_t i) {
    bits[i >> 3] &= (uint8_t)~(uint8_t)(1u << (i & 7u));
}

int main(void) {
    const size_t bytes = (LIMIT + 7u) >> 3;

    uint8_t *primes = malloc(bytes);
    if (!primes) {
        perror("malloc");
        return EXIT_FAILURE;
    }

    memset(primes, 0xFF, bytes);

    bit_clear(primes, 0);
    bit_clear(primes, 1);

    for (size_t p = 2; p <= LIMIT / p; ++p) {
        if (!bit_get(primes, p)) continue;

        for (size_t m = p * p; m < LIMIT; m += p) {
            bit_clear(primes, m);
        }
    }

    for (size_t i = 2; i < LIMIT; ++i) {
        if (bit_get(primes, i)) {
            printf("%zu\n", i);
        }
    }

    free(primes);
    return EXIT_SUCCESS;
}

Compile and run:

gcc -o sieve sieve.c -O3

The Sieve of Eratosthenes works as follows:

Assume all numbers ≥ 2 are prime.
Starting from the smallest prime p = 2, eliminate all multiples of p.
Move to the next number that is still marked prime.
Stop once p² > LIMIT.

The remaining marked numbers are primes.

Memory Representation: Bitset

Instead of storing a bool or char per number, this implementation uses one bit per number.

static const size_t LIMIT = (size_t)1000000000u;
const size_t bytes = (LIMIT + 7u) >> 3;

Each bit represents whether a number is prime.
LIMIT numbers require LIMIT / 8 bytes.
For LIMIT = 1,000,000,000, this is about 125 MB of memory.

This choice drastically reduces memory usage and improves cache locality compared to byte-based arrays.

Bit Access Helpers

Two inline functions abstract bit manipulation:

static inline bool bit_get(const uint8_t *bits, size_t i) {
    return (bits[i >> 3] >> (i & 7u)) & 1u;
}

static inline void bit_clear(uint8_t *bits, size_t i) {
    bits[i >> 3] &= (uint8_t)~(uint8_t)(1u << (i & 7u));
}

How this works:

i >> 3 selects the byte containing bit i (same as i / 8)
i & 7 selects the bit position within that byte (same as i % 8)
bit_get tests whether the bit is set
bit_clear marks a number as non-prime

Using static inline allows the compiler to eliminate function call overhead.

Initialization

uint8_t *primes = malloc(bytes);
memset(primes, 0xFF, bytes);

All bits are initially set to 1 (assume prime).

The numbers 0 and 1 are then explicitly cleared:

bit_clear(primes, 0);
bit_clear(primes, 1);

Sieving Loop

for (size_t p = 2; p <= LIMIT / p; ++p) {
    if (!bit_get(primes, p)) continue;

    for (size_t m = p * p; m < LIMIT; m += p) {
        bit_clear(primes, m);
    }
}

Loop bound of p <= LIMIT / p is equivalent to p ≤ √LIMIT, but avoids floating-point math and overflow. Start from p², all smaller multiples of p were already handled by smaller primes.

Output

for (size_t i = 2; i < LIMIT; ++i) {
    if (bit_get(primes, i)) {
        printf("%zu\n", i);
    }
}

After sieving completes, all remaining set bits correspond to prime numbers.

Note that printing primes will take orders of magnitude longer than computing them.