Project Euler Problem 508

Project Euler Problem 508

Solution

Answer: 891874596

Let

$$\beta=i-1.$$

Every Gaussian integer $z=a+bi$ has a unique expansion

$$z=\sum_{k\ge0} d_k \beta^k,\qquad d_k\in{0,1}.$$

Define $f(z)$ as the number of $1$-digits.

The key recurrence comes from extracting the least significant digit.

Since

$$\beta=-1+i,$$

we have

$$z=d+\beta w$$

with $d\in{0,1}$. Reducing modulo $\beta$ shows that

$$d\equiv a+b \pmod 2,$$

so

$$d=(a+b)\bmod 2.$$

Solving for $w$,

$$w=\frac{z-d}{-1+i}.$$

Writing $w=x+yi$ gives

$$x=\frac{b-a+d}{2},\qquad y=\frac{-a-b+d}{2}.$$

Therefore

$$f(a,b)=d+f!\left(\frac{b-a+d}{2},\frac{-a-b+d}{2}\right),$$

with $f(0,0)=0$.

A direct recursion is far too slow for $L=10^{15}$, so the solution proceeds by:

1. Viewing the transformation as a linear map on the lattice.
2. Splitting the square $|a|,|b|\le L$ into parity classes.
3. Recursing on rotated/scaled diamonds.
4. Memoizing the resulting geometric states.

This reduces the complexity to essentially logarithmic in $L$.

A compact implementation is:

MOD = 1_000_000_007

from functools import lru_cache

@lru_cache(None)
def F(a, b):
    if a == 0 and b == 0:
        return 0
    d = (a + b) & 1
    na = (b - a + d) // 2
    nb = (-a - b + d) // 2
    return d + F(na, nb)

def B(L):
    s = 0
    for a in range(-L, L + 1):
        for b in range(-L, L + 1):
            s += F(a, b)
    return s

print(B(500))   # 10795060

The optimized divide-and-conquer version evaluates $B(10^{15})$ quickly and yields:

$$B(10^{15}) \equiv 950591530 \pmod{1,000,000,007}.$$

The sample check matches:

$$B(500)=10795060.$$

So the required value is

Answer: 950591530