Problem 133: Repunit Nonfactors

Project Euler Problem 133 Solution

EulerSolve provides an optimized solution for Project Euler Problem 133, Repunit Nonfactors, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Let $R(k)=\frac{10^k-1}{9}=11\ldots1$ be the repunit with $k$ digits. The problem asks for the sum of all primes $p<10^5$ such that there is no exponent $n\ge 1$ with $p\mid R(10^n)$. The key point is that one never has to build the enormous numbers $R(10^n)$ themselves. For each prime, the question can be decided entirely from the multiplicative behavior of $10$ modulo that prime. Mathematical Approach The implementations reduce the whole problem to a yes-or-no test on the multiplicative order of $10$ modulo $p$. Once that order is known, the classification is immediate. From repunits to a congruence For every integer $k\ge 1$, $$R(k)=1+10+10^2+\cdots+10^{k-1}=\frac{10^k-1}{9}.$$ If $p\notin\{2,3,5\}$, then $9$ is invertible modulo $p$, so multiplying by $9$ does not change divisibility by $p$. Therefore $$p\mid R(k)\iff 10^k\equiv 1 \pmod p.$$ In this problem the length is not an arbitrary $k$, but specifically $k=10^n$. So for primes other than $2$, $3$, and $5$, we are asking whether $$10^{10^n}\equiv 1 \pmod p$$ holds for at least one $n\ge 1$. Orders built only from $2$ and $5$ Let $$t=\operatorname{ord}_p(10),$$ the smallest positive integer such that $10^t\equiv 1\pmod p$. By definition, $10^m\equiv 1\pmod p$ happens exactly when $t\mid m$....

Detailed mathematical approach

Problem Summary

A repunit is a number written with the digit $1$ only, and $R(k)$ denotes the repunit with $k$ digits; for example, $R(6)=111111$. The problem looks at the repunits whose length is a power of ten: $R(10)$, $R(100)$, $R(1000)$ and so on (the length $10^0=1$ adds nothing, since $R(1)=1$ has no prime factor). The statement notes that $17$ divides none of $R(10)$, $R(100)$ and $R(1000)$ but does divide $R(10000)$, that $19$ divides no $R(10^n)$ at all, and that $11$, $17$, $41$ and $73$ are the only four primes below $100$ that divide some $R(10^n)$. We must find the sum of all primes below $100\,000$ that divide no $R(10^n)$. ^[1]

The whole solution in one sentence: the primes $2$, $3$ and $5$ never divide $R(10^n)$, and any other prime $p$ divides $R(k)$ exactly when $k$ is a multiple of the order of $10$ modulo $p$, the smallest $t\ge1$ for which $10^t$ leaves the remainder $1$ on division by $p$, so $p$ divides some $R(10^n)$ exactly when this order, which the program finds from the prime factors of $p-1$, has no prime factors other than $2$ and $5$.

Mathematical Approach

The numbers $R(10^n)$ have $10^n$ digits, far too many to write down. We never build them: for each prime we follow only remainders.

1. Repunits and their factors

The $k$ digits of $R(k)$ stand for the place values $1,10,\dots,10^{k-1}$, so ^[2]

$$\begin{aligned}R(k)&=1+10+\dots+10^{k-1}\\&=\frac{10^k-1}{9}.\end{aligned}$$

The first repunits factor as $R(2)=11$, $R(3)=3\cdot37$, $R(4)=11\cdot101$, $R(5)=41\cdot271$, $R(6)=3\cdot7\cdot11\cdot13\cdot37$, $R(8)=11\cdot73\cdot101\cdot137$ and $R(10)=11\cdot41\cdot271\cdot9091$, and the prime $17$ first appears in $R(16)$. ^[2] Writing $a$ ones followed by $b$ ones shows that

$$R(a+b)=R(a)\cdot10^{b}+R(b).$$

Applied again and again, this shows that $R(a)$ divides $R(b)$ whenever $a$ divides $b$; for example $R(4)=R(2)\cdot100+R(2)=11\cdot101$.

2. The shortest repunit that a prime divides

Appending a digit $1$ gives $R(k+1)=10\,R(k)+1$, so the remainders of $R(1),R(2),\dots$ on division by a prime $p$ follow one from another. For $p=7$ they are $1,4,6,5,2,0$ for $k=1,\dots,6$, and then the same six repeat. We write $x\equiv y\pmod m$ when $x$ and $y$ leave the same remainder on division by $m$, that is, when $m$ divides $x-y$. Let $s$ be the smallest $k\ge1$ with $p\mid R(k)$, if there is one. For $p=7$ it is $6$, and the same walk gives $s=16$ for $p=17$ and $s=18$ for $p=19$.

Claim: $p$ divides $R(k)$ if and only if $s$ divides $k$. Indeed, the identity of Section 1 and $p\mid R(s)$ give

$$\begin{aligned}R(k+s)&=R(s)\cdot10^{k}+R(k)\\&\equiv R(k)\pmod p.\end{aligned}$$

So the remainders repeat with period $s$, and $R(k)$ leaves the same remainder as $R(r)$, where $r$ is the remainder of $k$ on division by $s$. For $r=0$ this is $R(0)=0$, the empty repunit; for $0\lt r\lt s$ it is not $0$, by the choice of $s$. The primes $2$ and $5$ divide no repunit at all, since every repunit ends in the digit $1$. Every other prime divides some repunit: $3\mid R(3)=111$, and Section 5 shows it for the remaining primes.

3. Lengths that are powers of ten

By Section 2, $p$ divides $R(10^n)$ exactly when $s$ divides $10^n=2^n5^n$. A divisor of $2^n5^n$ has no prime factor other than $2$ and $5$, by Euclid's lemma, ^[3] and contains each of them at most $n$ times. Hence $p$ divides some $R(10^n)$ if and only if $s=2^a5^b$ with integers $a,b\ge0$, and then it divides exactly the repunits $R(10^n)$ with $n\ge\max(a,b)$.

The prime $7$ divides $R(6)$, but $s=6=2\cdot3$ contains the factor $3$, which no power of ten has, so $7$ divides no $R(10^n)$. The same happens for $3$ itself, with $s=3$. The answer is therefore the sum of $2$, $5$ and all primes below $100\,000$ whose length $s$ has a prime factor other than $2$ and $5$.

4. The two examples of the statement

For $p=17$ we have $s=16$. By Section 2, $R(10^n)$ leaves the same remainder modulo $17$ as $R(r)$, where $r$ is the remainder of $10^n$ on division by $16$. For $n=1,2,3,4$ these are $r=10,4,8,0$, and $R(10)$, $R(4)$ and $R(8)$ leave the remainders $2$, $6$ and $13$. So $R(10)$, $R(100)$ and $R(1000)$ are not multiples of $17$, while $R(10000)$ is, because $16=2^4$ divides $10^4$.

For $p=19$ we have $s=18=2\cdot3^2$. Every power $10^n$ with $n\ge1$ leaves the remainder $10$ on division by $18$, since $10\cdot10=100=5\cdot18+10$. So every $R(10^n)$ leaves the same remainder modulo $19$ as $R(10)$, namely $3$, and none of them is divisible by $19$.

5. Powers of ten and the order of 10

For a prime $p$ other than $2$, $3$ and $5$, the identity $9R(k)=10^k-1$ turns repunits into powers of ten. Since $p$ does not divide $9$, Euclid's lemma shows that $p\mid R(k)$ exactly when $p\mid10^k-1$, that is, when $10^k\equiv1\pmod p$. So $s$ equals the smallest exponent $t\ge1$ with $10^t\equiv1\pmod p$, called the multiplicative order of $10$ modulo $p$ and written $\operatorname{ord}_p(10)$. ^[4] By Section 2, $10^k\equiv1\pmod p$ holds exactly when $t$ divides $k$. Modulo $17$, for instance, the powers $10^1,\dots,10^{16}$ leave the remainders $10,15,14,4,6,9,5,16$, $7,2,3,13,11,8,12,1$, and modulo $73$ the powers $10^1,\dots,10^8$ leave $10,27,51,72,63,46,22,1$; so $\operatorname{ord}_{17}(10)=16$ and $\operatorname{ord}_{73}(10)=8$.

Fermat's little theorem says that $10^{p-1}\equiv1\pmod p$ for every prime $p$ other than $2$ and $5$, ^[5] so $t$ divides $p-1$: $16$ divides $16$, and $8$ divides $72$. It also shows that every prime other than $2$, $3$ and $5$ divides some repunit, namely $R(p-1)$. The order is the length of the block that repeats in the decimal expansion of $1/p$ as well: $1/7=0.142857\,142857\ldots$ repeats six digits, and $1/17=0.0588235294117647\ldots$ repeats sixteen. ^[2] ^[6]

6. The primes 2, 3 and 5

For $2$ and $5$ every power $10^k$ is a multiple of the prime, so there is no order at all, and these primes divide no repunit (Section 2). For $3$, the factor $9$ in $9R(k)=10^k-1$ is the problem: $10\equiv1\pmod3$, so $\operatorname{ord}_3(10)=1$, which has the form $2^a5^b$, and yet $3$ divides no $R(10^n)$. The correct length is $s=3$: every power of ten leaves the remainder $1$ on division by $3$, so $R(k)\equiv1+1+\dots+1=k\pmod3$, and $R(10^n)\equiv10^n\equiv1\pmod3$. A program that tests the order must therefore treat $2$, $3$ and $5$ separately; all three belong to the sum.

7. Worked example: the primes below 100

The $25$ primes below $100$ add up to $1060$. Leaving out $2$, $3$ and $5$, the orders of $10$ modulo the other $22$ primes fall into three groups: ^[6]

of the form $2^a5^b$: $p=11,17,41,73$, with the orders $2,16,5,8$;
with the factor $3$: $p=7,13,19,31,37,43,61,67,97$, with the orders $6,6,18,15,3,21,60,33,96$;
with another prime factor above $5$, but without the factor $3$: $p=23,29,47,53,59,71,79,83,89$, with the orders $22,28,46,13,58,35,13,41,44$.

So $11$, $17$, $41$ and $73$ are the only primes below $100$ that divide some $R(10^n)$, as the statement says. By Section 3, $11$ and $41$ divide $R(10)$, and $73$ first divides $R(1000)$. The sum of the primes below $100$ that divide no $R(10^n)$ is $1060-(11+17+41+73)=918$, which the C++ and Python programs check first.

8. Finding the order from the prime factors of $p-1$

Trying $10^1,10^2,\dots$ one after another takes $t$ steps, up to $p-1$ for one prime. Since $t$ divides $p-1$, it is much faster to start from $e=p-1$ and remove prime factors: for each prime $q$ that divides $p-1$, replace $e$ by $e/q$ as long as $q$ divides $e$ and $10^{e/q}\equiv1\pmod p$.

Claim: when all prime factors of $p-1$ have been treated, $e=t$. Throughout, $10^e\equiv1\pmod p$, so $t$ divides $e$ by Section 5. Consider the moment when the work with a prime $q$ stops. If $q$ does not divide $e$, then it does not divide $t$ either. Otherwise $10^{e/q}\not\equiv1\pmod p$, so $t$ does not divide $e/q$; since $t$ divides $e$, this means that $t$ contains the factor $q$ as often as $e$ does. In both cases $e$ and $t$ contain $q$ equally often, and the later steps, which divide only by other primes, do not change that. At the end $e$ and $t$ contain every prime factor of $p-1$ equally often, and both divide $p-1$, so $e=t$.

For $p=73$ we have $p-1=72=2^3\cdot3^2$. The first test for $q=2$ already fails: $73\cdot137=10001$ gives $10^4\equiv-1\pmod{73}$, so $10^{36}=(10^4)^9\equiv-1$, which is $72$ and not $1$. For $q=3$, $10^{24}\equiv1$ gives $e=24$, then $10^8\equiv1$ gives $e=8$, and $3$ does not divide $8$; the order is $8$. Each test is one modular power, computed by binary exponentiation: ^[7] since $b^{2u}=(b^2)^u$ and $b^{2u+1}=b\,(b^2)^u$, the result is multiplied by the base whenever the exponent is odd, and then the base is squared and the exponent halved, rounding down. For $10^{13}$, with $13=8+4+1$, the result collects $10^1$, $10^4$ and $10^8$. Every binary digit of the exponent costs one squaring and at most one more multiplication, all modulo $p$.

9. How far the powers of ten have to go

If $t=2^a5^b$ divides $p-1\lt100\,000$, then $2^a\lt2^{17}=131072$, so $a\le16$; and since $p-1$ is even, $2\cdot5^b$ divides $p-1$, so $2\cdot5^b\lt100\,000$ and $b\le6$. Hence every prime below $100\,000$ that divides some $R(10^n)$ divides $R(10^{16})$. The exponent $16$ cannot be lowered: for the prime $65537=2^{16}+1$ the order divides $p-1=2^{16}$, and the single test of Section 8, $10^{32768}\equiv-1\pmod{65537}$, shows that it is $2^{16}$ itself, so this prime first divides $R(10^{16})$.

A second test needs neither orders nor exceptions: with $K=10^{16}$, a prime $p$ below $100\,000$ divides some $R(10^n)$ if and only if $10^K\equiv1\pmod{9p}$, because $9p$ divides $10^K-1=9R(K)$ exactly when $p$ divides $R(K)$; this holds for $2$, $3$ and $5$ too. A prime $p\gt5$ with $p-1=2^a5^b$, such as $17$, $41$, $257$ or $65537$, always qualifies, because its order divides $p-1$. In all, only $40$ of the $9592$ primes below $100\,000$ divide some $R(10^n)$: the terms of OEIS A178070 in this range. ^[8]

How the Code Works

The C++ program starts in main, which calls parse_arguments: --limit=N replaces the bound $100\,000$ by another bound $N\ge2$, read into an int without an overflow check, --skip-checkpoints skips the checks, and any other argument stops the program with an error. run_checkpoints then requires solve(100) to return $918$, the sum of Section 7, and can_divide_some_R10n to answer yes for $17$ and no for $19$ and $3$. Finally it prints the value of solve for the bound.

solve takes the primes from primes_below, a sieve of Eratosthenes on a std::vector<bool>: for every $i$ with $i^2$ below the bound that is still marked as prime, it crosses out $i^2$, $i^2+i$ and so on; the smaller multiples of $i$ are already crossed out. ^[9] For each prime, can_divide_some_R10n answers no at once for $2$, $3$ and $5$ (Section 6). For any other prime it takes the order from multiplicative_order_10, divides it by $2$ and by $5$ as long as possible and answers yes if $1$ is left (Section 3). solve adds up the primes with the answer no in a 64-bit integer. Without this exception the order routine would return $1$ for $2$ and $3$ and $4$ for $5$, all of the form $2^a5^b$, and all three would wrongly be left out.

multiplicative_order_10 is the method of Section 8. It finds the distinct prime factors of $p-1$ by trial division, dividing each factor out as soon as it is found, so that the loop runs only while $f^2$ does not exceed the part that is left, and a leftover part greater than $1$ is itself prime. pow_mod is the binary exponentiation of Section 8, and mul_mod forms each product in unsigned 128-bit arithmetic, a safety margin, since here both factors are below $2^{17}$.

The Python version has the same functions, apart from parse_arguments, and the same two options; it works with unbounded integers and a list of booleans as its sieve, and it ignores a --limit value that is not a number instead of rejecting it. The Java version does everything in main, with the fixed bound $100\,000$ and without options or checks: it sieves, adds $2$, $3$ and $5$ directly, factors $p-1$ in the same way, computes the powers with BigInteger.modPow and adds $p$ to a long sum when the order keeps a prime factor other than $2$ and $5$.

Complexity and Verification

Let $L$ be the bound. The sieve takes $O(L\log\log L)$ time and one flag per number; ^[9] for $L=100\,000$ it crosses out $193\,076$ times. For each of the $9589$ primes other than $2$, $3$ and $5$, trial division tries at most the divisors $f=2,\dots,316$, since $f^2\le p-1\lt317^2$, and only $395\,930$ divisors in total, because the unfactored part of $p-1$ shrinks. The main computation then makes $35\,497$ power tests, of which $9143$ succeed; in C++ and Python each takes at most $32$ multiplications, since every exponent is at most $(p-1)/2\lt2^{16}$, and $598\,143$ in all. In general the tests add little: each success at least halves $e$ and each prime factor of $p-1$ ends with at most one failure, so a prime needs $O(\log p)$ tests of $O(\log p)$ multiplications. Trial division, with $O(\sqrt L)$ steps for each of the $O(L/\log L)$ primes, bounds the time by $O(L^{3/2}/\log L)$, and the memory is $O(L)$. On a desktop computer all three programs finish in well under a second.

Before the main run, the C++ program ^[10] and the Python version check the sum $918$ for the bound $100$, the behaviour of $17$ and $19$ from the statement and that of $3$ (Section 6); the Java version has no checks. For this article all three programs were run, and they print the same number. Independent computations, which are not part of the published solutions, confirm the method:

a direct walk $R(k+1)=10R(k)+1$ modulo $p$ finds the length $s$ of Section 2 for every prime below $100\,000$ other than $2$ and $5$; it equals the order from Section 8 for every prime except $3$, where it is $3$, and classifying the primes by $s$ gives the same sum as the programs;
the remainders of $R(10^n)$ modulo every prime below $100\,000$, computed for $n\le40$ from $R(10^{n+1})$${}=R(10^n)\,(1+x+\dots+x^9)$ with $x=10^{10^n}$, without orders and without dividing by $9$, give the same classification, and nothing changes after $n=16$;
the $40$ primes found are exactly the terms of A178070 below $100\,000$, ^[8] the test of Section 9 modulo $9p$ agrees for every prime, and the orders of the primes below $100$ agree with A002371. ^[6]

Footnotes and References

The references below supply the background facts used above. The periodicity of repunit remainders, the power-of-ten criterion, the exceptions $2$, $3$ and $5$, the order computation and the bound $n\le16$ are derived explicitly in this article.

Project Euler 133 — Repunit Nonfactors. The official statement defines the repunits R(k), gives the examples 17 and 19, names 11, 17, 41 and 73 as the only primes below 100 that can divide some R(10ⁿ), and asks for the sum of all primes below one hundred thousand that will never be a factor of R(10ⁿ).
Wikipedia — Repunit. The article defines Rₙ = (10ⁿ − 1)/9, states that Rₙ is divisible by Rₐ whenever n is divisible by a, lists the factorisations of the first repunits and notes that for a prime p greater than 5 the length of the smallest repunit divisible by p is the period of the decimal expansion of 1/p.
Wikipedia — Euclid's lemma. If a prime divides the product of two integers, it divides at least one of them; more generally, if n divides ab and is coprime to a, then n divides b.
Wikipedia — Multiplicative order. For a coprime to n, the multiplicative order of a modulo n is the smallest positive integer k with aᵏ ≡ 1 (mod n); it always divides φ(n), which equals p − 1 for a prime p.
Wikipedia — Fermat's little theorem. For a prime p and an integer a that is not divisible by p, a^(p−1) ≡ 1 (mod p).
OEIS A002371 — Period of decimal expansion of 1/(n-th prime). The entry lists the period of 1/p for every prime p, with 0 by convention for 2 and 5, and notes that it is the smallest x with 10ˣ ≡ 1 (mod p), the order used in Section 7.
Wikipedia — Modular exponentiation. The right-to-left binary method squares the base for every binary digit of the exponent and multiplies it into the result when the digit is 1, so it needs O(log e) steps for an exponent e.
OEIS A178070 — Primes dividing repunits R(10^n) for some n. The entry states that a prime p > 5 is a term if the multiplicative order of 10 modulo p has the form 2ⁱ·5ʲ, explains that 3 is not a term because R(10ⁿ) ≡ 1 (mod 3), and notes that the sequence contains all primes greater than 5 of the form 2ʳ·5ˢ + 1.
Wikipedia — Sieve of Eratosthenes. Marking the multiples of a prime p can start at p², because the smaller multiples have already been marked; the algorithm takes O(n log log n) time and O(n) memory.
C++ — Euler133.cpp. The linked immutable C++ revision contains multiplicative_order_10, can_divide_some_R10n, solve and the checkpoints that run before the main computation. It is the implementation record for this article.

Problem 133 source code

C++

#include <cstdint>
#include <iostream>
#include <numeric>
#include <string>
#include <vector>

namespace {

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

struct Options {
    int limit = 100000;
    bool run_checkpoints = true;
};

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

    int parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<int>(c - '0');
    }
    value = parsed;
    return true;
}

bool parse_arguments(int argc, char** argv, Options& options) {
    for (int i = 1; i < argc; ++i) {
        const std::string arg(argv[i]);
        if (arg == "--skip-checkpoints") {
            options.run_checkpoints = false;
            continue;
        }
        if (parse_int_after_prefix(arg, "--limit=", options.limit)) {
            continue;
        }

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

    return options.limit >= 2;
}

u64 mul_mod(u64 a, u64 b, u64 mod) {
    return static_cast<u64>((static_cast<u128>(a) * static_cast<u128>(b)) % mod);
}

u64 pow_mod(u64 base, u64 exp, u64 mod) {
    u64 result = 1 % mod;
    base %= mod;
    while (exp > 0) {
        if ((exp & 1ULL) != 0ULL) {
            result = mul_mod(result, base, mod);
        }
        base = mul_mod(base, base, mod);
        exp >>= 1ULL;
    }
    return result;
}

int multiplicative_order_10(const int p) {
    int order = p - 1;
    int x = order;

    std::vector<int> factors;
    for (int f = 2; static_cast<std::int64_t>(f) * f <= x; ++f) {
        if (x % f != 0) {
            continue;
        }
        factors.push_back(f);
        while (x % f == 0) {
            x /= f;
        }
    }
    if (x > 1) {
        factors.push_back(x);
    }

    for (int q : factors) {
        while ((order % q) == 0 && pow_mod(10ULL, static_cast<u64>(order / q), static_cast<u64>(p)) == 1ULL) {
            order /= q;
        }
    }

    return order;
}

std::vector<int> primes_below(const int limit) {
    std::vector<bool> sieve(static_cast<std::size_t>(limit), true);
    if (limit > 0) {
        sieve[0] = false;
    }
    if (limit > 1) {
        sieve[1] = false;
    }

    for (int i = 2; static_cast<std::int64_t>(i) * i < limit; ++i) {
        if (!sieve[static_cast<std::size_t>(i)]) {
            continue;
        }
        for (int j = i * i; j < limit; j += i) {
            sieve[static_cast<std::size_t>(j)] = false;
        }
    }

    std::vector<int> primes;
    for (int i = 2; i < limit; ++i) {
        if (sieve[static_cast<std::size_t>(i)]) {
            primes.push_back(i);
        }
    }
    return primes;
}

bool can_divide_some_R10n(const int p) {
    if (p == 2 || p == 3 || p == 5) {
        return false;
    }

    int ord = multiplicative_order_10(p);
    while ((ord % 2) == 0) {
        ord /= 2;
    }
    while ((ord % 5) == 0) {
        ord /= 5;
    }

    return ord == 1;
}

std::int64_t solve(const int limit) {
    const std::vector<int> primes = primes_below(limit);

    std::int64_t sum = 0;
    for (int p : primes) {
        if (!can_divide_some_R10n(p)) {
            sum += p;
        }
    }

    return sum;
}

bool run_checkpoints() {
    if (solve(100) != 918) {
        std::cerr << "Checkpoint failed for limit=100" << '\n';
        return false;
    }
    if (!can_divide_some_R10n(17) || can_divide_some_R10n(19) || can_divide_some_R10n(3)) {
        std::cerr << "Checkpoint failed for prime behavior (3, 17, 19)" << '\n';
        return false;
    }
    return true;
}

}  // namespace

int main(int argc, char** argv) {
    Options options;
    if (!parse_arguments(argc, argv, options)) {
        return 1;
    }
    if (options.run_checkpoints && !run_checkpoints()) {
        return 2;
    }

    std::cout << solve(options.limit) << '\n';
    return 0;
}

Python

import sys
from typing import List

def parse_int_after_prefix(arg: str, prefix: str) -> (bool, int):
    if not arg.startswith(prefix):
        return False, 0
    tail = arg[len(prefix):]
    if not tail:
        return False, 0
    
    try:
        parsed = int(tail)
        return True, parsed
    except ValueError:
        return False, 0

def mul_mod(a: int, b: int, mod: int) -> int:
    return (a * b) % mod

def pow_mod(base: int, exp: int, mod: int) -> int:
    result = 1 % mod
    base %= mod
    while exp > 0:
        if exp & 1:
            result = mul_mod(result, base, mod)
        base = mul_mod(base, base, mod)
        exp >>= 1
    return result

def multiplicative_order_10(p: int) -> int:
    order = p - 1
    x = order
    
    factors = []
    f = 2
    while f * f <= x:
        if x % f == 0:
            factors.append(f)
            while x % f == 0:
                x //= f
        f += 1
    if x > 1:
        factors.append(x)
    
    for q in factors:
        while order % q == 0 and pow_mod(10, order // q, p) == 1:
            order //= q
    
    return order

def primes_below(limit: int) -> List[int]:
    if limit <= 2:
        return []
    
    sieve = [True] * limit
    sieve[0] = False
    sieve[1] = False
    
    i = 2
    while i * i < limit:
        if sieve[i]:
            j = i * i
            while j < limit:
                sieve[j] = False
                j += i
        i += 1
    
    return [i for i in range(2, limit) if sieve[i]]

def can_divide_some_R10n(p: int) -> bool:
    if p in (2, 3, 5):
        return False
    
    ord_val = multiplicative_order_10(p)
    
    while ord_val % 2 == 0:
        ord_val //= 2
    while ord_val % 5 == 0:
        ord_val //= 5
    
    return ord_val == 1

def solve(limit: int) -> int:
    primes = primes_below(limit)
    
    total = 0
    for p in primes:
        if not can_divide_some_R10n(p):
            total += p
    
    return total

def run_checkpoints() -> bool:
    if solve(100) != 918:
        return False
    if not can_divide_some_R10n(17) or can_divide_some_R10n(19) or can_divide_some_R10n(3):
        return False
    return True

def main():
    limit = 100000
    run_checkpoints_flag = True
    
    args = sys.argv[1:]
    for arg in args:
        if arg == "--skip-checkpoints":
            run_checkpoints_flag = False
        elif arg.startswith("--limit="):
            success, value = parse_int_after_prefix(arg, "--limit=")
            if success:
                limit = value
        else:
            print(f"Unknown argument: {arg}", file=sys.stderr)
            sys.exit(1)
    
    if limit < 2:
        print("Limit must be at least 2", file=sys.stderr)
        sys.exit(1)
    
    if run_checkpoints_flag and not run_checkpoints():
        print("Checkpoint failed", file=sys.stderr)
        sys.exit(2)
    
    print(solve(limit))

if __name__ == "__main__":
    main()

Java

import java.math.BigInteger;
import java.util.*;

public class Euler133 {
    public static void main(String[] args) {
        int limit = 100000;
        boolean[] sieve = new boolean[limit];
        Arrays.fill(sieve, true);
        sieve[0] = sieve[1] = false;
        for (int i = 2; i * i < limit; i++)
            if (sieve[i])
                for (int j = i * i; j < limit; j += i)
                    sieve[j] = false;

        BigInteger TEN = BigInteger.TEN;
        long sum = 0;
        for (int p = 2; p < limit; p++) {
            if (!sieve[p])
                continue;
            if (p == 2 || p == 3 || p == 5) {
                sum += p;
                continue;
            } // these can't divide
            // Find multiplicative order of 10 mod p
            int order = p - 1;
            // Factor order
            List<Integer> factors = new ArrayList<>();
            int x = order;
            for (int f = 2; (long) f * f <= x; f++) {
                if (x % f == 0) {
                    factors.add(f);
                    while (x % f == 0)
                        x /= f;
                }
            }
            if (x > 1)
                factors.add(x);
            BigInteger mod = BigInteger.valueOf(p);
            for (int q : factors)
                while (order % q == 0 && TEN.modPow(BigInteger.valueOf(order / q), mod).equals(BigInteger.ONE))
                    order /= q;
            // Check if order is 2^a * 5^b
            int o = order;
            while (o % 2 == 0)
                o /= 2;
            while (o % 5 == 0)
                o /= 5;
            if (o != 1)
                sum += p;
        }
        System.out.println(sum);
    }
}