Problem 910: Prime Constellation Counting

Mathematical Analysis

Twin Primes and the Sieve

Twin primes are pairs of primes differing by 2. The first few: $(3,5), (5,7), (11,13), (17,19), (29,31), (41,43), \ldots$

Observation. For $p \ge 5$, if $(p, p+2)$ are both prime, then $p \equiv 5 \pmod{6}$ (since one of $p, p+1, p+2$ must be divisible by 3, and it cannot be $p$ or $p+2$, so $p+1 \equiv 0 \pmod{3}$, giving $p \equiv 2 \pmod{3}$, combined with $p$ odd gives $p \equiv 5 \pmod{6}$).

The Twin Prime Conjecture

Conjecture (Twin Prime). There are infinitely many twin prime pairs.

This is one of the oldest unsolved problems in number theory. The best known result is Zhang’s 2013 breakthrough (later refined by Maynard and the Polymath project): there exist infinitely many pairs of primes differing by at most 246.

Hardy-Littlewood Conjecture

Conjecture (Hardy-Littlewood First). The number of twin prime pairs $(p, p+2)$ with $p \le N$ satisfies:

where $C_2 = \prod_{p \ge 3} \frac{p(p-2)}{(p-1)^2} = 0.6601618\ldots$ is the twin prime constant.

The product converges because the terms are $1 - 1/(p-1)^2 = 1 - O(1/p^2)$.

Brun’s Theorem

Theorem (Brun, 1919). The sum of reciprocals of all twin primes converges:

This is Brun’s constant. In contrast, $\sum_{p \text{ prime}} 1/p = \infty$ (Euler).

The Sieve of Eratosthenes

The standard algorithm for finding all primes up to $N$:

1. Create boolean array $\text{sieve}[0..N]$, initialized to true.
2. Set $\text{sieve}[0] = \text{sieve}[1] = \text{false}$.
3. For $i = 2, 3, \ldots, \lfloor\sqrt{N}\rfloor$: if $\text{sieve}[i]$ is true, mark all multiples $i^2, i^2+i, \ldots$ as composite.

Theorem. The sieve correctly identifies all primes $\le N$.

Proof. Every composite $n \le N$ has a prime factor $p \le \sqrt{N}$, so $n$ is marked composite when processing $p$. $\square$

Counting Twin Primes

After sieving, scan pairs $(p, p+2)$ for $p = 2, 3, \ldots, N-2$. The only even twin pair is $(2, 4)$, but 4 is not prime. So all twin pairs have $p$ odd, $p+2$ odd, and we only check odd $p$.

Concrete Values

The Hardy-Littlewood prediction improves with $N$ but overestimates the convergence rate.

Prime Gaps and Twin Primes

The prime gap $g_n = p_{n+1} - p_n$ equals 2 exactly for twin primes (plus $(2,3)$ where $g=1$). The average gap near $N$ is $\ln N$. Twin primes represent the minimum possible gap for $p > 3$.

Editorial

Count twin prime pairs (p, p+2) with p < 10^7. We sieve of Eratosthenes up to $N + 2$ (to check $p + 2$ for the largest $p$).

Pseudocode

Sieve of Eratosthenes up to $N + 2$ (to check $p + 2$ for the largest $p$)
Count pairs where both $\text{sieve}[p]$ and $\text{sieve}[p+2]$ are true

Proof of Correctness

Theorem. The algorithm correctly counts all twin prime pairs $(p, p+2)$ with $p \le N$.

Proof. The sieve identifies all primes $\le N + 2$ without error. For each $p \le N$, checking $\text{sieve}[p] \land \text{sieve}[p+2]$ is necessary and sufficient for $(p, p+2)$ to be a twin pair. $\square$

Complexity Analysis

• Sieve: $O(N \log \log N)$ time, $O(N)$ space.
• Counting: $O(N)$ single pass.
• Segmented sieve variant: $O(N)$ time, $O(\sqrt{N})$ space.

For $N = 10^7$: the sieve takes about 0.1 seconds in C++.

Answer

#include <bits/stdc++.h>
using namespace std;

/*
 * Problem 910: Prime Constellation Counting
 *
 * Count twin prime pairs (p, p+2) with p < 10^7.
 *
 * Method: Sieve of Eratosthenes + linear scan.
 *
 * Hardy-Littlewood conjecture: pi_2(N) ~ 2*C_2 * N / (ln N)^2
 * where C_2 = prod_{p>=3} p(p-2)/(p-1)^2 ≈ 0.6602.
 *
 * Brun's theorem: sum of 1/p over twin primes converges (B ≈ 1.9022).
 */

int main() {
    const int N = 10000000;

    // Sieve of Eratosthenes
    vector<bool> sieve(N + 3, true);
    sieve[0] = sieve[1] = false;
    for (int i = 2; (long long)i * i <= N + 2; i++) {
        if (sieve[i]) {
            for (int j = i * i; j <= N + 2; j += i) {
                sieve[j] = false;
            }
        }
    }

    // Count twin prime pairs
    int count = 0;
    for (int p = 2; p <= N; p++) {
        if (sieve[p] && sieve[p + 2]) {
            count++;
        }
    }

    // Verify small known values
    // Twin primes below 100: (3,5),(5,7),(11,13),(17,19),(29,31),(41,43),(59,61),(71,73) = 8
    int small_count = 0;
    for (int p = 2; p <= 100; p++)
        if (sieve[p] && sieve[p + 2])
            small_count++;
    assert(small_count == 8);

    // Below 1000: 35 pairs
    int med_count = 0;
    for (int p = 2; p <= 1000; p++)
        if (sieve[p] && sieve[p + 2])
            med_count++;
    assert(med_count == 35);

    cout << count << endl;
    return 0;
}

"""
Problem 910: Prime Constellation Counting

Count twin prime pairs (p, p+2) with p < 10^7.

Key results:
    - Hardy-Littlewood: pi_2(N) ~ 2*C_2 * N / (ln N)^2
    - Brun's theorem: sum of 1/p over twin primes converges
    - Twin prime constant C_2 = prod_{p>=3} p(p-2)/(p-1)^2

Methods:
    1. Sieve of Eratosthenes + pair counting
    2. Segmented sieve (memory-efficient)
    3. Verification via trial division (small N)
"""

from math import isqrt, log

def count_twin_primes_sieve(N: int):
    """Count twin prime pairs (p, p+2) with p <= N using Eratosthenes sieve."""
    sieve = bytearray(b'\x01') * (N + 3)
    sieve[0] = sieve[1] = 0
    for i in range(2, isqrt(N + 2) + 1):
        if sieve[i]:
            sieve[i * i::i] = bytearray(len(sieve[i * i::i]))
    count = 0
    for p in range(2, N + 1):
        if sieve[p] and sieve[p + 2]:
            count += 1
    return count

def is_prime_trial(n: int) -> bool:
    """Primality test by trial division."""
    if n < 2:
        return False
    if n < 4:
        return True
    if n % 2 == 0 or n % 3 == 0:
        return False
    i = 5
    while i * i <= n:
        if n % i == 0 or n % (i + 2) == 0:
            return False
        i += 6
    return True

def count_twin_primes_trial(N: int):
    """Count twin primes by trial division (slow, for verification)."""
    count = 0
    for p in range(2, N + 1):
        if is_prime_trial(p) and is_prime_trial(p + 2):
            count += 1
    return count

# Hardy-Littlewood prediction
def twin_prime_constant(num_primes: int = 100) -> float:
    """Compute twin prime constant C_2 = prod_{p>=3} p(p-2)/(p-1)^2."""
    # Get primes via sieve
    sieve = bytearray(b'\x01') * 1000
    sieve[0] = sieve[1] = 0
    for i in range(2, 32):
        if sieve[i]:
            sieve[i * i::i] = bytearray(len(sieve[i * i::i]))
    primes = [i for i in range(3, 1000) if sieve[i]]

    C2 = 1.0
    for p in primes[:num_primes]:
        C2 *= p * (p - 2) / (p - 1) ** 2
    return C2

def hardy_littlewood(N: int) -> float:
    """Hardy-Littlewood prediction for pi_2(N)."""
    C2 = twin_prime_constant()
    return 2 * C2 * N / log(N) ** 2

# Solve
N = 10**7
answer = count_twin_primes_sieve(N)

# Verify small cases
for test_n in [100, 1000, 5000]:
    sieve_ans = count_twin_primes_sieve(test_n)
    trial_ans = count_twin_primes_trial(test_n)
    assert sieve_ans == trial_ans, f"N={test_n}: sieve={sieve_ans}, trial={trial_ans}"

# Verify known values
assert count_twin_primes_sieve(100) == 8   # (3,5),(5,7),(11,13),(17,19),(29,31),(41,43),(59,61),(71,73)
assert count_twin_primes_sieve(1000) == 35

print(answer)