Problem 305: Reflexive Position

View on Project Euler

Project Euler Problem 305 Solution

EulerSolve provides an optimized solution for Project Euler Problem 305, Reflexive Position, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Consider the infinite concatenation $$C=12345678910111213\dots$$ For a given \(n\), let \(P=\text{str}(n)\). Define \(f(n)\) as the starting position of the \(n\)-th occurrence of \(P\) inside \(C\). The Project Euler sum is taken over \(n=3^k\). Mathematical Approach 1) Why this is a string-search problem on one continuous stream The pattern \(P\) can appear: 1. entirely inside one integer, 2. across the boundary between consecutive integers, 3. overlapping with itself. So we must search in the continuous digit stream \(123456789101112\dots\), not number by number. Resetting the match state at integer boundaries would miss valid occurrences. 2) KMP automaton for \(P=\text{str}(n)\) The code builds a Knuth-Morris-Pratt automaton for the decimal pattern \(P\). A state \(q\in\{0,\dots,L\}\), where \(L=|P|\), means: “the last processed digits currently match the first \(q\) digits of \(P\)”. For every state \(q\) and digit \(d\in\{0,\dots,9\}\), the automaton stores: $$q'=\delta(q,d),$$ and an output bit $$o(q,d)\in\{0,1\},$$ which tells whether reading digit \(d\) completes one new occurrence of \(P\). Because the KMP fallback is applied after a full match, overlapping occurrences are counted correctly. 3) Operator view of a digit block For any finite digit string \(w\), define an operator \(T_w\)....

Detailed mathematical approach

Problem Summary

Consider the infinite concatenation

$$C=12345678910111213\dots$$

For a given \(n\), let \(P=\text{str}(n)\). Define \(f(n)\) as the starting position of the \(n\)-th occurrence of \(P\) inside \(C\). The Project Euler sum is taken over \(n=3^k\).

Mathematical Approach

1) Why this is a string-search problem on one continuous stream

The pattern \(P\) can appear:

1. entirely inside one integer,

2. across the boundary between consecutive integers,

3. overlapping with itself.

So we must search in the continuous digit stream \(123456789101112\dots\), not number by number. Resetting the match state at integer boundaries would miss valid occurrences.

2) KMP automaton for \(P=\text{str}(n)\)

The code builds a Knuth-Morris-Pratt automaton for the decimal pattern \(P\). A state \(q\in\{0,\dots,L\}\), where \(L=|P|\), means: “the last processed digits currently match the first \(q\) digits of \(P\)”.

For every state \(q\) and digit \(d\in\{0,\dots,9\}\), the automaton stores:

$$q'=\delta(q,d),$$

and an output bit

$$o(q,d)\in\{0,1\},$$

which tells whether reading digit \(d\) completes one new occurrence of \(P\).

Because the KMP fallback is applied after a full match, overlapping occurrences are counted correctly.

3) Operator view of a digit block

For any finite digit string \(w\), define an operator \(T_w\). Starting from automaton state \(q\), this operator tells us two things after scanning \(w\):

$$T_w(q)=\bigl(q_{\text{end}}(q),c_w(q)\bigr),$$

where \(q_{\text{end}}(q)\) is the ending state and \(c_w(q)\) is the number of new matches found while reading \(w\).

This is exactly what the code stores as end_state and add_count.

4) Why operators compose so cleanly

If a block is split as \(w=uv\), then scanning \(u\) first and \(v\) second gives

$$T_{uv}(q)=\Bigl(q_v(q_u(q)),\,c_u(q)+c_v\bigl(q_u(q)\bigr)\Bigr).$$

So composition is associative. This is the central algebraic fact used by the solver: once we know the operator of a block, we never need to rescan its digits.

5) Full blocks and partial blocks

Let a decimal prefix already be fixed. Then:

1. full_block(prefix, rem_digits) means the concatenation of all numbers with that prefix followed by all \(10^{\text{rem\_digits}}\) possible suffixes;

2. partial_block(prefix, rem_digits, upper_suffix) means the same, but only up to a given suffix bound.

These operators are memoized. Since many queries reuse the same prefix-and-length structure, caching removes the exponential blowup.

6) Building the operator for \([1..N]\)

To process all integers from \(1\) to \(N\), the code does not concatenate them explicitly. Instead it composes operators in this order:

1. all complete digit-length blocks with fewer digits than \(N\);

2. within the top digit length, all complete leading-digit blocks smaller than the first digit of \(N\);

3. one final partial block for the prefix that leads exactly up to \(N\).

The result is an operator we may denote by

$$T_{[1..N]}.$$

Applying it to automaton state \(0\) yields the total number of occurrences of \(P\) in the finite prefix \(123\dots N\).

7) The monotone counting function \(A(N)\)

Define

$$A(N)=\text{number of occurrences of }P\text{ in }123\dots N.$$

Clearly \(A(N)\) is monotone nondecreasing in \(N\). Therefore we can binary-search the smallest \(N\) such that

$$A(N)\ge n.$$

That identifies the exact integer in which the \(n\)-th occurrence finishes.

8) Locating the occurrence inside the final number

After binary search finds this minimal \(N\), the code computes the automaton state after scanning \(123\dots(N-1)\). Then it scans the digits of \(N\) once more, counting only the matches produced inside that last number.

If the desired occurrence ends at digit index \(i\) inside \(\text{str}(N)\), and \(L=|P|\), then its global starting position is

$$\text{digits\_before}(N-1)+i-L+2.$$

Here

$$\text{digits\_before}(N)=\sum_{d=1}^{D-1}9\cdot 10^{d-1}\cdot d+\bigl(N-10^{D-1}+1\bigr)D,$$

where \(D\) is the number of digits of \(N\).

9) Checkpoints

The implementation verifies several known values:

$$f(1)=1,\qquad f(5)=81,\qquad f(12)=271,\qquad f(7780)=111111365.$$

It also brute-forces all cases \(1\le n\le 30\) against a direct stream simulation, which is a strong validation of the operator approach.

How the Code Works

The class PatternEngine does all heavy lifting. It builds the KMP automaton, interns operators so identical operators are reused, memoizes operator composition, and provides:

1. range_operator(N) for the whole prefix \(123\dots N\),

2. occ_and_state(N) for total match count and ending automaton state,

3. nth_occurrence_position(target) for the final binary search and local scan.

The outer solver simply evaluates this engine for patterns \(3,9,27,\dots,3^{13}\) and sums the returned positions.

Complexity Analysis

If \(L=|P|\), one operator composition costs \(O(L)\). Thanks to memoized full/partial blocks, evaluating \(A(N)\) is much closer to polylogarithmic in \(N\) than to linear in the total digit stream length. That is the whole reason the problem is feasible: the code never generates the massive prefix of \(C\) explicitly.

Further Reading

  1. Problem page: https://projecteuler.net/problem=305
  2. Knuth-Morris-Pratt algorithm: https://en.wikipedia.org/wiki/Knuth%E2%80%93Morris%E2%80%93Pratt_algorithm
  3. Champernowne constant: https://en.wikipedia.org/wiki/Champernowne_constant

Problem 305 source code

C++

#include <algorithm>
#include <array>
#include <cstdint>
#include <iostream>
#include <limits>
#include <numeric>
#include <string>
#include <unordered_map>
#include <utility>
#include <vector>

namespace {

using u64 = std::uint64_t;

struct Options {
    int k_max = 13;
    bool run_checkpoints = true;
    u64 single_n = 0;
};

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_u64_after_prefix(const std::string& arg, const std::string& prefix, u64& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }

    u64 parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<u64>(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, "--k-max=", options.k_max)) {
            continue;
        }
        if (parse_u64_after_prefix(arg, "--single=", options.single_n)) {
            continue;
        }
        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }

    return options.k_max >= 1;
}

struct Op {
    std::vector<std::uint8_t> end_state;
    std::vector<u64> add_count;

    bool operator==(const Op& other) const {
        return end_state == other.end_state && add_count == other.add_count;
    }
};

struct OpHash {
    std::size_t operator()(const Op& op) const {
        std::size_t h = 1469598103934665603ULL;
        for (std::uint8_t v : op.end_state) {
            h ^= static_cast<std::size_t>(v) + 0x9e3779b97f4a7c15ULL + (h << 6U) + (h >> 2U);
        }
        for (u64 v : op.add_count) {
            h ^= static_cast<std::size_t>(v ^ (v >> 32U)) + 0x9e3779b97f4a7c15ULL + (h << 6U) + (h >> 2U);
        }
        return h;
    }
};

class PatternEngine {
public:
    explicit PatternEngine(const std::string& pattern)
        : pattern_(pattern), L_(static_cast<int>(pattern.size())) {
        build_automaton();
        build_digit_operators();
        build_pow10(20);
    }

    u64 nth_occurrence_position(const u64 target_occurrence) {
        u64 hi = 1;
        while (occ_and_state(hi).first < target_occurrence) {
            if (hi > std::numeric_limits<u64>::max() / 2ULL) {
                hi = std::numeric_limits<u64>::max();
                break;
            }
            hi *= 2ULL;
        }

        u64 lo = 1;
        while (lo < hi) {
            const u64 mid = lo + (hi - lo) / 2ULL;
            if (occ_and_state(mid).first >= target_occurrence) {
                hi = mid;
            } else {
                lo = mid + 1ULL;
            }
        }

        const u64 number = lo;
        const auto [before_count, before_state] = occ_and_state(number - 1ULL);
        const u64 need_inside = target_occurrence - before_count;

        u64 prefix_digits = digits_before(number - 1ULL);
        u64 found_inside = 0;
        int state = before_state;

        const std::string num = std::to_string(number);
        for (std::size_t i = 0; i < num.size(); ++i) {
            const int d = num[i] - '0';
            if (out_[state][d] != 0) {
                ++found_inside;
                if (found_inside == need_inside) {
                    const u64 end_pos = prefix_digits + static_cast<u64>(i + 1);
                    return end_pos - static_cast<u64>(L_) + 1ULL;
                }
            }
            state = next_[state][d];
        }

        return 0;
    }

private:
    std::string pattern_;
    int L_ = 0;

    std::vector<std::vector<int>> next_;
    std::vector<std::vector<std::uint8_t>> out_;

    std::vector<Op> ops_;
    std::unordered_map<Op, int, OpHash> op_to_id_;

    int id_identity_ = -1;
    std::array<int, 10> digit_op_id_{};

    std::unordered_map<u64, int> compose_cache_;
    std::unordered_map<u64, int> append_cache_;
    std::unordered_map<u64, int> full_block_cache_;

    std::vector<u64> pow10_;

    static u64 make_key(const int a, const int b) {
        return (static_cast<u64>(a) << 32U) | static_cast<std::uint32_t>(b);
    }

    int intern_op(const Op& op) {
        const auto it = op_to_id_.find(op);
        if (it != op_to_id_.end()) {
            return it->second;
        }
        const int id = static_cast<int>(ops_.size());
        ops_.push_back(op);
        op_to_id_.emplace(ops_[static_cast<std::size_t>(id)], id);
        return id;
    }

    void build_automaton() {
        std::vector<int> pi(static_cast<std::size_t>(L_), 0);
        for (int i = 1; i < L_; ++i) {
            int j = pi[static_cast<std::size_t>(i - 1)];
            while (j > 0 && pattern_[static_cast<std::size_t>(i)] != pattern_[static_cast<std::size_t>(j)]) {
                j = pi[static_cast<std::size_t>(j - 1)];
            }
            if (pattern_[static_cast<std::size_t>(i)] == pattern_[static_cast<std::size_t>(j)]) {
                ++j;
            }
            pi[static_cast<std::size_t>(i)] = j;
        }

        next_.assign(static_cast<std::size_t>(L_ + 1), std::vector<int>(10, 0));
        out_.assign(static_cast<std::size_t>(L_ + 1), std::vector<std::uint8_t>(10, 0));

        for (int st = 0; st <= L_; ++st) {
            for (int d = 0; d <= 9; ++d) {
                const char c = static_cast<char>('0' + d);
                int j = st;
                while (j > 0 && (j == L_ || pattern_[static_cast<std::size_t>(j)] != c)) {
                    j = pi[static_cast<std::size_t>(j - 1)];
                }
                if (j < L_ && pattern_[static_cast<std::size_t>(j)] == c) {
                    ++j;
                }

                if (j == L_) {
                    out_[static_cast<std::size_t>(st)][static_cast<std::size_t>(d)] = 1;
                    j = pi[static_cast<std::size_t>(L_ - 1)];
                }
                next_[static_cast<std::size_t>(st)][static_cast<std::size_t>(d)] = j;
            }
        }
    }

    void build_digit_operators() {
        Op identity;
        identity.end_state.resize(static_cast<std::size_t>(L_ + 1));
        identity.add_count.assign(static_cast<std::size_t>(L_ + 1), 0);
        for (int s = 0; s <= L_; ++s) {
            identity.end_state[static_cast<std::size_t>(s)] = static_cast<std::uint8_t>(s);
        }
        id_identity_ = intern_op(identity);

        for (int d = 0; d <= 9; ++d) {
            Op op;
            op.end_state.resize(static_cast<std::size_t>(L_ + 1));
            op.add_count.resize(static_cast<std::size_t>(L_ + 1));
            for (int s = 0; s <= L_; ++s) {
                op.end_state[static_cast<std::size_t>(s)] =
                    static_cast<std::uint8_t>(next_[static_cast<std::size_t>(s)][static_cast<std::size_t>(d)]);
                op.add_count[static_cast<std::size_t>(s)] =
                    static_cast<u64>(out_[static_cast<std::size_t>(s)][static_cast<std::size_t>(d)]);
            }
            digit_op_id_[static_cast<std::size_t>(d)] = intern_op(op);
        }
    }

    void build_pow10(const int max_power) {
        pow10_.assign(static_cast<std::size_t>(max_power + 1), 1ULL);
        for (int i = 1; i <= max_power; ++i) {
            pow10_[static_cast<std::size_t>(i)] = pow10_[static_cast<std::size_t>(i - 1)] * 10ULL;
        }
    }

    int compose_ids(const int id_a, const int id_b) {
        const u64 key = make_key(id_a, id_b);
        const auto it = compose_cache_.find(key);
        if (it != compose_cache_.end()) {
            return it->second;
        }

        const Op& A = ops_[static_cast<std::size_t>(id_a)];
        const Op& B = ops_[static_cast<std::size_t>(id_b)];

        Op C;
        C.end_state.resize(static_cast<std::size_t>(L_ + 1));
        C.add_count.resize(static_cast<std::size_t>(L_ + 1));

        for (int s = 0; s <= L_; ++s) {
            const int mid = A.end_state[static_cast<std::size_t>(s)];
            C.end_state[static_cast<std::size_t>(s)] = B.end_state[static_cast<std::size_t>(mid)];
            C.add_count[static_cast<std::size_t>(s)] =
                A.add_count[static_cast<std::size_t>(s)] + B.add_count[static_cast<std::size_t>(mid)];
        }

        const int id = intern_op(C);
        compose_cache_.emplace(key, id);
        return id;
    }

    int append_digit(const int op_id, const int d) {
        const u64 key = make_key(op_id, d);
        const auto it = append_cache_.find(key);
        if (it != append_cache_.end()) {
            return it->second;
        }
        const int id = compose_ids(op_id, digit_op_id_[static_cast<std::size_t>(d)]);
        append_cache_.emplace(key, id);
        return id;
    }

    int full_block(const int prefix_op, const int rem_digits) {
        const u64 key = make_key(prefix_op, rem_digits);
        const auto it = full_block_cache_.find(key);
        if (it != full_block_cache_.end()) {
            return it->second;
        }

        int id;
        if (rem_digits == 0) {
            id = prefix_op;
        } else {
            int result = id_identity_;
            for (int d = 0; d <= 9; ++d) {
                const int child = full_block(append_digit(prefix_op, d), rem_digits - 1);
                result = compose_ids(result, child);
            }
            id = result;
        }

        full_block_cache_.emplace(key, id);
        return id;
    }

    int partial_block(const int prefix_op, const int rem_digits, const u64 upper_suffix) {
        if (rem_digits == 0) {
            return prefix_op;
        }

        const u64 base = pow10_[static_cast<std::size_t>(rem_digits - 1)];
        const int first = static_cast<int>(upper_suffix / base);
        const u64 rest = upper_suffix % base;

        int result = id_identity_;
        for (int d = 0; d < first; ++d) {
            const int child = full_block(append_digit(prefix_op, d), rem_digits - 1);
            result = compose_ids(result, child);
        }

        const int last = partial_block(append_digit(prefix_op, first), rem_digits - 1, rest);
        result = compose_ids(result, last);
        return result;
    }

    int range_operator(const u64 N) {
        if (N == 0) {
            return id_identity_;
        }

        const std::string s = std::to_string(N);
        const int D = static_cast<int>(s.size());

        int result = id_identity_;

        for (int d = 1; d < D; ++d) {
            const int rem = d - 1;
            for (int a = 1; a <= 9; ++a) {
                const int block = full_block(digit_op_id_[static_cast<std::size_t>(a)], rem);
                result = compose_ids(result, block);
            }
        }

        const int first = s[0] - '0';
        const int rem = D - 1;

        for (int a = 1; a < first; ++a) {
            const int block = full_block(digit_op_id_[static_cast<std::size_t>(a)], rem);
            result = compose_ids(result, block);
        }

        if (rem == 0) {
            result = compose_ids(result, digit_op_id_[static_cast<std::size_t>(first)]);
        } else {
            u64 suffix = 0;
            for (int i = 1; i < D; ++i) {
                suffix = suffix * 10ULL + static_cast<u64>(s[static_cast<std::size_t>(i)] - '0');
            }
            const int part = partial_block(digit_op_id_[static_cast<std::size_t>(first)], rem, suffix);
            result = compose_ids(result, part);
        }

        return result;
    }

    std::pair<u64, int> occ_and_state(const u64 N) {
        const int op_id = range_operator(N);
        const Op& op = ops_[static_cast<std::size_t>(op_id)];
        return {op.add_count[0], static_cast<int>(op.end_state[0])};
    }

    u64 digits_before(const u64 N) const {
        if (N == 0) {
            return 0;
        }

        const std::string s = std::to_string(N);
        const int D = static_cast<int>(s.size());

        u64 total = 0;
        for (int d = 1; d < D; ++d) {
            total += 9ULL * pow10_[static_cast<std::size_t>(d - 1)] * static_cast<u64>(d);
        }

        total += (N - pow10_[static_cast<std::size_t>(D - 1)] + 1ULL) * static_cast<u64>(D);
        return total;
    }
};

u64 brute_nth_position_small(const int n) {
    const std::string token = std::to_string(n);
    const int L = static_cast<int>(token.size());

    std::string stream;
    stream.reserve(10000);
    std::vector<u64> starts;
    starts.reserve(static_cast<std::size_t>(n));

    u64 checked_until = 0;

    for (int x = 1; static_cast<int>(starts.size()) < n; ++x) {
        const std::string s = std::to_string(x);
        const u64 old_size = static_cast<u64>(stream.size());
        stream += s;

        const u64 begin = (old_size >= static_cast<u64>(L - 1)) ? (old_size - static_cast<u64>(L - 1)) : 0ULL;
        for (u64 i = begin; i + static_cast<u64>(L) <= static_cast<u64>(stream.size()); ++i) {
            bool match = true;
            for (int j = 0; j < L; ++j) {
                if (stream[static_cast<std::size_t>(i + static_cast<u64>(j))] != token[static_cast<std::size_t>(j)]) {
                    match = false;
                    break;
                }
            }
            if (match) {
                if (starts.empty() || starts.back() != i + 1ULL) {
                    starts.push_back(i + 1ULL);
                }
            }
        }
        checked_until = old_size;
        (void)checked_until;
    }

    return starts[static_cast<std::size_t>(n - 1)];
}

bool run_checkpoints() {
    struct Sample {
        u64 n;
        u64 expected;
    };

    const std::vector<Sample> samples = {
        {1ULL, 1ULL},
        {5ULL, 81ULL},
        {12ULL, 271ULL},
        {7780ULL, 111111365ULL},
    };

    for (const auto& s : samples) {
        PatternEngine engine(std::to_string(s.n));
        const u64 got = engine.nth_occurrence_position(s.n);
        if (got != s.expected) {
            std::cerr << "Sample checkpoint failed for n=" << s.n << ": got " << got
                      << ", expected " << s.expected << '\n';
            return false;
        }
    }

    for (int n = 1; n <= 30; ++n) {
        PatternEngine engine(std::to_string(n));
        const u64 fast = engine.nth_occurrence_position(static_cast<u64>(n));
        const u64 brute = brute_nth_position_small(n);
        if (fast != brute) {
            std::cerr << "Brute checkpoint failed for n=" << n << ": fast=" << fast
                      << ", brute=" << brute << '\n';
            return false;
        }
    }

    return true;
}

u64 solve_sum_k(const int k_max) {
    u64 sum = 0;
    u64 n = 1;
    for (int k = 1; k <= k_max; ++k) {
        n *= 3ULL;
        PatternEngine engine(std::to_string(n));
        sum += engine.nth_occurrence_position(n);
    }
    return sum;
}

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

    if (options.single_n > 0) {
        PatternEngine engine(std::to_string(options.single_n));
        std::cout << engine.nth_occurrence_position(options.single_n) << '\n';
        return 0;
    }

    std::cout << solve_sum_k(options.k_max) << '\n';
    return 0;
}

Python

import sys

class PatternEngine:
    def __init__(self, pattern):
        self.pattern = pattern
        self.L = len(pattern)
        self.next_state = [[0] * 10 for _ in range(self.L + 1)]
        self.out = [[0] * 10 for _ in range(self.L + 1)]
        self.build_automaton()
        
        self.ops = []
        self.op_to_id = {}
        
        self.id_identity = -1
        self.digit_op_id = [0] * 10
        self.build_digit_operators()
        
        self.compose_cache = {}
        self.append_cache = {}
        self.full_block_cache = {}
        
        self.pow10 = [1] * 21
        for i in range(1, 21):
            self.pow10[i] = self.pow10[i - 1] * 10

    def build_automaton(self):
        pi = [0] * self.L
        for i in range(1, self.L):
            j = pi[i - 1]
            while j > 0 and self.pattern[i] != self.pattern[j]:
                j = pi[j - 1]
            if self.pattern[i] == self.pattern[j]:
                j += 1
            pi[i] = j
            
        for st in range(self.L + 1):
            for d in range(10):
                c = str(d)
                j = st
                while j > 0 and (j == self.L or self.pattern[j] != c):
                    j = pi[j - 1]
                if j < self.L and self.pattern[j] == c:
                    j += 1
                if j == self.L:
                    self.out[st][d] = 1
                    j = pi[self.L - 1]
                self.next_state[st][d] = j

    def intern_op(self, end_state, add_count):
        key = (tuple(end_state), tuple(add_count))
        if key in self.op_to_id:
            return self.op_to_id[key]
        op_id = len(self.ops)
        self.ops.append((end_state, add_count))
        self.op_to_id[key] = op_id
        return op_id

    def build_digit_operators(self):
        end_state = list(range(self.L + 1))
        add_count = [0] * (self.L + 1)
        self.id_identity = self.intern_op(end_state, add_count)
        
        for d in range(10):
            es = [0] * (self.L + 1)
            ac = [0] * (self.L + 1)
            for s in range(self.L + 1):
                es[s] = self.next_state[s][d]
                ac[s] = self.out[s][d]
            self.digit_op_id[d] = self.intern_op(es, ac)

    def compose_ids(self, id_a, id_b):
        key = (id_a, id_b)
        if key in self.compose_cache:
            return self.compose_cache[key]
            
        A_es, A_ac = self.ops[id_a]
        B_es, B_ac = self.ops[id_b]
        
        C_es = [0] * (self.L + 1)
        C_ac = [0] * (self.L + 1)
        for s in range(self.L + 1):
            mid = A_es[s]
            C_es[s] = B_es[mid]
            C_ac[s] = A_ac[s] + B_ac[mid]
            
        res = self.intern_op(C_es, C_ac)
        self.compose_cache[key] = res
        return res

    def append_digit(self, op_id, d):
        key = (op_id, d)
        if key in self.append_cache:
            return self.append_cache[key]
        res = self.compose_ids(op_id, self.digit_op_id[d])
        self.append_cache[key] = res
        return res

    def full_block(self, prefix_op, rem_digits):
        key = (prefix_op, rem_digits)
        if key in self.full_block_cache:
            return self.full_block_cache[key]
            
        if rem_digits == 0:
            res = prefix_op
        else:
            res = self.id_identity
            for d in range(10):
                child = self.full_block(self.append_digit(prefix_op, d), rem_digits - 1)
                res = self.compose_ids(res, child)
                
        self.full_block_cache[key] = res
        return res

    def partial_block(self, prefix_op, rem_digits, upper_suffix):
        if rem_digits == 0:
            return prefix_op
            
        base = self.pow10[rem_digits - 1]
        first = int(upper_suffix // base)
        rest = int(upper_suffix % base)
        
        res = self.id_identity
        for d in range(first):
            child = self.full_block(self.append_digit(prefix_op, d), rem_digits - 1)
            res = self.compose_ids(res, child)
            
        last = self.partial_block(self.append_digit(prefix_op, first), rem_digits - 1, rest)
        res = self.compose_ids(res, last)
        return res

    def range_operator(self, N):
        if N == 0:
            return self.id_identity
            
        s = str(N)
        D = len(s)
        res = self.id_identity
        
        for d in range(1, D):
            rem = d - 1
            for a in range(1, 10):
                block = self.full_block(self.digit_op_id[a], rem)
                res = self.compose_ids(res, block)
                
        first = int(s[0])
        rem = D - 1
        for a in range(1, first):
            block = self.full_block(self.digit_op_id[a], rem)
            res = self.compose_ids(res, block)
            
        if rem == 0:
            res = self.compose_ids(res, self.digit_op_id[first])
        else:
            suffix = int(s[1:])
            part = self.partial_block(self.digit_op_id[first], rem, suffix)
            res = self.compose_ids(res, part)
            
        return res

    def occ_and_state(self, N):
        op_id = self.range_operator(N)
        es, ac = self.ops[op_id]
        return ac[0], es[0]

    def digits_before(self, N):
        if N == 0:
            return 0
            
        s = str(N)
        D = len(s)
        total = 0
        for d in range(1, D):
            total += 9 * self.pow10[d - 1] * d
            
        total += (N - self.pow10[D - 1] + 1) * D
        return total

    def nth_occurrence_position(self, target_occurrence):
        hi = 1
        while self.occ_and_state(hi)[0] < target_occurrence:
            hi *= 2
            
        lo = 1
        while lo < hi:
            mid = (lo + hi) // 2
            if self.occ_and_state(mid)[0] >= target_occurrence:
                hi = mid
            else:
                lo = mid + 1
                
        number = lo
        before_count, before_state = self.occ_and_state(number - 1)
        need_inside = target_occurrence - before_count
        
        prefix_digits = self.digits_before(number - 1)
        found_inside = 0
        state = before_state
        
        num_str = str(number)
        for i, ch in enumerate(num_str):
            d = int(ch)
            if self.out[state][d]:
                found_inside += 1
                if found_inside == need_inside:
                    end_pos = prefix_digits + i + 1
                    return end_pos - self.L + 1
            state = self.next_state[state][d]
            
        return 0

def solve(k_max=13):
    total = 0
    n = 1
    for k in range(1, k_max + 1):
        n *= 3
        engine = PatternEngine(str(n))
        total += engine.nth_occurrence_position(n)
        
    return str(total)

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

Java

import java.util.*;

public class Euler305 {

    static class Op {
        byte[] endState;
        long[] addCount;

        Op(byte[] endState, long[] addCount) {
            this.endState = endState;
            this.addCount = addCount;
        }

        @Override
        public boolean equals(Object obj) {
            Op other = (Op) obj;
            return Arrays.equals(endState, other.endState) && Arrays.equals(addCount, other.addCount);
        }

        @Override
        public int hashCode() {
            int h = 1;
            for (byte b : endState)
                h = 31 * h + b;
            for (long c : addCount)
                h = 31 * h + (int) (c ^ (c >>> 32));
            return h;
        }
    }

    static class PatternEngine {
        String pattern;
        int L;

        int[][] nextState;
        byte[][] out;

        List<Op> ops = new ArrayList<>();
        Map<Op, Integer> opToId = new HashMap<>();

        int idIdentity = -1;
        int[] digitOpId = new int[10];

        Map<Long, Integer> composeCache = new HashMap<>();
        Map<Long, Integer> appendCache = new HashMap<>();
        Map<Long, Integer> fullBlockCache = new HashMap<>();

        long[] pow10 = new long[21];

        PatternEngine(String pattern) {
            this.pattern = pattern;
            this.L = pattern.length();
            nextState = new int[L + 1][10];
            out = new byte[L + 1][10];

            buildAutomaton();
            buildDigitOperators();

            pow10[0] = 1;
            for (int i = 1; i <= 20; ++i) {
                pow10[i] = pow10[i - 1] * 10L;
            }
        }

        long makeKey(int a, int b) {
            return (((long) a) << 32) | (b & 0xFFFFFFFFL);
        }

        int internOp(Op op) {
            Integer id = opToId.get(op);
            if (id != null)
                return id;
            int newId = ops.size();
            ops.add(op);
            opToId.put(op, newId);
            return newId;
        }

        void buildAutomaton() {
            int[] pi = new int[L];
            for (int i = 1; i < L; ++i) {
                int j = pi[i - 1];
                while (j > 0 && pattern.charAt(i) != pattern.charAt(j)) {
                    j = pi[j - 1];
                }
                if (pattern.charAt(i) == pattern.charAt(j)) {
                    j++;
                }
                pi[i] = j;
            }

            for (int st = 0; st <= L; ++st) {
                for (int d = 0; d < 10; ++d) {
                    char c = (char) ('0' + d);
                    int j = st;
                    while (j > 0 && (j == L || pattern.charAt(j) != c)) {
                        j = pi[j - 1];
                    }
                    if (j < L && pattern.charAt(j) == c) {
                        j++;
                    }
                    if (j == L) {
                        out[st][d] = 1;
                        j = pi[L - 1];
                    }
                    nextState[st][d] = j;
                }
            }
        }

        void buildDigitOperators() {
            byte[] es = new byte[L + 1];
            long[] ac = new long[L + 1];
            for (int s = 0; s <= L; ++s)
                es[s] = (byte) s;
            idIdentity = internOp(new Op(es, ac));

            for (int d = 0; d < 10; ++d) {
                byte[] esd = new byte[L + 1];
                long[] acd = new long[L + 1];
                for (int s = 0; s <= L; ++s) {
                    esd[s] = (byte) nextState[s][d];
                    acd[s] = out[s][d];
                }
                digitOpId[d] = internOp(new Op(esd, acd));
            }
        }

        int composeIds(int idA, int idB) {
            long key = makeKey(idA, idB);
            Integer cached = composeCache.get(key);
            if (cached != null)
                return cached;

            Op A = ops.get(idA);
            Op B = ops.get(idB);

            byte[] CEs = new byte[L + 1];
            long[] CAc = new long[L + 1];

            for (int s = 0; s <= L; ++s) {
                int mid = A.endState[s];
                CEs[s] = B.endState[mid];
                CAc[s] = A.addCount[s] + B.addCount[mid];
            }

            int id = internOp(new Op(CEs, CAc));
            composeCache.put(key, id);
            return id;
        }

        int appendDigit(int opId, int d) {
            long key = makeKey(opId, d);
            Integer cached = appendCache.get(key);
            if (cached != null)
                return cached;
            int id = composeIds(opId, digitOpId[d]);
            appendCache.put(key, id);
            return id;
        }

        int fullBlock(int prefixOp, int remDigits) {
            long key = makeKey(prefixOp, remDigits);
            Integer cached = fullBlockCache.get(key);
            if (cached != null)
                return cached;

            int id;
            if (remDigits == 0) {
                id = prefixOp;
            } else {
                int res = idIdentity;
                for (int d = 0; d < 10; ++d) {
                    int child = fullBlock(appendDigit(prefixOp, d), remDigits - 1);
                    res = composeIds(res, child);
                }
                id = res;
            }
            fullBlockCache.put(key, id);
            return id;
        }

        int partialBlock(int prefixOp, int remDigits, long upperSuffix) {
            if (remDigits == 0)
                return prefixOp;

            long base = pow10[remDigits - 1];
            int first = (int) (upperSuffix / base);
            long rest = upperSuffix % base;

            int res = idIdentity;
            for (int d = 0; d < first; ++d) {
                int child = fullBlock(appendDigit(prefixOp, d), remDigits - 1);
                res = composeIds(res, child);
            }

            int last = partialBlock(appendDigit(prefixOp, first), remDigits - 1, rest);
            res = composeIds(res, last);
            return res;
        }

        int rangeOperator(long N) {
            if (N == 0)
                return idIdentity;

            String s = String.valueOf(N);
            int D = s.length();
            int res = idIdentity;

            for (int d = 1; d < D; ++d) {
                int rem = d - 1;
                for (int a = 1; a < 10; ++a) {
                    int block = fullBlock(digitOpId[a], rem);
                    res = composeIds(res, block);
                }
            }

            int first = s.charAt(0) - '0';
            int rem = D - 1;
            for (int a = 1; a < first; ++a) {
                int block = fullBlock(digitOpId[a], rem);
                res = composeIds(res, block);
            }

            if (rem == 0) {
                res = composeIds(res, digitOpId[first]);
            } else {
                long suffix = Long.parseLong(s.substring(1));
                int part = partialBlock(digitOpId[first], rem, suffix);
                res = composeIds(res, part);
            }

            return res;
        }

        long[] occAndState(long N) {
            int opId = rangeOperator(N);
            Op op = ops.get(opId);
            return new long[] { op.addCount[0], op.endState[0] };
        }

        long digitsBefore(long N) {
            if (N == 0)
                return 0;

            String s = String.valueOf(N);
            int D = s.length();

            long total = 0;
            for (int d = 1; d < D; ++d) {
                total += 9L * pow10[d - 1] * d;
            }

            total += (N - pow10[D - 1] + 1) * D;
            return total;
        }

        long nthOccurrencePosition(long targetOccurrence) {
            long hi = 1;
            while (occAndState(hi)[0] < targetOccurrence) {
                hi *= 2;
            }

            long lo = 1;
            while (lo < hi) {
                long mid = lo + (hi - lo) / 2;
                if (occAndState(mid)[0] >= targetOccurrence) {
                    hi = mid;
                } else {
                    lo = mid + 1;
                }
            }

            long number = lo;
            long[] before = occAndState(number - 1);
            long beforeCount = before[0];
            int beforeState = (int) before[1];

            long needInside = targetOccurrence - beforeCount;
            long prefixDigits = digitsBefore(number - 1);

            long foundInside = 0;
            int state = beforeState;
            String numStr = String.valueOf(number);

            for (int i = 0; i < numStr.length(); ++i) {
                int d = numStr.charAt(i) - '0';
                if (out[state][d] != 0) {
                    foundInside++;
                    if (foundInside == needInside) {
                        long endPos = prefixDigits + i + 1;
                        return endPos - L + 1;
                    }
                }
                state = nextState[state][d];
            }

            return 0;
        }
    }

    public static String solve() {
        int kMax = 13;
        long sum = 0;
        long n = 1;
        for (int k = 1; k <= kMax; ++k) {
            n *= 3;
            PatternEngine engine = new PatternEngine(String.valueOf(n));
            sum += engine.nthOccurrencePosition(n);
        }
        return String.valueOf(sum);
    }

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