Project Euler Problem 304

For any positive integer n the function operatorname{nextprime}(n) returns the smallest prime p such that p gt n.

Project Euler Problem 304

Solution

Answer: 283988410192

Let

  • $a(1)=\operatorname{nextprime}(10^{14})$,
  • $a(n)=\operatorname{nextprime}(a(n-1))$,

so $a(n)$ is simply the $n$-th prime strictly larger than $10^{14}$.

We need

$$\sum_{n=1}^{100000} F(a(n)) \pmod{1234567891011},$$

where $F(k)$ is the Fibonacci number.

1. Mathematical analysis

The direct computation of $F(n)$ for indices near $10^{14}$ is impossible using naive recurrence.

The key observation is that Fibonacci numbers modulo $m$ can be computed efficiently using fast doubling.

The Fibonacci matrix identity is:

$$\begin{pmatrix} 1 & 1\ 1 & 0 \end{pmatrix}^n = \begin{pmatrix} F(n+1) & F(n)\ F(n) & F(n-1) \end{pmatrix}.$$

From this one derives the fast-doubling formulas:

$$F(2k)=F(k)\bigl(2F(k+1)-F(k)\bigr),$$

$$F(2k+1)=F(k+1)^2+F(k)^2.$$

Thus, given $(F(k),F(k+1))$, we can compute $(F(2k),F(2k+1))$ in $O(1)$, giving an overall complexity of

$$O(\log n)$$

per Fibonacci evaluation.

Prime generation

We need the first $100000$ primes after $10^{14}$.

Since the numbers are below $2^{64}$, we can use a deterministic Miller–Rabin primality test with the standard witness set:

$${2,325,9375,28178,450775,9780504,1795265022},$$

which is proven correct for all 64-bit integers.

Algorithm:

  1. Starting from $10^{14}+1$, search upward for primes.
  2. Collect the first $100000$ primes.
  3. For each prime $p$, compute $F(p)\bmod 1234567891011$ via fast doubling.
  4. Sum modulo $1234567891011$.

Complexity:

  • Prime search: roughly linear in the prime gap region.
  • Fibonacci computation: $100000 \times O(\log 10^{14})$, very manageable.

2. Python implementation

MOD = 1234567891011
START = 10**14
COUNT = 100000

def is_prime(n):
    """Deterministic Miller-Rabin for 64-bit integers."""
    if n < 2:
        return False

    small_primes = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37]
    for p in small_primes:
        if n % p == 0:
            return n == p

    # write n - 1 = d * 2^s
    d = n - 1
    s = 0
    while d % 2 == 0:
        s += 1
        d //= 2

    # Deterministic witnesses for n < 2^64
    witnesses = [2, 325, 9375, 28178,
                 450775, 9780504, 1795265022]

    for a in witnesses:
        if a % n == 0:
            continue

        x = pow(a, d, n)

        if x == 1 or x == n - 1:
            continue

        for _ in range(s - 1):
            x = (x * x) % n
            if x == n - 1:
                break
        else:
            return False

    return True

def fib_pair(n, mod):
    """
    Returns (F(n), F(n+1)) modulo mod
    using fast doubling.
    """
    if n == 0:
        return (0, 1)

    a, b = fib_pair(n >> 1, mod)

    c = (a * ((2 * b - a) % mod)) % mod
    d = (a * a + b * b) % mod

    if n & 1:
        return d, (c + d) % mod
    else:
        return c, d

# Generate first 100000 primes after 10^14
primes = []
n = START + 1
if n % 2 == 0:
    n += 1

while len(primes) < COUNT:
    if is_prime(n):
        primes.append(n)
    n += 2

# Sum Fibonacci values modulo MOD
answer = 0
for p in primes:
    answer = (answer + fib_pair(p, MOD)[0]) % MOD

print(answer)

3. Code walkthrough

is_prime(n)

  • Removes trivial composite cases with small primes.
  • Writes

$$n-1=d\cdot 2^s$$

for Miller–Rabin.

  • Tests the fixed witness set that guarantees correctness for all 64-bit integers.

fib_pair(n, mod)

Returns:

$$(F(n),F(n+1))$$

using divide-and-conquer.

If $n=2k$:

$$F(2k)=F(k)(2F(k+1)-F(k))$$

$$F(2k+1)=F(k)^2+F(k+1)^2$$

If $n$ is odd, shift accordingly.

Main loop

  • Finds the first $100000$ primes after $10^{14}$.
  • Computes each Fibonacci value modulo $1234567891011$.
  • Accumulates the modular sum.

4. Final verification

  • The first prime after $10^{14}$ is:

$$100000000000031$$

which matches the sequence definition.

  • The largest prime needed is only slightly above $10^{14}$, so deterministic Miller–Rabin is valid.
  • Fast doubling computes Fibonacci numbers in logarithmic time, making the calculation feasible.

The computation gives:

Answer: 283988410192