Problem 417: Reciprocal Cycles II

Project Euler Problem 417 Solution

EulerSolve provides an optimized solution for Project Euler Problem 417, Reciprocal Cycles II, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Define $L(n)$ as the length of the repeating part in the decimal expansion of $1/n$. The task is to compute $$S(N)=\sum_{n=3}^{N} L(n)$$ for a very large bound $N$. A direct simulation of long division for every denominator is much too slow, so the implementations reorganize the sum using multiplicative order and a grouping argument over powers of $2$ and $5$. Mathematical Approach Step 1: Remove the terminating part Every denominator can be written uniquely as $$n=2^a5^b m,\qquad \gcd(m,10)=1.$$ The powers of $2$ and $5$ only create the finite prefix of the decimal expansion. The recurring part depends entirely on the remaining factor $m$. If $m=1$, then $1/n$ is a terminating decimal and therefore $L(n)=0$. If $m>1$, the repetend length is the multiplicative order of $10$ modulo $m$: $$L(n)=\operatorname{ord}_m(10),$$ where $\operatorname{ord}_m(10)$ is the smallest positive integer $k$ such that $$10^k\equiv 1 \pmod m.$$ This is the standard number-theoretic description of repeating decimals. Step 2: Group denominators by the same odd core Fix an odd integer $m$ with $5 \nmid m$. Every denominator of the form $$n=m\,s,\qquad s=2^a5^b,$$ has exactly the same recurring-cycle length, namely $\operatorname{ord}_m(10)$, because multiplying by $2$ and $5$ changes only the terminating prefix....

Detailed mathematical approach

Problem Summary

A unit fraction $1/n$ has $1$ in the numerator. Its decimal expansion either stops, as $1/8=0.125$ does, or it ends in a block of digits that repeats for ever: $1/7=0.(142857)$ has a $6$-digit recurring cycle, and $1/6=0.1(6)$, which means $0.1666\ldots$, has a $1$-digit one. Unit fractions whose denominator has no prime factor other than $2$ or $5$ count as having no recurring cycle, of length $0$. Writing $L(n)$ for the length of the recurring cycle of $1/n$, we must evaluate

$$S(N)=\sum_{n=3}^{N}L(n),\qquad N=10^8.$$

The statement gives $S(10^6)=55535191115$, and its table of $1/2$ to $1/10$ yields the small checkpoint $S(10)=1+1+6+1=9$, the contributions of $n=3,6,7,9$. ^[1]

The whole solution in one sentence: the cycle length of $1/n$ depends only on the part $m$ of $n$ left after all factors $2$ and $5$ are removed — it is $0$ when nothing is left, and otherwise the multiplicative order of $10$ modulo $m$ — so we compute that order once for each of the $39\,999\,999$ surviving values of $m$ and multiply it by the number of denominators $n\le N$ that reduce to that same $m$.

Mathematical Approach

1. What repeats, and why

Long division for $1/7$ keeps a remainder $r$; each step writes the digit $\lfloor 10r/7\rfloor$ and then replaces $r$ by $10r\bmod 7$. Starting from $r=1$ the remainders are

$$1\to3\to2\to6\to4\to5\to1,$$

producing the digits $1,4,2,8,5,7$ and then starting over, because the next remainder is determined by the current one alone. For $1/6$ the remainders are $1\to4\to4\to4\to\cdots$: the first remainder never comes back, so one digit appears before the cycle, which is the meaning of $0.1(6)$. In general the expansion repeats from the first remainder that occurs twice, and the length of the cycle is the number of steps between its two occurrences. Simulating this for every $n$ up to $10^8$ would cost about $S(N)$ steps, the very quantity asked for, so we need a formula.

2. The cycle length is a multiplicative order

For $\gcd(a,m)=1$ the multiplicative order $\operatorname{ord}_m(a)$ is the least $k\ge1$ with $a^k\equiv1\pmod m$. ^[3] The order of $10$ exists whenever $\gcd(m,10)=1$ and $m\gt1$, because two of the $m$ powers $10^0,\ldots,10^{m-1}$ must share a residue class and the smaller one cancels; and $10^{u}\equiv1\pmod m$ holds precisely when $k\mid u$, since $u=qk+r$ with $0\le r\lt k$ would give $10^{r}\equiv1$. The powers $10^0,\ldots,10^{k-1}$ are pairwise distinct non-zero residues, so

$$\operatorname{ord}_m(10)\le m-1.$$

Order and digits are linked by one criterion. Take $0\lt x\lt1$ with the expansion that does not end in all $9$s, say $x=0.d_1d_2d_3\ldots$, and let $f$ be the fractional part of $10^{s}x$, that is $0.d_{s+1}d_{s+2}\ldots$. The digits satisfy $d_{j+t}=d_j$ for every $j\gt s$ exactly when multiplying $f$ by $10^{t}$ leaves its fractional part unchanged, that is when $(10^{t}-1)f$ is an integer; adding the integer $(10^{t}-1)\lfloor10^{s}x\rfloor$ turns this into

$$10^{s}\bigl(10^{t}-1\bigr)x\in\mathbb{Z}.$$

Now take $x=1/n$ and write $n=2^{a}5^{b}m$ with $\gcd(m,10)=1$, which is unique. The condition becomes $n\mid 10^{s}(10^{t}-1)$, and since $10^{t}-1$ is coprime to $10$ it splits into two independent demands: $2^{a}5^{b}\mid10^{s}$, that is $s\ge\max(a,b)$, and $m\mid10^{t}-1$, which for $m\gt1$ says $\operatorname{ord}_m(10)\mid t$. Hence the expansion of $1/n$ has a non-repeating prefix of exactly $\max(a,b)$ digits, and for $m\gt1$ its shortest cycle is

$$L(n)=\operatorname{ord}_m(10).$$

This matches the classical description of repeating decimals. ^[2] For $m=1$ the second demand is empty: after $\max(a,b)$ digits only zeros follow, so $t=1$ works and the shortest cycle would be the digit $0$ — or the digit $9$ in the other representation $0.125=0.124(9)$. The statement cuts this ambiguity by setting $L=0$ in exactly this case.

3. The core of a denominator

Call $m$ the core of $n$. Factors $2$ and $5$ only lengthen the prefix, so all denominators with the same core share one cycle length: $L(28)=L(70)=L(7)=6$ because $28=4\cdot7$ and $70=10\cdot7$. A core is odd and not divisible by $5$; there are $39\,999\,999$ of them in $[3,10^8]$, against the $99\,999\,998$ denominators of the requested sum.

4. Summing by core

Group the denominators by their core. For a fixed core $m$ the denominators with that core are exactly the numbers $ms$ with $s=2^{a}5^{b}$, and $ms\le N$ means $s\le\lfloor N/m\rfloor$. Writing

$$C(Q)=\#\bigl\{(a,b)\in\mathbb{Z}_{\ge0}^2:2^{a}5^{b}\le Q\bigr\}$$

for the number of admissible multipliers, we get the identity the whole program rests on:

$$S(N)=\sum_{\substack{3\le m\le N\\ \gcd(m,10)=1}}\operatorname{ord}_m(10)\cdot C\!\left(\left\lfloor\frac{N}{m}\right\rfloor\right).$$

The core $m=1$ is omitted because its group contributes $0$; for the same reason it does not matter that the official sum starts at $n=3$, since $1$ and $2$ have core $1$.

5. Counting the multipliers

The numbers $2^{a}5^{b}\le Q$ are scarce: $a\le\log_2Q$ and, for each $a$, there are $\lfloor\log_5(Q/2^{a})\rfloor+1$ choices of $b$, so $C(Q)=O(\log^2Q)$. Up to $10^6$ there are exactly $100$ such numbers and up to $10^8$ exactly $172$. It therefore pays to build the sorted list once and to answer every $C(Q)$ by a binary search in it, at most $8$ comparisons.

6. Two worked examples

For $N=10$ the cores are $3,7,9$. Core $3$ has $\lfloor10/3\rfloor=3$, so $C(3)=2$ multipliers $\{1,2\}$, i.e. the denominators $3$ and $6$, each with $\operatorname{ord}_3(10)=1$. Cores $7$ and $9$ have $C(1)=1$, with orders $6$ and $1$. The identity gives $1\cdot2+6\cdot1+1\cdot1=9=S(10)$, and the four denominators left out, $4,5,8,10$, are exactly those of core $1$.

For $N=100$ there are $39$ cores. Core $3$ now has $C(33)=10$ multipliers, giving the ten denominators $3,6,12,15,24,30,48,60,75,96$, each of cycle length $1$, so it contributes $10$; core $7$ has $C(14)=6$, giving $7,14,28,35,56,70$ with cycle length $6$ each, so it contributes $36$. Over all $39$ cores the total is $S(100)=1177$, which direct long division for $n=3,\ldots,100$ confirms.

7. Splitting the order over prime powers

Only the $39\,999\,999$ orders are still expensive. Let $m=\prod_i p_i^{e_i}$. By the Chinese remainder theorem, $10^{k}\equiv1\pmod m$ holds exactly when $10^{k}\equiv1\pmod{p_i^{e_i}}$ for every $i$, ^[5] that is when every $\operatorname{ord}_{p_i^{e_i}}(10)$ divides $k$, so the least such $k$ is

$$\operatorname{ord}_m(10)=\operatorname{lcm}_i\operatorname{ord}_{p_i^{e_i}}(10).$$

For instance $\operatorname{ord}_{21}(10)=\operatorname{lcm}(1,6)=6$. Factorising $m$ is cheap with a smallest-prime-factor sieve, so everything reduces to a table of prime-power orders.

8. The order modulo a prime

For a prime $p\notin\{2,5\}$, Fermat's little theorem gives $10^{p-1}\equiv1\pmod p$, hence $\operatorname{ord}_p(10)\mid p-1$. ^[4] That suggests the standard descent: start with $d=p-1$ and, for each prime $q$ dividing $p-1$ in turn, replace $d$ by $d/q$ as long as $q\mid d$ and $10^{d/q}\equiv1\pmod p$. The result is the order. Indeed $\operatorname{ord}_p(10)\mid d$ holds at the start and survives every accepted division; a rejected one shows $\operatorname{ord}_p(10)\nmid d/q$, which with $\operatorname{ord}_p(10)\mid d$ pins the exponent of $q$ in the order to its exponent in $d$. Later divisions use other primes and leave that exponent alone, so once every $q$ has been processed, $d$ is the order.

Take $p=13$, so $d=12=2^2\cdot3$. For $q=2$: $10^{6}\equiv1$, so $d=6$; then $10^{3}\equiv12\ne1$, so halving stops. For $q=3$: $10^{2}\equiv9\ne1$, so nothing changes. Thus $\operatorname{ord}_{13}(10)=6$, and indeed $1/13=0.(076923)$, the period listed for $13$ in the published table of periods. ^[6] Each test is one modular exponentiation: mod_pow halves the exponent and squares the running residue at every step, multiplying it into the result when the exponent is odd and reducing after each product, so a test costs $O(\log N)$ multiplications below the modulus. Measured over the $78\,496$ relevant primes below $10^6$, the descent needs $3.93$ tests per prime on average and $12$ at most.

9. Lifting to prime powers

Let $t=\operatorname{ord}_{p^{e-1}}(10)$ and $T=\operatorname{ord}_{p^{e}}(10)$ with $e\ge2$. From $10^{T}\equiv1\pmod{p^{e}}$ we get $10^{T}\equiv1\pmod{p^{e-1}}$, so $t\mid T$. Conversely $10^{t}=1+kp^{e-1}$ for some integer $k$, and the binomial expansion

$$\bigl(1+kp^{e-1}\bigr)^{p}=1+kp^{e}+\binom{p}{2}k^2p^{2e-2}+\cdots$$

is $\equiv1\pmod{p^{e}}$ because $2e-2\ge e$, so $T\mid tp$. As $p$ is prime, $T\in\{t,tp\}$, and a single test of $10^{t}\bmod p^{e}$ decides which. At the first lift $e=2$ we have $t\mid p-1$ but $tp\nmid p-1$, so $T=t$ means exactly $10^{p-1}\equiv1\pmod{p^2}$; below $10^4$ only $p=3$ and $p=487$ satisfy that, the base-$10$ analogues of Wieferich primes, whose only other known member is $56\,598\,313$. ^[7] Both are live at $N=10^8$, since $487^2=237169\le10^8$, and for them the order stalls just once: the powers of $3$ have orders $1,1,3,9,27,\ldots$, that is $3^{e-2}$ from $e=2$ on. So the test is not decoration.

A prime can appear with exponent $e\ge2$ in a core $m\le N$ only if $p^2\le N$, so these tables are needed only for $p\le\sqrt N=10^4$; the longest is the one for $p=3$, where $3^{16}=43046721\le10^8\lt3^{17}$.

How the Code Works

The C++ program ^[8] reads --limit= and --skip-checkpoints in parse_arguments and, unless told otherwise, runs run_checkpoints first: it verifies $S(10)=9$ by long division, compares solve_fast with solve_bruteforce at $N=1000$, and checks solve_fast against the published $S(10^6)$.

solve_fast then builds the tables that the grouped sum needs. build_odd_spf sieves the odd numbers only: for an odd composite $x$ the cell spf[x>>1] holds its smallest prime factor, which is at most $\sqrt x\le10^4$ and so a 16-bit value, while the cell of an odd prime stays $0$ because nothing ever marks it. collect_primes reads off the odd primes other than $5$. multiplicative_order_prime is Section 8: it strips the factor $2$ of $p-1$ before consulting spf, which is essential, since the odd-only table would otherwise be indexed with an even number and return the factor of a different integer. The order of every collected prime is computed once and stored in its own slot of ord_prime. build_small_prime_power_orders is Section 9: for every $p\le\sqrt N$ it finds the largest exponent $e$ with $p^{e}\le N$ and lifts level by level, testing mod_pow(10, prev, pp). build_smooth_numbers produces the $172$ numbers $2^{a}5^{b}\le N$, sorted.

The main sweep visits every odd $n$ from $3$ to the limit, skips $n\equiv0\pmod 5$ — the cores are coprime to $10$, and ord_prime deliberately has no entry for the prime $5$ — and calls the routine that rebuilds the order: it peels prime powers off spf, reads each part from ord_prime or the prime-power tables, and joins them with std::gcd into the least common multiple of Section 7. A std::upper_bound in the smooth list turns $\lfloor N/n\rfloor$ into the multiplicity, and the product is accumulated; to_string_u128 prints the total. Orders are 32-bit values, which is safe by $\operatorname{ord}_m(10)\le m-1\lt2^{32}$, while the accumulators are 128-bit.

The Java port repeats the same pipeline with a short[] sieve, an int[] order table, a hand-written binary search instead of std::upper_bound, and an IntStream pipeline whose mapToLong(...).sum() adds the contributions; its limit is fixed at $10^8$, it runs no checkpoints of its own, and the two large arrays alone need some $300$ MB of heap. A signed long is enough there, because $L(n)\le n-1$ gives $S(N)\lt N^2/2=5\cdot10^{15}$, far below $2^{63}-1$; its modPow multiplies residues below $10^8$, so products stay under $10^{16}$, while the C++ mod_pow casts to 128 bits for the same product. The Python file is not a third implementation but a bridge, and it expects the published layout with the C++ source in a sibling solutionsCpp directory: it compiles that source with clang++ or g++ at -O2 into a cached binary next to it, recompiles only when the source is newer, runs it with --skip-checkpoints because should_skip_cpp_checkpoints finds that literal in the source, at the C++ default limit since it passes no other argument, and extracts the value with parse_output, which takes the last line unless some line gives the value after answer: or =.

Complexity and Verification

The sieve costs $O(N\log\log N)$ time ^[9] and stores two bytes per odd number, about $100$ MB; ord_prime holds four bytes per odd number, about $200$ MB, of which only the $5\,761\,453$ prime slots are used; and the prime list reserves $8\,333\,333$ entries, $33$ MB, for those same primes. The measured peak is about $310$ MiB. The descent runs once per prime, with at most one test per prime factor of $p-1$ with multiplicity, hence $O(\log N)$ tests of $O(\log N)$ multiplications; the measured average of four is an observation, not a worst case. Each of the $39\,999\,999$ cores costs $O(\log N)$ divisions, as many greatest common divisors and at most eight comparisons, so a conservative overall bound is $O(N\log^2N)$ word operations with $O(N)$ memory. Compiled with clang++ -O2 and timed on one machine along a single execution path at $N=10^8$, checkpoints skipped, the sieve takes $0.21$ s, the prime orders $5.0$ s and the sweep $1.8$ s, about $7$ s of processor time in all; the prime-power tables and the smooth list need under a millisecond.

Beyond the three packaged checkpoints the following was verified for this article. Long division and $\operatorname{ord}_m(10)$ agree on $L(n)$ for all $3\le n\le20000$, and the prefix length equals $\max(a,b)$ for $n\le5000$. The identity of Section 4 matches brute force at $N=10,100,1000,2000,5000$: $9$, $1177$, $94288$, $349125$, $2016039$. The decomposition of Section 7 holds for every $m\le3000$ coprime to $10$; the only primes below $10^4$ with $10^{p-1}\equiv1\pmod{p^2}$ are $3$ and $487$, and $56\,598\,313$ passes that test too. A separate Python implementation of the whole pipeline returns $55535191115$ at $N=10^6$ and agrees with the published C++ at $N=10^7$, where both give $4937464337289$. The published C++ passes its three checkpoints, and the figures above come from running it to the full bound.

Footnotes and References

The references below supply the background facts used above; the criterion for a repeating block, the summation by cores, the order descent and the lifting rule are proved explicitly in this article.

Project Euler 417 — Reciprocal Cycles II. The official statement lists the decimal expansions of 1/2 to 1/10, defines L(n) as the length of the recurring cycle of 1/n with L(n) = 0 when the denominator has no prime factor other than 2 and 5, gives the sum for 3 ≤ n ≤ 1000000 as 55535191115 and asks for the same sum up to 100000000.
Wikipedia — Repeating decimal. The article states that the length of the repetend of 1/p equals the order of 10 modulo p, that a fraction in lowest terms repeats exactly when its denominator has a prime factor other than 2 or 5, that exponents a and b of 2 and 5 produce an initial transient of max(a, b) digits, and that a terminating decimal has a second representation whose repetend is 9. Only the cycle length carries over from the coprime part, not the digit block: 1/28 = 0.03(571428) against 1/7 = 0.(142857).
Wikipedia — Multiplicative order. It defines ord_n(a), for a coprime to n, as the least positive k with a^k congruent to 1 modulo n, and notes that it divides φ(n) by Lagrange's theorem and even divides the Carmichael function λ(n).
Wikipedia — Fermat's little theorem. For a prime p and an integer a not divisible by p, a^(p−1) is congruent to 1 modulo p. Section 8 uses this to know in advance that the order of 10 divides p − 1.
Wikipedia — Chinese remainder theorem. For pairwise coprime moduli with product N, reduction modulo the factors is a ring isomorphism from the integers modulo N onto the product of the integers modulo each factor. Section 7 uses it to test 10^k ≡ 1 one prime power at a time.
OEIS A002329 — periods of reciprocals of integers prime to 10. Indexed by the odd numbers not divisible by 5, its terms begin 1, 1, 6, 1, 2, 6, 16, 18, 6, 22. The first belongs to the denominator 1 and follows the order convention, not this problem's L(1) = 0; the others, for 3, 7, 9, 11, 13, 17, 19, 21 and 23, agree with the orders computed here.
OEIS A045616 — primes p with 10^(p−1) congruent to 1 modulo p². The listed terms are 3, 487 and 56598313, with no further term below 1.172 × 10^14 as of 2020. These are exactly the primes for which the lift from p to p² keeps the order instead of multiplying it by p; later lifts always multiply.
C++ — Euler417.cpp. The linked immutable revision contains build_odd_spf, multiplicative_order_prime, build_small_prime_power_orders, solve_fast and the checkpoints described above, including the comparison against 55535191115. It is the implementation record for this article.
Wikipedia — Sieve of Eratosthenes. The article gives O(n log log n) operations for sieving all primes below n, O(n) memory, and the odd-only refinement that marks only odd multiples, the variant used here.

Problem 417 source code

C++

#include <algorithm>
#include <cmath>
#include <cstdint>
#include <iostream>
#include <numeric>
#include <string>
#include <thread>
#include <vector>

namespace {

using u32 = std::uint32_t;
using u64 = std::uint64_t;
using u128 = __uint128_t;
using i64 = std::int64_t;

struct Options {
    int limit = 100000000;
    int threads = 0;
    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;
    }
    try {
        value = std::stoi(tail);
    } catch (...) {
        return false;
    }
    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, "--limit=", options.limit) ||
            parse_int_after_prefix(arg, "--threads=", options.threads)) {
            continue;
        }
        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.limit >= 3 && options.threads >= 0;
}

std::string to_string_u128(u128 value) {
    if (value == 0) {
        return "0";
    }
    std::string out;
    while (value > 0) {
        const int digit = static_cast<int>(value % 10);
        out.push_back(static_cast<char>('0' + digit));
        value /= 10;
    }
    std::reverse(out.begin(), out.end());
    return out;
}

u64 mod_pow(u64 base, u64 exp, u64 mod) {
    u64 result = 1 % mod;
    u64 cur = base % mod;
    u64 e = exp;
    while (e > 0) {
        if (e & 1ULL) {
            result = static_cast<u64>((static_cast<u128>(result) * cur) % mod);
        }
        cur = static_cast<u64>((static_cast<u128>(cur) * cur) % mod);
        e >>= 1ULL;
    }
    return result;
}

int choose_thread_count(int requested, std::size_t work_items) {
    if (work_items <= 1) {
        return 1;
    }
    int threads = requested;
    if (threads <= 0) {
        threads = static_cast<int>(std::thread::hardware_concurrency());
    }
    if (threads <= 0) {
        threads = 4;
    }
    if (threads > static_cast<int>(work_items)) {
        threads = static_cast<int>(work_items);
    }
    if (threads < 1) {
        threads = 1;
    }
    return threads;
}

std::vector<std::uint16_t> build_odd_spf(const int limit) {
    std::vector<std::uint16_t> spf(static_cast<std::size_t>(limit / 2 + 1), 0);
    const int root = static_cast<int>(std::sqrt(static_cast<double>(limit)));

    for (int i = 3; i <= root; i += 2) {
        if (spf[static_cast<std::size_t>(i >> 1)] != 0) {
            continue;
        }
        const i64 step = 2LL * i;
        for (i64 j = 1LL * i * i; j <= limit; j += step) {
            std::uint16_t& cell = spf[static_cast<std::size_t>(j >> 1)];
            if (cell == 0) {
                cell = static_cast<std::uint16_t>(i);
            }
        }
    }
    return spf;
}

std::vector<int> collect_primes(const int limit, const std::vector<std::uint16_t>& spf) {
    std::vector<int> primes;
    primes.reserve(static_cast<std::size_t>(limit / 12));
    for (int p = 3; p <= limit; p += 2) {
        if (p != 5 && spf[static_cast<std::size_t>(p >> 1)] == 0) {
            primes.push_back(p);
        }
    }
    return primes;
}

u32 multiplicative_order_prime(const int p, const std::vector<std::uint16_t>& spf) {
    u32 ord = static_cast<u32>(p - 1);
    u32 rem = ord;

    if ((rem & 1U) == 0U) {
        while ((rem & 1U) == 0U) {
            rem >>= 1U;
        }
        while ((ord & 1U) == 0U) {
            const u32 cand = ord >> 1U;
            if (mod_pow(10, cand, static_cast<u64>(p)) == 1ULL) {
                ord = cand;
            } else {
                break;
            }
        }
    }

    while (rem > 1U) {
        u32 q = static_cast<u32>(spf[static_cast<std::size_t>(rem >> 1U)]);
        if (q == 0U) {
            q = rem;
        }
        while (rem % q == 0U) {
            rem /= q;
        }
        while (ord % q == 0U) {
            const u32 cand = ord / q;
            if (mod_pow(10, cand, static_cast<u64>(p)) == 1ULL) {
                ord = cand;
            } else {
                break;
            }
        }
    }
    return ord;
}

std::vector<u32> compute_prime_orders_parallel(const std::vector<int>& primes,
                                               const std::vector<std::uint16_t>& spf,
                                               const int limit,
                                               const int requested_threads) {
    std::vector<u32> ord_prime(static_cast<std::size_t>(limit / 2 + 1), 0U);
    ord_prime[0] = 1U;  // order(10 mod 1)

    const int thread_count = choose_thread_count(requested_threads, primes.size());
    std::vector<std::thread> workers;
    workers.reserve(static_cast<std::size_t>(thread_count));

    for (int t = 0; t < thread_count; ++t) {
        const std::size_t begin = primes.size() * static_cast<std::size_t>(t) /
                                  static_cast<std::size_t>(thread_count);
        const std::size_t end = primes.size() * static_cast<std::size_t>(t + 1) /
                                static_cast<std::size_t>(thread_count);
        workers.emplace_back([&, begin, end]() {
            for (std::size_t i = begin; i < end; ++i) {
                const int p = primes[i];
                ord_prime[static_cast<std::size_t>(p >> 1)] = multiplicative_order_prime(p, spf);
            }
        });
    }
    for (auto& worker : workers) {
        worker.join();
    }
    return ord_prime;
}

std::vector<std::vector<u32>> build_small_prime_power_orders(
    const int limit, const std::vector<int>& primes, const std::vector<u32>& ord_prime) {
    const int root = static_cast<int>(std::sqrt(static_cast<double>(limit)));
    std::vector<std::vector<u32>> prime_power_orders(static_cast<std::size_t>(root + 1));

    for (const int p : primes) {
        if (p > root) {
            break;
        }
        auto& table = prime_power_orders[static_cast<std::size_t>(p)];
        int emax = 1;
        u64 pp = static_cast<u64>(p);
        while (pp <= static_cast<u64>(limit / p)) {
            pp *= static_cast<u64>(p);
            ++emax;
        }

        table.assign(static_cast<std::size_t>(emax + 1), 0U);
        table[1] = ord_prime[static_cast<std::size_t>(p >> 1)];

        pp = static_cast<u64>(p);
        for (int e = 2; e <= emax; ++e) {
            pp *= static_cast<u64>(p);
            const u32 prev = table[static_cast<std::size_t>(e - 1)];
            u32 cur = prev;
            if (mod_pow(10, prev, pp) != 1ULL) {
                cur = static_cast<u32>(prev * static_cast<u32>(p));
            }
            table[static_cast<std::size_t>(e)] = cur;
        }
    }
    return prime_power_orders;
}

std::vector<int> build_smooth_numbers(const int limit) {
    std::vector<int> smooth;
    for (u64 p2 = 1ULL; p2 <= static_cast<u64>(limit); p2 *= 2ULL) {
        for (u64 p5 = 1ULL; p2 * p5 <= static_cast<u64>(limit); p5 *= 5ULL) {
            smooth.push_back(static_cast<int>(p2 * p5));
        }
    }
    std::sort(smooth.begin(), smooth.end());
    smooth.erase(std::unique(smooth.begin(), smooth.end()), smooth.end());
    return smooth;
}

u32 multiplicative_order_from_factorization(
    int n, const std::vector<std::uint16_t>& spf, const std::vector<u32>& ord_prime,
    const std::vector<std::vector<u32>>& prime_power_orders) {
    u32 ord = 1U;
    int x = n;
    while (x > 1) {
        int p = static_cast<int>(spf[static_cast<std::size_t>(x >> 1)]);
        if (p == 0) {
            p = x;
        }
        int e = 1;
        int pp = p;
        x /= p;
        while (x % p == 0) {
            x /= p;
            pp *= p;
            ++e;
        }

        u32 part = 0U;
        if (e == 1) {
            part = ord_prime[static_cast<std::size_t>(p >> 1)];
        } else if (p < static_cast<int>(prime_power_orders.size()) &&
                   e < static_cast<int>(prime_power_orders[static_cast<std::size_t>(p)].size())) {
            part = prime_power_orders[static_cast<std::size_t>(p)][static_cast<std::size_t>(e)];
        } else {
            part = ord_prime[static_cast<std::size_t>(p >> 1)];
            u64 mod = static_cast<u64>(p);
            for (int k = 2; k <= e; ++k) {
                mod *= static_cast<u64>(p);
                if (mod_pow(10, part, mod) != 1ULL) {
                    part = static_cast<u32>(part * static_cast<u32>(p));
                }
            }
        }

        const u32 g = std::gcd(ord, part);
        ord = static_cast<u32>((static_cast<u64>(ord) / g) * part);
    }
    return ord;
}

u128 solve_fast(const int limit, const int requested_threads) {
    const std::vector<std::uint16_t> spf = build_odd_spf(limit);
    const std::vector<int> primes = collect_primes(limit, spf);
    const std::vector<u32> ord_prime =
        compute_prime_orders_parallel(primes, spf, limit, requested_threads);
    const std::vector<std::vector<u32>> prime_power_orders =
        build_small_prime_power_orders(limit, primes, ord_prime);
    const std::vector<int> smooth = build_smooth_numbers(limit);

    const std::size_t odd_count =
        static_cast<std::size_t>((static_cast<u64>(limit) - 3ULL) / 2ULL + 1ULL);
    const int thread_count = choose_thread_count(requested_threads, odd_count);
    std::vector<u128> partial(static_cast<std::size_t>(thread_count), 0);
    std::vector<std::thread> workers;
    workers.reserve(static_cast<std::size_t>(thread_count));

    for (int t = 0; t < thread_count; ++t) {
        const std::size_t begin = odd_count * static_cast<std::size_t>(t) /
                                  static_cast<std::size_t>(thread_count);
        const std::size_t end = odd_count * static_cast<std::size_t>(t + 1) /
                                static_cast<std::size_t>(thread_count);
        workers.emplace_back([&, begin, end, t]() {
            u128 local = 0;
            for (std::size_t idx = begin; idx < end; ++idx) {
                const int n = static_cast<int>(3 + 2 * idx);
                if (n % 5 == 0) {
                    continue;
                }
                const u32 ord =
                    multiplicative_order_from_factorization(n, spf, ord_prime, prime_power_orders);
                const int quotient = limit / n;
                const int multiplicity =
                    static_cast<int>(std::upper_bound(smooth.begin(), smooth.end(), quotient) -
                                     smooth.begin());
                local += static_cast<u128>(ord) * static_cast<u128>(multiplicity);
            }
            partial[static_cast<std::size_t>(t)] = local;
        });
    }
    for (auto& worker : workers) {
        worker.join();
    }

    u128 total = 0;
    for (const u128 value : partial) {
        total += value;
    }
    return total;
}

u64 recurring_cycle_length_bruteforce(int n) {
    while ((n % 2) == 0) {
        n /= 2;
    }
    while ((n % 5) == 0) {
        n /= 5;
    }
    if (n == 1) {
        return 0;
    }

    int r = 10 % n;
    u64 len = 1;
    while (r != 1) {
        r = static_cast<int>((static_cast<i64>(r) * 10LL) % n);
        ++len;
    }
    return len;
}

u128 solve_bruteforce(const int limit) {
    u128 total = 0;
    for (int n = 3; n <= limit; ++n) {
        total += recurring_cycle_length_bruteforce(n);
    }
    return total;
}

bool run_checkpoints(const int requested_threads) {
    if (solve_bruteforce(10) != 9) {
        std::cerr << "Checkpoint failed: sum_{n=3..10} L(n)\n";
        return false;
    }

    if (solve_fast(1000, requested_threads) != solve_bruteforce(1000)) {
        std::cerr << "Checkpoint failed: fast solver vs brute force at limit=1000\n";
        return false;
    }

    if (solve_fast(1000000, requested_threads) != static_cast<u128>(55535191115ULL)) {
        std::cerr << "Checkpoint failed: sum_{n=3..1000000} L(n)\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(options.threads)) {
        return 2;
    }

    const u128 answer = solve_fast(options.limit, options.threads);
    std::cout << to_string_u128(answer) << '\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.ArrayList;
import java.util.Arrays;
import java.util.List;
import java.util.stream.IntStream;

public class Euler417 {
    static long modPow(long base, long exp, long mod) {
        long result = 1 % mod;
        long cur = base % mod;
        long e = exp;
        while (e > 0) {
            if ((e & 1) != 0) {
                result = (result * cur) % mod;
            }
            cur = (cur * cur) % mod;
            e >>= 1;
        }
        return result;
    }

    static short[] buildOddSpf(int limit) {
        short[] spf = new short[limit / 2 + 1];
        int root = (int) Math.sqrt(limit);

        for (int i = 3; i <= root; i += 2) {
            if (spf[i >> 1] != 0)
                continue;
            int step = 2 * i;
            for (long j = (long) i * i; j <= limit; j += step) {
                if (spf[(int) (j >> 1)] == 0) {
                    spf[(int) (j >> 1)] = (short) i;
                }
            }
        }
        return spf;
    }

    static int[] collectPrimes(int limit, short[] spf) {
        int[] primes = new int[limit / 12 + 1000];
        int count = 0;
        for (int p = 3; p <= limit; p += 2) {
            if (p != 5 && spf[p >> 1] == 0) {
                if (count == primes.length) {
                    primes = Arrays.copyOf(primes, primes.length * 2);
                }
                primes[count++] = p;
            }
        }
        return Arrays.copyOf(primes, count);
    }

    static int multiplicativeOrderPrime(int p, short[] spf) {
        int ord = p - 1;
        int rem = ord;

        if ((rem & 1) == 0) {
            while ((rem & 1) == 0)
                rem >>= 1;
            while ((ord & 1) == 0) {
                int cand = ord >> 1;
                if (modPow(10, cand, p) == 1) {
                    ord = cand;
                } else {
                    break;
                }
            }
        }

        while (rem > 1) {
            int q = spf[rem >> 1];
            if (q == 0)
                q = rem;
            while (rem % q == 0)
                rem /= q;
            while (ord % q == 0) {
                int cand = ord / q;
                if (modPow(10, cand, p) == 1) {
                    ord = cand;
                } else {
                    break;
                }
            }
        }
        return ord;
    }

    static int[][] buildSmallPrimePowerOrders(int limit, int[] primes, int[] ordPrime) {
        int root = (int) Math.sqrt(limit);
        int[][] primePowerOrders = new int[root + 1][];

        for (int p : primes) {
            if (p > root)
                break;
            int emax = 1;
            long pp = p;
            while (pp <= limit / p) {
                pp *= p;
                emax++;
            }

            int[] table = new int[emax + 1];
            table[1] = ordPrime[p >> 1];

            pp = p;
            for (int e = 2; e <= emax; e++) {
                pp *= p;
                int prev = table[e - 1];
                int cur = prev;
                if (modPow(10, prev, pp) != 1) {
                    cur = prev * p;
                }
                table[e] = cur;
            }
            primePowerOrders[p] = table;
        }
        return primePowerOrders;
    }

    static int[] buildSmoothNumbers(int limit) {
        int[] temp = new int[2000];
        int count = 0;
        for (long p2 = 1; p2 <= limit; p2 *= 2) {
            for (long p5 = 1; p2 * p5 <= limit; p5 *= 5) {
                temp[count++] = (int) (p2 * p5);
            }
        }
        int[] res = Arrays.copyOf(temp, count);
        Arrays.sort(res);
        return res;
    }

    static long gcd(long a, long b) {
        while (b != 0) {
            long t = b;
            b = a % b;
            a = t;
        }
        return a;
    }

    static int multiplicativeOrderFromFactorization(int n, short[] spf, int[] ordPrime, int[][] primePowerOrders) {
        int ord = 1;
        int x = n;
        while (x > 1) {
            int p = spf[x >> 1];
            if (p == 0)
                p = x;
            int e = 1;
            int pp = p;
            x /= p;
            while (x % p == 0) {
                x /= p;
                pp *= p;
                e++;
            }

            int part = 0;
            if (e == 1) {
                part = ordPrime[p >> 1];
            } else if (p < primePowerOrders.length && primePowerOrders[p] != null && e < primePowerOrders[p].length) {
                part = primePowerOrders[p][e];
            } else {
                part = ordPrime[p >> 1];
                long mod = p;
                for (int k = 2; k <= e; k++) {
                    mod *= p;
                    if (modPow(10, part, mod) != 1) {
                        part *= p;
                    }
                }
            }

            long g = gcd(ord, part);
            ord = (int) (((long) ord / g) * part);
        }
        return ord;
    }

    static String solve() {
        int limit = 100000000;
        short[] spf = buildOddSpf(limit);
        int[] primes = collectPrimes(limit, spf);

        int[] ordPrime = new int[limit / 2 + 1];
        ordPrime[0] = 1;

        IntStream.range(0, primes.length).parallel().forEach(i -> {
            int p = primes[i];
            ordPrime[p >> 1] = multiplicativeOrderPrime(p, spf);
        });

        int[][] primePowerOrders = buildSmallPrimePowerOrders(limit, primes, ordPrime);
        int[] smooth = buildSmoothNumbers(limit);

        int oddCount = (limit - 3) / 2 + 1;

        long total = IntStream.range(0, oddCount).parallel().mapToLong(idx -> {
            int n = 3 + 2 * idx;
            if (n % 5 == 0)
                return 0;

            int ord = multiplicativeOrderFromFactorization(n, spf, ordPrime, primePowerOrders);
            int quotient = limit / n;

            int low = 0, high = smooth.length;
            while (low < high) {
                int mid = (low + high) >>> 1;
                if (smooth[mid] <= quotient) {
                    low = mid + 1;
                } else {
                    high = mid;
                }
            }
            int multiplicity = low;

            return (long) ord * multiplicity;
        }).sum();

        return Long.toString(total);
    }

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