Problem 529: $10$-substrings

View on Project Euler

Project Euler Problem 529 Solution

EulerSolve provides an optimized solution for Project Euler Problem 529, $10$-substrings, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary A decimal string \(d_1d_2\dots d_\ell\) with \(d_1\ne 0\) is valid if every position belongs to at least one contiguous substring whose digit sum is \(10\). Let \(a_\ell\) be the number of valid strings of length \(\ell\), and define $$T(n)=\sum_{\ell=1}^{n} a_\ell \pmod{10^9+7}.$$ The target input is enormous, so the implementation cannot iterate lengths up to \(n\) directly. The solution turns the coverage condition into a finite automaton, extracts a linear recurrence for \(T(n)\), and evaluates that recurrence at very large indices. Mathematical Approach Write the prefix sums of the digits as \(P_0=0\) and \(P_t=d_1+\cdots+d_t\). A substring \(d_i\dots d_j\) has digit sum \(10\) exactly when $$P_j-P_{i-1}=10.$$ The whole method is built around scanning the string from left to right while tracking only the information that can still matter for future sum-\(10\) substrings. Step 1: Track the Leftmost Unresolved Position Suppose we have processed digits up to position \(j\), and let \(c\) be the first position that is not yet guaranteed to be covered by a sum-\(10\) substring. Then the unresolved window is \(d_c\dots d_j\), whose total digit sum is $$\Delta=P_j-P_{c-1}.$$ All digits are nonnegative, so \(\Delta\) can only increase as we append more digits....

Detailed mathematical approach

Problem Summary

A decimal string \(d_1d_2\dots d_\ell\) with \(d_1\ne 0\) is valid if every position belongs to at least one contiguous substring whose digit sum is \(10\). Let \(a_\ell\) be the number of valid strings of length \(\ell\), and define

$$T(n)=\sum_{\ell=1}^{n} a_\ell \pmod{10^9+7}.$$

The target input is enormous, so the implementation cannot iterate lengths up to \(n\) directly. The solution turns the coverage condition into a finite automaton, extracts a linear recurrence for \(T(n)\), and evaluates that recurrence at very large indices.

Mathematical Approach

Write the prefix sums of the digits as \(P_0=0\) and \(P_t=d_1+\cdots+d_t\). A substring \(d_i\dots d_j\) has digit sum \(10\) exactly when

$$P_j-P_{i-1}=10.$$

The whole method is built around scanning the string from left to right while tracking only the information that can still matter for future sum-\(10\) substrings.

Step 1: Track the Leftmost Unresolved Position

Suppose we have processed digits up to position \(j\), and let \(c\) be the first position that is not yet guaranteed to be covered by a sum-\(10\) substring. Then the unresolved window is \(d_c\dots d_j\), whose total digit sum is

$$\Delta=P_j-P_{c-1}.$$

All digits are nonnegative, so \(\Delta\) can only increase as we append more digits. Therefore, if \(\Delta>10\), no future extension can create a sum-\(10\) substring containing position \(c\), and that branch is impossible. This gives the crucial bound

$$0\le \Delta \le 10.$$

Step 2: Record Only Finite Data

Inside the current unresolved window we keep the set of partial sums

$$S=\{P_t-P_{c-1}: c-1\le t\le j\}\subseteq\{0,1,\dots,\Delta\}.$$

These values describe every internal boundary of the unresolved window. We also keep a second set

$$H\subseteq\{0,1,\dots,10\},$$

whose elements are offsets from earlier prefix positions that can still serve as left boundaries of a future sum-\(10\) substring. When a new digit \(x\) is appended, the unresolved total becomes

$$\Delta'=\Delta+x.$$

If \(\Delta'>10\), the transition is rejected. Otherwise we enlarge the partial-sum set to

$$S'=S\cup\{\Delta'\}.$$

Step 3: Detect When the Unresolved Window Closes

The new right endpoint closes the whole unresolved window exactly when there is a complementary historical offset:

$$10-\Delta' \in H.$$

In that case some earlier prefix position \(u\) satisfies \(P_{c-1}-P_u=10-\Delta'\), so

$$P_{j+1}-P_u=10.$$

Hence the substring from \(u+1\) to \(j+1\) has digit sum \(10\) and covers every position from \(c\) through \(j+1\). After this closure, the unresolved window becomes empty again, so the window data resets to

$$S_{\text{new}}=\{0\}.$$

The historical offsets relative to the new boundary are updated by translation and by every internal cut of the closed window:

$$H_{\text{new}}=\{h+\Delta' : h\in H,\ h+\Delta'\le 10\}\cup\{\Delta'-s : s\in S'\}.$$

Step 4: Finite Automaton and Acceptance

Because both \(H\) and \(S\) live inside \(\{0,1,\dots,10\}\), only finitely many states can ever occur. Each appended digit gives a deterministic state transition. A state is accepting exactly when no unresolved window remains, which in this encoding means

$$S=\{0\}.$$

So counting valid strings of length \(\ell\) becomes counting length-\(\ell\) paths in a deterministic finite automaton. The first digit uses the alphabet \(\{1,\dots,9\}\), while all later digits use \(\{0,\dots,9\}\).

Step 5: From the Automaton to a Linear Recurrence

Let \(a_\ell\) be the number of accepting paths of length \(\ell\). Dynamic programming on the minimized automaton gives the initial values of \(a_\ell\) and therefore of

$$T(n)=\sum_{\ell=1}^{n} a_\ell.$$

Any fixed finite automaton induces repeated matrix multiplication modulo \(10^9+7\), so the cumulative sequence \(T(0),T(1),T(2),\dots\) satisfies a linear recurrence. The implementation samples enough initial terms, recovers the shortest recurrence with Berlekamp-Massey, and then evaluates the \(n\)-th term with polynomial binary exponentiation.

Worked Example

Consider the string \(352\). Its prefix sums are \(0,3,8,10\). While reading it, the unresolved total moves from \(0\) to \(3\), then to \(8\), then to \(10\), never exceeding the allowed bound. At the last digit we have \(10-\Delta'=0\), and the initial historical set already contains \(0\), so the whole string closes as one sum-\(10\) block. Therefore \(352\) is valid.

As a small counting check, no length-\(1\) string can be valid, so \(a_1=0\). For length \(2\), the only valid strings are \(19,28,37,46,55,64,73,82,91\), giving

$$a_2=9,\qquad T(2)=9.$$

How the Code Works

The C++, Python, and Java implementations use the same mathematical pipeline. The compiled solver first enumerates every reachable automaton state from the start configuration by repeatedly applying the digit-transition rule described above. After that, it minimizes the raw automaton by partition refinement so that equivalent states are merged without changing the accepted language.

Next, the implementation runs dynamic programming on the minimized automaton. The first step allows digits \(1\) through \(9\); each later step allows \(0\) through \(9\). Summing the counts in accepting states yields \(a_\ell\), and cumulative sums yield \(T(\ell)\). The code generates \(2s+5\) initial values when the minimized automaton has \(s\) states, which is enough to recover the recurrence used in practice.

Finally, Berlekamp-Massey extracts the shortest linear recurrence for the cumulative sequence modulo \(10^9+7\). The solver then evaluates that recurrence at index \(n\) using Kitamasa-style polynomial reduction with binary exponentiation. The Python entry point delegates to the compiled solver, so its numerical behavior matches the compiled versions.

Complexity Analysis

Let \(s_{\text{raw}}\) be the number of reachable raw states, \(s\) the number of minimized states, and \(r\) the recovered recurrence order. Building the reachable raw automaton requires \(O(10s_{\text{raw}})\) transitions. Minimization is near-linear in the transition graph, typically written as \(O(10s_{\text{raw}}\log s_{\text{raw}})\).

The initial dynamic programming phase computes \(2s+5\) terms, so it costs \(O(10s(2s+5))=O(s^2)\) time. Once the recurrence is known, the large-index evaluation costs

$$O(r^2\log n)$$

time and \(O(r)\) extra memory, on top of the automaton tables.

Footnotes and References

  1. Project Euler Problem 529
  2. Wikipedia - Prefix sum
  3. Wikipedia - Finite-state machine
  4. Wikipedia - Berlekamp-Massey algorithm
  5. Wikipedia - Linear recurrence with constant coefficients

Problem 529 source code

C++

#include <algorithm>
#include <array>
#include <cstdint>
#include <cstdlib>
#include <deque>
#include <iostream>
#include <thread>
#include <unordered_map>
#include <vector>

using int64 = long long;

namespace {

constexpr int kTargetSum = 10;
constexpr int kDigits = 10;
constexpr int64 kMod = 1000000007LL;

struct State {
    uint16_t history = 0;  // bits for reachable prefix sums in [S_c-10, S_c]
    uint16_t segment = 0;  // bits for segment sums from c to current position
};

struct Automaton {
    int start = 0;
    std::vector<std::array<int, kDigits>> trans;
    std::vector<char> accept;
};

constexpr uint32_t kHistoryMask = (1u << (kTargetSum + 1)) - 1u;

uint32_t encode_state(uint16_t history, uint16_t segment) {
    return static_cast<uint32_t>(history) | (static_cast<uint32_t>(segment) << 11);
}

int highest_bit(uint16_t mask) {
    return 31 - __builtin_clz(static_cast<unsigned>(mask));
}

bool next_state(uint16_t history, uint16_t segment, int digit, State& out) {
    int delta = highest_bit(segment);
    int new_delta = delta + digit;
    if (new_delta > kTargetSum) return false;
    uint16_t new_segment = static_cast<uint16_t>(segment | (1u << new_delta));
    if ((history >> (kTargetSum - new_delta)) & 1u) {
        uint32_t shifted = (static_cast<uint32_t>(history) << new_delta) & kHistoryMask;
        uint16_t new_history = static_cast<uint16_t>(shifted);
        for (int k = 0; k <= new_delta; ++k) {
            if ((new_segment >> k) & 1u) {
                new_history |= static_cast<uint16_t>(1u << (new_delta - k));
            }
        }
        out.history = new_history;
        out.segment = 1u;
    } else {
        out.history = history;
        out.segment = new_segment;
    }
    return true;
}

struct RawAutomaton {
    std::vector<State> states;
    std::vector<std::array<int, kDigits>> trans;
    std::vector<char> accept;
    int start = 0;
};

RawAutomaton build_raw() {
    RawAutomaton raw;
    raw.states.reserve(20000);
    std::unordered_map<uint32_t, int> index;
    index.reserve(20000);

    State start{1u, 1u};
    raw.states.push_back(start);
    index[encode_state(start.history, start.segment)] = 0;

    for (size_t i = 0; i < raw.states.size(); ++i) {
        const State& cur = raw.states[i];
        for (int d = 0; d < kDigits; ++d) {
            State nxt;
            if (!next_state(cur.history, cur.segment, d, nxt)) continue;
            uint32_t key = encode_state(nxt.history, nxt.segment);
            auto it = index.find(key);
            if (it == index.end()) {
                int id = static_cast<int>(raw.states.size());
                index[key] = id;
                raw.states.push_back(nxt);
            }
        }
    }

    int n = static_cast<int>(raw.states.size());
    raw.trans.assign(n, {});
    raw.accept.assign(n, 0);
    for (int i = 0; i < n; ++i) {
        const State& cur = raw.states[i];
        raw.accept[i] = (cur.segment == 1u) ? 1 : 0;
        for (int d = 0; d < kDigits; ++d) {
            State nxt;
            if (!next_state(cur.history, cur.segment, d, nxt)) {
                raw.trans[i][d] = -1;
            } else {
                uint32_t key = encode_state(nxt.history, nxt.segment);
                raw.trans[i][d] = index[key];
            }
        }
    }
    raw.start = 0;
    return raw;
}

Automaton minimize_automaton(const RawAutomaton& raw) {
    int n = static_cast<int>(raw.states.size());
    int dead = n;
    int total = n + 1;

    std::vector<std::array<int, kDigits>> trans(total);
    std::vector<char> accept(total, 0);
    for (int i = 0; i < n; ++i) {
        accept[i] = raw.accept[i];
        for (int d = 0; d < kDigits; ++d) {
            int t = raw.trans[i][d];
            trans[i][d] = (t == -1) ? dead : t;
        }
    }
    accept[dead] = 0;
    trans[dead].fill(dead);

    std::vector<std::vector<int>> inv[kDigits];
    for (int d = 0; d < kDigits; ++d) {
        inv[d].assign(total, {});
    }
    for (int s = 0; s < total; ++s) {
        for (int d = 0; d < kDigits; ++d) {
            inv[d][trans[s][d]].push_back(s);
        }
    }

    std::vector<std::vector<int>> classes;
    classes.reserve(total);
    classes.push_back({});
    classes.push_back({});
    for (int s = 0; s < total; ++s) {
        classes[accept[s] ? 0 : 1].push_back(s);
    }
    if (classes[0].empty()) classes.erase(classes.begin());
    if (classes.size() == 2 && classes[1].empty()) classes.pop_back();

    std::vector<int> class_of(total, -1);
    for (int i = 0; i < static_cast<int>(classes.size()); ++i) {
        for (int s : classes[i]) class_of[s] = i;
    }

    std::deque<int> work;
    std::vector<char> in_work(classes.size(), 0);
    for (int i = 0; i < static_cast<int>(classes.size()); ++i) {
        work.push_back(i);
        in_work[i] = 1;
    }

    std::vector<char> marked(total, 0);
    std::vector<std::vector<int>> bucket(classes.size());

    while (!work.empty()) {
        int A = work.front();
        work.pop_front();
        in_work[A] = 0;
        for (int d = 0; d < kDigits; ++d) {
            std::vector<int> X;
            X.reserve(total / 4);
            for (int t : classes[A]) {
                for (int s : inv[d][t]) {
                    if (!marked[s]) {
                        marked[s] = 1;
                        X.push_back(s);
                    }
                }
            }
            if (X.empty()) continue;
            std::vector<int> touched;
            for (int s : X) {
                int c = class_of[s];
                if (bucket[c].empty()) touched.push_back(c);
                bucket[c].push_back(s);
            }
            for (int c : touched) {
                if (bucket[c].size() < classes[c].size()) {
                    std::vector<int> inX;
                    inX.swap(bucket[c]);
                    std::vector<int> outX;
                    outX.reserve(classes[c].size() - inX.size());
                    for (int s : classes[c]) {
                        if (!marked[s]) outX.push_back(s);
                    }
                    classes[c].swap(inX);
                    int new_id = static_cast<int>(classes.size());
                    classes.push_back(std::move(outX));
                    if (bucket.size() < classes.size()) bucket.resize(classes.size());
                    if (in_work.size() < classes.size()) in_work.resize(classes.size(), 0);
                    for (int s : classes[new_id]) class_of[s] = new_id;
                    if (in_work[c]) {
                        if (!in_work[new_id]) {
                            work.push_back(new_id);
                            in_work[new_id] = 1;
                        }
                    } else {
                        int pick = (classes[c].size() <= classes[new_id].size()) ? c : new_id;
                        if (!in_work[pick]) {
                            work.push_back(pick);
                            in_work[pick] = 1;
                        }
                    }
                } else {
                    bucket[c].clear();
                }
            }
            for (int s : X) marked[s] = 0;
        }
    }

    Automaton min;
    min.trans.assign(classes.size(), {});
    min.accept.assign(classes.size(), 0);
    min.start = class_of[raw.start];
    for (int i = 0; i < static_cast<int>(classes.size()); ++i) {
        int rep = classes[i][0];
        min.accept[i] = accept[rep];
        for (int d = 0; d < kDigits; ++d) {
            min.trans[i][d] = class_of[trans[rep][d]];
        }
    }
    return min;
}

struct Sequences {
    std::vector<int64> counts;
    std::vector<int64> prefix;
};

Sequences build_sequences(const Automaton& aut, int terms) {
    int n = static_cast<int>(aut.trans.size());
    std::vector<int> accept_list;
    accept_list.reserve(n);
    for (int i = 0; i < n; ++i) {
        if (aut.accept[i]) accept_list.push_back(i);
    }

    std::vector<int64> counts(terms + 1, 0);
    std::vector<int64> prefix(terms + 1, 0);
    std::vector<int64> cur(n, 0), next(n, 0);
    cur[aut.start] = 1;

    auto apply_digits = [&](int d_start, int d_end) {
        std::fill(next.begin(), next.end(), 0);
        for (int i = 0; i < n; ++i) {
            int64 val = cur[i];
            if (!val) continue;
            for (int d = d_start; d <= d_end; ++d) {
                int j = aut.trans[i][d];
                int64& slot = next[j];
                slot += val;
                if (slot >= kMod) slot -= kMod;
            }
        }
        cur.swap(next);
    };

    apply_digits(1, 9);
    int64 sum = 0;
    for (int idx : accept_list) {
        sum += cur[idx];
        if (sum >= kMod) sum -= kMod;
    }
    counts[1] = sum;
    prefix[1] = sum;

    for (int len = 2; len <= terms; ++len) {
        apply_digits(0, 9);
        sum = 0;
        for (int idx : accept_list) {
            sum += cur[idx];
            if (sum >= kMod) sum -= kMod;
        }
        counts[len] = sum;
        prefix[len] = prefix[len - 1] + sum;
        if (prefix[len] >= kMod) prefix[len] -= kMod;
    }
    return {counts, prefix};
}

int64 mod_pow(int64 base, uint64_t exp) {
    int64 result = 1 % kMod;
    int64 cur = base % kMod;
    while (exp > 0) {
        if (exp & 1) result = (result * cur) % kMod;
        cur = (cur * cur) % kMod;
        exp >>= 1;
    }
    return result;
}

int64 mod_inv(int64 x) {
    return mod_pow(x, kMod - 2);
}

std::vector<int64> berlekamp_massey(const std::vector<int64>& seq) {
    std::vector<int64> C(1, 1), B(1, 1);
    int L = 0;
    int m = 1;
    int64 b = 1;
    for (int n = 0; n < static_cast<int>(seq.size()); ++n) {
        int64 d = 0;
        for (int i = 0; i <= L; ++i) {
            d = (d + C[i] * seq[n - i]) % kMod;
        }
        if (d == 0) {
            ++m;
            continue;
        }
        std::vector<int64> T = C;
        int64 coef = d * mod_inv(b) % kMod;
        if (C.size() < B.size() + m) C.resize(B.size() + m, 0);
        for (size_t i = 0; i < B.size(); ++i) {
            int64 val = (coef * B[i]) % kMod;
            C[i + m] -= val;
            if (C[i + m] < 0) C[i + m] += kMod;
        }
        if (2 * L <= n) {
            L = n + 1 - L;
            B.swap(T);
            b = d;
            m = 1;
        } else {
            ++m;
        }
    }
    C.erase(C.begin());
    for (auto& x : C) x = (kMod - x) % kMod;
    return C;
}

std::vector<int64> combine_poly(const std::vector<int64>& a, const std::vector<int64>& b,
                                const std::vector<int64>& coef, int threads) {
    int k = static_cast<int>(coef.size());
    std::vector<int64> res(2 * k, 0);
    if (threads <= 1 || k < 512) {
        for (int i = 0; i < k; ++i) {
            if (a[i] == 0) continue;
            for (int j = 0; j < k; ++j) {
                if (b[j] == 0) continue;
                int64 add = (a[i] * b[j]) % kMod;
                int64& slot = res[i + j];
                slot += add;
                if (slot >= kMod) slot -= kMod;
            }
        }
    } else {
        int tcount = std::min(threads, k);
        std::vector<std::vector<int64>> partial(tcount, std::vector<int64>(2 * k, 0));
        std::vector<std::thread> pool;
        pool.reserve(tcount);
        for (int t = 0; t < tcount; ++t) {
            int start = t * k / tcount;
            int end = (t + 1) * k / tcount;
            pool.emplace_back([&, t, start, end]() {
                auto& loc = partial[t];
                for (int i = start; i < end; ++i) {
                    if (a[i] == 0) continue;
                    for (int j = 0; j < k; ++j) {
                        if (b[j] == 0) continue;
                        int64 add = (a[i] * b[j]) % kMod;
                        int64& slot = loc[i + j];
                        slot += add;
                        if (slot >= kMod) slot -= kMod;
                    }
                }
            });
        }
        for (auto& th : pool) th.join();
        for (int i = 0; i < 2 * k - 1; ++i) {
            int64 sum = 0;
            for (int t = 0; t < tcount; ++t) {
                sum += partial[t][i];
            }
            res[i] = sum % kMod;
        }
    }

    for (int i = 2 * k - 2; i >= k; --i) {
        int64 val = res[i];
        if (val == 0) continue;
        for (int j = 0; j < k; ++j) {
            int64 add = (val * coef[j]) % kMod;
            int64& slot = res[i - 1 - j];
            slot += add;
            if (slot >= kMod) slot -= kMod;
        }
    }
    res.resize(k);
    return res;
}

int64 linear_recurrence(const std::vector<int64>& coef, const std::vector<int64>& init,
                        uint64_t n, int threads) {
    int k = static_cast<int>(coef.size());
    if (k == 0) return 0;
    if (n < init.size()) return init[n];
    if (k == 1) {
        return init[0] * mod_pow(coef[0], n) % kMod;
    }
    std::vector<int64> pol(k, 0), e(k, 0);
    pol[0] = 1;
    e[1] = 1;
    while (n > 0) {
        if (n & 1) pol = combine_poly(pol, e, coef, threads);
        e = combine_poly(e, e, coef, threads);
        n >>= 1;
    }
    int64 res = 0;
    for (int i = 0; i < k; ++i) {
        res = (res + pol[i] * init[i]) % kMod;
    }
    return res;
}

bool run_validations(const Sequences& seq) {
    struct Check {
        int n;
        int64 expected;
    };
    const std::vector<Check> checks = {
        {1, 0},
        {2, 9},
        {5, 3492},
    };

    for (const auto& check : checks) {
        if (check.n >= static_cast<int>(seq.prefix.size())) {
            std::cerr << "Validation skipped for T(" << check.n << "): term not computed.\n";
            continue;
        }
        int64 got = seq.prefix[check.n];
        if (got != check.expected) {
            std::cerr << "Validation failed for T(" << check.n << "): got=" << got
                      << " expected=" << check.expected << "\n";
            return false;
        }
    }
    if (seq.counts.size() > 2 && seq.counts[2] != 9) {
        std::cerr << "Validation failed for length 2: got=" << seq.counts[2] << " expected=9\n";
        return false;
    }
    return true;
}

}  // namespace

int main(int argc, char** argv) {
    std::ios::sync_with_stdio(false);
    std::cin.tie(nullptr);

    uint64_t n = 1000000000000000000ULL;
    if (argc > 1) n = static_cast<uint64_t>(std::strtoull(argv[1], nullptr, 10));

    int threads = static_cast<int>(std::thread::hardware_concurrency());
    if (threads == 0) threads = 4;
    if (argc > 2) threads = std::max(1, std::atoi(argv[2]));

    RawAutomaton raw = build_raw();
    Automaton aut = minimize_automaton(raw);

    int terms = 2 * static_cast<int>(aut.trans.size()) + 5;
    Sequences seq = build_sequences(aut, terms);
    if (!run_validations(seq)) return 1;

    std::vector<int64> prefix(seq.prefix.begin(), seq.prefix.end());
    std::vector<int64> coef = berlekamp_massey(prefix);
    std::vector<int64> init(prefix.begin(), prefix.begin() + static_cast<int>(coef.size()));
    int64 answer = linear_recurrence(coef, init, n, threads);
    std::cout << answer % kMod << "\n";
    return 0;
}

Python

from __future__ import annotations

import re
import shutil
import subprocess
from pathlib import Path

ANSWER_RE = re.compile(r"answer\s*:\s*(.+)$", re.IGNORECASE)
EQUAL_RE = re.compile(r"=\s*(.+)$")


def parse_output(stdout: str) -> str:
    lines = [line.strip() for line in stdout.splitlines() if line.strip()]
    if not lines:
        return ""
    answers = []
    equals = []
    for line in lines:
        m1 = ANSWER_RE.search(line)
        if m1:
            answers.append(m1.group(1).strip())
        m2 = EQUAL_RE.search(line)
        if m2:
            equals.append(m2.group(1).strip())
    if answers:
        return answers[-1]
    if equals:
        return equals[-1]
    return lines[-1]


def should_skip_cpp_checkpoints(src: Path) -> bool:
    try:
        text = src.read_text(encoding="utf-8", errors="ignore")
    except OSError:
        return False
    return "--skip-checkpoints" in text


def run_cpp(binary: Path, src: Path, root: Path) -> str:
    cmd = [str(binary)]
    if should_skip_cpp_checkpoints(src):
        cmd.append("--skip-checkpoints")

    try:
        return subprocess.check_output(cmd, text=True, cwd=root)
    except subprocess.CalledProcessError:
        return subprocess.check_output(cmd, text=True, cwd=src.parent)


def solve() -> str:
    problem_id = __file__.split("Euler")[-1].split(".")[0]
    root = Path(__file__).resolve().parent.parent
    src = root / "solutionsCpp" / f"Euler{problem_id}.cpp"
    binary = root / "solutionsCpp" / f".euler{problem_id}_py_bridge"

    if not binary.exists() or src.stat().st_mtime > binary.stat().st_mtime:
        compiler = shutil.which("clang++") or shutil.which("g++")
        if not compiler:
            raise RuntimeError("No C++ compiler found (clang++/g++).")
        subprocess.check_call([compiler, "-std=c++17", "-O2", str(src), "-o", str(binary)])

    output = run_cpp(binary=binary, src=src, root=root)
    parsed = parse_output(output)
    if not parsed:
        raise RuntimeError(f"Euler{problem_id} bridge produced empty output.")
    return parsed


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

Java

import java.util.*;

public class Euler529 {

    static final int kTargetSum = 10;
    static final int kDigits = 10;
    static final long kMod = 1000000007L;
    static final int kHistoryMask = (1 << (kTargetSum + 1)) - 1;

    static int encodeState(int history, int segment) {
        return history | (segment << 11);
    }

    static int highestBit(int mask) {
        return 31 - Integer.numberOfLeadingZeros(mask);
    }

    static class State {
        int history, segment;

        State(int h, int s) {
            history = h;
            segment = s;
        }
    }

    static State nextState(int history, int segment, int digit) {
        int delta = highestBit(segment);
        int newDelta = delta + digit;
        if (newDelta > kTargetSum)
            return null;
        int newSegment = segment | (1 << newDelta);
        if (((history >> (kTargetSum - newDelta)) & 1) == 1) {
            int shifted = (history << newDelta) & kHistoryMask;
            int newHistory = shifted;
            for (int k = 0; k <= newDelta; k++) {
                if (((newSegment >> k) & 1) == 1) {
                    newHistory |= (1 << (newDelta - k));
                }
            }
            return new State(newHistory, 1);
        } else {
            return new State(history, newSegment);
        }
    }

    static class RawAutomaton {
        int[][] trans;
        int[] accept;
        int start;
    }

    static RawAutomaton buildRaw() {
        List<State> states = new ArrayList<>();
        Map<Integer, Integer> index = new HashMap<>();

        State start = new State(1, 1);
        states.add(start);
        index.put(encodeState(1, 1), 0);

        for (int i = 0; i < states.size(); i++) {
            State cur = states.get(i);
            for (int d = 0; d < kDigits; d++) {
                State nxt = nextState(cur.history, cur.segment, d);
                if (nxt == null)
                    continue;
                int key = encodeState(nxt.history, nxt.segment);
                if (!index.containsKey(key)) {
                    index.put(key, states.size());
                    states.add(nxt);
                }
            }
        }

        int n = states.size();
        RawAutomaton raw = new RawAutomaton();
        raw.trans = new int[n][kDigits];
        raw.accept = new int[n];
        raw.start = 0;

        for (int i = 0; i < n; i++) {
            State cur = states.get(i);
            raw.accept[i] = (cur.segment == 1) ? 1 : 0;
            for (int d = 0; d < kDigits; d++) {
                State nxt = nextState(cur.history, cur.segment, d);
                if (nxt == null) {
                    raw.trans[i][d] = -1;
                } else {
                    raw.trans[i][d] = index.get(encodeState(nxt.history, nxt.segment));
                }
            }
        }
        return raw;
    }

    static class Automaton {
        int[][] trans;
        int[] accept;
        int start;
    }

    static Automaton minimizeAutomaton(RawAutomaton raw) {
        int n = raw.trans.length;
        int dead = n;
        int total = n + 1;

        int[][] trans = new int[total][kDigits];
        int[] accept = new int[total];
        for (int i = 0; i < n; i++) {
            accept[i] = raw.accept[i];
            for (int d = 0; d < kDigits; d++) {
                int t = raw.trans[i][d];
                trans[i][d] = (t == -1) ? dead : t;
            }
        }
        accept[dead] = 0;
        Arrays.fill(trans[dead], dead);

        List<List<Integer>>[] inv = new List[kDigits];
        for (int d = 0; d < kDigits; d++) {
            inv[d] = new ArrayList<>(total);
            for (int i = 0; i < total; i++)
                inv[d].add(new ArrayList<>());
        }
        for (int s = 0; s < total; s++) {
            for (int d = 0; d < kDigits; d++) {
                inv[d].get(trans[s][d]).add(s);
            }
        }

        List<List<Integer>> classes = new ArrayList<>();
        classes.add(new ArrayList<>());
        classes.add(new ArrayList<>());
        for (int s = 0; s < total; s++) {
            classes.get(accept[s] != 0 ? 0 : 1).add(s);
        }
        if (classes.get(0).isEmpty())
            classes.remove(0);
        else if (classes.size() == 2 && classes.get(1).isEmpty())
            classes.remove(1);

        int[] classOf = new int[total];
        for (int i = 0; i < classes.size(); i++) {
            for (int s : classes.get(i))
                classOf[s] = i;
        }

        Deque<Integer> work = new ArrayDeque<>();
        boolean[] inWork = new boolean[10000];
        for (int i = 0; i < classes.size(); i++) {
            work.add(i);
            inWork[i] = true;
        }

        boolean[] marked = new boolean[total];
        List<List<Integer>> bucket = new ArrayList<>();
        for (int i = 0; i < 10000; i++)
            bucket.add(new ArrayList<>());

        while (!work.isEmpty()) {
            int A = work.poll();
            inWork[A] = false;
            for (int d = 0; d < kDigits; d++) {
                List<Integer> X = new ArrayList<>();
                for (int t : classes.get(A)) {
                    for (int s : inv[d].get(t)) {
                        if (!marked[s]) {
                            marked[s] = true;
                            X.add(s);
                        }
                    }
                }
                if (X.isEmpty())
                    continue;

                List<Integer> touched = new ArrayList<>();
                for (int s : X) {
                    int c = classOf[s];
                    if (bucket.get(c).isEmpty())
                        touched.add(c);
                    bucket.get(c).add(s);
                }

                for (int c : touched) {
                    if (bucket.get(c).size() < classes.get(c).size()) {
                        List<Integer> inX = new ArrayList<>(bucket.get(c));
                        List<Integer> outX = new ArrayList<>();
                        for (int s : classes.get(c)) {
                            if (!marked[s])
                                outX.add(s);
                        }
                        classes.set(c, inX);
                        int newId = classes.size();
                        classes.add(outX);

                        for (int s : outX)
                            classOf[s] = newId;

                        if (inWork[c]) {
                            if (!inWork[newId]) {
                                work.add(newId);
                                inWork[newId] = true;
                            }
                        } else {
                            int pick = (classes.get(c).size() <= classes.get(newId).size()) ? c : newId;
                            if (!inWork[pick]) {
                                work.add(pick);
                                inWork[pick] = true;
                            }
                        }
                    }
                    bucket.get(c).clear();
                }
                for (int s : X)
                    marked[s] = false;
            }
        }

        Automaton min = new Automaton();
        min.trans = new int[classes.size()][kDigits];
        min.accept = new int[classes.size()];
        min.start = classOf[raw.start];

        for (int i = 0; i < classes.size(); i++) {
            int rep = classes.get(i).get(0);
            min.accept[i] = accept[rep];
            for (int d = 0; d < kDigits; d++) {
                min.trans[i][d] = classOf[trans[rep][d]];
            }
        }
        return min;
    }

    static long[] buildSequences(Automaton aut, int terms) {
        int n = aut.trans.length;
        List<Integer> acceptList = new ArrayList<>();
        for (int i = 0; i < n; i++)
            if (aut.accept[i] == 1)
                acceptList.add(i);

        long[] prefix = new long[terms + 1];
        long[] cur = new long[n];
        cur[aut.start] = 1;
        long[] nxt = new long[n];

        for (int i = 0; i < n; i++)
            nxt[i] = 0;
        for (int i = 0; i < n; i++) {
            if (cur[i] == 0)
                continue;
            for (int d = 1; d <= 9; d++) {
                int j = aut.trans[i][d];
                nxt[j] = (nxt[j] + cur[i]) % kMod;
            }
        }
        for (int i = 0; i < n; i++)
            cur[i] = nxt[i];

        long sum = 0;
        for (int idx : acceptList)
            sum = (sum + cur[idx]) % kMod;
        prefix[1] = sum;

        for (int len = 2; len <= terms; len++) {
            for (int i = 0; i < n; i++)
                nxt[i] = 0;
            for (int i = 0; i < n; i++) {
                if (cur[i] == 0)
                    continue;
                for (int d = 0; d <= 9; d++) {
                    int j = aut.trans[i][d];
                    nxt[j] = (nxt[j] + cur[i]) % kMod;
                }
            }
            for (int i = 0; i < n; i++)
                cur[i] = nxt[i];

            sum = 0;
            for (int idx : acceptList)
                sum = (sum + cur[idx]) % kMod;
            prefix[len] = (prefix[len - 1] + sum) % kMod;
        }
        return prefix;
    }

    static long modPow(long base, long exp) {
        long res = 1;
        long cur = base % kMod;
        while (exp > 0) {
            if ((exp & 1) == 1)
                res = (res * cur) % kMod;
            cur = (cur * cur) % kMod;
            exp >>= 1;
        }
        return res;
    }

    static List<Long> berlekampMassey(long[] seq) {
        List<Long> C = new ArrayList<>();
        C.add(1L);
        List<Long> B = new ArrayList<>();
        B.add(1L);
        int L = 0, m = 1;
        long b = 1;

        for (int n = 0; n < seq.length; n++) {
            long d = 0;
            for (int i = 0; i <= L; i++) {
                d = (d + C.get(i) * seq[n - i]) % kMod;
            }
            if (d == 0) {
                m++;
                continue;
            }
            List<Long> T = new ArrayList<>(C);
            long coef = (d * modPow(b, kMod - 2)) % kMod;
            while (C.size() < B.size() + m)
                C.add(0L);
            for (int i = 0; i < B.size(); i++) {
                long val = (C.get(i + m) - coef * B.get(i)) % kMod;
                if (val < 0)
                    val += kMod;
                C.set(i + m, val);
            }
            if (2 * L <= n) {
                L = n + 1 - L;
                B = T;
                b = d;
                m = 1;
            } else {
                m++;
            }
        }
        C.remove(0);
        for (int i = 0; i < C.size(); i++)
            C.set(i, (kMod - C.get(i)) % kMod);
        return C;
    }

    static long[] combinePoly(long[] a, long[] b, long[] coef) {
        int k = coef.length;
        long[] res = new long[2 * k];
        for (int i = 0; i < k; i++) {
            if (a[i] == 0)
                continue;
            for (int j = 0; j < k; j++) {
                if (b[j] == 0)
                    continue;
                res[i + j] = (res[i + j] + a[i] * b[j]) % kMod;
            }
        }
        for (int i = 2 * k - 2; i >= k; i--) {
            long val = res[i];
            if (val == 0)
                continue;
            for (int j = 0; j < k; j++) {
                res[i - 1 - j] = (res[i - 1 - j] + val * coef[j]) % kMod;
            }
        }
        return Arrays.copyOf(res, k);
    }

    static long linearRecurrence(long[] coef, long[] init, long n) {
        int k = coef.length;
        if (k == 0)
            return 0;
        if (n < init.length)
            return init[(int) n];
        if (k == 1)
            return (init[0] * modPow(coef[0], n)) % kMod;

        long[] pol = new long[k];
        pol[0] = 1;
        long[] e = new long[k];
        e[1] = 1;

        while (n > 0) {
            if ((n & 1) == 1)
                pol = combinePoly(pol, e, coef);
            e = combinePoly(e, e, coef);
            n >>= 1;
        }

        long res = 0;
        for (int i = 0; i < k; i++)
            res = (res + pol[i] * init[i]) % kMod;
        return res;
    }

    public static void main(String[] args) {
        RawAutomaton raw = buildRaw();
        Automaton aut = minimizeAutomaton(raw);
        int terms = 2 * aut.trans.length + 5;
        long[] prefix = buildSequences(aut, terms);

        List<Long> coefList = berlekampMassey(prefix);
        long[] coef = new long[coefList.size()];
        for (int i = 0; i < coef.length; i++)
            coef[i] = coefList.get(i);

        long[] init = Arrays.copyOf(prefix, coef.length);
        long answer = linearRecurrence(coef, init, 1000000000000000000L);
        System.out.println(answer);
    }
}