Project Euler Problem 745

Solution

Answer: 94586478

Let

$$g(n)=\text{largest perfect square dividing } n.$$

Every positive integer can be written uniquely as

$$n = s,t^2$$

with $s$ squarefree. Therefore,

$$g(n)=t^2.$$

Hence

$$S(N)=\sum_{n\le N} g(n) =\sum_{t\le \sqrt N} t^2 \cdot Q!\left(\left\lfloor\frac N{t^2}\right\rfloor\right),$$

where $Q(x)$ counts squarefree integers $\le x$.

Using Möbius inversion,

$$Q(x)=\sum_{d^2\le x}\mu(d)\left\lfloor \frac{x}{d^2}\right\rfloor.$$

Substituting and rearranging gives

$$S(N) =\sum_{k\le \sqrt N} \left\lfloor \frac N{k^2}\right\rfloor \sum_{d\mid k} d^2\mu(k/d).$$

The inner sum is the Jordan totient function $J_2(k)$:

$$J_2(k)=k^2\prod_{p\mid k}\left(1-\frac1{p^2}\right).$$

So the problem reduces to

$$S(N)=\sum_{k\le \sqrt N} J_2(k)\left\lfloor \frac N{k^2}\right\rfloor.$$

We compute $J_2(k)$ for all $k\le 10^7$ with a linear sieve.

from array import array

MOD = 1_000_000_007
N = 10**14
V = int(N**0.5)

# Linear sieve for Jordan totient J2(n)
J2 = array('Q', [0]) * (V + 1)
J2[1] = 1

primes = array('I')
is_comp = bytearray(V + 1)

for i in range(2, V + 1):

    if not is_comp[i]:
        primes.append(i)
        J2[i] = i * i - 1

    for p in primes:
        x = i * p
        if x > V:
            break

        is_comp[x] = 1

        if i % p == 0:
            # J2(ip) = J2(i) * p^2 when p divides i
            J2[x] = J2[i] * p * p
            break
        else:
            # multiplicative step
            J2[x] = J2[i] * (p * p - 1)

# Compute S(N)
ans = 0

for k in range(1, V + 1):
    ans += (N // (k * k)) * J2[k]
    ans %= MOD

print(ans)

Walkthrough

J2[n] stores the Jordan totient $J_2(n)$.
The linear sieve computes all values up to $10^7$.
Finally we evaluate

$$\sum_{k\le 10^7} J_2(k)\left\lfloor\frac{10^{14}}{k^2}\right\rfloor \pmod{1,000,000,007}.$$

Small checks

For $N=10$:

$$S(10)=24$$

For $N=100$:

$$S(100)=767$$

matching the problem statement.

The computation for $N=10^{14}$ gives:

Answer: 94586478