One of the earliest algorithms ever invented. Eratosthenes found a quick way to find all the primes up to a given number n.

Start with the numbers 2,3,...,n. Repeatedly remove all the multiples of the first number in the list, except that number. The first time round, that's 2, so strike 4,6,8... Then the first number is 3, so you strike 6,9,12,... (Note that 6 is already out, from the first stage -- that's OK). Stop when the square of the first number in the list is greater than n. That's it. You can code this pretty efficiently with various tricks of the trade, but the principle is the same.

Obviously, we only delete numbers from the list which are not prime (they're all non-trivial multiples of some number). On the other hand, any composite number up to n is necessarily deleted at some point: such a number may be written xy, and either x2<=n or y2<=n (otherwise xy>n). So the list we end up with is exactly the primes up to n.

There are two equations which are an algebraic expression of this sieve:

k+1 is prime iff for all natural numbers u,v:

(1)k <> u*v + u + v

Equation (1) is read "k does not equal u times v plus u plus v". The second equation yields only the odd primes:

2*k+1 is prime iff for all natural numbers u,v:

(2)k <> 2*u*v + u + v

One thing to notice. k=2*n satisfies equation (1) for some n iff k=n satisfies equation (2). There are similar equations describing primes in terms of b*k+/-c ("+/-" is to be read "plus or minus") for all naturals b, and each one has an accompanying equation 2*b*k+/-c which describes the primes in the same manner. As b increases, however, more and more work must be done with modulo operations to ensure that the numbers are indeed prime, and that all primes are covered.


The mathematics for the above equations are as follows:

A number n+1 is composite for n in the natural numbers iff there exists a k>1 in the natural numbers such that k divides n and k<n.

A number k divides a number n iff there exists an x in the integers such that k*x = n.

So test n with k which are greater than one. If k is in the set of naturals, then use k+1 for tests. Similarly, use x+1 for x in the naturals.

Thus, a number n+1 is composite iff there exist an x and a k in the naturals such that n+1 = (k+1)*(x+1).

The negation of the above is to say that a number n+1 is prime iff for all x and k in the naturals, n+1 != (k+1)*(x+1), which is the same as saying n+1 != k*x+k+x+1, which is the same as n != k*x+k+x. Similar logic applies for 2*n+1.

Here is a pretty clear explanation of this concept.

Goal: Find all primes which are less than 30.

Step 1: List 2-N

2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

Step 2: Take 2 as prime, and remove all multiples of 2 up to N.

2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

Step 3: Take 3 as prime, and remove all multiples of 3 up to N.

2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 2122 23 24 25 26 27 28 29 30

Step 4: Take 5 as prime, and remove all multiples of 5 up to N.

2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 2122 23 24 25 26 27 28 29 30

Since 5 is the largest integer that is smaller than the square root of N, we can stop here. All number that have not been removed are primes. This leaves us with the following set...

{2, 3, 5, 7, 11, 13, 17, 19, 23, 29}

This algorithm is considered good for small values of N (small in computer terms, so up to 10,000,000 should be ok). You can further enhance the algorithm by skipping all multiples between the current prime and its square. So for instance, in this example when we took 5 as prime, we could skip straight to 25, since we can be sure all other multiples between 5 and 25 (10, 15, 20) would have already been removed. This may not seem like much now, but it can save a bunch of time when your range spans millions.

Here is a simple implementation of the Eratosthenes' Sieve in C. It finds all the primes from 2 to PRIME_MAX - 1 and prints them out.

This implementation takes the smallest known prime that hasn't yet been sieved "currprime" and removes all multiples less than PRIME_MAX from the array. Once the square of currprime is greater than or equal to PRIME_MAX, it stops. As noted above, we can start with the square of currprime, since all lesser multiples will have already been removed.

The program uses an array of chars to hold the primeness of each number from 0 to PRIME_MAX - 1. Each char holds 1 if the index is a prime, or zero otherwise. Since only a bit of storage is needed, storing the primeness in a char is wasteful, but simplifies the code considerably.

#include <stdio.h>

#define PRIME_MAX 100

char sieve[PRIME_MAX];

int nextPrime(int lastprime) {
    int i;

    for (i = lastprime + 1; i < PRIME_MAX; i++)
        if (sieve[i] == 1)
            return i;

    return 0;
}

int main() {
    int i, currprime;

    sieve[0] = 0;
    sieve[1] = 0;
    for (i = 2; i < PRIME_MAX; i++)
        sieve[i] = 1;

    currprime = nextPrime(0);
    while ((currprime * currprime) < PRIME_MAX) {

        for (i = (currprime * currprime); i < PRIME_MAX; i += currprime)
            sieve[i] = 0;

        currprime = nextPrime(currprime);
    }

    for (i = 0; i < PRIME_MAX; i++)
        if (sieve[i] == 1)
            printf("%d\n", i);

    return 0;
}

Log in or register to write something here or to contact authors.