Sieves & Bulk Precomputation

Intuition, Vocabulary and Foundations

A sieve treats a range of integers like a checklist. Begin with every integer from $2$ through $n$ as a possible prime, then cross out each number once a known prime proves it composite. For example, sieving through $30$ leaves

$$ 2,3,5,7,11,13,17,19,23,29. $$

This resembles eliminating impossible candidates from a search: each discovered prime supplies a rule that rejects all of its larger multiples. A prime has exactly two positive divisors, while a composite integer can be expressed as a product of smaller positive integers.

In the Sieve of Eratosthenes, a composite $m$ is marked when a prime $p$ divides it. Marking may begin at $p^2$, because every smaller multiple $kp$ with $2\le k\lt p$ already has a prime factor smaller than $p$.

The smallest prime factor, written $\operatorname{spf}(m)$, is the least prime dividing $m$. Repeatedly dividing by this value recovers the full factorization. For example,

$$ 84\longrightarrow42\longrightarrow21\longrightarrow7\longrightarrow1 $$

extracts $2,2,3,7$, so $84=2^2\cdot3\cdot7$.

A segmented sieve marks only an interval $[L,R]$. It first computes all primes through $\lfloor\sqrt R\rfloor$, because every composite not exceeding $R$ has a prime factor no greater than its square root. Those base primes can then process one memory-sized block at a time.

Bulk precomputation can also produce arithmetic-function tables. Euler's totient $\varphi(n)$ counts the integers in $\{1,\ldots,n\}$ coprime to $n$. The Möbius function $\mu(n)$ is $0$ when a squared prime divides $n$; otherwise it is $(-1)^k$ when $n$ is a product of $k$ distinct primes.

Step-by-Step Derivations and Algorithms

For the ordinary sieve, allocate flags for $0$ through $n$, mark $0$ and $1$ non-prime, and inspect candidates through $\lfloor\sqrt n\rfloor$. Whenever $p$ remains unmarked, mark

$$ p^2,\ p^2+p,\ p^2+2p,\ldots $$

as composite. If some composite $m\le n$ remained unmarked, its least prime factor $p$ would satisfy $p\le\sqrt m\le\sqrt n$, so the iteration for $p$ would mark $m$. This contradiction proves completeness. The total marking work is

$$ O\!\left(n\sum_{p\le n}\frac1p\right)=O(n\log\log n), $$

with $O(n)$ storage. An odd-only representation removes even candidates after handling $2$, improving constants without changing the asymptotic bounds.

For a block $[L,R]$, each base prime $p\le\sqrt R$ begins at

$$ \max\!\left(p^2,\left\lceil\frac Lp\right\rceil p\right). $$

The ceiling term selects the first multiple inside the block, while $p^2$ prevents the prime $p$ itself from being marked. A block of width $B$ needs $O(B)$ working storage in addition to the base-prime table.

The linear sieve maintains a prime list and an SPF table. For each $i$, it generates $ip$ in increasing prime order and stops after processing the first prime $p$ dividing $i$. Thus every composite is assigned to its least prime factor exactly once. SPF factorization then repeatedly replaces $x$ by $x/\operatorname{spf}(x)$, requiring at most $\lfloor\log_2x\rfloor$ divisions.

flowchart TD
    A["Choose range and required outputs"] --> B{"Dense prefix fits in memory?"}
    B -->|Yes| C{"Need factors or functions?"}
    B -->|No| D["Generate base primes through sqrt R"]
    C -->|No| E["Run Eratosthenes or odd-only sieve"]
    C -->|Yes| F["Run linear sieve with SPF, phi, and mu"]
    D --> G["Mark one segment at a time"]
    E --> H["Answer prime queries"]
    F --> I["Answer factor and function queries"]
    G --> H

Implementation, Complexity and Engineering

A robust implementation rejects invalid bounds, represents $0$ and $1$ explicitly as non-prime, and uses integer square roots. In a fixed-width language, the loop condition $p\le n/p$ avoids overflowing while evaluating $p^2\le n$. Packed bitsets minimize memory, whereas byte arrays often offer simpler indexing and good cache behavior.

Segment width is an engineering parameter. Blocks that are too small repeat setup work; blocks that are too large may leave cache. Benchmark representative ranges and query batches instead of selecting a size from asymptotic notation alone.

from math import isqrt


def eratosthenes(limit: int) -> list[int]:
    """Return every prime p satisfying 2 <= p <= limit."""
    if limit < 0:
        raise ValueError("limit must be non-negative")

    flags = bytearray(b"\x01") * (limit + 1)
    flags[0] = 0
    if limit >= 1:
        flags[1] = 0

    for prime in range(2, isqrt(limit) + 1):
        if flags[prime]:
            start = prime * prime
            count = (limit - start) // prime + 1
            flags[start : limit + 1 : prime] = b"\x00" * count

    return [value for value, flag in enumerate(flags) if flag]


def linear_sieve(limit: int):
    """Return primes, SPF, Euler phi, and Moebius mu through limit."""
    if limit < 0:
        raise ValueError("limit must be non-negative")

    primes: list[int] = []
    spf = [0] * (limit + 1)
    phi = [0] * (limit + 1)
    mu = [0] * (limit + 1)
    if limit >= 1:
        phi[1] = 1
        mu[1] = 1

    for value in range(2, limit + 1):
        if spf[value] == 0:
            spf[value] = value
            primes.append(value)
            phi[value] = value - 1
            mu[value] = -1

        for prime in primes:
            product = value * prime
            if product > limit:
                break
            spf[product] = prime
            if value % prime == 0:
                phi[product] = phi[value] * prime
                mu[product] = 0
                break
            phi[product] = phi[value] * (prime - 1)
            mu[product] = -mu[value]

    return primes, spf, phi, mu


def factor_with_spf(value: int, spf: list[int]) -> dict[int, int]:
    if value < 1 or value >= len(spf):
        raise ValueError("value must be positive and covered by the table")
    factors: dict[int, int] = {}
    while value > 1:
        prime = spf[value]
        if prime == 0:
            raise ValueError("SPF table is incomplete")
        factors[prime] = factors.get(prime, 0) + 1
        value //= prime
    return factors


def segmented_primes(left: int, right: int) -> list[int]:
    if left > right:
        raise ValueError("left must not exceed right")
    if right < 2:
        return []

    low = max(left, 2)
    flags = bytearray(b"\x01") * (right - low + 1)
    for prime in eratosthenes(isqrt(right)):
        first = max(prime * prime, ((low + prime - 1) // prime) * prime)
        if first <= right:
            count = (right - first) // prime + 1
            offset = first - low
            flags[offset : right - low + 1 : prime] = b"\x00" * count

    return [low + offset for offset, flag in enumerate(flags) if flag]


assert eratosthenes(1) == []
assert eratosthenes(30) == [2, 3, 5, 7, 11, 13, 17, 19, 23, 29]

primes, spf, phi, mu = linear_sieve(100)
assert primes == eratosthenes(100)
assert factor_with_spf(1, spf) == {}
assert factor_with_spf(84, spf) == {2: 2, 3: 1, 7: 1}
assert (phi[1], phi[10], phi[36]) == (1, 4, 12)
assert (mu[1], mu[6], mu[12], mu[30]) == (1, 1, 0, -1)
assert segmented_primes(100, 130) == [101, 103, 107, 109, 113, 127]
assert segmented_primes(-10, 3) == [2, 3]
assert segmented_primes(14, 16) == []

for invalid in (
    lambda: eratosthenes(-1),
    lambda: linear_sieve(-1),
    lambda: segmented_primes(10, 5),
    lambda: factor_with_spf(101, spf),
):
    try:
        invalid()
        raise AssertionError("expected ValueError")
    except ValueError:
        pass

print("all sieve tests passed")

Applications, Case Studies and Exercises

Dense sieves are ideal when many queries share a moderate upper bound. One preprocessing pass can support primality checks, factorizations, divisor statistics, totients, and Möbius values. Segmented sieves target a different workload: the right endpoint may be enormous, but only a manageable interval needs to reside in memory.

A strong test oracle compares optimized output with trial division on small random ranges. Useful invariants include a strictly increasing prime list, $\operatorname{spf}(x)\mid x$, primality of every nonzero SPF value, and reconstruction of each input from its factorization.