Sieve of Eratosthenes (C Plus Plus)

From LiteratePrograms

Other implementations: Alice ML | Bash | C++ | Forth | Haskell | Java | Python | Python, arrays | Scala

The Sieve of Eratosthenes is an algorithm for rapidly locating all the prime numbers in a certain range of integers. It operates by marking as composite all nontrivial multiples of each prime in sequence until only primes remain unmarked. It is most well-suited to implementation using arrays, which are compact and feature random access.

Contents 1 Design 2 Main implementation 3 Files 4 Test driver 5 Performance

Design

We have three main portable choices for representing the sieve supplied by the Standard Template Library:

• bitset: provides a fixed size bit array whose size is fixed at compile-time. This is not suitable for a general algorithm because the integer range being searched must be fixed at compile time.
• : each byte is set to zero or one to indicate prime or composite, respectively. On a byte-addressed machine no address computation or bit manipulation is required. However, this uses 8 bits per integer.
• : a template specialization that uses only one bit per entry. This requires slightly more computation per access for address computation and bit manipulation, but is much more compact, and utilizes the data cache more efficiently. For larger ranges of integers, we expect memory access time to dominate computation time, making this the ideal solution.

Main implementation

To save additional space and time, at the cost of some additional complexity, we choose not to represent even integers in our sieve. Instead, the element at index k will indicate the primarily of the number 2k + 3. We abstract this complexity away by wrapping vector<bool> in a class exposing special indexer methods:

<<Sieve class>>=
class Sieve
{
private:
    std::vector<bool> sieve;
public:
    Sieve(size_t size);
    bool is_composite(size_t k);
    void set_composite(size_t k);
};
inline Sieve::Sieve(size_t size)
  : sieve((size+1)/2)
{ }
inline bool Sieve::is_composite(size_t k)
{
    assert(k >= 3 && (k % 2) == 1);
    return sieve[(k-3)/2];
}
inline void Sieve::set_composite(size_t k)
{
    assert(k >= 3 && (k % 2) == 1);
    sieve[(k-3)/2] = true;
}

We rely on inlining and compiler optimization to eliminate the additional computation in this method. Standard headers are required for vector, size_t, and assert:

<<header files>>=
#include <vector>
#include <cstdlib> // size_t
#include <cassert>

The sieve itself begins by creating a sieve array and then searching for the first zero entry. All odd multiples of this entry are set to true to indicate that they are composite. We start with the square i*i of i, because all multiples of i smaller than that will have already been set to true while dealing with smaller primes. Note that by default a vector<bool> is initialized to all zeros, as the algorithm requires.

<<sieve of Eratosthenes main body>>=
Sieve sieve(max + 1 /* +1 to include max itself*/);
size_t i;
for(i=3; i*i <= max; i += 2)
{
    if (sieve.is_composite(i))
    {
        continue;
    }
    // We increment by 2*i to skip even multiples of i
    for(size_t multiple_i=i*i; multiple_i <= max; multiple_i += 2*i)
    {
        sieve.set_composite(multiple_i);
    }
}

When we're done, we iterate through the sieve, forming a std::vector of elements still set to false, which we return:

<<sieve of Eratosthenes prototype>>=
std::vector<size_t> sieve_of_eratosthenes(size_t max);
<<sieve of Eratosthenes>>=
std::vector<size_t> sieve_of_eratosthenes(size_t max)
{
    sieve of Eratosthenes main body
    std::vector<size_t> primes;
    primes.push_back(2);
    for(size_t i=3; i <= max; i += 2)
    {
        if (!sieve.is_composite(i))
        {
            primes.push_back(i);
        }
    }
    return primes;
}

Files

We place the sieve class and the prototype for the sieve of Eratosthenes in a header file:

<<Sieve.h>>=
#ifndef _SIEVE_H_
#define _SIEVE_H_
header files
Sieve class
sieve of Eratosthenes prototype
#endif // #ifndef _SIEVE_H_

We place the implementation of the sieving in the corresponding source file:

<<Sieve.cpp>>=
#include "Sieve.h"
sieve of Eratosthenes

Test driver

We write a simple test driver that takes a maximum value on the command-line, determines the time required for the computation in total and per integer, and writes out the resulting list of primes. We offer an optional second parameter indicating the number of times to repeat the sieve, for better precision in the time measurement.

<<test.cpp>>=
#include <cstdlib>  // For atoi()
#include <ctime>    // For clock(), clock_t
#include <iostream>
#include "Sieve.h"
using std::cout;
using std::endl;
int main(int argc, char* argv[])
{
    size_t max = (size_t)atoi(argv[1]);
    int num_times = 1;
    if (argc > 2)
    {
        num_times = atoi(argv[2]);
    }
    std::vector<size_t> result;
    clock_t start_time = clock();
    for(int i=0; i < num_times; i++)
    {
        result = sieve_of_eratosthenes(max);
    }
    double time_in_ms = ((double)(clock() - start_time)) / CLOCKS_PER_SEC / num_times * 1000;
    double time_per_integer_ns = time_in_ms / max * 1000000;
    cout << "Sieved over integers 1 to " << max
         << " in " << time_in_ms << " ms (" << time_per_integer_ns << " ns per integer)"
         << endl << endl;
    for(size_t i=0; i < result.size(); i++)
    {
        cout << result[i] << endl;
    }
    return 0;
}

Performance

We compiled the implementation with gcc with the "-O3" flag and ran it for various input sizes on a machine with a Pentium 4 2.40 GhZ CPU with 512KB of cache and a 2GB of RAM. In order to determine times accurately for small integers, a larger number of times was supplied. Times follow:

Input size	Sieving time (ms)	Sieving time per integer (ns)
10	0.000887	88.7
100	0.00252	25.2
1 000	0.0130	13.0
10 000	0.102	10.2
100 000	0.98	9.8
1 000 000	12.15	12.2
10 000 000	216.5	21.7
100 000 000	3900	39.0

Initially, the time per element is relatively high due to the overhead of setting up the sieve and producing the result list. Optimal performance is reached around 100 000 values. After this it begins to go up again due to nonlinear algorithmic complexity and a larger number of cache misses. A more efficient approach to the larger sieves might be a segmenting scheme that performs all the necessary sieving over one block of numbers before moving on to the next, making better use of the cache.

Download code