Sieve of Sundaram in Python
List the primes up to n by crossing out each x = i + j + 2ij with 1 <= i <= j. The numbers 2x + 1 left are the odd primes. It takes Θ(n log n) word operations and about n/2 flags.
- File
sieve_of_sundaram.py, 22 lines- Level
- I0 not checked yet
- Toolchain
- Python 3.14, pytest 9
- Test framework
- pytest
- Run the tests
python3 -m pytest
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/sieve-of-sundaram/. The authors are listed in prov/KQBYRB.jsonl.
def sieve_of_sundaram(n: int) -> list[int]:
"""Return the increasing list of all primes p <= n.
Precondition: n is an integer (n < 2 gives the empty list).
Cost: Theta(n log n) word operations and about n/2 bytes of flags.
"""
if n < 2: # @line L1
return []
k = (n - 1) // 2 # @line L2
marked = bytearray(k + 1) # @line L3: index 0 is unused
i = 1 # @line L4
while 2 * i + 2 * i * i <= k: # @line L5
j = i # @line L6
while i + j + 2 * i * j <= k: # @line L7
marked[i + j + 2 * i * j] = 1 # @line L8
j += 1 # @line L9
i += 1 # @line L10
return [2] + [2 * x + 1 for x in range(1, k + 1) if not marked[x]] # @line L11
test_sieve_of_sundaram.py 55 lines
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/sieve-of-sundaram/. The authors are listed in prov/KQBYRB.jsonl.
import random
from sieve_of_sundaram import sieve_of_sundaram
def is_prime_by_trial_division(m: int) -> bool:
if m < 2:
return False
d = 2
while d * d <= m:
if m % d == 0:
return False
d += 1
return True
LIMIT = 30000
ORACLE = [m for m in range(LIMIT + 1) if is_prime_by_trial_division(m)]
def oracle_primes(n: int) -> list[int]:
return [p for p in ORACLE if p <= n]
def test_example() -> None:
assert sieve_of_sundaram(30) == [2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
def test_edges() -> None:
assert sieve_of_sundaram(-3) == []
assert sieve_of_sundaram(0) == []
assert sieve_of_sundaram(1) == []
assert sieve_of_sundaram(2) == [2]
assert sieve_of_sundaram(3) == [2, 3]
# even n with n + 1 prime: the bound k = (n - 1) // 2 must not reach n + 1
for n in (4, 6, 10, 12, 16, 18, 22, 28, 30, 100):
assert sieve_of_sundaram(n) == oracle_primes(n)
# n equal to a square of a prime
for n in (9, 25, 49, 121, 169, 289, 361, 529, 841, 961):
assert sieve_of_sundaram(n) == oracle_primes(n)
def test_every_n_up_to_2000() -> None:
for n in range(2001):
assert sieve_of_sundaram(n) == oracle_primes(n)
def test_random_larger_n() -> None:
rng = random.Random(1934)
for _ in range(15):
n = rng.randint(2001, LIMIT)
assert sieve_of_sundaram(n) == oracle_primes(n)
What the next level needs
- the tests did not run in this build