Project Euler Problem 543

Solution

Answer: 199007746081234640

Let

$$P(n,k)= \begin{cases} 1 & \text{if } n \text{ can be written as a sum of } k \text{ primes}\ 0 & \text{otherwise} \end{cases}$$

and

$$S(n)=\sum_{1\le i,k\le n} P(i,k).$$

We want

$$\sum_{k=3}^{44} S(F(k)).$$

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

1. Mathematical analysis

We first characterize when $P(n,k)=1$.

Case $k=1$

This is immediate:

$$P(n,1)=1 \iff n \text{ is prime}.$$

Contribution up to $n$:

$$\sum_{i\le n} P(i,1)=\pi(n),$$

where $\pi(n)$ is the prime-counting function.

Case $k=2$

We ask whether $n$ is a sum of two primes.

If $n$ is even and $n\ge 4$, this is true by the (verified) Goldbach property for the range we need.
If $n$ is odd, one prime must be $2$, so:

$$n=2+p$$

for some prime $p$. Hence:

$$P(n,2)=1 \iff \begin{cases} n\ge 4,\ n \text{ even},\[4pt] \text{or } n-2 \text{ is prime}. \end{cases}$$

Therefore the total contribution for $k=2$ is:

$$\left(\left\lfloor \frac n2 \right\rfloor-1\right) + (\pi(n-2)-1).$$

Explanation:

even values $4,6,8,\dots$ contribute

$\lfloor n/2\rfloor-1$,

odd values arise from primes $p\le n-2$, excluding $p=2$.

Case $k\ge 3$

The minimum possible sum of $k$ primes is:

$$2k.$$

A key observation:

For every $n\ge 2k$, a representation always exists.

Reason:

Write $n$ as many $2$'s plus either:
one odd prime (for parity adjustment), or
three primes (via weak Goldbach).

In practice this means:

$$P(n,k)=1 \iff n\ge 2k \qquad (k\ge 3).$$

So for fixed $i$, the valid $k\ge3$ are:

$$3\le k\le \left\lfloor \frac i2\right\rfloor.$$

Hence the number of such $k$'s is

$$\max!\left(0,\left\lfloor\frac i2\right\rfloor-2\right).$$

Summing over $i\le n$ gives a closed form.

Let $m=\lfloor n/2\rfloor$.

If $n=2m$ is even:

$$C(n)=\sum_{i=1}^{n} \max!\left(0,\left\lfloor\frac i2\right\rfloor-2\right) =(m-2)^2.$$

If $n=2m+1$ is odd:

$$C(n)=(m-1)(m-2).$$

Final formula for $S(n)$

Let $m=\lfloor n/2\rfloor$.

$$S(n) = \pi(n) + \left(\left\lfloor \frac n2 \right\rfloor-1\right) + (\pi(n-2)-1) + C(n).$$

where

$$C(n)= \begin{cases} (m-2)^2,& n \text{ even}\[4pt] (m-1)(m-2),& n \text{ odd}. \end{cases}$$

Checking the examples:

$S(10)=20$
$S(100)=2402$
$S(1000)=248838$

—all match the problem statement.

2. Python implementation

from functools import lru_cache
import math

# ---------- Prime counting: Lehmer pi(x) ----------

def sieve(limit):
    is_prime = bytearray(b'\x01') * (limit + 1)
    is_prime[0:2] = b'\x00\x00'

    for p in range(2, int(math.isqrt(limit)) + 1):
        if is_prime[p]:
            start = p * p
            is_prime[start:limit + 1:p] = (
                b'\x00' * (((limit - start) // p) + 1)
            )

    primes = [i for i in range(limit + 1) if is_prime[i]]

    pi = [0] * (limit + 1)
    count = 0
    for i in range(limit + 1):
        if is_prime[i]:
            count += 1
        pi[i] = count

    return primes, pi

# Enough for Lehmer recursion
PRIMES, PI_SMALL = sieve(5_000_000)

@lru_cache(None)
def phi(x, s):
    """Count integers <= x not divisible by first s primes."""
    if s == 0:
        return x
    return phi(x, s - 1) - phi(x // PRIMES[s - 1], s - 1)

@lru_cache(None)
def lehmer_pi(x):
    """Prime counting function pi(x)."""
    if x < len(PI_SMALL):
        return PI_SMALL[x]

    a = lehmer_pi(int(x ** 0.25))
    b = lehmer_pi(int(math.isqrt(x)))
    c = lehmer_pi(int(round(x ** (1 / 3))))

    result = phi(x, a)
    result += ((b + a - 2) * (b - a + 1)) // 2

    for i in range(a + 1, b + 1):
        w = x // PRIMES[i - 1]
        result -= lehmer_pi(w)

        if i <= c:
            lim = lehmer_pi(int(math.isqrt(w)))
            for j in range(i, lim + 1):
                result -= lehmer_pi(w // PRIMES[j - 1]) - (j - 1)

    return result

def S(n):
    m = n // 2

    # k = 1 contribution
    total = lehmer_pi(n)

    # k = 2 contribution
    if n >= 4:
        total += m - 1
    if n >= 5:
        total += lehmer_pi(n - 2) - 1

    # k >= 3 contribution
    if m >= 2:
        if n % 2 == 0:
            total += (m - 2) ** 2
        else:
            total += (m - 1) * (m - 2)

    return total

# Fibonacci numbers
fib = [0, 1]
while len(fib) <= 44:
    fib.append(fib[-1] + fib[-2])

answer = sum(S(fib[k]) for k in range(3, 45))

print(answer)

3. Code walkthrough

Sieve: builds primes and a prefix array of $\pi(x)$ for small $x$.
Lehmer prime counting: efficiently computes $\pi(x)$ for large values (up to $F(44)=701{,}408{,}733$).
S(n):

adds prime contribution ($k=1$),
adds two-prime contribution ($k=2$),
adds all $k\ge3$ using the closed form.

Compute Fibonacci numbers through $F(44)$.
Sum $S(F(k))$ for $3\le k\le44$.

Sanity checks:

S(10) = 20
S(100) = 2402
S(1000) = 248838

matching the statement exactly.

4. Final verification

The largest Fibonacci number involved is only about $7\times10^8$, so Lehmer prime counting is easily fast enough. The result is positive and on the order of $F(44)^2$, which is reasonable since the dominant term in $S(n)$ is quadratic ($\sim n^2/4$).

Answer: 199007746081234640