The versatile sieve of Eratosthenes
Introduction
The sieve of Eratosthenes is an efficient algorithm to mark whether each integer between 0 to n is prime. It turns out that the algorithm is useful for calculating other interesting values that are also based on the prime factorization of a number.
This page will show how the basic sieve of Eratosthenes can be modified to calculate other numbertheoretic functions. A prior understanding of the basic algorithm will be assumed (please browse other resources for an introduction).
Designing and applying these modified forms of the sieve of Eratosthenes has been very helpful in a number of my Project Euler solution implementations.
Technical notes
The sieve of Eratosthenes (and its variants) takes an upper limit integer n, and computes a table of some particular numbertheoretic function f(k) for every value of k from 0 to n. Hence if you wish to evaluate f(k) for many values of k between 0 and n (i.e. the usage is dense), then it helps to run a sieve first to precompute all the values in that domain. If the queries of k are sparse however (e.g. only testing if a couple of numbers are prime), then sieving is not efficient in time or memory.
All code examples are given in the Python programming language, for clarity and to avoid issues with integer overflow. These algorithms can certainly be implemented in languages with fixedwidth integer types (such as C/C++/C#/Java/JavaScript), but extreme care must be taken to prove that every calculation does not overflow the chosen integer type’s range.
Each example sieving function takes a nonnegative integer named limit
and returns a list/array of limit+1
elements (usually integers). Then, for any query number k
in the range 0 <= k <= limit
, the expression result[k]
gives the value of the desired numbertheoretic function.
Most numbertheoretic functions don’t have a meaningful or consistent definition at 0 or 1. But because Python arrays start at index 0, I made an earnest attempt to give sensible values to indexes 0 and 1. For the most part, you can ignore these first two values, and only consider values starting at index 2.
Source code
 Python (plain, compatible with 2 and 3)
 Python (NumPy, compatible with 2 and 3)
 Java

 EratosthenesSieves.java (SE 1.2+)
 EratosthenesSievesTest.java (SE 8+)
 C#
 C++ (C++11 and above)
 C (C99 and above)
License: Project Nayuki hereby places all code on this page regarding the sieve of Eratosthenes in the public domain. Retaining the credit notice containing the author and URL is encouraged but not required.
Primeness (the basic sieve)
def sieve_primeness(limit): result = [True] * (limit + 1) result[0] = False result[1] = False for i in range(2, len(result)): if result[i]: for j in range(i * i, len(result), i): result[j] = False return result
isPrime(k) is true if k is a prime number; otherwise it’s false if k is 0 or 1 or composite.
Table of example values:
k  isPrime(k) 

0  false 
1  false 
2  true 
3  true 
4  false 
5  true 
6  false 
7  true 
8  false 
9  false 
10  false 
11  true 
12  false 
13  true 
14  false 
15  false 
16  false 
17  true 
18  false 
19  true 
20  false 
21  false 
22  false 
23  true 
24  false 
25  false 
26  false 
27  false 
28  false 
29  true 
30  false 
Smallest prime factor
def sieve_smallest_prime_factor(limit): result = [0] * (limit + 1) result[1] = 1 for i in range(2, len(result)): if result[i] == 0: result[i] = i for j in range(i * i, len(result), i): if result[j] == 0: result[j] = i return result
smallestPrimeFactor(k) is the smallest integer p such that p is greater than 1 and p divides k. The smallest proper divisor is necessarily a prime number.
When given a composite number k, we can repeatedly divide out its smallest prime factor in order to obtain its prime factorization.
Table of example values:
k  smallestPrimeFactor(k) 

0  0 
1  1 
2  2 
3  3 
4  2 
5  5 
6  2 
7  7 
8  2 
9  3 
10  2 
11  11 
12  2 
13  13 
14  2 
15  3 
16  2 
17  17 
18  2 
19  19 
20  2 
21  3 
22  2 
23  23 
24  2 
25  5 
26  2 
27  3 
28  2 
29  29 
30  2 
Notice that smallestPrimeFactor(k) = k if and only if k is a prime number. Hence, this information can be viewed as a richer version of the plain sieve of Eratosthenes (which only gives a Boolean value per number).
Totient function
def sieve_totient(limit): result = list(range(limit + 1)) for i in range(2, len(result)): if result[i] == i: for j in range(i, len(result), i): result[j] = result[j] // i return result
totient(k) (or φ(k)) is defined as the number of integers in the range [1, k] that are coprime with k. This value is useful for modular arithmetic.
Table of example values:
k  totient(k) 

0  0 
1  1 
2  1 
3  2 
4  2 
5  4 
6  2 
7  6 
8  4 
9  6 
10  4 
11  10 
12  4 
13  12 
14  6 
15  8 
16  8 
17  16 
18  6 
19  18 
20  8 
21  12 
22  10 
23  22 
24  8 
25  20 
26  12 
27  18 
28  12 
29  28 
30  8 
Notice that totient(k) = k−1 if and only if k is a prime number. Hence, this information can be viewed as a richer version of the plain sieve of Eratosthenes (which only gives a Boolean value per number).
Omega function
def sieve_omega(limit): result = [0] * (limit + 1) for i in range(2, len(result)): if result[i] == 0: for j in range(i, len(result), i): result[j] += 1 return result
omega(k) (or ω(k)) is defined as the number of distinct prime factors that k has.
Table of example values:
k  omega(k) 

0  0 
1  0 
2  1 
3  1 
4  1 
5  1 
6  2 
7  1 
8  1 
9  1 
10  2 
11  1 
12  2 
13  1 
14  2 
15  2 
16  1 
17  1 
18  2 
19  1 
20  2 
21  2 
22  2 
23  1 
24  2 
25  1 
26  2 
27  1 
28  2 
29  1 
30  3 
Note that for a given m, the smallest k such that omega(k) = m is given by k = primorial(m). Furthermore the function grows rather slowly, because omega(k) ≤ log_{2}(k).
Radical function
def sieve_radical(limit): result = [1] * (limit + 1) result[0] = 0 for i in range(2, len(result)): if result[i] == 1: for j in range(i, len(result), i): result[j] *= i return result
radical(k) is defined as the product of the unique prime factors of k.
Table of example values:
k  radical(k) 

0  0 
1  1 
2  2 
3  3 
4  2 
5  5 
6  6 
7  7 
8  2 
9  3 
10  10 
11  11 
12  6 
13  13 
14  14 
15  15 
16  2 
17  17 
18  6 
19  19 
20  10 
21  21 
22  22 
23  23 
24  6 
25  5 
26  26 
27  3 
28  14 
29  29 
30  30 
Notice that radical(k) = k if and only if k is squarefree.