Problem 125: Palindromic Sums

Project Euler Problem 125 Solution

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

Problem Summary We seek every palindromic integer $n<10^8$ that can be written as a sum of consecutive positive squares with at least two terms: $$n=a^2+(a+1)^2+\cdots+b^2,\qquad 1\le a<b.$$ The final answer is the sum of those palindromic values, not the number of representations. That distinction matters because the same palindrome may arise from more than one interval of consecutive squares, so valid values must be counted only once. Mathematical Approach The implementations do not search blindly. They convert every consecutive-square sum into a difference of square-prefix totals, use monotonic growth to stop each search branch as soon as it becomes too large, and store successful palindromes in a set so repeated representations do not change the final total. Interval sums from square-prefix totals Define the prefix sum of squares $$P_n=\sum_{k=1}^{n} k^2=\frac{n(n+1)(2n+1)}{6}.$$ Then every block of consecutive squares from $a$ through $b$ can be written as $$S(a,b)=\sum_{k=a}^{b} k^2=P_b-P_{a-1}=\frac{b(b+1)(2b+1)-(a-1)a(2a-1)}{6}.$$ This is the central algebraic object of the problem. Once the prefix totals are known, the sum for any interval $[a,b]$ is obtained in constant time from one subtraction instead of recomputing every square inside that interval....

Detailed mathematical approach

Problem Summary

A palindrome is a positive integer whose decimal digits read the same in both directions, such as $5$, $77$ or $595$. We must find the sum of all palindromes $n\lt10^8$ that can be written as a sum of two or more consecutive squares of positive integers,

$$\begin{aligned}n&=a^2+(a+1)^2+\dots+b^2,\\1&\le a\lt b,\end{aligned}$$

where a palindrome with several such representations is added only once. The statement gives the example $595=6^2+7^2+\dots+12^2$, says that exactly eleven palindromes below $1000$ have the property and that their sum is $4164$, and excludes $1=0^2+1^2$, because only squares of positive integers are allowed. ^[1] Section 1 shows that this count of eleven also rules out a single square on its own.

The whole solution in one sentence: obtain every sum of two or more consecutive squares below $10^8$ as a difference of two running totals of squares, ending each run of longer and longer sums as soon as it reaches $10^8$, and add up the palindromic sums, each value only once.

Mathematical Approach

The work lies in organising the search so that it is finite, misses nothing and counts every palindrome once; Section 7 proves that the method does all three.

1. Small cases by hand

Start with $1^2$ and keep adding the next square: $5$, $14$, $30$, $55$, $91$, $140$, $204$, $285$, $385$, $506$, $650$, $819$, and the next sum, $1015$, is above the limit $1000$. Starting with $2^2$ gives $13$, $29$, $54$, $90$, $139$, $203$, $284$, $384$, $505$, $649$, $818$. Doing the same for every start and keeping the palindromes yields exactly eleven numbers below $1000$:

$5=1^2+2^2$, $55=1^2+\dots+5^2$, $77=4^2+5^2+6^2$, $181=9^2+10^2$;
$313=12^2+13^2$, $434=11^2+12^2+13^2$, $505=2^2+\dots+11^2$, $545=16^2+17^2$;
$595=6^2+\dots+12^2$, $636=4^2+\dots+12^2$, $818=2^2+\dots+13^2$.

Their sum is $4164$, as the statement says. The palindromic squares $1$, $4$, $9$, $121$, $484$ and $676$ are missing: none of them is a sum of two or more consecutive squares, and if a single square counted, there would be $17$ palindromes below $1000$ instead of eleven. The same eleven numbers begin the OEIS sequence A216446. ^[5]

2. Running totals of squares

Let $P_m=1^2+2^2+\dots+m^2$, with $P_0=0$, be the running total of the squares. A run of squares from $a^2$ to $b^2$ is what remains of $P_b$ after its first $a-1$ terms are removed:

$$\begin{aligned}S(a,b)&=a^2+(a+1)^2+\dots+b^2\\&=P_b-P_{a-1}.\end{aligned}$$

With a table of $P_0,P_1,\dots,P_m$, a so-called prefix-sum table, ^[3] every such sum costs one subtraction; for the example of the statement, $S(6,12)=P_{12}-P_5$${}=650-55=595$. The running totals are the square pyramidal numbers, with the closed form ^[2]

$$P_m=\frac{m(m+1)(2m+1)}{6}.$$

It holds for $m=0$, and from $m-1$ to $m$ the right-hand side grows by exactly $m^2$: after the common factor $m/6$ is taken out of the difference, what remains is

$$\begin{aligned}&(m+1)(2m+1)\\&\quad{}-(m-1)(2m-1)=6m.\end{aligned}$$

By induction the formula holds for every $m$. The programs build their table by adding one square at a time, but the bounds below need the formula.

3. How far the search must go

Fix the limit $L=10^8$ and call a sum $S(a,b)$ with $a\lt b$ useful if it is below $L$. A useful sum contains its last two terms, so $(b-1)^2+b^2\lt L$. Since $7071^2+7072^2=100012225$ is already too large, every useful end satisfies $b\le7071$ and every useful start $a\le7070$; indeed $7070^2+7071^2=99983941\lt L$. The number of terms $k=b-a+1$ is bounded too: the $j$-th term $(a+j-1)^2$ is at least $j^2$, so $S(a,b)\ge P_k$, and $P_{668}=99582434$ while $P_{669}=100029995$. No useful sum has more than $668$ terms, so the search is finite.

4. Stopping early

For a fixed start $a$, each longer run adds one more positive square,

$$\begin{aligned}S(a,b+1)&=S(a,b)+(b+1)^2\\&\gt S(a,b),\end{aligned}$$

so the sums from one start increase strictly. Once one of them reaches $L$, every longer run from the same start does too, and the search can go on with the next start. For $L=1000$ and $a=2$, the run $13,29,\dots,818$ of Section 1 is followed by $818+14^2=1014$: after eleven useful sums, one sum overshoots and the run ends.

5. Recognising palindromes

To test a number, write its decimal digits as a string, reverse the string and compare: $595$ stays $595$, while $149=6^2+7^2+8^2$ becomes $941$. Palindromes are rare. One with $d$ digits is fixed by its first $\lceil d/2\rceil$ digits, the first of them not $0$, so there are $9\cdot10^{\lceil d/2\rceil-1}$ of them: $9$, $9$, $90$, $90$, $900$, $900$, $9000$ and $9000$ for $d=1,\dots,8$, together $19998$ palindromes below $10^8$. ^[4] One could also start from them and test each one; the complexity section compares the two routes.

6. Counting each palindrome once

A palindrome can be a sum of consecutive squares in more than one way. The smallest example is

$$\begin{aligned}554455&=9^2+10^2+\dots+118^2\\&=331^2+332^2+\dots+335^2,\end{aligned}$$

and below $10^8$ there is exactly one more, ^[6]

$$\begin{aligned}&9343439\\&\quad=102^2+103^2+\dots+307^2\\&\quad=657^2+658^2+\dots+677^2.\end{aligned}$$

No palindrome below $10^8$ has three representations. The task asks for a sum of numbers, not of representations, so a running total that adds every palindromic sum as soon as it is found would be wrong: for the limit $10^6$ it would give $15447478$ instead of the correct $14893023$, exactly $554455$ too much. The programs therefore collect the palindromic sums in a set, which keeps one copy of each value, and add up the set at the end.

7. Why the search is complete and exact

Claim. Let $M$ be the least positive integer with $M^2\ge L$. For every start $a=1,2,\dots,M$, try the ends $b=a+1,a+2,\dots,M$ in increasing order until $S(a,b)\ge L$, and put every palindromic sum below $L$ into a set. Then the set contains exactly the palindromes below $L$ that are sums of two or more consecutive squares of positive integers.

Proof. Every value put into the set is a palindrome below $L$ of the form $S(a,b)$ with $1\le a\lt b$. Conversely, let $n\lt L$ be such a palindrome, $n=S(a,b)$. Then $b^2\lt n\lt L\le M^2$, so $b\lt M$, and the table reaches $P_b$. In the pass of the start $a$, every earlier end $e\lt b$ gives $S(a,e)\lt S(a,b)=n\lt L$ by Section 4, so the pass does not stop before $e=b$, where it finds $n$ and puts it into the set. A set ignores a second insertion of a value it already holds, so every such palindrome is in the set exactly once, and the sum of the set is the requested total. $\square$

How the Code Works

In the C++ program, ^[9] parse_arguments accepts --skip-checkpoints and --limit= followed by digits, which parse_i64_after_prefix reads; any other argument, or a limit below $2$, gives exit code $1$. Unless it is skipped, run_checkpoints requires solve(1000) to return the statement's $4164$, with exit code $2$ otherwise, and then main prints solve(limit), by default for $L=10^8$.

solve carries out Sections 2 to 7. It raises max_base until its square reaches the limit, which gives $M=10000$, and fills the std::vector<i64> named prefix with $P_0,\dots,P_M$. The outer loop runs start from $1$ to max_base, the inner loop end from start + 1; the sum is prefix[end] - prefix[start - 1], the inner loop breaks at the first sum that reaches the limit, and every palindromic sum goes into the std::set<i64> named values, whose elements are added up at the end. is_palindrome converts the number with std::to_string, reverses a copy with std::reverse and compares. The table is longer than Section 3 requires, so each start from $7071$ to $9999$ makes one useless step, which changes nothing.

The Python and Java programs run the same loops with the limit fixed, without arguments or checkpoint. Python reverses the string with ss[::-1] and collects the values in a built-in set; Java's isPalin uses StringBuilder.reverse, and the values go into a HashSet<Long>. The C++ std::set is a sorted tree with logarithmic insertion, ^[7] the other two sets are hash tables, ^[8] and all three keep one copy of each value. C++ and Java store the table as 64-bit integers, i64 and long: its last entry $P_{10000}=333383335000$ exceeds $2^{31}-1$, although every kept sum is below $10^8$. Python's integers have no fixed size. The C++ parser does not check the range of the digits after --limit=, and for limits above $3024616^2\approx9.1\cdot10^{12}$ the table would exceed the signed 64-bit range; both cases are far from the default run.

Complexity and Verification

With $M=\lceil\sqrt L\,\rceil$, the table takes $M$ additions; for $L=10^8$ the C++ and Java tables hold $10001$ entries of $8$ bytes, about $80$ kB. The loops compute one sum for every useful pair $(a,b)$ and one overshooting sum for every start below $M$. The number of useful pairs is at most of order $L^{2/3}$: a run of $k$ terms from the start $a$ is at least $k\,a^2$ and at least $P_k\ge k^3/3$, so a start has fewer than $\min\bigl(L/a^2,(3L)^{1/3}\bigr)$ useful ends. With $T=\lceil L^{1/3}\rceil$, the starts $a\le T$ therefore give fewer than $T(3L)^{1/3}$, about $3^{1/3}L^{2/3}$, pairs, and the starts $a\gt T$ fewer than $L/T\le L^{2/3}$, because $1/a^2\lt1/(a-1)-1/a$ and these differences add up to $1/T$.

For $L=10^8$ the loops compute $371234$ useful sums and $9999$ overshoots, $381233$ in all; $168$ sums are palindromes, and the set ends with $166$ values. The set never holds more entries than there are palindromes below $L$, $O(\sqrt L)$ by Section 5, and an insertion costs at most $O(\log L)$. Counting arithmetic on the stored numbers as one step, the time is therefore $O(L^{2/3}\log L)$, as a palindrome test reads $O(\log L)$ digits, and the memory $O(\sqrt L)$; on a desktop computer each program finishes in less than a tenth of a second. The palindrome route of Section 5 would need $7124319$ integer square roots, even after skipping the lengths $k$ with $P_k\gt n$.

The C++ program checks the total $4164$ for $L=1000$, where $M=32$ and the loops compute $106$ useful sums and $31$ overshoots. For this article all three programs were run and print the same number. An independent program took the palindrome route: for each palindrome $n\lt10^8$ and each number of terms $k$ from $2$ to $668$, it solved $S(a,a+k-1)=n$, that is

$$k\,a^2+k(k-1)\,a+P_{k-1}=n,$$

for a positive integer $a$. Such an $a$ exists exactly when the discriminant $D=k^2(k-1)^2-4k(P_{k-1}-n)$ is a perfect square $q^2$ and $q-k(k-1)$ is a positive multiple of $2k$; for $n=595$ and $k=7$ this gives $D=15876=126^2$ and $a=(126-42)/14=6$. It found the same $166$ palindromes, the two of Section 6 with two representations among them, and the same total, and its list agrees term by term with the table of A216446 below $10^8$. ^[5] A sliding window over the squares, applied to every palindrome below $10^6$, gives the totals $4164$, $20268$, $667526$ and $14893023$ for the limits $10^3$ to $10^6$, the same as the C++ program with --limit=. These checks are not part of the published solutions.

Footnotes and References

The references below supply the background facts used above. The search bounds, the stopping rule, the completeness proof and the operation counts are derived in this article.

Project Euler 125 — Palindromic Sums. The official statement writes 595 as 6² + 7² + ⋯ + 12², says that exactly eleven palindromes below one thousand are such sums and add up to 4164, excludes 1 = 0² + 1², and asks for the sum of all palindromes below 10⁸ that are sums of consecutive squares.
Wikipedia — Square pyramidal number. The numbers 1, 5, 14, 30, 55, 91, … are the sums 1² + 2² + ⋯ + n² = n(n + 1)(2n + 1)/6; the article notes that the formula may be proved by induction.
Wikipedia — Prefix sum. Prefix sums are the running totals of a sequence; summed-area tables use them to compute sums of rectangular blocks.
Wikipedia — Palindromic number. A palindromic number stays the same when its digits are reversed; the article's table counts 19999 palindromic natural numbers below 10⁸, including 0.
OEIS A216446. Palindromic numbers which can be written as the sum of two or more consecutive squares: 5, 55, 77, 181, 313, 434, 505, 545, 595, 636, 818, 1001, …; the linked table of terms contains the 166 terms below 10⁸.
OEIS A267600. Palindromes that are sums of consecutive squares in more than one way; the entry gives the two representations of 554455 and of 9343439, and the next term is 923222222329.
cppreference — std::set. A sorted set of unique objects with logarithmic search, insertion and removal, usually implemented as a red–black tree.
Python — set types: an unordered collection of distinct hashable objects, commonly used to remove duplicates. Java — HashSet: backed by a hash table; adding an element that is already present leaves the set unchanged.
C++ — Euler125.cpp. The linked immutable C++ revision contains parse_arguments, is_palindrome, solve and run_checkpoints. It is the implementation record for this article.

Problem 125 source code

C++

#include <algorithm>
#include <cstdint>
#include <iostream>
#include <set>
#include <string>
#include <vector>

namespace {

using i64 = std::int64_t;

struct Options {
    i64 limit = 100000000LL;
    bool run_checkpoints = true;
};

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

    i64 parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<i64>(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_i64_after_prefix(arg, "--limit=", options.limit)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.limit >= 2;
}

bool is_palindrome(i64 n) {
    std::string s = std::to_string(n);
    std::string t = s;
    std::reverse(t.begin(), t.end());
    return s == t;
}

i64 solve(const i64 limit) {
    int max_base = 1;
    while (static_cast<i64>(max_base) * max_base < limit) {
        ++max_base;
    }

    std::vector<i64> prefix(static_cast<std::size_t>(max_base + 1), 0);
    for (int i = 1; i <= max_base; ++i) {
        prefix[static_cast<std::size_t>(i)] = prefix[static_cast<std::size_t>(i - 1)] +
                                              static_cast<i64>(i) * static_cast<i64>(i);
    }

    std::set<i64> values;

    for (int start = 1; start <= max_base; ++start) {
        for (int end = start + 1; end <= max_base; ++end) {
            const i64 sum = prefix[static_cast<std::size_t>(end)] - prefix[static_cast<std::size_t>(start - 1)];
            if (sum >= limit) {
                break;
            }
            if (is_palindrome(sum)) {
                values.insert(sum);
            }
        }
    }

    i64 total = 0;
    for (i64 v : values) {
        total += v;
    }
    return total;
}

bool run_checkpoints() {
    if (solve(1000) != 4164) {
        std::cerr << "Checkpoint failed for limit=1000" << '\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.limit) << '\n';
    return 0;
}

Python

# Problem 125: Palindromic sums
# Sum of palindromes below 10^8 that are sums of consecutive squares.

def solve():
    limit = 10**8
    max_base = 1
    while max_base * max_base < limit:
        max_base += 1
    
    prefix = [0] * (max_base + 1)
    for i in range(1, max_base + 1):
        prefix[i] = prefix[i-1] + i * i
    
    values = set()
    for start in range(1, max_base + 1):
        for end in range(start + 1, max_base + 1):
            s = prefix[end] - prefix[start - 1]
            if s >= limit: break
            ss = str(s)
            if ss == ss[::-1]:
                values.add(s)
    
    print(sum(values))

solve()

Java

import java.util.*;

public class Euler125 {
    static boolean isPalin(long n) {
        String s = Long.toString(n);
        return s.equals(new StringBuilder(s).reverse().toString());
    }

    public static void main(String[] args) {
        long limit = 100000000L;
        int maxBase = 1;
        while ((long) maxBase * maxBase < limit)
            maxBase++;
        long[] prefix = new long[maxBase + 1];
        for (int i = 1; i <= maxBase; i++)
            prefix[i] = prefix[i - 1] + (long) i * i;
        Set<Long> vals = new HashSet<>();
        for (int s = 1; s <= maxBase; s++)
            for (int e = s + 1; e <= maxBase; e++) {
                long sum = prefix[e] - prefix[s - 1];
                if (sum >= limit)
                    break;
                if (isPalin(sum))
                    vals.add(sum);
            }
        long total = 0;
        for (long v : vals)
            total += v;
        System.out.println(total);
    }
}