Problem 512: Sums of Totients of Powers

View on Project Euler

Project Euler Problem 512 Solution

EulerSolve provides an optimized solution for Project Euler Problem 512, Sums of Totients of Powers, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary For this problem, the required total can be rewritten as the odd totient prefix sum $$G(n)=\sum_{\substack{1\le m\le n\\ m\text{ odd}}}\varphi(m),$$ where \(\varphi\) is Euler's totient function. At the target scale \(n=500000000\), directly evaluating every odd \(\varphi(m)\) is too expensive. The implementation therefore replaces the raw summation by a divisor identity, a closed formula for an auxiliary prefix sum, and a memoized divide-and-conquer recurrence. Mathematical Approach We keep the notation $$G(n)=\sum_{\substack{1\le m\le n\\ m\text{ odd}}}\varphi(m)$$ and define the odd part of a positive integer by $$\operatorname{odd}(x)=\frac{x}{2^{v_2(x)}}.$$ The auxiliary summatory function used by the implementation is $$A(n)=\sum_{x=1}^{n}\operatorname{odd}(x).$$ Step 1: Relate Odd Divisors to the Odd Part Every odd divisor of \(x\) is exactly a divisor of \(\operatorname{odd}(x)\). The classical identity $$\sum_{d\mid m}\varphi(d)=m$$ therefore gives $$\sum_{\substack{d\mid x\\ d\text{ odd}}}\varphi(d)=\sum_{d\mid \operatorname{odd}(x)}\varphi(d)=\operatorname{odd}(x).$$ This is the key arithmetic fact behind the whole algorithm: instead of summing totients directly, we first sum the odd parts of all integers up to \(n\)....

Detailed mathematical approach

Problem Summary

For this problem, the required total can be rewritten as the odd totient prefix sum

$$G(n)=\sum_{\substack{1\le m\le n\\ m\text{ odd}}}\varphi(m),$$

where \(\varphi\) is Euler's totient function. At the target scale \(n=500000000\), directly evaluating every odd \(\varphi(m)\) is too expensive. The implementation therefore replaces the raw summation by a divisor identity, a closed formula for an auxiliary prefix sum, and a memoized divide-and-conquer recurrence.

Mathematical Approach

We keep the notation

$$G(n)=\sum_{\substack{1\le m\le n\\ m\text{ odd}}}\varphi(m)$$

and define the odd part of a positive integer by

$$\operatorname{odd}(x)=\frac{x}{2^{v_2(x)}}.$$

The auxiliary summatory function used by the implementation is

$$A(n)=\sum_{x=1}^{n}\operatorname{odd}(x).$$

Step 1: Relate Odd Divisors to the Odd Part

Every odd divisor of \(x\) is exactly a divisor of \(\operatorname{odd}(x)\). The classical identity

$$\sum_{d\mid m}\varphi(d)=m$$

therefore gives

$$\sum_{\substack{d\mid x\\ d\text{ odd}}}\varphi(d)=\sum_{d\mid \operatorname{odd}(x)}\varphi(d)=\operatorname{odd}(x).$$

This is the key arithmetic fact behind the whole algorithm: instead of summing totients directly, we first sum the odd parts of all integers up to \(n\).

Step 2: Convert the Divisor Sum into a Summatory Recurrence

Summing the previous identity over \(1\le x\le n\) and exchanging the order of summation yields

$$A(n)=\sum_{x=1}^{n}\sum_{\substack{d\mid x\\ d\text{ odd}}}\varphi(d)=\sum_{\substack{d\le n\\ d\text{ odd}}}\varphi(d)\left\lfloor\frac{n}{d}\right\rfloor.$$

Now use the standard floor-sum rearrangement

$$\sum_{d\le n} a(d)\left\lfloor\frac{n}{d}\right\rfloor=\sum_{t=1}^{n}\sum_{d\le \lfloor n/t\rfloor}a(d).$$

With \(a(d)=\varphi(d)\) for odd \(d\) and \(a(d)=0\) for even \(d\), this becomes

$$A(n)=\sum_{t=1}^{n}G\left(\left\lfloor\frac{n}{t}\right\rfloor\right).$$

The term for \(t=1\) is exactly \(G(n)\), so we isolate it and obtain the recurrence

$$G(n)=A(n)-\sum_{t=2}^{n}G\left(\left\lfloor\frac{n}{t}\right\rfloor\right).$$

Step 3: Compute \(A(n)\) in Closed Form

Write every positive integer uniquely as \(x=2^j u\) with \(u\) odd. Then \(\operatorname{odd}(x)=u\), so grouping by the exponent \(j\) gives

$$A(n)=\sum_{j\ge 0}\ \sum_{\substack{u\le \lfloor n/2^j\rfloor\\ u\text{ odd}}}u.$$

If \(r=\left\lfloor\frac{u+1}{2}\right\rfloor\), then the odd numbers up to \(u\) are \(1,3,\dots,2r-1\), and

$$1+3+\cdots+(2r-1)=r^2.$$

Therefore the inner sum has the closed form

$$\sum_{\substack{v\le u\\ v\text{ odd}}}v=\left\lfloor\frac{u+1}{2}\right\rfloor^2,$$

and hence

$$A(n)=\sum_{j\ge 0}\left(\left\lfloor\frac{\lfloor n/2^j\rfloor+1}{2}\right\rfloor\right)^2.$$

Only \(O(\log n)\) terms are present, because \(\lfloor n/2^j\rfloor=0\) once \(2^j>n\).

Step 4: Group Equal Quotients

The recurrence still seems to sum over all \(t\), but the quotient \(\left\lfloor n/t\right\rfloor\) repeats on long intervals. If

$$q=\left\lfloor\frac{n}{\ell}\right\rfloor,$$

then the same quotient persists for every

$$t\in[\ell,r],\qquad r=\left\lfloor\frac{n}{q}\right\rfloor.$$

So one whole block contributes

$$\sum_{t=\ell}^{r}G\left(\left\lfloor\frac{n}{t}\right\rfloor\right)=(r-\ell+1)\,G(q).$$

This reduces the cost of one recursive evaluation from linear in \(n\) to the number of distinct quotients, which is \(O(\sqrt n)\).

Step 5: Memoize the Recurrence

The same quotient values appear again and again inside different blocks, so the recurrence is evaluated only once for each distinct argument. After a value \(G(q)\) has been computed, later requests reuse the cached result instead of recomputing the entire subtree.

This turns the formula above into a practical divide-and-conquer algorithm for very large \(n\).

Worked Example: \(n=10\)

First compute \(A(10)\) from the closed form:

$$A(10)=\left\lfloor\frac{10+1}{2}\right\rfloor^2+\left\lfloor\frac{5+1}{2}\right\rfloor^2+\left\lfloor\frac{2+1}{2}\right\rfloor^2+\left\lfloor\frac{1+1}{2}\right\rfloor^2=25+9+1+1=36.$$

Now use quotient blocks for \(t\ge 2\):

$$\left\lfloor\frac{10}{t}\right\rfloor=5,3,2,2,1,1,1,1,1 \qquad (t=2,\dots,10).$$

Hence

$$G(10)=36-\bigl(G(5)+G(3)+2G(2)+5G(1)\bigr).$$

The smaller values are

$$G(1)=1,\qquad G(2)=1,\qquad G(3)=1+\varphi(3)=3,\qquad G(5)=1+\varphi(3)+\varphi(5)=7.$$

Therefore

$$G(10)=36-(7+3+2+5)=19.$$

A direct check agrees:

$$\varphi(1)+\varphi(3)+\varphi(5)+\varphi(7)+\varphi(9)=1+2+4+6+6=19.$$

How the Code Works

The C++, Python, and Java implementations all follow the same structure. For a requested \(n\), the implementation first computes \(A(n)\) by repeatedly halving the current bound and adding the square of the number of odd integers up to that bound. It then evaluates the recurrence for \(G(n)\) by scanning quotient blocks \([\ell,r]\), recursively requesting the smaller value \(G(q)\), and subtracting \((r-\ell+1)G(q)\) from the running total.

A cache keyed by the argument stores every completed subproblem, so each distinct value is solved only once. The three language versions differ only in surface syntax and data-structure choice; mathematically they implement the same recurrence, the same block decomposition, and the same exact integer arithmetic.

Complexity Analysis

Computing \(A(n)\) costs \(O(\log n)\) time because the argument is halved at each step. A single recursive evaluation on argument \(m\) processes one interval for each distinct quotient \(\lfloor m/t\rfloor\), so that call uses \(O(\sqrt m)\) block iterations. If \(M\) denotes the set of distinct cached arguments generated from the original input, the total running time is

$$O\left(\sum_{m\in M}\sqrt{m}\right),$$

and the memory usage is \(O(|M|)\). In particular, the top-level quotient list already has only \(O(\sqrt n)\) distinct values, which is why this approach is dramatically faster than summing \(\varphi(m)\) for every odd \(m\le n\).

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=512
  2. Euler's totient function and the divisor identity \(\sum_{d\mid m}\varphi(d)=m\).
  3. The odd part of an integer, \(\operatorname{odd}(x)=x/2^{v_2(x)}\).
  4. Dirichlet-style divisor reordering for summatory functions.
  5. Floor-division block decomposition in recursive arithmetic summation.

Problem 512 source code

C++

#include <cstdint>
#include <iostream>
#include <string>
#include <unordered_map>

namespace {

using u64 = std::uint64_t;
using i128 = __int128_t;

struct Options {
    u64 n = 500'000'000ULL;
    bool run_checkpoints = true;
};

bool parse_u64_after_prefix(const std::string& arg, const std::string& prefix, u64& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }
    u64 parsed = 0ULL;
    for (char ch : tail) {
        if (ch < '0' || ch > '9') {
            return false;
        }
        parsed = parsed * 10ULL + static_cast<u64>(ch - '0');
    }
    value = parsed;
    return true;
}

bool parse_arguments(int argc, char** argv, Options& options) {
    for (int i = 1; i < argc; ++i) {
        std::string arg(argv[i]);
        if (arg == "--skip-checkpoints") {
            options.run_checkpoints = false;
            continue;
        }
        if (parse_u64_after_prefix(arg, "--n=", options.n)) {
            continue;
        }
        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.n >= 1ULL;
}

class OddTotientSummatory {
public:
    u64 sum(const u64 n) {
        if (n == 0ULL) {
            return 0ULL;
        }
        auto it = memo_.find(n);
        if (it != memo_.end()) {
            return it->second;
        }

        i128 result = static_cast<i128>(odd_part_prefix_sum(n));
        for (u64 l = 2ULL; l <= n;) {
            const u64 q = n / l;
            const u64 r = n / q;
            result -= static_cast<i128>(r - l + 1ULL) * static_cast<i128>(sum(q));
            l = r + 1ULL;
        }

        const u64 out = static_cast<u64>(result);
        memo_.emplace(n, out);
        return out;
    }

private:
    static u64 odd_numbers_sum(const u64 x) {
        const u64 m = (x + 1ULL) / 2ULL;
        return m * m;
    }

    static u64 odd_part_prefix_sum(u64 n) {
        u64 total = 0ULL;
        while (n > 0ULL) {
            total += odd_numbers_sum(n);
            n >>= 1ULL;
        }
        return total;
    }

    std::unordered_map<u64, u64> memo_;
};

u64 solve(const u64 n) {
    OddTotientSummatory solver;
    return solver.sum(n);
}

u64 brute_phi_odd_sum(const u64 n) {
    auto phi = [](u64 x) {
        u64 result = x;
        u64 t = x;
        for (u64 p = 2ULL; p * p <= t; ++p) {
            if (t % p != 0ULL) {
                continue;
            }
            while (t % p == 0ULL) {
                t /= p;
            }
            result -= result / p;
        }
        if (t > 1ULL) {
            result -= result / t;
        }
        return result;
    };
    u64 s = 0ULL;
    for (u64 i = 1ULL; i <= n; i += 2ULL) {
        s += phi(i);
    }
    return s;
}

bool run_checkpoints() {
    if (solve(100ULL) != 2007ULL) {
        std::cerr << "Checkpoint failed: g(100)=2007" << '\n';
        return false;
    }
    if (solve(10'000ULL) != brute_phi_odd_sum(10'000ULL)) {
        std::cerr << "Checkpoint failed: brute-force odd-phi sum for n=10000" << '\n';
        return false;
    }
    return true;
}

}  // namespace

int main(int argc, char** argv) {
    Options options;
    if (!parse_arguments(argc, argv, options)) {
        return 1;
    }
    if (options.run_checkpoints && !run_checkpoints()) {
        return 2;
    }
    std::cout << solve(options.n) << '\n';
    return 0;
}

Python

import sys

# Increase recursion depth just in case, though it only goes to O(log N) depth.
sys.setrecursionlimit(2000)

memo = {}

def odd_numbers_sum(x):
    m = (x + 1) // 2
    return m * m

def odd_part_prefix_sum(n):
    total = 0
    x = n
    while x > 0:
        total += odd_numbers_sum(x)
        x //= 2
    return total

def get_sum(n):
    if n == 0:
        return 0
    if n in memo:
        return memo[n]
    
    result = odd_part_prefix_sum(n)
    l = 2
    while l <= n:
        q = n // l
        r = n // q
        result -= (r - l + 1) * get_sum(q)
        l = r + 1
        
    memo[n] = result
    return result

def solve():
    n = 500000000
    ans = get_sum(n)
    return str(ans)

if __name__ == '__main__':
    print(solve())

Java

import java.util.HashMap;
import java.util.Map;

public class Euler512 {

    static Map<Integer, Long> memo = new HashMap<>();

    static long oddNumbersSum(long x) {
        long m = (x + 1) / 2;
        return m * m;
    }

    static long oddPartPrefixSum(long n) {
        long total = 0;
        long x = n;
        while (x > 0) {
            total += oddNumbersSum(x);
            x /= 2;
        }
        return total;
    }

    static long getSum(int n) {
        if (n == 0) {
            return 0;
        }
        if (memo.containsKey(n)) {
            return memo.get(n);
        }

        long result = oddPartPrefixSum(n);
        int l = 2;
        while (l <= n) {
            int q = n / l;
            int r = n / q;
            result -= (long) (r - l + 1) * getSum(q);
            l = r + 1;
        }

        memo.put(n, result);
        return result;
    }

    public static void main(String[] args) {
        int n = 500000000;
        long ans = getSum(n);
        System.out.println(ans);
    }
}