Problem 119: Digit Power Sum

View on Project Euler

Project Euler Problem 119 Solution

EulerSolve provides an optimized solution for Project Euler Problem 119, Digit Power Sum, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Let \(\operatorname{digitSum}(N)\) denote the sum of the decimal digits of \(N\). We want the increasing sequence of integers \(N\) with at least two digits such that $$N=\operatorname{digitSum}(N)^k,\qquad k\ge 2,\qquad N\ge 10.$$ The task is to determine the 30th term of that sequence. Typical examples are \(81=9^2\), \(512=8^3\), \(2401=7^4\), and \(4913=17^3\). The important feature is that the base is not chosen freely: it is forced to be the digit sum of the final number. A brute-force scan over all integers is the wrong viewpoint, because almost every number fails immediately. The implementations instead reverse the logic: they iterate over plausible digit sums \(s\), generate the power chain \(s^2,s^3,s^4,\dots\), and keep only the values whose digit sum comes back to \(s\). Mathematical Approach The whole solution is a bounded search over digit sums. Once the correct search space is identified, the rest is just a complete enumeration of short power chains. The Base Is Determined by the Digit Sum If a number \(N\) belongs to the sequence, then its defining base is $$s=\operatorname{digitSum}(N).$$ So every valid term has the form \(N=s^k\) for some \(k\ge 2\). That means the real objects of interest are not arbitrary integers, but the chains of powers attached to possible digit sums....

Detailed mathematical approach

Problem Summary

Write \(S(N)\) for the digit sum of a whole number \(N\), the sum of its decimal digits. The number \(512\) is a power of its own digit sum: \(S(512)=5+1+2=8\) and \(8^3=512\). So is \(614656=28^4\), whose digits add up to \(28\). Let \(a_1\lt a_2\lt a_3\lt\dots\) be all numbers with this property, in increasing order, where, as the statement insists, a number needs at least two digits to have a sum: the terms are the numbers \(N\ge10\) with \(N=S(N)^k\) for some whole number \(k\ge1\). We are given \(a_2=512\) and \(a_{10}=614656\) and must find \(a_{30}\). [1]

The whole solution in one sentence: the base of the power has to be the digit sum itself, and a number with \(d\) digits has a digit sum of at most \(9d\), so instead of examining every number we go through every possible digit sum \(s\), test its powers \(s^2,s^3,\dots\) up to a ceiling such as \(10^{18}\), keep the powers whose digits add up to \(s\) again, and sort what is left.

Mathematical Approach

We start with small numbers, turn the question around, prove that the resulting short search misses no term, and then see why the terms are so rare.

1. Small numbers by hand

Among the two-digit numbers only \(81\) qualifies: \(S(81)=8+1=9\) and \(9^2=81\). A power alone is not enough: \(64=8^2\) has \(S(64)=10\), and \(64\) is no power of \(10\). The digits are added only once: \(614656\) has the digit sum \(28\), and adding again, \(2+8=10\) and \(1+0=1\), would give the digital root, a different quantity. [2]

The most direct method tests one number after another: compute its digit sum \(s\). A digit sum of \(1\), as for \(10\) or \(100\), never qualifies, since every power of \(1\) is \(1\), and must be rejected at once, because dividing by \(1\) would never stop. Otherwise divide the number by \(s\) as long as no remainder is left; the number is a term exactly when this ends at \(1\). For \(2401\) the digit sum is \(7\), and \(2401\to343\to49\to7\to1\), so \(2401=7^4\). This test gives the first ten terms, each written as a power of its digit sum:

\(81=9^2\), \(512=8^3\), \(2401=7^4\), \(4913=17^3\), \(5832=18^3\),
\(17576=26^3\), \(19683=27^3\), \(234256=22^4\),
\(390625=25^4\), \(614656=28^4\).

The given \(a_2\) and \(a_{10}\) are in their places. The terms thin out quickly: five lie below \(10^4\) and only twelve more between \(10^4\) and \(10^9\), so testing number after number is correct but slow.

2. The base and the exponent are at least \(2\)

Let \(d_i\) be the digit of \(N\) in the place worth \(10^i\), so that \(N=\sum_i d_i10^i\) and \(S(N)=\sum_i d_i\). Subtracting gives

$$N-S(N)=\sum_{i\ge1}d_i\,(10^i-1),$$

because the units place, \(i=0\), contributes \(d_0-d_0=0\). No term of this sum is negative, and a number with at least two digits has a nonzero digit in some place \(i\ge1\), which contributes at least \(9\). So such a number is larger than its digit sum, and \(N=S(N)^1\) is impossible. The exponent \(0\) gives only \(1\), and so does the base \(1\). Every term is therefore \(s^k\) with \(s\ge2\) and \(k\ge2\); the rule about two digits only removes the one-digit numbers, each of which is trivially its own digit sum to the power \(1\).

Moreover, each \(10^i-1\) is one of \(9,99,999,\dots\), a multiple of \(9\), so \(N-S(N)\) is a multiple of \(9\): a number and its digit sum leave the same remainder after division by \(9\). This is the old rule of casting out nines. [3] For example, \(512-8=504=9\cdot56\). Section 7 turns this into a test.

3. Turning the question around

If \(N\) is a term, the base of its power is not free: it must be \(s=S(N)\). So instead of asking for every number whether it is a power of its digit sum, we ask for every base \(s\) which of its powers \(s^2,s^3,\dots\) have the digit sum \(s\); each power is the previous one times \(s\). For \(s=8\) this gives

\(S(64)=10\), \(S(512)=8\), \(S(4096)=19\), \(S(32768)=26\),
\(S(262144)=19\), \(S(2097152)=26\), \(S(16777216)=37\),

and only \(512\) is kept. A base can also give several terms. For \(s=28\) the powers \(784\), \(21952\), \(614656\), \(17210368\) and \(481890304\) have the digit sums \(19\), \(19\), \(28\), \(28\) and \(37\), so \(28^4\) and \(28^5\) are both terms, and below \(10^{18}\) the base \(18\) even gives three: \(18^3=5832\), \(18^6=34012224\) and \(18^7=612220032\). So every chain must be followed to the end.

Two facts make this search exact. Every accepted power is a term: if \(S(s^k)=s\) with \(s,k\ge2\), then \(s^k\) is a power of its digit sum, and it has at least two digits, because a one-digit number is its own digit sum while \(s^k\gt s\). And no term is found twice. A candidate can lie on several chains, as \(64=2^6=4^3=8^2\) does, but a power accepted for the base \(s\) has the digit sum \(s\), so no other base accepts it, and the powers of one base are all different.

4. How large the digit sum can be

A number with \(d\) digits has a digit sum of at most \(9d\), and only \(99\dots9\) reaches it. [2] Take a ceiling \(10^D\). Every number below it has at most \(D\) digits, and \(10^D\) itself has the digit sum \(1\) and is no term. Hence every term \(N\le10^D\) has a base \(s=S(N)\) with \(2\le s\le9D\), and its exponent \(k\) satisfies \(s^k\le10^D\). Together with Section 3 this proves the claim on which the search rests:

Claim. The terms up to \(10^D\) are exactly the powers \(s^k\le10^D\) with \(2\le s\le9D\), \(k\ge2\) and \(S(s^k)=s\).

Indeed, every such power is a term by Section 3, and every term up to \(10^D\) is such a power by the bound above. For a fixed base \(s\) the exponents run from \(2\) to \(\lfloor\log_s10^D\rfloor\), the largest \(k\) with \(s^k\le10^D\), where \(\lfloor x\rfloor\) is the largest whole number not above \(x\). So the number of powers to test is

$$\sum_{s=2}^{9D}\bigl(\lfloor\log_s10^D\rfloor-1\bigr).$$

For \(D=18\), the ceiling of the Python and Java programs, the base runs from \(2\) to \(162\); the base \(2\) has the exponents \(2\) to \(59\), since \(2^{59}\le10^{18}\lt2^{60}\), the base \(162\) only \(2\) to \(8\), and the sum is \(1533\).

5. Worked example: all terms below \(10^4\)

For the ceiling \(10^4\) the bound gives \(s\le36\). The base \(2\) has the twelve candidates \(2^2,\dots,2^{13}\), because \(2^{13}=8192\le10^4\lt2^{14}\); the base \(3\) has seven, up to \(3^8=6561\); the bases \(4\) to \(10\) have \(5,4,4,3,3,3,3\); the bases \(11\) to \(21\) have their squares and cubes, up to \(21^3=9261\); and the bases \(22\) to \(36\) only their squares, since \(22^3=10648\gt10^4\). That makes \(81\) candidates. The cubes of \(11\) to \(21\), for example, have the digit sums \(8, 18, 19, 17, 18, 19,\) \(17, 18, 28, 8, 18\), so only \(17^3=4913\) and \(18^3=5832\) hit their base. In all, five candidates survive: \(9^2=81\), \(8^3=512\), \(7^4=2401\), \(17^3\) and \(18^3\). Sorted, they are \(a_1,\dots,a_5\), in agreement with Section 1.

6. Sorting the hits

The search visits the bases in increasing order, but the sequence is ordered by size: the bases \(7\), \(8\) and \(9\) give \(7^4=2401\), \(8^3=512\) and \(9^2=81\), exactly in reverse order, and the base \(18\) gives \(a_5\), \(a_{13}\) and \(a_{17}\). So all accepted powers are collected first and sorted afterwards. If the collection has at least \(n\) elements, its \(n\)-th smallest is \(a_n\), because by the claim of Section 4 it holds exactly the terms up to the ceiling, the first terms of the sequence. The ceiling \(10^{12}\) holds \(24\) terms, too few for \(a_{30}\); the ceiling \(10^{18}\) holds \(34\), enough. The ceiling, not a cap on the base, makes the list complete: the bases up to \(30\) alone, however far their powers are followed, miss \(36^4=1679616=a_{11}\).

7. A test modulo \(9\)

Section 2 showed that \(N\) and \(S(N)\) leave the same remainder after division by \(9\). In the notation of modular arithmetic, \(a\equiv b\pmod m\) means that \(a-b\) is a multiple of \(m\), and such congruences may be multiplied, so \(a\equiv b\) gives \(a^k\equiv b^k\). [4] For a term \(N=s^k\) with \(S(N)=s\), the congruence \(N\equiv S(N)\) becomes the necessary condition

$$s^k\equiv s\pmod 9.$$

It depends only on the remainder \(r\) of \(s\) after division by \(9\):

  • \(r=0\) or \(r=1\): every exponent passes;
  • \(r=3\) or \(r=6\): no exponent passes, because \(9\) divides \(s^2\) but not \(s\);
  • \(r=8\): the powers leave \(8,1,8,1,\dots\), so the odd exponents pass;
  • \(r=4\) or \(r=7\): the powers leave \(4,7,1,\dots\) or \(7,4,1,\dots\), so \(k\equiv1\pmod 3\) is needed;
  • \(r=2\) or \(r=5\): the powers leave \(2,4,8,7,5,1,\dots\) or \(5,7,8,4,2,1,\dots\), so \(k\equiv1\pmod 6\) is needed.

For example, the digit sums of \(7^2,\dots,7^{10}\) are \(13, 10, 7, 22, 28,\) \(25, 31, 28, 43\), with the remainders \(4,1,7\) repeating, so only \(k=4,7,10\) can give the digit sum \(7\). The exponent \(4\) does; \(7\) and \(10\) pass the test but fail, so the test is necessary, not sufficient. Bases such as \(12\) or \(15\) never give a term. Of the \(1533\) candidates up to \(10^{18}\) only \(534\) pass, and of the \(81\) up to \(10^4\) only \(19\). The published programs do not use the test, since a digit sum is cheap anyway; it explains the gaps among the bases, and all \(58\) terms below \(10^{30}\) pass it.

8. Each exponent gives only finitely many terms

For a fixed exponent \(k\), apply the digit bound of Section 4 to the power itself. The number \(s^k\) has at most \(k\log_{10}s+1\) digits, where \(\log_{10}\) is the logarithm to base \(10\), so

$$s=S(s^k)\le9\,(k\log_{10}s+1).$$

The left side grows faster. To see this, let \(s\) have \(m\) digits. Then \(s\ge10^{m-1}\), while \(s^k\) has at most \(km\) digits, so a base that works satisfies \(10^{m-1}\le s\le9km\). When \(m\) grows by one, the left end grows tenfold and \(9km\) at most doubles, so once this inequality fails, it fails for all larger \(m\): for each \(k\) only finitely many bases qualify. For \(k\le6\) it fails from \(m=4\) on, since \(10^3\gt9\cdot6\cdot4\), so these exponents leave only bases with at most three digits. For squares the bound is easy to make exact. Here \(m=3\) already fails, since \(10^2\gt9\cdot2\cdot3\), and for \(37\le s\le99\) the square has at most four digits, since \(99^2=9801\), so its digit sum is at most \(36\lt s\). So only the bases \(2\) to \(36\) remain; their squares have digit sums of at most \(19\), and only \(9^2=81\) hits its base. Thus \(81\) is the only term with \(N=S(N)^2\). Other terms can still be perfect squares, such as \(2401=49^2\); its base must be \(S(2401)=7\), so its exponent is \(4\).

For cubes the same reasoning leaves the bases up to \(54\), and exactly five of them work: \(8\), \(17\), \(18\), \(26\) and \(27\), with the cubes \(512\), \(4913\), \(5832\), \(17576\) and \(19683\). Together with \(1=1^3\) these are the Dudeney numbers, named after Henry Dudeney, who noted them in his puzzle Root Extraction. [5] For the exponents \(4\), \(5\) and \(6\), the condition that \(s\) is at most nine times the number of digits of \(s^k\) leaves only the bases up to \(72\), \(90\) and \(117\), and the terms come from the bases \(7,22,25,28,36\), then \(28,35,36,46\) and \(18,45,54,64\). So exactly \(19\) terms have an exponent of at most \(6\), the largest of them \(64^6=68719476736=a_{23}\), and every other term has an exponent of at least \(7\).

9. Rarity, and other number bases

Up to \(10^6\) there are \(10\) terms, up to \(10^{12}\) there are \(24\), up to \(10^{18}\) \(34\), up to \(10^{24}\) \(44\) and up to \(10^{30}\) \(58\). A heuristic reason, not a proof: if the digits of a power behaved like random digits, they would average \(4.5\), so a power with \(D\) digits would usually have a digit sum near \(4.5D\), and only rarely would it land exactly on its base. The terms fit this picture: below \(10^{18}\) the base of a term is on average \(3.9\) and at most \(5.82\) times its number of digits, far below the factor \(9\) of the bound. Digit sums also drop at every carry, as in \(S(99)=18\) but \(S(100)=1\), so along a chain of powers they jump up and down, as for the base \(8\) in Section 3, and the search has to test every candidate. The method works in every base \(b\ge2\): digits are found by dividing by \(b\), a number with \(D\) digits has a digit sum of at most \((b-1)D\), and the argument of Section 2 gives the remainder test modulo \(b-1\).

How the Code Works

The C++ program [12] starts in main with the default index \(30\) and the checkpoints switched on. parse_arguments accepts --index=N, read by parse_int_after_prefix, and --skip-checkpoints; any other argument, or an index below \(1\), gives exit code \(1\). Unless the checkpoints are skipped, run_checkpoints checks that solve(2) returns \(a_2=512\) and exits with code \(2\) otherwise; \(a_{10}\) is not checked. solve(index) starts with the ceiling limit \(=10^6\) and calls generate_terms; while fewer than index terms come back, it multiplies the ceiling by \(10\) and searches again, so it never relies on a guessed ceiling, and it returns terms[index - 1].

generate_terms(limit) implements Sections 3 and 4. Nine times digits_count(limit) bounds the base; the ceiling \(10^D\) has \(D+1\) digits, so this allows one digit too many, \(63\) instead of \(54\) for \(10^6\). For every base \(s\) the running power is an unsigned __int128, a 128-bit type that GCC and Clang provide as an extension; [6] it starts at \(s\cdot s\) and is multiplied by \(s\) while it does not exceed the ceiling, so at every test it equals \(s^k\) for the current exponent \(k\), as the claim of Section 4 requires. The ceiling is a 64-bit number with at most \(20\) digits, so the base is at most \(180\), every product stays below \(180\cdot2^{64}\lt2^{128}\), and the comparison with the ceiling is exact, unlike an exponent from floating-point logarithms, which could round the wrong way at an exact power; only a power that passed it is converted back to 64 bits. A power is stored if it has at least two digits and digit_sum returns \(s\); this function adds the remainders of repeated division by \(10\), for \(512\) the remainders \(2\), \(1\) and \(5\). By Section 3 the two-digit test is redundant. The powers go into a std::set, which keeps its elements sorted and without repeats, [7] and copying it into a std::vector gives the terms in increasing order.

Two limits concern only other indices. The ceiling is an unsigned 64-bit integer, which wraps around modulo \(2^{64}\approx1.8\cdot10^{19}\): after \(10^{19}\) comes \(10^{20}-5\cdot2^{64}\)\({}\approx7.8\cdot10^{18}\). So only the \(35\) terms below \(2^{64}\) can ever be returned; \(a_{35}\), which lies above \(10^{19}\), is found only because a later wrapped ceiling, \(10^{26}\bmod2^{64}\), happens to exceed it, and for --index=36 or more the loop never ends. And parse_int_after_prefix collects the digits in a signed int without a range check, so an index above \(2^{31}-1\) overflows, which C++ leaves undefined. Neither matters for the default index.

The Python program uses the fixed ceiling \(10^{18}\) and the base bound max_ds \(=9\cdot18=162\), exactly that of Section 4. The loop is as in C++; Python integers have unlimited precision, [8] so val *= s cannot overflow, and the last power of each base, up to \(162\) times the ceiling, just ends the loop. The terms go into a set, which has no order, sorted puts them in order, and result[29] is \(a_{30}\) because the ceiling holds \(34\) terms.

The Java program also uses the ceiling \(10^{18}\), but lets the base run up to \(200\), more than the \(162\) needed, so it tests \(1776\) genuine powers instead of \(1533\). A TreeSet<Long> keeps the terms in increasing order, [9] toArray copies them in the order of the set's iterator, [10] and the element with index \(29\) is printed. The product val*=s can exceed the largest long, \(2^{63}-1\approx9.2\cdot10^{18}\), when val is close to the ceiling. Java then keeps only the low 64 bits of the product, so the sign of the result may differ from that of the true product, but need not. [11] The program ends a base when the product becomes negative, which is not a complete overflow test. Of the \(95\) bases that overflow, \(58\) give a negative value and \(36\) a positive one above the ceiling, which also ends the loop; but for \(s=139\) the true value \(139^9\approx1.94\cdot10^{19}\) wraps around to \(923415668714480043\), below the ceiling, and is tested as if it were the next power of \(139\), which makes \(1777\) tests in all. Its digit sum is \(75\), so it is rejected, and the next product is negative. The stored terms and the output are therefore not affected, but the guard is incomplete: a sound loop would stop as soon as val > limit / s instead of multiplying.

Complexity and Verification

For the ceiling \(10^D\) there are at most \(9D\) bases, each with at most \(\log_2 10^D\approx3.3D\) exponents, and a digit sum costs at most \(D+1\) loop steps, so the work is at most of the order of \(D^3\), counting each loop step as one operation, fair while the numbers fit into machine words; the direct test of Section 1 needs at least \(10^D\) steps. The Python program tests \(1533\) powers, none with more than \(19\) digits; their digit sums take \(16471\) loop steps in all, each a remainder, a division and an addition, and the program stores at most \(34\) numbers. The C++ program repeats its search for the ceilings \(10^6,10^7,\dots\) until the list is long enough, at most thirteen times for the default index, since \(10^{18}\) is enough; even thirteen rounds together would test only \(10536\) powers. On the machine used for this article, each program finishes within a few hundredths of a second, including the start of the Python interpreter or the Java virtual machine.

The C++ program checks \(a_2=512\) before it searches; the Python and Java programs have no checks. For this article, all three programs were run and print the same number, and the C++ program prints \(a_{10}=614656\) for --index=10. A separate search with Python's unbounded integers agrees with the C++ program for every index from \(1\) to \(35\), and a small harness around generate_terms confirms that the ceiling itself is included: it finds no term up to \(80\), one up to \(81\) and two up to \(512\). Agreement alone cannot show that no term is missing; the claim of Section 4 does. The direct test of Section 1, written in C for every number from \(10\) to \(10^9\), with the digit sum updated from one number to the next, finds exactly the \(17\) terms that the search finds below \(10^9\); it takes about \(1.4\) seconds, so at the same rate a scan up to \(10^{18}\) would take more than \(40\) years. A simulation of Java's 64-bit arithmetic and an instrumented copy of the Java loop, compared with exact powers, reproduce the overflow cases above, and a C++ build with a detector for undefined behaviour runs the default calculation cleanly and confirms the overflow of the index parser. All other counts and lists in this article were computed in the same way, those of Section 8 for every base up to \(2000\).

Footnotes and References

The references supply the background facts; the claim of Section 4, the test of Section 7, the bounds of Section 8 and the Java overflow analysis are derived here.

  1. Project Euler 119 — Digit Power Sum. The official statement: 512 = 8³ and 614656 = 28⁴ equal a power of their digit sums 5 + 1 + 2 = 8 and 28; aₙ is the n-th such number, a number needs at least two digits to have a sum, a₂ = 512 and a₁₀ = 614656 are given, and a₃₀ is asked for.
  2. Wikipedia — Digit sum. The digit sum is the sum of all digits of a number; applied repeatedly until one digit remains, it gives the digital root. For a given number of digits it is largest when every digit equals b − 1, which in base 10 gives the bound 9d.
  3. Wikipedia — Casting out nines. Replacing a number by its digit sum leaves the same remainder after division by nine, because the difference is a sum of digit multiples of the numbers 10ⁱ − 1, which are all divisible by 9.
  4. Wikipedia — Modular arithmetic. Two integers are congruent modulo m if their difference is a multiple of m; congruence is compatible with addition, multiplication and exponentiation, so a ≡ b (mod m) gives aᵏ ≡ bᵏ (mod m).
  5. Wikipedia — Dudeney number. In base 10 exactly six positive integers equal the cube of their digit sum: 1, 512, 4913, 5832, 17576 and 19683. The name comes from Henry Dudeney's puzzle Root Extraction, and a bound on the roots leaves only finitely many for each power and base.
  6. GCC — 128-bit Integers. As an extension, __int128 is supported on targets with an integer mode wide enough to hold 128 bits; unsigned __int128 is the unsigned 128-bit type.
  7. cppreference — std::set. std::set is an associative container that contains a sorted set of unique objects; search, removal and insertion take logarithmic time.
  8. Python — Numeric types. Integers have unlimited precision.
  9. Java — TreeSet. The elements are ordered by their natural ordering, the iterator returns them in ascending order, and add, remove and contains take guaranteed logarithmic time.
  10. Java — AbstractCollection. Its toArray methods return the elements in the order of the collection's iterator whenever the collection guarantees such an order.
  11. The Java Language Specification, §15.17.1 — Multiplication Operator *. An overflowing integer multiplication yields the low-order bits of the mathematical product, so the sign may differ from that of the true product; no exception is thrown.
  12. C++ — Euler119.cpp. The linked immutable revision contains generate_terms, solve and run_checkpoints; it is the implementation record for this article.

Problem 119 source code

C++

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

namespace {

using u64 = std::uint64_t;
using u128 = unsigned __int128;

struct Options {
    int index = 30;
    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, "--index=", options.index)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }

    return options.index >= 1;
}

int digit_sum(u64 n) {
    int sum = 0;
    while (n > 0) {
        sum += static_cast<int>(n % 10ULL);
        n /= 10ULL;
    }
    return sum;
}

int digits_count(u64 n) {
    int count = 0;
    while (n > 0) {
        ++count;
        n /= 10ULL;
    }
    return std::max(1, count);
}

std::vector<u64> generate_terms(const u64 limit) {
    const int max_digit_sum = 9 * digits_count(limit);
    std::set<u64> values;

    for (int s = 2; s <= max_digit_sum; ++s) {
        u128 value = static_cast<u128>(s) * static_cast<u128>(s);
        while (value <= limit) {
            const u64 v = static_cast<u64>(value);
            if (v >= 10ULL && digit_sum(v) == s) {
                values.insert(v);
            }
            value *= static_cast<u128>(s);
        }
    }

    return std::vector<u64>(values.begin(), values.end());
}

u64 solve(const int index) {
    u64 limit = 1000000ULL;

    for (;;) {
        std::vector<u64> terms = generate_terms(limit);
        if (static_cast<int>(terms.size()) >= index) {
            return terms[static_cast<std::size_t>(index - 1)];
        }
        limit *= 10ULL;
    }
}

bool run_checkpoints() {
    if (solve(2) != 512ULL) {
        std::cerr << "Checkpoint failed for second term" << '\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.index) << '\n';
    return 0;
}

Python

# Problem 119: Digit power sum
# Find the 30th term where n = (digit_sum(n))^k for some k.

def digit_sum(n):
    s = 0
    while n > 0:
        s += n % 10
        n //= 10
    return s

def solve():
    limit = 10**18
    values = set()
    max_ds = 9 * 18  # max digit sum for 18-digit numbers
    for s in range(2, max_ds + 1):
        val = s * s
        while val <= limit:
            if val >= 10 and digit_sum(val) == s:
                values.add(val)
            val *= s
    
    result = sorted(values)
    print(result[29])  # 30th term (0-indexed)

solve()

Java

import java.util.*;

public class Euler119 {
    static int digitSum(long n){int s=0;while(n>0){s+=(int)(n%10);n/=10;}return s;}
    public static void main(String[] args) {
        long limit=1000000000000000000L; TreeSet<Long> vals=new TreeSet<>();
        for(int s=2;s<=200;s++){long val=(long)s*s;while(val<=limit&&val>0){if(val>=10&&digitSum(val)==s)vals.add(val);val*=s;if(val<0)break;}}
        Long[] arr=vals.toArray(new Long[0]);
        System.out.println(arr[29]);
    }
}