Problem 284: Steady Squares

Project Euler Problem 284 Solution

EulerSolve provides an optimized solution for Project Euler Problem 284, Steady Squares, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary In base $14$, a positive integer $x$ is a steady square of length $n$ if its square ends in the same $n$ base-$14$ digits, i.e. $$x^2 \equiv x \pmod{14^n}.$$ The code sums the base-$14$ digit sums of all positive steady squares with at most $n_{\max}$ digits and prints the final total in base $14$. The Project Euler final value is intentionally omitted here. Mathematical Approach 1) Idempotents modulo $14^n$. The congruence $x^2 \equiv x \pmod{14^n}$ is equivalent to $$x(x-1)\equiv 0 \pmod{2^n7^n}.$$ Because $\gcd(x,x-1)=1$, for the factor $2^n$ one of $x$ and $x-1$ must be divisible by $2^n$, so $$x \equiv 0 \text{ or } 1 \pmod{2^n}.$$ The same argument modulo $7^n$ gives $$x \equiv 0 \text{ or } 1 \pmod{7^n}.$$ By the Chinese Remainder Theorem there are exactly four idempotent residue classes modulo $14^n$: $$ (0,0),\quad (1,1),\quad (0,1),\quad (1,0) $$ with respect to $(\bmod\,2^n,\bmod\,7^n)$. For $n=1$ these are exactly $$0,\ 1,\ 7,\ 8 \pmod{14}.$$ The branches starting at $0$ and $1$ are trivial; the two nontrivial infinite branches start at $7$ and $8$. That is why the implementation only tracks those two sequences, while adding the one-digit steady square $1$ separately. 2) One-step lifting in base $14$....

Detailed mathematical approach

Problem Summary

In base $14$, a positive integer $x$ is a steady square of length $n$ if its square ends in the same $n$ base-$14$ digits, i.e.

$$x^2 \equiv x \pmod{14^n}.$$

The code sums the base-$14$ digit sums of all positive steady squares with at most $n_{\max}$ digits and prints the final total in base $14$. The Project Euler final value is intentionally omitted here.

Mathematical Approach

1) Idempotents modulo $14^n$. The congruence $x^2 \equiv x \pmod{14^n}$ is equivalent to

$$x(x-1)\equiv 0 \pmod{2^n7^n}.$$

Because $\gcd(x,x-1)=1$, for the factor $2^n$ one of $x$ and $x-1$ must be divisible by $2^n$, so

$$x \equiv 0 \text{ or } 1 \pmod{2^n}.$$

The same argument modulo $7^n$ gives

$$x \equiv 0 \text{ or } 1 \pmod{7^n}.$$

By the Chinese Remainder Theorem there are exactly four idempotent residue classes modulo $14^n$:

$$ (0,0),\quad (1,1),\quad (0,1),\quad (1,0) $$

with respect to $(\bmod\,2^n,\bmod\,7^n)$. For $n=1$ these are exactly

$$0,\ 1,\ 7,\ 8 \pmod{14}.$$

The branches starting at $0$ and $1$ are trivial; the two nontrivial infinite branches start at $7$ and $8$. That is why the implementation only tracks those two sequences, while adding the one-digit steady square $1$ separately.

2) One-step lifting in base $14$. Suppose $x$ is already an idempotent modulo

$$M=14^k,$$

so that

$$x^2-x=fM$$

for some integer $f$. Any lift to modulus $14M=14^{k+1}$ has the form

$$x'=x+tM,\qquad t\in\{0,1,\dots,13\}.$$

Expanding gives

$$x'^2-x'=(x^2-x)+(2x-1)tM+t^2M^2=M\bigl(f+(2x-1)t+t^2M\bigr).$$

For divisibility by $14M$, the bracket must vanish modulo $14$. Since $M$ is already a multiple of $14$, the term $t^2M$ disappears modulo $14$, so we only need

$$f+(2x-1)t\equiv 0 \pmod{14}.$$

3) Why the next digit is unique. For an idempotent, $x\equiv 0$ or $1\pmod 2$ and also $x\equiv 0$ or $1\pmod 7$. Therefore

$$2x-1\equiv \pm 1 \pmod 2,\qquad 2x-1\equiv \pm 1 \pmod 7,$$

so $\gcd(2x-1,14)=1$. Hence $2x-1$ has a unique inverse modulo $14$, and the congruence above determines one and only one digit $t\in\{0,\dots,13\}$. This is the whole reason the algorithm has only one successor per branch at each step.

4) The carry-like state update. After choosing $t$, divide the bracket by $14$:

$$f'=\frac{f+(2x-1)t+t^2M}{14}.$$

Then the lifted value satisfies

$$x'^2-x'=f'(14M).$$

So the state for the next step is again just $(x',f',14M)$. This is why the implementation can advance each branch using constant-time arithmetic per length.

5) Why $t$ is the new leading digit. Because $0\le x<M=14^k$, the number $x$ occupies at most $k$ base-$14$ digits. Adding $tM=t\cdot 14^k$ inserts one new digit to the left. Therefore the lift

$$x'=x+t14^k$$

has base-$14$ representation obtained by prefixing the old $k$-digit block with the new digit $t$.

If $t=0$, then the lifted idempotent is still valid modulo $14^{k+1}$, but it is not a genuine $(k+1)$-digit positive number. That is exactly why the code stores leading_digit and ignores states with leading digit $0$ when counting steady squares of a fixed length.

6) Concrete branch examples. Starting from the nontrivial roots modulo $14$:

$$7 \to 37 \to \mathrm{c37} \to 0\mathrm{c37} \to \mathrm{a0c37} \to \cdots$$

and

$$8 \to \mathrm{a8} \to 1\mathrm{a8} \to \mathrm{d1a8} \to 3\mathrm{d1a8} \to \cdots$$

The value $0\mathrm{c37}$ is a valid idempotent modulo $14^4$, but its new leading digit is $0$, so it is not counted as a $4$-digit steady square. The branch rooted at $1$ always lifts with $t=0$, so it only contributes the one-digit number $1$.

7) Digit-sum accumulation. If a branch currently has digit sum $\sigma$ and the next leading digit is $t$, then the new digit sum is simply

$$\sigma'=\sigma+t.$$

So the code never recomputes digit sums from scratch; it just keeps one running total for the $7$-branch and one for the $8$-branch. The overall answer starts from

$$1$$

to count the steady square $1$, then adds the digit sum of each nonzero-leading branch at every length.

Worked Checkpoints

For the first few lengths, the positive steady squares are:

$$ n=1:\ 1,7,8;\qquad n=2:\ 37,\mathrm{a8};\qquad n=3:\ \mathrm{c37},1\mathrm{a8};\qquad n=4:\ \mathrm{d1a8}. $$

The cumulative digit-sum totals are therefore

$$16,\ 44,\ 85,\ 117,\dots$$

and the source file checks that

$$S(9)=582,$$

which is

$$582_{10}=2\mathrm{d}8_{14}.$$

It also verifies directly during the first 20 lifts that each tracked branch still satisfies $x^2\equiv x\pmod{14^n}$.

Complexity Analysis

There are only two nontrivial active branches, and each step performs a constant amount of modular arithmetic and big-integer updates. Therefore the time complexity is

$$O(n_{\max}),$$

and the extra memory is

$$O(1).$$

The values of $x$ themselves become enormous, which is why the code uses cpp_int, but the state space never branches.

Problem 284 source code

C++

#include <array>
#include <cstdint>
#include <iostream>
#include <string>

#include <boost/multiprecision/cpp_int.hpp>

namespace {

using u64 = std::uint64_t;
using boost::multiprecision::cpp_int;

constexpr int kBase = 14;
constexpr int kDefaultN = 10000;

struct Options {
    int n_max = kDefaultN;
    bool run_checkpoints = true;
};

bool parse_int_after_prefix(const std::string& arg, const std::string& prefix, int& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }

    int parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<int>(c - '0');
    }
    value = parsed;
    return true;
}

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

int mod14(const cpp_int& v) {
    cpp_int r = v % kBase;
    if (r < 0) {
        r += kBase;
    }
    return static_cast<int>(r);
}

std::string to_base14(u64 value) {
    constexpr std::array<char, 14> digits = {
        '0', '1', '2', '3', '4', '5', '6', '7', '8', '9', 'a', 'b', 'c', 'd'};
    if (value == 0) {
        return "0";
    }
    std::string out;
    while (value > 0) {
        out.push_back(digits[static_cast<std::size_t>(value % 14)]);
        value /= 14;
    }
    std::reverse(out.begin(), out.end());
    return out;
}

struct SequenceState {
    cpp_int x;
    cpp_int f;
    cpp_int mod;
    int leading_digit;
    u64 digit_sum;

    void step(const std::array<int, kBase>& inv_mod14) {
        const cpp_int x_old = x;
        const cpp_int mod_old = mod;

        const int a = mod14(2 * x_old - 1);
        const int inv = inv_mod14[static_cast<std::size_t>(a)];
        const int f_mod = mod14(f);
        const int t = (kBase - f_mod) % kBase * inv % kBase;

        x = x_old + static_cast<cpp_int>(t) * mod_old;
        f = (f + (2 * x_old - 1) * t + static_cast<cpp_int>(t) * t * mod_old) / kBase;
        mod = mod_old * kBase;

        leading_digit = t;
        digit_sum += static_cast<u64>(t);
    }
};

SequenceState make_initial_state(int root) {
    SequenceState s;
    s.x = root;
    s.mod = kBase;
    s.f = (s.x * s.x - s.x) / kBase;
    s.leading_digit = root;
    s.digit_sum = static_cast<u64>(root);
    return s;
}

u64 solve_sum_of_digit_sums(const int n_max) {
    std::array<int, kBase> inv_mod14{};
    inv_mod14.fill(-1);
    for (int a = 1; a < kBase; ++a) {
        for (int b = 1; b < kBase; ++b) {
            if ((a * b) % kBase == 1) {
                inv_mod14[static_cast<std::size_t>(a)] = b;
                break;
            }
        }
    }

    SequenceState seq7 = make_initial_state(7);
    SequenceState seq8 = make_initial_state(8);

    u64 total = 1;  // steady square "1" for n = 1

    for (int n = 1; n <= n_max; ++n) {
        if (seq7.leading_digit != 0) {
            total += seq7.digit_sum;
        }
        if (seq8.leading_digit != 0) {
            total += seq8.digit_sum;
        }

        if (n == n_max) {
            break;
        }

        seq7.step(inv_mod14);
        seq8.step(inv_mod14);
    }

    return total;
}

bool run_checkpoints() {
    if (solve_sum_of_digit_sums(9) != 582ULL) {
        std::cerr << "Checkpoint failed for n<=9" << '\n';
        return false;
    }

    // Verify steady-square condition during the first few lifts.
    SequenceState seq7 = make_initial_state(7);
    SequenceState seq8 = make_initial_state(8);
    std::array<int, kBase> inv_mod14{};
    inv_mod14.fill(-1);
    for (int a = 1; a < kBase; ++a) {
        for (int b = 1; b < kBase; ++b) {
            if ((a * b) % kBase == 1) {
                inv_mod14[static_cast<std::size_t>(a)] = b;
                break;
            }
        }
    }

    for (int n = 1; n <= 20; ++n) {
        if ((seq7.x * seq7.x - seq7.x) % seq7.mod != 0) {
            std::cerr << "Idempotence checkpoint failed for root 7 at n=" << n << '\n';
            return false;
        }
        if ((seq8.x * seq8.x - seq8.x) % seq8.mod != 0) {
            std::cerr << "Idempotence checkpoint failed for root 8 at n=" << n << '\n';
            return false;
        }
        seq7.step(inv_mod14);
        seq8.step(inv_mod14);
    }

    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;
    }

    const u64 answer = solve_sum_of_digit_sums(options.n_max);
    std::cout << to_base14(answer) << '\n';
    return 0;
}

Python

def mod14(v):
    return v % 14

def to_base14(value):
    digits = "0123456789abcd"
    if value == 0:
        return "0"
    out = []
    while value > 0:
        out.append(digits[value % 14])
        value //= 14
    return "".join(reversed(out))

class SequenceState:
    def __init__(self, root):
        self.x = root
        self.mod = 14
        self.f = (self.x * self.x - self.x) // 14
        self.leading_digit = root
        self.digit_sum = root

    def step(self, inv_mod14):
        x_old = self.x
        mod_old = self.mod
        
        a = mod14(2 * x_old - 1)
        inv = inv_mod14[a]
        f_mod = mod14(self.f)
        t = ((14 - f_mod) % 14 * inv) % 14
        
        self.x = x_old + t * mod_old
        self.f = (self.f + (2 * x_old - 1) * t + t * t * mod_old) // 14
        self.mod = mod_old * 14
        
        self.leading_digit = t
        self.digit_sum += t

def solve_sum_of_digit_sums(n_max):
    inv_mod14 = [-1] * 14
    for a in range(1, 14):
        for b in range(1, 14):
            if (a * b) % 14 == 1:
                inv_mod14[a] = b
                break
                
    seq7 = SequenceState(7)
    seq8 = SequenceState(8)
    
    total = 1
    
    for n in range(1, n_max + 1):
        if seq7.leading_digit != 0:
            total += seq7.digit_sum
        if seq8.leading_digit != 0:
            total += seq8.digit_sum
            
        if n == n_max:
            break
            
        seq7.step(inv_mod14)
        seq8.step(inv_mod14)
        
    return total

def solve():
    ans = solve_sum_of_digit_sums(10000)
    return to_base14(ans)

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

Java

import java.math.BigInteger;

public class Euler284 {

    static final int kBase = 14;

    static String toBase14(long value) {
        char[] digits = "0123456789abcd".toCharArray();
        if (value == 0)
            return "0";
        StringBuilder sb = new StringBuilder();
        while (value > 0) {
            sb.append(digits[(int) (value % 14)]);
            value /= 14;
        }
        return sb.reverse().toString();
    }

    static int mod14(BigInteger v) {
        BigInteger r = v.remainder(BigInteger.valueOf(kBase));
        if (r.signum() < 0) {
            r = r.add(BigInteger.valueOf(kBase));
        }
        return r.intValue();
    }

    static class SequenceState {
        BigInteger x;
        BigInteger f;
        BigInteger mod;
        int leadingDigit;
        long digitSum;

        SequenceState(int root) {
            this.x = BigInteger.valueOf(root);
            this.mod = BigInteger.valueOf(kBase);
            this.f = this.x.multiply(this.x).subtract(this.x).divide(BigInteger.valueOf(kBase));
            this.leadingDigit = root;
            this.digitSum = root;
        }

        void step(int[] invMod14) {
            BigInteger xOld = x;
            BigInteger modOld = mod;

            int a = mod14(xOld.multiply(BigInteger.TWO).subtract(BigInteger.ONE));
            int inv = invMod14[a];
            int fMod = mod14(f);
            int t = ((kBase - fMod) % kBase * inv) % kBase;

            x = xOld.add(modOld.multiply(BigInteger.valueOf(t)));
            BigInteger term2 = xOld.multiply(BigInteger.TWO).subtract(BigInteger.ONE).multiply(BigInteger.valueOf(t));
            BigInteger term3 = modOld.multiply(BigInteger.valueOf((long) t * t));
            f = f.add(term2).add(term3).divide(BigInteger.valueOf(kBase));
            mod = modOld.multiply(BigInteger.valueOf(kBase));

            leadingDigit = t;
            digitSum += t;
        }
    }

    static long solveSumOfDigitSums(int nMax) {
        int[] invMod14 = new int[kBase];
        for (int i = 0; i < kBase; i++)
            invMod14[i] = -1;
        for (int a = 1; a < kBase; a++) {
            for (int b = 1; b < kBase; b++) {
                if ((a * b) % kBase == 1) {
                    invMod14[a] = b;
                    break;
                }
            }
        }

        SequenceState seq7 = new SequenceState(7);
        SequenceState seq8 = new SequenceState(8);

        long total = 1;

        for (int n = 1; n <= nMax; n++) {
            if (seq7.leadingDigit != 0) {
                total += seq7.digitSum;
            }
            if (seq8.leadingDigit != 0) {
                total += seq8.digitSum;
            }

            if (n == nMax)
                break;

            seq7.step(invMod14);
            seq8.step(invMod14);
        }

        return total;
    }

    public static void main(String[] args) {
        System.out.println(toBase14(solveSumOfDigitSums(10000)));
    }
}