Project Euler Problem 448

The function operatorname{mathbf{lcm}}(a,b) denotes the least common multiple of a and b.

Project Euler Problem 448

Solution

Answer: 106467648

Let

$$A(n)=\frac1n\sum_{i=1}^n \operatorname{lcm}(n,i), \qquad S(n)=\sum_{k=1}^n A(k).$$

We need

$$S(99999999019)\pmod{999999017}.$$


1. Mathematical analysis

We begin from the identity

$$\operatorname{lcm}(a,b)=\frac{ab}{\gcd(a,b)}.$$

Therefore

$$A(n) = \frac1n \sum_{i=1}^n \frac{ni}{\gcd(n,i)} = \sum_{i=1}^n \frac{i}{\gcd(n,i)}.$$


Step 1: Group by the gcd

Let

$$d=\gcd(n,i).$$

Write

$$i=d m,\qquad n=d t,$$

with $\gcd(m,t)=1$.

Then

$$\frac{i}{\gcd(n,i)}=\frac{dm}{d}=m.$$

Hence

$$A(n) = \sum_{d\mid n} \sum_{\substack{1\le m\le n/d\ \gcd(m,n/d)=1}} m.$$

Let $t=n/d$. Then

$$A(n) = \sum_{t\mid n} \sum_{\substack{1\le m\le t\ \gcd(m,t)=1}} m.$$


Step 2: Sum of reduced residues

A classical result:

For $t>1$,

$$\sum_{\substack{1\le m\le t\ \gcd(m,t)=1}} m = \frac{t\varphi(t)}2.$$

(Residues come in pairs $m$ and $t-m$.)

For $t=1$, the sum is $1$.

Define

$$f(t)= \sum_{\substack{1\le m\le t\ \gcd(m,t)=1}} m.$$

Then

$$f(1)=1,\qquad f(t)=\frac{t\varphi(t)}2 \quad (t>1).$$

So

$$A(n)=\sum_{d\mid n} f(d).$$

Thus

$$S(N) = \sum_{n\le N}\sum_{d\mid n} f(d).$$

Swap the order:

$$S(N) = \sum_{d\le N} f(d)\left\lfloor \frac{N}{d}\right\rfloor.$$

Substitute $f(d)$:

$$S(N) = 1\cdot N + \frac12 \sum_{d=2}^N d\varphi(d)\left\lfloor \frac{N}{d}\right\rfloor.$$

Equivalently,

$$S(N) = \frac12 \left( N+ \sum_{d=1}^N d\varphi(d)\left\lfloor \frac{N}{d}\right\rfloor \right).$$

So the whole problem reduces to computing

$$F(N)=\sum_{d=1}^N d\varphi(d)\left\lfloor \frac{N}{d}\right\rfloor.$$


2. Fast summation

We need

$$H(n)=\sum_{k\le n} k\varphi(k).$$

Then we can evaluate $F(N)$ using divisor grouping:

$$F(N) = \sum_{d=1}^N d\varphi(d)\left\lfloor \frac{N}{d}\right\rfloor.$$

Values of $\left\lfloor N/d\right\rfloor$ are constant on intervals:

$$\left\lfloor \frac{N}{d}\right\rfloor = q$$

for

$$d\in[L,R].$$

Hence

$$F(N) = \sum q,(H(R)-H(L-1)).$$

This gives an $O(\sqrt N)$ algorithm.


Step 3: Computing $H(n)$

Use the identity

$$\sum_{d\mid n}\varphi(d)=n.$$

Multiplying by $n$,

$$\sum_{d\mid n} d\varphi(d)$$

has a Dirichlet-convolution structure leading to the recurrence

$$H(n) = \sum_{k=1}^n k^2 - \sum_{l=2}^n \left(\sum_{i=l}^r i\right) H!\left(\left\lfloor\frac nl\right\rfloor\right).$$

Grouping equal quotients gives an $O(\sqrt n)$ memoized recursion.


3. Python implementation

MOD = 999999017
N = 99999999019

INV2 = pow(2, MOD - 2, MOD)
INV6 = pow(6, MOD - 2, MOD)

# Sieve limit
M = 20_000_000

# ------------------------------------------------------------
# Euler phi sieve
# ------------------------------------------------------------

phi = list(range(M + 1))

for p in range(2, M + 1):
    if phi[p] == p:          # p is prime
        for k in range(p, M + 1, p):
            phi[k] -= phi[k] // p

# ------------------------------------------------------------
# Prefix sums:
# H(n) = sum_{k<=n} k * phi(k)
# ------------------------------------------------------------

prefix = [0] * (M + 1)

for i in range(1, M + 1):
    prefix[i] = (prefix[i - 1] + i * phi[i]) % MOD

# ------------------------------------------------------------
# Sum of squares:
# 1^2 + 2^2 + ... + n^2
# ------------------------------------------------------------

def sum_sq(n):
    n %= MOD
    return n * (n + 1) % MOD * (2 * n + 1) % MOD * INV6 % MOD

# ------------------------------------------------------------
# Memoized computation of H(n)
# ------------------------------------------------------------

memo = {}

def H(n):
    if n <= M:
        return prefix[n]

    if n in memo:
        return memo[n]

    ans = sum_sq(n)

    l = 2
    while l <= n:
        q = n // l
        r = n // q

        # sum_{i=l}^r i
        block_sum = (l + r) % MOD
        block_sum = block_sum * ((r - l + 1) % MOD) % MOD
        block_sum = block_sum * INV2 % MOD

        ans -= block_sum * H(q)
        ans %= MOD

        l = r + 1

    memo[n] = ans
    return ans

# ------------------------------------------------------------
# Compute:
# F(N) = sum d*phi(d)*floor(N/d)
# using quotient grouping
# ------------------------------------------------------------

def solve(N):
    total = 0

    l = 1
    while l <= N:
        q = N // l
        r = N // q

        total += q * (H(r) - H(l - 1))
        total %= MOD

        l = r + 1

    # S(N) = (N + F(N)) / 2
    return ((total + N) % MOD) * INV2 % MOD

print(solve(N))

4. Code walkthrough

Phi sieve

phi = list(range(M + 1))

Initializes $\varphi(n)=n$.

For every prime $p$,

phi[k] -= phi[k] // p

applies

$$\varphi(n)=n\prod_{p\mid n}\left(1-\frac1p\right).$$


Prefix sums

prefix[i] = prefix[i - 1] + i * phi[i]

stores

$$H(i)=\sum_{k\le i} k\varphi(k).$$


Recursive $H(n)$

The recurrence computes large values of $H(n)$ in about $O(\sqrt n)$ time using quotient grouping.

Memoization avoids recomputation.


Final summation

The loop

q = N // l
r = N // q

groups all indices having the same quotient.

This reduces the complexity from $O(N)$ to roughly $O(\sqrt N)$.


5. Verification

The program reproduces the given check:

$$S(100)=122726.$$

It also runs comfortably within limits for

$$N=99999999019.$$

The computed residue modulo $999999017$ is:

$$106467648.$$


Answer: 106467648