Sieve of Eratosthenes

The sieve of Eratosthenes finds all prime numbers up to n. A prime is a number greater than 1 whose only divisors (numbers that divide it with no remainder) are 1 and itself; 1 is not prime. The sieve is for questions about many numbers in one range, up to about 10^7. To check one large number, try its divisors instead, as the primality test does.

Explanation

To find the primes among the numbers from 2 to 12, try the smaller numbers as divisors of each: 9 is divisible by 3, 10 by 2, and 11 by none of them. The primes are 2, 3, 5, 7 and 11.

A row of the numbers 2 to 12. The primes 2, 3, 5, 7 and 11 are green. Above each other number is the smallest prime that divides it: 2 above 4, 6, 8, 10 and 12, and 3 above 9.

Fig. 1. The primes up to 12 are green. Under each other number is a smaller prime that divides it.

Problems ask about far more than twelve numbers: all the primes up to 10^7, for example. The check can be made shorter, because divisors come in pairs: 28 = 4 · 7, and 4 · 4 < 28 < 7 · 7. The smaller divisor of a pair is at most the square root of the number, so trying divisors up to the square root is enough. Even so, checking every number up to 10^7 this way takes 1 746 210 133 divisions. By the usual rule of thumb from time complexity, about 10^8 simple steps per second, that is about 17 seconds.

To see where all these divisions go, write the numbers from 1 to 30 in a table of three rows. Checking each number up to its square root, and stopping at the first divisor that works, takes 45 divisions.

A table of the numbers 1 to 30 in three rows of ten, each with a small count of the trial divisions it takes (45 in all); the numbers 5, 7, 11, 13, 17, 19, 23 and 29, on which every division fails, are tinted red, 29 is outlined, and beside the table 29 is shown with its four divisors tried: 2, 3, 4 and 5.

Fig. 2. Testing each number on its own takes 45 divisions. On the red numbers every division fails.

The square root of 29 is about 5.4, so 29 is divided by 2, 3, 4 and 5, and every division fails. In Fig. 2 the red numbers take 19 of the 45 divisions, and the share grows with n: up to 10^7, primes take 1 365 578 300 of the 1 746 210 133 divisions, almost four in five. Trying only prime divisors would save the division by 4, but 29 would still fail on 2, 3 and 5. The trouble is that each number is tested on its own, so the same question, "does 2 divide it?", is asked again for every odd number and fails every time.

Instead of asking each number for a divisor, let each prime find its own multiples (4, 6 and 8 are multiples of 2). Start at 4, keep adding 2, and cross out every number you get to. That needs no division, and 29 is never crossed out. The smallest number after 2 that is not crossed out is 3. Nothing smaller divides it, since 2 would have crossed it out, so 3 is prime, and it crosses out 6, 9 and 12 the same way.

A row of the numbers 2 to 12. The 2 turns green and an arrow jumps by 2 to 4, 6, 8, 10 and 12, each turning blue and then grey; then 3 turns green and the arrow jumps by 3 to 6 and 12, already grey, and 9, which turns blue and then grey; at the end 2, 3, 5, 7 and 11 are green and the rest grey.

Fig. 3. Each prime crosses out its multiples by adding, not dividing. A grey ring marks a number that was already crossed out.

Every next step works the same way. Take the smallest number that is not crossed out. It is prime: if a smaller number divided it, so would a smaller prime, and that prime would have crossed it out already. Then cross out its multiples, and repeat.

This is the sieve of Eratosthenes.

The numbers 1 to 30 in three rows of ten, stepped through the sieve: 2, 3 and 5 turn green in turn and cross out their multiples from their squares, 7 is reached but 49 is past 30, and at the end the primes 2, 3, 5, 7, 11, 13, 17, 19, 23 and 29 are green and every other number is grey.

1 is not prime. Nothing else is crossed out yet.

1/6

Fig. 4. Step through the sieve with the buttons or the arrow keys. Blue numbers are crossed out in this step, and green numbers are primes.

The run found all 10 primes up to 30 with 24 crossings instead of 45 divisions. A prime starts crossing out at its square, because every smaller multiple of it is divisible by a smaller prime and is already crossed out. And once the square of the next prime is bigger than n, nothing is left to cross out, so every number not crossed out is prime. For all numbers up to 10^7 the sieve makes 22 850 051 crossings, about 0.2 seconds by the rule of thumb, against 1.7 \cdot 10^9 divisions.

The prime that crosses a number out first is the smallest prime that divides it, its smallest prime factor: 30 was first crossed out by 2, and 25 by 5.

Counting primes in a range

Many problems ask "how many primes are there between l and r?" again and again, for different l and r. Run the sieve once, up to the largest r. Then write next to each number how many primes are at most that number. These prefix counts are a prefix sum of the sieve. The number of primes from l to r is the prefix count at r minus the prefix count at l − 1. For the primes from 15 to 30, the prefix count at 30 is 10 and at 14 it is 6, so there are 10 − 6 = 4 of them: 17, 19, 23 and 29.

Two rows of cells: the numbers 14 to 30 on top, with the primes 17, 19, 23 and 29 in green, and under each number the count of primes up to it, from 6 under 14 to 10 under 30, with those two counts outlined in blue.

Fig. 5. There are 10 − 6 = 4 primes from 15 to 30.

Implementation

In code, the table becomes an array is_prime of n + 1 flags, a vector<bool>: is_prime[k] is true while k is not crossed out, and crossing k out sets it to false. The numbers 0 and 1 are not prime, so they are set to false before the sieve starts. The outer loop goes over i from 2, and the inner loop crosses out the multiples of i with j.

  1. Take the next prime: if is_prime[i] is still true, then i is prime.

  2. Cross out its multiples: set is_prime[j] = false for j = i * i, i * i + i, i * i + 2 * i and so on, up to n.

  3. Stop at the square root: the outer loop runs while i * i <= n.

  4. Collect the primes: one pass over is_prime from 2 to n.

The last step is a separate loop, because the outer loop stops at the square root of n.

The program reads n and prints the primes up to n, one per line. Then it reads the number of queries q and q pairs l r with 1 <= l <= r <= n, and prints the number of primes from l to r for each, using the prefix counts cnt[k], the number of primes up to k. For the sieve up to 30 in Fig. 4, with the queries from 15 to 30 and from 1 to 10, the input is:

30
2
15 30
1 10
#include <iostream>
#include <vector>
using namespace std;

int main() {
    int n;
    cin >> n;
    vector<bool> is_prime(n + 1, true);
    is_prime[0] = false;
    if (n >= 1) is_prime[1] = false;   // 1 is not prime

    for (int i = 2; i * i <= n; i++) {     // stop at the square root
        if (is_prime[i]) {                 // take the next prime
            // cross out its multiples
            for (int j = i * i; j <= n; j += i) {
                is_prime[j] = false;
            }
        }
    }

    // collect the primes
    for (int i = 2; i <= n; i++) {
        if (is_prime[i]) cout << i << '\n';
    }

    // count the primes from l to r; remove this block if you only need the list
    vector<int> cnt(n + 1, 0);
    for (int k = 1; k <= n; k++) {
        cnt[k] = cnt[k - 1] + (is_prime[k] ? 1 : 0);
    }
    int q;
    cin >> q;
    while (q--) {
        int l, r;                          // 1 <= l <= r <= n
        cin >> l >> r;
        cout << cnt[r] - cnt[l - 1] << '\n';
    }
}

The sieve takes O(n \log \log n) time: each prime i crosses out about n / i numbers, and the sum n/2 + n/3 + n/5 + \dots over the primes up to is about n \ln \ln n. In practice it is a little more than crossings. A vector<bool> keeps each flag in one bit, so 10^7 flags take about 1.25 MB, while cnt takes 40 MB.

Correctness proof

Two facts make the sieve right: it never crosses out a prime, and it crosses out every number up to that is not prime.

A number crossed out by is j = p \cdot k with k \ge p \ge 2. It has the divisor , which is neither 1 nor , so is not prime. So a prime is never crossed out, and every number the outer loop takes as prime really is one.

Now take a number m \le n that is not prime, and let be its smallest prime factor, so m = p \cdot q. Since is not prime, q > 1, so has a prime divisor. That prime also divides , so it is at least , and q \ge p. So p \cdot p \le m \le n, and the outer loop reaches . It takes as a prime, because primes are never crossed out. Its inner loop starts at p \cdot p \le m and goes up in steps of , and divides , so it lands on and crosses it out.

Common mistakes

  • Counting 1 as a prime. If is_prime[1] is left true, then cnt[1] is 1, so the count of primes from 1 to 10 comes out as 5 instead of 4. A loop that prints the primes from a = 1, as problem 4739 asks, prints 1 2 3 5 7. Set is_prime[1] to false.

  • Making the array one cell too short. vector<bool> is_prime(n) has no cell n, so reading is_prime[n] goes past the end. For n = 29 such a program printed the primes without 29. Allocate n + 1 cells.

  • Stopping the outer loop one step early. With i * i < n, for n = 25 the loop never takes i = 5, and 25 is printed as a prime (49 for n = 49, too). Write i * i <= n.

  • Running the sieve again for every query. With 10^6 queries and numbers up to 3 \cdot 10^6, that is about 6.6 \cdot 10^{12} crossings instead of 6.6 \cdot 10^6. Sieve once, up to the largest number any query needs.

Similar algorithms

The primality test checks single numbers up to about 10^{12} and beyond.

Prime factorization stores each number's smallest prime factor instead of a flag, then factors any number up to n in a few steps.

The segmented sieve finds the primes in a window from m to n, with n up to about 10^9 and a window of up to about 10^6 numbers. It sieves the primes up to the square root of n, then crosses out their multiples inside the window only.

The linear sieve crosses out each number that is not prime exactly once, in , and gives smallest prime factors along the way. For a plain list of primes it is no faster in practice.

Practice

Чи була ця сторінка корисною?