Problem 358: Cyclic Numbers

View on Project Euler

Project Euler Problem 358 Solution

EulerSolve provides an optimized solution for Project Euler Problem 358, Cyclic Numbers, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary We seek the unique prime \(p\) such that the decimal expansion of \(1/p\) begins with the digits \(00000000137\) and the full repetend defines a cyclic number whose last five digits are \(56789\). Once that prime is identified, the required Project Euler quantity is the sum of the digits of the cyclic number. Mathematical Approach Step 1: Write the repetend as an integer For a prime \(p \neq 2,5\), the decimal expansion of \(1/p\) is purely repeating. If the period is maximal, then its length is \(p-1\), and the repeating block can be written as $$N=\frac{10^{p-1}-1}{p}.$$ The number \(N\) is the integer formed by the repeating digits of \(1/p\). When \(p\) is a full reptend prime in base \(10\), the multiples \(N,2N,\dots,(p-1)N\) are cyclic rotations of the same digit string, which is exactly why the problem is about cyclic numbers. Step 2: Convert the decimal prefix into an interval for \(p\) The first eleven digits after the decimal point are \(00000000137\)....

Detailed mathematical approach

Problem Summary

A cyclic number with \(n\) digits is a number whose products with \(1,2,3,\dots,n\) all consist of exactly the same digits in the same order, rotated circularly. The smallest is the six-digit \(142857\), with multiples \(285714\), \(428571\), \(571428\), \(714285\), \(857142\); the next is \(0588235294117647\), whose sixteenth multiple is \(9411764705882352\). Leading zeros count: a cyclic number is a digit string of fixed length, read as a positive integer, and \(n\ge 2\) throughout, since for \(n=1\) the condition is vacuous and every digit would qualify. Exactly one cyclic number has the eleven leftmost digits \(00000000137\) and the five rightmost digits \(56789\), with an unknown number of digits in between, and we must find the sum of all its digits. [1]

The whole solution in one sentence: every cyclic number is the repeating block of \(1/p\) for a prime \(p\) of longest possible decimal period, the leading digits trap \(p\) in a short interval and the trailing digits in one arithmetic progression, and the digit sum of such a block is \(9(p-1)/2\), so the answer follows from \(p\) alone.

Mathematical Approach

1. Rotating digits is multiplying by a power of ten

Write \(N\) as a string of exactly \(n\) digits, padded with leading zeros, and put \(M=10^{n}-1\). Split it after the first \(j\) digits: the front is \(A=\lfloor N/10^{n-j}\rfloor\), the back is \(B=N\bmod 10^{n-j}\), so \(N=A\cdot 10^{n-j}+B\). Rotating left by \(j\) places moves the front to the end and gives \(R_j=B\cdot 10^{j}+A\). Since \(10^{j}N=A\cdot 10^{n}+B\cdot 10^{j}\),

$$10^{j}N=A\,M+R_j,\qquad A=\lfloor N/10^{n-j}\rfloor,$$

so \(R_j\equiv 10^{j}N\pmod M\), meaning that \(R_j-10^{j}N\) is a multiple of \(M\). For \(N=142857\) and \(j=1\) this reads \(1428570=1\cdot 999999+428571\). As \(0\le R_j\le M\), a rotation is pinned down by its residue modulo \(M\) once we know it is not \(M\) itself, which happens only when every digit is a nine.

2. Long-period primes produce cyclic numbers

Let \(p\) be a prime other than \(2\) and \(5\), \(n=p-1\), \(M=10^{n}-1\). Fermat's little theorem gives \(10^{n}\equiv 1\pmod p\), so

$$N=\frac{10^{p-1}-1}{p}$$

is an integer with \(pN=M\), so \(N\lt M\) and the string of \(N\) is not all nines. Take \(k\) with \(1\le k\le p-1\) and suppose \(k\equiv 10^{j}\pmod p\); then \(p\) divides \(k-10^{j}\) and

$$kN-10^{j}N=\frac{k-10^{j}}{p}\cdot M$$

is a multiple of \(M\), while \(kN\le (p-1)N\lt M\). By Section 1 the rotation \(R_j\) is the one value in \(\{0,\dots,M-1\}\) congruent to \(10^{j}N\), so \(kN=R_j\) and the product is a rotation.

The multiplicative order \(\operatorname{ord}_p(10)\) is the least \(t\ge 1\) with \(10^{t}\equiv 1\pmod p\), and \(10^{0},\dots,10^{t-1}\) are then distinct, since an earlier repetition could be cancelled to give a smaller such \(t\); modulo \(7\) they leave \(1,3,2,6,4,5\), every nonzero remainder. [5] So a suitable \(j\) exists for every \(k\) exactly when \(\operatorname{ord}_p(10)=p-1\), that is when \(10\) is a primitive root modulo \(p\). Those are the full reptend primes, beginning \(7\), \(17\), \(19\), \(23\), \(29\), \(47\), \(59\), \(61\), \(97\), \(109\). [3] [6] Then \(k\mapsto j\) is a bijection, so the \(p-1\) products are the \(p-1\) rotations and \(N\) is cyclic with \(p-1\) digits. From \(pN=M\) we get \(1/p=N/M\), which is \(N\cdot(10^{-n}+10^{-2n}+\cdots)\), so the string is the repeating block of \(1/p\): \(1/7=0.\overline{142857}\), \(1/17=0.\overline{0588235294117647}\). [4]

Without a primitive root this collapses. For \(p=13\) the order is \(6\), and \((10^{12}-1)/13=076923076923\) is the block \(076923\) twice: \(6\) distinct rotations instead of \(12\), and \(2\cdot 076923076923=153846153846\) is none of them. [2]

3. Nothing else is a cyclic number

The converse is what the search rests on, and it is usually quoted without proof. [2] Let \(N\ge 1\) be cyclic with \(n\ge 2\) digits, \(M=10^{n}-1\). The products \(N,2N,\dots,nN\) are \(n\) distinct rotations, and there are at most \(n\) rotations, so they are all of them. Hence \(nN\le M\) and, by Section 1, there are exponents \(j_k\) with

$$kN\equiv 10^{j_k}N\pmod M\quad(1\le k\le n),$$

where \(k\mapsto j_k\) is a bijection onto \(\{0,\dots,n-1\}\). Put \(d=\gcd(N,M)\), \(N=dN_1\), \(M=dM_1\), so \(\gcd(N_1,M_1)=1\); dividing by \(d\) and cancelling \(N_1\) leaves

$$k\equiv 10^{j_k}\pmod{M_1}\quad(1\le k\le n).$$

  • \(M_1\gt n\). From \(nN\le M\), \(ndN_1\le dM_1\), so \(M_1\ge nN_1\ge n\ge 2\). As \(M_1\mid 10^{n}-1\), powers of \(10\) are invertible modulo \(M_1\), whereas \(M_1=n\) would make \(k=n\) both \(\equiv 0\) and invertible. So \(1,2,\dots,n\) are \(n\) distinct residues, each a power of \(10\).
  • \(\operatorname{ord}_{M_1}(10)=n\). Call it \(t\). The powers \(10^{0},\dots,10^{n-1}\) repeat with period \(t\) and so take \(\min(n,t)\) distinct values, and we just found \(n\) of them, so \(t\ge n\); and \(10^{n}\equiv 1\) with \(n=at+r\), \(0\le r\lt t\), gives \(10^{r}\equiv 1\), so \(r=0\) and \(t\mid n\). Hence \(t=n\) and the powers of \(10\) are precisely \(1,2,\dots,n\).
  • Every prime factor of \(M_1\) exceeds \(n\), since each \(k\le n\) is congruent to an invertible element, so \(\gcd(k,M_1)=1\).

Now count the digits twice, with \(s\) the digit sum. Across the \(n\) rotations each position receives each of the \(n\) digits once, so every column sums to \(s\) and the rotations total \(sM/9\), while as multiples they total \(N\cdot n(n+1)/2\). Substituting \(N=dN_1\), \(M=dM_1\) gives \(9N_1n(n+1)/2=sM_1\), so

$$M_1\ \Bigm|\ 9\,\frac{n(n+1)}{2}.$$

Let \(n\ge 3\) and let \(q\) be a prime factor of \(M_1\), so \(q\gt n\ge 3\). Then \(q\nmid 2\), \(q\nmid 9\), \(q\nmid n\), so \(q\mid n+1\), which with \(q\gt n\) forces \(q=n+1\). Thus \(M_1\) is a power of the prime \(n+1\), which occurs in \(9n(n+1)/2\) to the first power only, so \(M_1=n+1\). Finally \(N_1\le M_1/n\lt 2\), so \(N_1=1\) and, with \(p=n+1\),

$$\begin{aligned}N&=\frac{M}{M_1}=\frac{10^{p-1}-1}{p},\\ \operatorname{ord}_p(10)&=p-1.\end{aligned}$$

For \(n=2\) the divisibility gives \(M_1\mid 27\) with \(M_1\gt 2\), so \(M_1\in\{3,9,27\}\), whose orders of \(10\) are \(1\), \(1\), \(3\) and never \(2\): no two-digit cyclic number exists. With Section 2, the cyclic numbers of \(n\ge 2\) digits are exactly the blocks \((10^{p-1}-1)/p\) for full reptend primes \(p=n+1\), so hunting cyclic numbers is hunting such primes.

4. Testing whether ten is a primitive root

By Lagrange's theorem \(t=\operatorname{ord}_p(10)\) divides \(p-1\). [5] If \(t\lt p-1\) then \((p-1)/t\) has a prime factor \(q\) and \(t\mid (p-1)/q\), so \(10^{(p-1)/q}\equiv 1\pmod p\); conversely that congruence forces \(t\mid (p-1)/q\lt p-1\). Hence

$$\begin{aligned}&\operatorname{ord}_p(10)=p-1\\ &\iff 10^{(p-1)/q}\not\equiv 1\pmod p\end{aligned}$$

for every prime \(q\) dividing \(p-1\). [5] Only the distinct prime factors matter, and one modular exponentiation settles each. Most primes fail: the order of \(10\) is \(1\) for \(p=3\), \(2\) for \(11\), \(6\) for \(13\), \(15\) for \(31\), \(3\) for \(37\) and \(5\) for \(41\), against maxima \(2\), \(10\), \(12\), \(30\), \(36\), \(40\).

5. The digit sum without the digits

Let \(p\) be full reptend, \(n=p-1\), \(N=(10^{n}-1)/p\) and \(s\) the digit sum. The same double count, now with the products known to be the rotations, gives

$$\begin{aligned}\sum_{j=0}^{n-1}R_j&=\frac{sM}{9},\\ \sum_{k=1}^{p-1}kN&=N\,\frac{p(p-1)}{2}.\end{aligned}$$

Both sides count the same numbers and \(M=pN\), so \(spN/9=Np(p-1)/2\). Cancelling the nonzero \(pN/9\) leaves

$$s=\frac{9(p-1)}{2},$$

an integer because \(p\) is odd. It checks out on every small case: \(p=7\) gives \(27=1+4+2+8+5+7\), and \(p=17\), \(19\), \(23\), \(29\), \(47\) give \(72\), \(81\), \(99\), \(126\), \(207\), which is what adding the digits of their blocks returns. This is the lever of the problem: the sum is a one-line function of \(p\), although the number has hundreds of millions of digits.

6. The eleven leading digits trap the prime in an interval

The first eleven digits of the \(n\)-digit string, \(n=p-1\), read as a number are \(\lfloor N/10^{\,n-11}\rfloor\), and

$$\begin{aligned}\frac{N}{10^{\,n-11}}&=\frac{10^{n}-1}{p\cdot 10^{\,n-11}}\\ &=\frac{10^{11}-10^{\,11-n}}{p}.\end{aligned}$$

For \(n\ge 12\) the subtracted term is below \(1\), and \(p\) divides neither \(2\) nor \(5\), so \(10^{11}\) leaves a remainder of at least \(1\) on division by \(p\) and the floor does not move:

$$\bigl\lfloor N/10^{\,n-11}\bigr\rfloor=\bigl\lfloor 10^{11}/p\bigr\rfloor.$$

Demanding \(00000000137\) means \(137\le 10^{11}/p\lt 138\), that is \(10^{11}/138\lt p\le 10^{11}/137\). Numerically \(10^{11}/138=724637681.159\ldots\) and \(10^{11}/137=729927007.299\ldots\), so

$$724637682\le p\le 729927007,$$

a window of \(5289326\) integers. Equivalently: \(1/p\) must begin \(0.00000000137\).

7. The five trailing digits trap it in a progression

The last five digits say \(N\equiv 56789\pmod{10^{5}}\). The number has at least \(16\) digits, so \(p-1\ge 5\) and \(10^{p-1}\equiv 0\pmod{10^{5}}\); then \(pN=10^{p-1}-1\) reduces to \(pN\equiv -1\) and

$$56789\,p\equiv -1\pmod{10^{5}}.$$

As \(56789\) is odd and not a multiple of \(5\), it is invertible modulo \(10^{5}\), with inverse \(90109\), so

$$p\equiv -90109\equiv 9891\pmod{10^{5}},$$

and indeed \(56789\cdot 9891\equiv 99999\). Inside the window this progression starts at \(724709891\) and has \(53\) members. The argument also runs backwards: for every prime of the progression the block \((10^{p-1}-1)/p\) ends in \(56789\), full reptend or not, so the suffix alone says nothing about being cyclic.

8. A free test the suffix happens to pass

Because \(10^{5}=2^{5}\cdot 5^{5}\) is a multiple of \(40\), every member of the progression has the same residue modulo \(40\) — here \(9891\equiv 11\), so \(p\equiv 3\pmod 8\) and \(p\equiv 1\pmod 5\). By the supplementary laws of quadratic reciprocity, \(2\) is the remainder of a square modulo \(p\) exactly for \(p\equiv\pm 1\pmod 8\) and \(5\) exactly for \(p\equiv\pm 1\pmod 5\); here \(5\) is such a remainder and \(2\) is not, and a square times a non-square is a non-square, [7] so \(10\) is a non-square modulo every prime of the progression.

That is necessary for a primitive root, because Euler's criterion gives \(a^{(p-1)/2}\equiv 1\) for a square \(a\), whose order then divides \((p-1)/2\). [7] Had these five digits led to a class in which \(10\) is a square, no cyclic number could have ended in them at all. The test is the case \(q=2\) of Section 4, which the programs run anyway, so this explains the search without shortening it.

9. The whole method on a miniature

Ask the same question with leading digits \(001\) and trailing digits \(21\), small enough to check by hand.

  • Prefix: \(\lfloor 10^{3}/p\rfloor=1\) means \(10^{3}/2\lt p\le 10^{3}\), so \(501\le p\le 1000\).
  • Suffix: \(21p\equiv -1\pmod{100}\), and \(21^{-1}=81\), so \(p\equiv -81\equiv 19\pmod{100}\); indeed \(21\cdot 19=399\equiv 99\).
  • Candidates \(519\), \(619\), \(719\), \(819\), \(919\), of which \(519=3\cdot 173\) and \(819=3^{2}\cdot 7\cdot 13\) are composite.
  • Order test: \(618=2\cdot 3\cdot 103\), so Section 4 asks for \(10^{309}\), \(10^{206}\), \(10^{6}\) modulo \(619\), and repeated squaring gives \(618\), \(252\), \(315\) — none of them \(1\), so \(\operatorname{ord}_{619}(10)=618\). The other two primes fail already at \(q=2\): \(10^{359}\equiv 1\pmod{719}\) and \(10^{459}\equiv 1\pmod{919}\), so their orders are half of \(p-1\).

So \(p=619\), and the cyclic number is \((10^{618}-1)/619\), a \(618\)-digit string beginning \(001615\) and ending \(273021\), as demanded. Section 5 predicts the digit sum \(9\cdot 618/2=2781\), and adding the \(618\) digits one at a time gives \(2781\) too.

10. The real search

The actual instance has \(53\) candidates instead of five. Three are prime, and of those exactly one has \(10\) as a primitive root — the uniqueness the statement promises; the other two have order \((p-1)/5\), and by Artin's conjecture only some \(37\) percent of primes qualify at all. [8] That prime goes into \(s=9(p-1)/2\). Its cyclic number has between \(7.24\cdot 10^{8}\) and \(7.30\cdot 10^{8}\) digits, some \(700\) megabytes at one digit per byte, none of which is ever needed.

How the Code Works

All three programs run the same search: C++ in find_prime_p, called from main after the optional run_checkpoints, Python and Java directly inside solve. The names below are the C++ ones; Python shares them but writes is_full_reptend, and Java camel-cases them as modPow, isPrime, primeFactorsDistinct, isFullReptendPrimeBase10. The linked C++ revision is the implementation record. [10]

mod_pow(base, exp, mod) is binary modular exponentiation: at each step it multiplies the result by the current base when the lowest bit of exp is set, then squares the base and shifts exp right, so each exponent bit costs at most two modular multiplications and a test from Section 4 costs \(O(\log p)\) of them. The languages differ only in how they form those products. Two numbers near \(7.3\cdot 10^{8}\) multiply to near \(5.3\cdot 10^{17}\), still within \(64\) bits, yet C++ widens the product to __int128, Java routes it through BigInteger.multiply and remainder and back with longValue, [9] and Python just multiplies, its integers being unbounded.

is_prime(n) is trial division: handle \(2\), then odd d while d*d <= n. prime_factors_distinct(n) returns each prime factor once, dividing it out before moving on, into a std::vector, an ArrayList or a Python list, one candidate at a time. is_full_reptend_prime_base10(p) is Section 4 verbatim: refuse \(p\) unless it is prime and neither \(2\) nor \(5\), factor \(p-1\), and return false as soon as mod_pow(10, phi/q, p) gives \(1\).

The search is Sections 6 and 7: integer division gives low = 100000000000/138 + 1 and high = 100000000000/137, the start is shifted onto the class \(9891\) modulo mod = 100000, and the loop steps by mod. The in-loop checks that 100000000000/p is \(137\) and that (p*56789 + 1) % mod is zero follow from how the loop was built, so they guard a mis-set bound rather than filter. Here the versions part company: C++ keeps scanning and reports an error if a second candidate survives, establishing uniqueness rather than assuming it, Java stops at the first survivor with break, Python keeps the last. The result becomes (9*(p-1))/2, Section 5, which fits in std::uint64_t and runs in long arithmetic in Java, so the factor \(9\) cannot overflow an int.

Unless --skip-checkpoints is passed, the C++ run_checkpoints — the only self-check in the three files — verifies that \(7\) is full reptend in base ten, that the search finds a unique prime, that it is the expected one and that the digit-sum formula reproduces the expected sum; it calls find_prime_p itself, so a default run searches twice.

Complexity and Verification

The window holds \(5289326\) integers and the progression cuts it to \(C=53\) candidates, so all cost is per candidate. Trial division dominates, with \(\lfloor\sqrt{729927007}\rfloor=27017\) and so at most \(13508\) odd divisors each. In one full sweep the odd-divisor loops run \(46552\) iterations for the primality tests and \(1172\) for factoring \(p-1\) at the three primes, a composite with a small factor exiting at once — \(38\) of the \(50\) composites have a prime factor below \(100\). The order tests stop at the first failing factor, so they perform ten modular exponentiations and \(390\) modular multiplications. Python makes one sweep, C++ two by default, Java stops at the survivor after \(52\) candidates. Counting word arithmetic as one operation the method is \(O(C\sqrt{p})\), and it holds a few words plus the distinct prime factors of \(p-1\), at most \(\log_2 p\) of them. All three finish in milliseconds here, with outputs that agree byte for byte.

Verified for this article, beyond those checkpoints:

  • The statement's data: \((10^{6}-1)/7=142857\) with its six products and \((10^{16}-1)/17=0588235294117647\) with its sixteenth, \(9411764705882352\); and the identity of Section 1 for every rotation of the blocks of \(7\), \(17\), \(19\).
  • Section 3 by brute force: for each \(n\) from \(2\) to \(8\), every \(N\) up to \((10^{n}-1)/n\) — a necessary bound, since \(nN\) must still have \(n\) digits — was tested for \(kN\) being a rotation at all \(k\le n\). The only cyclic number is \(142857\) at \(n=6\), and \(7\) is the only full reptend prime with \(p-1\le 8\).
  • The digit-sum formula against direct summation for \(p=7\), \(17\), \(19\), \(23\), \(29\), \(47\) — \(27\), \(72\), \(81\), \(99\), \(126\), \(207\) — and for the \(618\)-digit block of Section 9, \(2781\).
  • The prefix and suffix on the surviving prime, modularly instead of by building the string: \(\lfloor 10^{11}/p\rfloor\) is \(137\), and \(-p^{-1}\bmod 10^{5}\), which is \(N\bmod 10^{5}\) because \(pN\equiv -1\), is \(56789\).
  • Section 8, and the specialness of the uniqueness, against a sample: among \(300\) randomly drawn residue classes modulo \(10^{5}\) coprime to \(10\), each supplying \(52\) or \(53\) candidates from the same window, all \(142\) classes in which \(10\) is a square contain no full reptend prime, \(143\) contain none at all, and only \(6\) contain exactly one.

These extra checks are not packaged in the published programs.

Footnotes and References

The rotation identity, the converse of Section 3, the digit-sum formula, the interval and the progression are derived in full above; the references supply the background facts.

  1. Project Euler 358 — Cyclic Numbers. The official statement defines a cyclic number with n digits by its products with 1, 2, …, n being rotations of the same digit string, gives 142857 and 0588235294117647 as the first two, notes that leading zeros matter, and asks for the digit sum of the unique one of the form 00000000137⋯56789.
  2. Wikipedia — Cyclic number. States that cyclic numbers have the form (b^(p−1) − 1)/p for a prime p not dividing the base, and gives 076923076923 for p = 13 as a case where that fails. It only says the form "can be shown", hence Section 3.
  3. Wikipedia — Full reptend prime. An odd prime is full reptend in base b when the order of b is p − 1, that is when b is a primitive root, and only then does p have a cyclic number of p − 1 digits.
  4. Wikipedia — Repeating decimal. The repetend length of 1/p is the order of 10 modulo p, hence p − 1 for a primitive root and a proper divisor otherwise; in the first case the repetend read as an integer is a cyclic number, as with 142857/999999.
  5. Wikipedia — Multiplicative order defines the order as the least positive k with a^k ≡ 1 (mod n) and states, from Lagrange's theorem, that it divides φ(n), which is p − 1 for a prime. Wikipedia — Primitive root modulo n gives the code's test: for g coprime to n, g is a primitive root exactly when no residue g^(φ(n)/qᵢ) mod n is 1, over the distinct prime factors qᵢ of φ(n).
  6. OEIS A001913 — Full reptend primes: primes with primitive root 10. Terms 7, 17, 19, 23, 29, 47, 59, 61, 97, 109, …, also called long period primes and linked to the cyclic numbers of A180340.
  7. Wikipedia — Quadratic reciprocity gives the supplementary laws that x² ≡ 2 (mod p) is solvable exactly for p ≡ ±1 (mod 8) and that 5 is a quadratic residue modulo p exactly for p ≡ ±1 (mod 5), and the product rules by which a residue times a non-residue is a non-residue. Wikipedia — Euler's criterion: a^((p−1)/2) is 1 modulo an odd prime for a residue, −1 otherwise.
  8. Wikipedia — Artin's conjecture on primitive roots. An integer neither a perfect square nor −1 is conjectured to be a primitive root modulo infinitely many primes, with density Artin's constant 0.3739558136… or a rational multiple of it; the page records that this is still unresolved, so the 37 percent above is heuristic.
  9. GCC manual — 128-bit integers documents __int128 as an extension on targets with a 128-bit integer mode. Java — BigInteger documents arbitrary-precision integers with multiply, remainder and a longValue that keeps only the low 64 bits.
  10. C++ — Euler358.cpp. The immutable revision with mod_pow, is_prime, prime_factors_distinct, is_full_reptend_prime_base10, find_prime_p and run_checkpoints.

Problem 358 source code

C++

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

namespace {

using u64 = std::uint64_t;

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

bool is_prime(u64 n) {
    if (n < 2) {
        return false;
    }
    if ((n % 2ULL) == 0ULL) {
        return n == 2;
    }
    for (u64 d = 3; d * d <= n; d += 2) {
        if ((n % d) == 0ULL) {
            return false;
        }
    }
    return true;
}

std::vector<u64> prime_factors_distinct(u64 n) {
    std::vector<u64> factors;
    if ((n % 2ULL) == 0ULL) {
        factors.push_back(2ULL);
        while ((n % 2ULL) == 0ULL) {
            n /= 2ULL;
        }
    }
    for (u64 d = 3; d * d <= n; d += 2) {
        if ((n % d) == 0ULL) {
            factors.push_back(d);
            while ((n % d) == 0ULL) {
                n /= d;
            }
        }
    }
    if (n > 1ULL) {
        factors.push_back(n);
    }
    return factors;
}

bool is_full_reptend_prime_base10(const u64 p) {
    if (!is_prime(p) || p == 2ULL || p == 5ULL) {
        return false;
    }
    const u64 phi = p - 1ULL;
    const std::vector<u64> factors = prime_factors_distinct(phi);
    for (const u64 q : factors) {
        if (mod_pow(10ULL, phi / q, p) == 1ULL) {
            return false;
        }
    }
    return true;
}

u64 find_prime_p() {
    // First 11 digits of repetend of 1/p are 00000000137.
    const u64 low = 100000000000ULL / 138ULL + 1ULL;
    const u64 high = 100000000000ULL / 137ULL;

    // Last 5 digits of repetend are 56789.
    // If N is repetend integer, p*N = 10^(p-1)-1 => p*N == -1 (mod 100000).
    // Therefore p*56789 == -1 (mod 100000), i.e. p == 9891 (mod 100000).
    const u64 mod = 100000ULL;
    const u64 residue = 9891ULL;

    u64 start = low;
    const u64 rem = start % mod;
    if (rem <= residue) {
        start += residue - rem;
    } else {
        start += mod - (rem - residue);
    }

    u64 found = 0;
    for (u64 p = start; p <= high; p += mod) {
        if (100000000000ULL / p != 137ULL) {
            continue;
        }
        if (!is_full_reptend_prime_base10(p)) {
            continue;
        }
        if ((p * 56789ULL + 1ULL) % mod != 0ULL) {
            continue;
        }

        if (found != 0ULL) {
            std::cerr << "More than one candidate prime found\n";
            return 0ULL;
        }
        found = p;
    }

    return found;
}

bool run_checkpoints() {
    if (is_full_reptend_prime_base10(7ULL) == false) {
        std::cerr << "Checkpoint failed: p=7 should be full reptend in base 10\n";
        return false;
    }

    const u64 p = find_prime_p();
    if (p == 0ULL) {
        std::cerr << "Checkpoint failed: no unique prime found\n";
        return false;
    }

    if (p != 729809891ULL) {
        std::cerr << "Checkpoint failed: unexpected prime value\n";
        return false;
    }

    if ((9ULL * (p - 1ULL)) / 2ULL != 3284144505ULL) {
        std::cerr << "Checkpoint failed: digit sum formula mismatch\n";
        return false;
    }

    return true;
}

}  // namespace

int main(int argc, char** argv) {
    bool skip_checkpoints = false;
    for (int i = 1; i < argc; ++i) {
        const std::string arg(argv[i]);
        if (arg == "--skip-checkpoints") {
            skip_checkpoints = true;
        } else {
            std::cerr << "Unknown argument: " << arg << '\n';
            return 1;
        }
    }

    if (!skip_checkpoints && !run_checkpoints()) {
        return 2;
    }

    const u64 p = find_prime_p();
    const u64 answer = (9ULL * (p - 1ULL)) / 2ULL;
    std::cout << answer << '\n';
    return 0;
}

Python

def solve():
    def mod_pow(base, exp, mod):
        result = 1 % mod
        base %= mod
        while exp > 0:
            if exp & 1:
                result = (result * base) % mod
            base = (base * base) % mod
            exp >>= 1
        return result

    def is_prime(n):
        if n < 2:
            return False
        if n % 2 == 0:
            return n == 2
        d = 3
        while d * d <= n:
            if n % d == 0:
                return False
            d += 2
        return True

    def prime_factors_distinct(n):
        factors = []
        if n % 2 == 0:
            factors.append(2)
            while n % 2 == 0:
                n //= 2
        d = 3
        while d * d <= n:
            if n % d == 0:
                factors.append(d)
                while n % d == 0:
                    n //= d
            d += 2
        if n > 1:
            factors.append(n)
        return factors

    def is_full_reptend(p):
        if not is_prime(p) or p == 2 or p == 5:
            return False
        phi = p - 1
        factors = prime_factors_distinct(phi)
        for q in factors:
            if mod_pow(10, phi // q, p) == 1:
                return False
        return True

    # Find prime p with specific constraints
    low = 100000000000 // 138 + 1
    high = 100000000000 // 137
    mod = 100000
    residue = 9891

    start = low
    rem = start % mod
    if rem <= residue:
        start += residue - rem
    else:
        start += mod - (rem - residue)

    found = 0
    p = start
    while p <= high:
        if 100000000000 // p == 137:
            if is_full_reptend(p):
                if (p * 56789 + 1) % mod == 0:
                    found = p
        p += mod

    answer = (9 * (found - 1)) // 2
    return str(answer)

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

Java

import java.util.*;

public class Euler358 {
    static long modPow(long base, long exp, long mod) {
        long result = 1 % mod;
        base %= mod;
        while (exp > 0) {
            if ((exp & 1) != 0) {
                result = (long) (((java.math.BigInteger.valueOf(result)).multiply(java.math.BigInteger.valueOf(base)))
                        .remainder(java.math.BigInteger.valueOf(mod)).longValue());
            }
            base = (long) (((java.math.BigInteger.valueOf(base)).multiply(java.math.BigInteger.valueOf(base)))
                    .remainder(java.math.BigInteger.valueOf(mod)).longValue());
            exp >>= 1;
        }
        return result;
    }

    static boolean isPrime(long n) {
        if (n < 2)
            return false;
        if (n % 2 == 0)
            return n == 2;
        for (long d = 3; d * d <= n; d += 2) {
            if (n % d == 0)
                return false;
        }
        return true;
    }

    static List<Long> primeFactorsDistinct(long n) {
        List<Long> factors = new ArrayList<>();
        if (n % 2 == 0) {
            factors.add(2L);
            while (n % 2 == 0)
                n /= 2;
        }
        for (long d = 3; d * d <= n; d += 2) {
            if (n % d == 0) {
                factors.add(d);
                while (n % d == 0)
                    n /= d;
            }
        }
        if (n > 1)
            factors.add(n);
        return factors;
    }

    static boolean isFullReptendPrimeBase10(long p) {
        if (!isPrime(p) || p == 2 || p == 5)
            return false;
        long phi = p - 1;
        for (long q : primeFactorsDistinct(phi)) {
            if (modPow(10, phi / q, p) == 1)
                return false;
        }
        return true;
    }

    public static String solve() {
        long low = 100000000000L / 138 + 1;
        long high = 100000000000L / 137;
        long mod = 100000L;
        long residue = 9891L;

        long start = low;
        long rem = start % mod;
        if (rem <= residue) {
            start += residue - rem;
        } else {
            start += mod - (rem - residue);
        }

        long found = 0;
        for (long p = start; p <= high; p += mod) {
            if (100000000000L / p != 137)
                continue;
            if (!isFullReptendPrimeBase10(p))
                continue;
            if ((p * 56789L + 1) % mod != 0)
                continue;
            found = p;
            break;
        }

        return String.valueOf((9 * (found - 1)) / 2);
    }

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