Detailed mathematical approach

Problem Summary

Let

$M=11^8=214358881.$

Problem 478 asks for a counting function \(E(n)\) modulo \(M\). The C++, Python, and Java implementations do not enumerate mixtures directly. Instead, they rewrite the problem in terms of primitive lattice data inside the box \(\{0,\dots,n\}^3\), together with a structured correction built from reduced two-dimensional boundary profiles. The checkpoints embedded in the implementations are

$E(1)=103,\qquad E(2)=520447,\qquad E(10)=82608406,\qquad E(500)=13801403.$

Mathematical Approach

The solution is organized around four quantities: a primitive-triple count \(P(n)\), a totient-weighted boundary term \(B(n)\), an inner coprime count \(H_n(L)\), and a final correction sum \(T(n)\). Once these are known, \(E(n)\) follows from one modular identity.

Step 1: Count Primitive Triples with Möbius Inversion

Define

$\mathcal{P}_n=\left\{(a,b,c)\in \{0,\dots,n\}^3\setminus\{(0,0,0)\}:\gcd(a,b,c)=1\right\}.$

Let \(P(n)=|\mathcal{P}_n|\). For a fixed divisor \(d\), the number of triples whose three coordinates are all divisible by \(d\) is

$\left(\left\lfloor\frac{n}{d}\right\rfloor+1\right)^3-1,$

because each coordinate has \(\left\lfloor n/d\right\rfloor+1\) multiples of \(d\) between \(0\) and \(n\), and the all-zero triple must be excluded. Möbius inversion therefore gives

$P(n)=\sum_{d=1}^{n}\mu(d)\left(\left(\left\lfloor\frac{n}{d}\right\rfloor+1\right)^3-1\right).$

This is the large exponent that ultimately drives the term \(2^{P(n)}\).

Step 2: Explain the Totient Weight \(6\varphi(L)\)

For each \(L\ge 1\), the number of coprime positive pairs \((u,v)\) satisfying

$u+v=L$

is exactly \(\varphi(L)\). Indeed, choosing \(u\) with \(1\le u\le L\) and \(\gcd(u,L)=1\) uniquely determines \(v=L-u\), and then \(\gcd(u,v)=1\). The final correction formula treats each reduced pair with a sixfold symmetry factor, which explains the coefficient \(6\varphi(L)\). Hence the global boundary weight is

$B(n)=6\sum_{L=1}^{n}\varphi(L).$

Step 3: Derive the Inner Count \(H_n(L)\)

Fix \(L\). The inner quantity used by the implementations is

$H_n(L)=1+\#\left\{(x,y)\in \mathbb{Z}_{>0}^2:\gcd(x,y)=1,\ x+Ly\le n\right\}.$

The initial \(1\) is explicit. To count the coprime pairs, apply Möbius inversion again. If \(d\mid x\) and \(d\mid y\), write \(x=da\), \(y=db\). Then

$a+Lb\le \left\lfloor\frac{n}{d}\right\rfloor.$

Set

$q_{d,L}=\left\lfloor\frac{n}{Ld}\right\rfloor.$

For each \(b=1,\dots,q_{d,L}\), the variable \(a\) has \(\left\lfloor n/d\right\rfloor-Lb\) positive choices. Summing those choices yields

$\sum_{b=1}^{q_{d,L}}\left(\left\lfloor\frac{n}{d}\right\rfloor-Lb\right)=q_{d,L}\left\lfloor\frac{n}{d}\right\rfloor-L\frac{q_{d,L}(q_{d,L}+1)}{2}.$

Therefore

$H_n(L)=1+\sum_{d=1}^{\lfloor n/L\rfloor}\mu(d)\left(q_{d,L}\left\lfloor\frac{n}{d}\right\rfloor-L\frac{q_{d,L}(q_{d,L}+1)}{2}\right),\qquad q_{d,L}=\left\lfloor\frac{n}{Ld}\right\rfloor.$

This matches the inner loop term for term.

Step 4: Reduce the Huge Exponents Safely

Let

$\Phi=\varphi(M)=\varphi(11^8)=10\cdot 11^7=194871710.$

Since \(\gcd(2,M)=1\), Euler's theorem tells us that every power of \(2\) modulo \(M\) depends only on the exponent modulo \(\Phi\). However, the correction terms also need

$\frac{P(n)-1}{2},$

so parity must be preserved before dividing by \(2\). For that reason the implementations first compute \(P(n)\) modulo \(2\Phi\), subtract \(1\), verify that the result is even, and only then reduce the half-exponent modulo \(\Phi\).

Step 5: Assemble the Final Formula

With the reduced half-exponent understood modulo \(\Phi\), define

$T(n)\equiv \sum_{L=1}^{n}6\varphi(L)\,2^{\left(\frac{P(n)-1}{2}-H_n(L)\right)\bmod \Phi}\pmod{M}.$

Then the shared correction term is

$F(n)\equiv 1+2^{\left(\frac{P(n)-1}{2}\right)\bmod \Phi}B(n)-T(n)\pmod{M},$

and the required value is

$\boxed{E(n)\equiv 2^{P(n)\bmod \Phi}-F(n)\pmod{M}.}$

This is exactly the algebra used by the C++, Python, and Java implementations.

Worked Example: \(n=2\)

Here \(\mu(1)=1\), \(\mu(2)=-1\), and \(\varphi(1)=\varphi(2)=1\). First,

$P(2)=\left((2+1)^3-1\right)-\left((1+1)^3-1\right)=26-7=19.$

Hence

$\frac{P(2)-1}{2}=9,\qquad B(2)=6(\varphi(1)+\varphi(2))=12.$

For \(L=1\),

$H_2(1)=1+\left(2\cdot 2-\frac{2\cdot 3}{2}\right)-\left(1\cdot 1-\frac{1\cdot 2}{2}\right)=2.$

For \(L=2\),

$H_2(2)=1+\left(1\cdot 2-2\cdot\frac{1\cdot 2}{2}\right)=1.$

Therefore

$T(2)=6\cdot 2^{9-2}+6\cdot 2^{9-1}=768+1536=2304,$

$F(2)=1+2^9\cdot 12-2304=3841,$

and

$E(2)=2^{19}-3841=520447.$

This agrees with the checkpoint used by the implementations.

How the Code Works

The implementation begins with a linear sieve that computes \(\mu(d)\) and \(\varphi(d)\) for every \(d\le n\). It also stores the quotient table \(\left\lfloor n/d\right\rfloor\), because both \(P(n)\) and \(H_n(L)\) reuse those values constantly.

Next, it evaluates \(P(n)\) modulo \(2\Phi\), accumulates \(B(n)=6\sum_{L\le n}\varphi(L)\) modulo \(M\), and precomputes modular powers of \(2\) by fast exponentiation. Then it loops over \(L=1,\dots,n\), computes \(H_n(L)\) from the Möbius-weighted inner sum, forms the exponent residue \(\left(\frac{P(n)-1}{2}-H_n(L)\right)\bmod \Phi\), and adds the contribution \(6\varphi(L)\cdot 2^{((P(n)-1)/2-H_n(L))\bmod\Phi}\) to \(T(n)\).

Finally, it combines \(2^{P(n)\bmod\Phi}\), \(B(n)\), and \(T(n)\) through the boxed formula. The C++ and Java implementations parallelize the outer loop over \(L\), while the Python implementation evaluates the same mathematics sequentially.

Complexity Analysis

The sieve and the quotient table each take \(O(n)\) time and \(O(n)\) memory. The dominant cost is the double summation for \(H_n(L)\): for each \(L\), the inner loop runs to \(\left\lfloor n/L\right\rfloor\), so the total number of iterations is

$\sum_{L=1}^{n}\left\lfloor\frac{n}{L}\right\rfloor=n\log n+O(n).$

Thus the overall running time is \(O(n\log n)\) arithmetic operations with \(O(n)\) memory. Parallel execution improves wall-clock time but does not change the asymptotic complexity.

Footnotes and References

1. Problem page: https://projecteuler.net/problem=478
2. Möbius inversion formula: Wikipedia — Möbius inversion formula
3. Euler's totient function: Wikipedia — Euler's totient function
4. Inclusion-exclusion principle: Wikipedia — Inclusion-exclusion principle
5. Greatest common divisor: Wikipedia — Greatest common divisor