MathAlgorithmaka prime sieve, linear sieve, SPF sieve, smallest prime factor sieve

Sieve of Eratosthenes

Find every prime up to n by crossing out multiples of each prime starting from its square, in O(n log log n).

▶ VisualizePattern: Math & Number TheoryPractice (1)
Progress

Overview

The Sieve of Eratosthenes marks composites instead of testing primes. Start with every number 2..n unmarked. For each i from 2 upward, if i is still unmarked it is prime — cross out i², i² + i, i² + 2i, …. Stop the outer loop at √n; whatever remains unmarked is prime. Total work is O(n log log n), nearly linear.

Two upgrades are common. The smallest-prime-factor (SPF) sieve records, for every x ≤ n, its smallest prime divisor; with it, any x factors in O(log x) by repeatedly dividing by spf[x] — see Prime Factorization. The linear sieve (Euler's sieve) marks every composite exactly once via its smallest prime factor, achieving O(n) and producing the SPF table as a by-product.

Memory is a bit or byte per number: n = 10^7 fits in 10 MB of bytes or ~1.2 MB of bits. Beyond ~10^8, use a segmented sieve that processes windows of size √n with the small primes.

primessievenumber theoryO(n log log n)precomputation

Intuition

A mental model before the formal terms.

Write the numbers 2 to 100 on a grid. Circle 2 and strike out every second number after 4. The next unstruck number, 3, is prime — strike every third number from 9. Then 5 from 25, 7 from 49. After 7 (7² = 49 ≤ 100 < 11² = 121) nothing more needs to be struck; the survivors are the 25 primes below 100. Every composite was hit by its smallest prime factor, which is at most √100.

How it works

  1. Allocate isPrime[0..n] set to true; mark 0 and 1 false.
  2. For i = 2 while i·i ≤ n: if isPrime[i], for j = i·i; j ≤ n; j += i: isPrime[j] = false.
  3. Collect all i with isPrime[i] true.
  4. SPF variant: allocate spf[0..n] = 0. For each i with spf[i] == 0 (prime), set spf[i] = i and for multiples j = i·i .. n step i, set spf[j] = i only if spf[j] == 0.
  5. Linear sieve: keep a list primes. For i = 2..n: if spf[i] == 0 then spf[i] = i and append. Then for each prime p in primes with p ≤ spf[i] and i·p ≤ n: spf[i·p] = p. Each composite c is set exactly once, when i = c / spf[c].

Why it works

Every composite c ≤ n has a prime factor p ≤ √c ≤ √n, so it is crossed out during the pass for p (or earlier). Every prime is never a multiple of a smaller prime, so it survives.

Starting at is safe because i·k for k < i was already crossed out by a prime factor of k, which is smaller than i.

Cost: the pass for prime p touches n/p cells. Summing n/p over primes p ≤ n gives n · Σ 1/p ≈ n · ln ln n (Mertens), so O(n log log n).

Linear sieve: the loop over primes stops at spf[i], so composite i·p is only written with p ≤ spf[i], i.e. p is its smallest prime factor. Each composite has one smallest prime factor, so it is written exactly once.

Recognition

How to tell a problem wants this.

  • "Count primes ≤ n", "list primes", "is each of these q numbers prime" with many queries and n ≤ ~10^7.
  • Repeated factorization of many numbers — SPF sieve.
  • Multiplicative functions over a range (Euler's φ, number of divisors, Möbius μ) — the linear sieve computes them alongside.

Interactive visualization

Play, step, change the input. ← → and space work too.

2
0
3
1
4
2
5
3
6
4
7
5
8
6
9
7
10
8
11
9
12
10
13
11
14
12
15
13
16
14
17
15
18
16
19
17
20
18
21
19
22
20
23
21
24
22
25
23
26
24
27
25
28
26
29
27
30
28
1/30List every number from 2 to 30, all assumed prime. We only need to sieve with p up to sqrt(30) ≈ 5, because any composite ≤ n has a factor that small.
Current prime pMultiple being crossed outComposite (crossed out)Prime
1isPrime = [true] * (n + 1)
2for p in 2 .. floor(sqrt(n)):
3 if isPrime[p]:
4 for m in p*p, p*p+p, .. n:
5 isPrime[m] = false
6primes = [p for p in 2..n if isPrime[p]]
Variables
n30
limit5
Complexity
best O(n)
avg O(n log log n)
worst O(n log log n)
space O(n)
Speed

Pseudocode

1isPrime[0..n] = true; isPrime[0] = isPrime[1] = false
2for i = 2; i * i <= n; i++:
3 if isPrime[i]:
4 for j = i * i; j <= n; j += i:
5 isPrime[j] = false
6return [i for i in 2..n if isPrime[i]]

Implementations

1import math
2
3
4def sieve(n: int) -> list[int]:
5 """Cross out every multiple of every prime. Finds all primes below n in
6 O(n log log n) — effectively linear for any practical n."""
7 if n < 2:
8 return []
9
101 · One flag per number; 0 and 1 are not prime by definition
11 is_composite = bytearray(n)
12 primes: list[int] = []
13
14 for p in range(2, n):
15 if is_composite[p]:
16 continue
17 primes.append(p)
18
192 · Start at p*p — smaller multiples already have a smaller factor
20 if p * p >= n:
21 continue
22 is_composite[p * p :: p] = b"\x01" * ((n - p * p + p - 1) // p)
23 return primes
24
25
263 · Smallest-prime-factor sieve: same cost, but factorises in O(log n)
27def smallest_prime_factor(n: int) -> list[int]:
28 spf = [0] * n
29 for i in range(2, n):
30 if spf[i] != 0:
31 continue
32 for m in range(i, n, i):
33 if spf[m] == 0:
34 spf[m] = i
35 return spf
36
37
384 · Segmented sieve: primes in [lo, hi) without allocating hi flags
39def segmented_sieve(lo: int, hi: int) -> list[int]:
40 limit = math.isqrt(hi) + 1
41 base = sieve(limit + 1)
42
43 composite = bytearray(hi - lo)
44 for p in base:
45 start = max(p * p, -(-lo // p) * p) # -(-a // b) is ceiling division
46 for m in range(start, hi, p):
47 composite[m - lo] = 1
48
495 · Collect what survived, skipping 0 and 1 if the range includes them
50 return [v for v in range(max(lo, 2), hi) if not composite[v - lo]]
Walkthrough
  1. bytearray(n) is the dense flag array: one byte per entry, zero-initialised, mutable.
  2. is_composite[p * p :: p] = b"\x01" * count is the crucial Python idiom — slice assignment with a step crosses out every multiple in one C-level operation instead of a Python loop.
  3. The count (n - p*p + p - 1) // p is the number of slots the extended slice covers; slice assignment requires the right-hand side to match that length exactly.
  4. math.isqrt(hi) gives the exact integer square root with no float rounding, which matters when hi exceeds 2^53.
  5. -(-lo // p) * p is ceiling division written with the negation trick, since Python // floors rather than truncates.
Complexity (this implementation)
time O(n log log n) for the classic sieve; O((hi - lo) log log hi + sqrt(hi)) segmented · space O(n) bytes; O(hi - lo + sqrt(hi)) segmented

The slice-assignment trick moves the inner loop into C and is typically 5-10x faster than an explicit for m in range(...).

Language notes
  • Extended slice assignment on a bytearray is the fastest pure-Python sieve idiom; numpy with arr[p*p::p] = 1 is faster still and reads the same.
  • math.isqrt is exact for arbitrarily large integers, unlike int(math.sqrt(n)) which loses precision past 2^53.
  • // is floor division, so -(-a // b) is the standard ceiling-division idiom; math.ceil(a / b) goes through a float and is inexact for huge values.
  • sympy.primerange and sympy.factorint are the library answers when correctness matters more than owning the code.
Common mistakes in this language
  • Mis-computing the slice-assignment length, which raises ValueError: attempt to assign bytes of size X to extended slice of size Y.
  • Using int(math.sqrt(n)) on a large n and getting a bound one too small, so a prime factor is missed.
  • Writing math.ceil(lo / p) for huge lo, where the float division has already lost the low bits.
Language differences that matter here
  • Dense flags: C++ std::vector<char>, JS/TS Uint8Array, Python bytearray — and in Python the crossing-out loop can be replaced entirely by extended slice assignment, which has no equivalent in the other three.
  • Exact integer square root exists only in Python (math.isqrt); C++ and JS/TS go through a double sqrt, which is exact only below 2^53.
  • Overflow of p * p is a real hazard only in C++ with 32-bit int; JavaScript doubles are exact to 2^53 and Python integers are unbounded.
  • Ceiling division: C++ (a + b - 1) / b, JS/TS Math.ceil(a / b), Python -(-a // b) — three idioms for one operation, and only the Python one stays exact at arbitrary magnitude.

Complexity

Best
O(n)
Average
O(n log log n)
Worst
O(n log log n)
Space
O(n)

Linear sieve is O(n) time. Segmented sieve reduces memory to O(sqrt n) plus one window. Trial division of a single number is O(sqrt n).

Compare growth rates in the Complexity Explorer →

When to use — and when not to

Use it when
  • All primes up to n ≤ ~10^7–10^8, or many primality queries in that range.
  • Factorizing many numbers quickly — precompute SPF once, then O(log x) per number.
  • Computing multiplicative functions (φ, μ, divisor counts) for every number up to n.
Avoid it when
  • A single large number (10^12 or 10^18): use trial division to √n, Miller–Rabin, or Pollard rho — a sieve cannot allocate that range.
  • Only a few small queries: O(√n) trial division per query is simpler and uses no memory.
  • A range [L, R] with huge L — use a segmented sieve over that window rather than sieving from 2.

Alternatives

Common mistakes

  • Starting the inner loop at 2·i instead of i·i (correct but ~2× slower), or bounding the outer loop at n instead of √n.
  • i * i overflowing 32-bit int when n is near 2^31 — cast to 64-bit or loop with i <= n / i.
  • Marking with a List<Boolean> or Array of boxed values in Java/JS — use primitive arrays or typed arrays; the sieve is memory-bound.
  • Forgetting to mark 0 and 1 as non-prime.
  • Linear sieve: omitting the p > spf[i] break, which makes it O(n log n) and marks composites multiple times.

Interview patterns

  • Count Primes (LeetCode 204) — the sieve is the intended solution; trial division per number times out.
  • Prime pairs with a target sum, or "closest prime numbers in range": sieve then scan.
  • Number of distinct prime factors / divisors for all x ≤ n via SPF.
  • Sum of Euler's totient over 1..n with the linear sieve.

Example problems