Project Euler Problem 423

Let n be a positive integer.

Project Euler Problem 423

Solution

Answer: 653972374

Let

$$C(n)=#{\text{length-}n\text{ die sequences with }c\le \pi(n)},$$

where $c$ is the number of adjacent equal pairs.

We must compute

$$S(50,000,000)=\sum_{n=1}^{50,000,000} C(n) \pmod{1,000,000,007}.$$


1. Mathematical analysis

Counting sequences with exactly $k$ equal adjacencies

A sequence of length $n$ has $n-1$ adjacent gaps.

Suppose exactly $k$ of those gaps are equal.

  • Choose the $k$ equal gaps:

$$\binom{n-1}{k}.$$

  • Choose the first die value:

$$6.$$

  • Every non-equal gap forces a change, giving $5$ choices each.

There are $n-1-k$ such gaps.

Hence:

$$\boxed{ N(n,k)=6\binom{n-1}{k}5^{,n-1-k} }$$

Therefore

$$\boxed{ C(n)=6\sum_{k=0}^{\pi(n)} \binom{n-1}{k}5^{,n-1-k} }$$

This immediately verifies:

$$C(4)=6\left(5^3+3\cdot 5^2+3\cdot 5\right)=1290.$$


2. A fast recurrence

Define

$$D(n)=\frac{C(n)}{6} = \sum_{k=0}^{\pi(n)} \binom{n-1}{k}5^{,n-1-k}.$$

Also define

$$p=\pi(n), \qquad B_n=\binom{n-1}{p}5^{,n-1-p}.$$


Composite transition

If $n+1$ is composite, then $\pi(n+1)=p$.

Using Pascal:

$$\binom{n}{k}=\binom{n-1}{k}+\binom{n-1}{k-1}.$$

Then

$$D(n+1) = 6D(n)-B_n.$$

Thus

$$\boxed{ C(n+1)=6C(n)-6B_n }$$


Prime transition

If $n+1$ is prime, then $\pi(n+1)=p+1$.

The new upper summation limit contributes one extra term:

$$D(n+1) = 6D(n) + \binom{n-1}{p+1}5^{,n-1-p}.$$

Using

$$\binom{n-1}{p+1} = \binom{n-1}{p}\frac{n-1-p}{p+1},$$

we obtain

$$\boxed{ C(n+1) = 6C(n) + 6B_n\frac{n-1-p}{p+1} }$$


3. Updating $B_n$

We also need a recurrence for $B_n$.


If $n+1$ is composite

Then $p$ stays fixed:

$$B_{n+1} = \binom{n}{p}5^{n-p} = B_n\cdot \frac{5n}{n-p}.$$

So

$$\boxed{ B_{n+1}=B_n\frac{5n}{n-p} }$$


If $n+1$ is prime

Then $p\to p+1$:

$$B_{n+1} = \binom{n}{p+1}5^{,n-1-p} = B_n\frac{n}{p+1}.$$

Hence

$$\boxed{ B_{n+1}=B_n\frac{n}{p+1} }$$


4. Complexity

Everything becomes an $O(L)$ loop with only constant work per step.

We sieve primes up to $50,000,000$, precompute modular inverses, and iterate once.

This easily runs within a few seconds in optimized C++.


5. Python implementation

MOD = 1_000_000_007
L = 50_000_000

# ---------- sieve ----------
is_prime = bytearray(b"\x01") * (L + 1)
is_prime[0] = is_prime[1] = 0

p = 2
while p * p <= L:
    if is_prime[p]:
        step = p
        start = p * p
        is_prime[start:L + 1:step] = b"\x00" * (((L - start) // step) + 1)
    p += 1

# ---------- modular inverses ----------
inv = [0] * (L + 1)
inv[1] = 1

for i in range(2, L + 1):
    inv[i] = MOD - (MOD // i) * inv[MOD % i] % MOD

# ---------- recurrence ----------
# C(1) = 6
C = 6

# B_1 = binom(0,0) * 5^0 = 1
B = 1

S = 6
pi_n = 0

for n in range(1, L):

    if is_prime[n + 1]:
        # prime transition

        extra = (
            6
            * B
            % MOD
            * (n - 1 - pi_n)
            % MOD
            * inv[pi_n + 1]
        ) % MOD

        C = (6 * C + extra) % MOD

        pi_n += 1

        B = B * n % MOD * inv[pi_n] % MOD

    else:
        # composite transition

        C = (6 * C - 6 * B) % MOD

        B = (
            B
            * 5
            % MOD
            * n
            % MOD
            * inv[n - pi_n]
        ) % MOD

    S += C
    S %= MOD

print(S)

6. Walkthrough of the code

Prime sieve

We compute primality for all integers up to $50,000,000$.

This lets us know whether $\pi(n)$ increases at each step.


Modular inverses

The recurrence divides by numbers like $p+1$ and $n-p$.

Using:

$$\text{inv}[i] = -(M//i)\cdot \text{inv}[M\bmod i] \pmod M$$

we precompute all inverses in linear time.


Main loop

We maintain:

  • C = C(n)
  • B = B_n
  • pi_n = π(n)

At each step:

  • if $n+1$ is prime:

use the prime recurrence,

  • otherwise:

use the composite recurrence.

We accumulate into S.


7. Verification on the given examples

The formulas reproduce:

$$C(3)=216$$

$$C(4)=1290$$

$$C(11)=361912500$$

$$C(24)=4727547363281250000$$

and

$$S(50)\bmod 1,000,000,007 = 832833871.$$

All match the problem statement.


8. Final computation

Running the program for

$$L=50,000,000$$

gives:

$$\boxed{653972374}$$

Answer: 653972374