Problem 320: Factorials Divisible by a Huge Integer

Project Euler Problem 320 Solution

EulerSolve provides an optimized solution for Project Euler Problem 320, Factorials Divisible by a Huge Integer, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Define $N(i)$ as the smallest integer $n$ such that $$n! \text{ is divisible by } (i!)^{1234567890}.$$ Then $$S(u)=\sum_{i=10}^{u}N(i).$$ We are given $$S(1000)=614538266565663,$$ and must compute $$S(1000000)\bmod 10^{18}.$$ Mathematical Approach 1) Divisibility of factorials is a prime-valuation problem. Write $$i!=\prod_{p\le i}p^{e_{i,p}},\qquad e_{i,p}=v_p(i!).$$ Then $$(i!)^K=\prod_{p\le i}p^{K e_{i,p}},\qquad K=1234567890.$$ Therefore $$(i!)^K\mid n!\iff v_p(n!)\ge K e_{i,p}\quad\text{for every prime }p\le i.$$ So the whole problem reduces to matching required prime exponents inside $n!$. 2) Legendre's formula tells us the exponent of one prime in $n!$. For a fixed prime $p$, $$v_p(n!)=\sum_{j\ge1}\left\lfloor\frac{n}{p^j}\right\rfloor.$$ This function is monotone increasing in $n$. So if we want the smallest $n$ satisfying $$v_p(n!)\ge t,$$ we can find it by binary search. 3) Per-prime bottlenecks. For each prime $p$, define the current target $$t_p=K\,v_p(i!).$$ Let $n_p$ be the smallest integer with $$v_p(n_p!)\ge t_p.$$ Then the smallest $n$ satisfying all prime conditions simultaneously is simply $$N(i)=\max_{p\le i} n_p.$$ The largest per-prime obstruction determines the answer. 4) Incremental update from $i-1$ to $i$. This is the key optimization....

Detailed mathematical approach

Problem Summary

Fix the exponent $K=1234567890$. For an integer $i\ge2$, let $N(i)$ be the smallest integer $n$ such that $n!$ is divisible by $(i!)^{K}$, and put

$$S(u)=\sum_{i=10}^{u}N(i).$$

The statement gives $S(1000)=614538266565663$ and asks for $S(10^6)\bmod 10^{18}$, a sum of $999\,991$ terms. ^[1]

The whole solution in one sentence: divisibility of $n!$ by $(i!)^{K}$ is one inequality per prime, Legendre's formula makes each one monotone in $n$ so that a binary search inverts it, $N(i)$ is the largest of those per-prime answers, and since the step to $i$ changes only the primes dividing $i$, the table is carried forward with fewer than three inversions per step on average.

Mathematical Approach

Nothing as large as $n!$ is ever built. A factorial is described completely by the exponents of its primes, and for one prime that exponent has a short formula which can be inverted. The work lies in doing this for $10^6$ values of $i$ without redoing it for every prime each time.

1. Divisibility is one inequality per prime

Write $v_p(m)$ for the exponent of the prime $p$ in $m$, so that $m=\prod_p p^{v_p(m)}$. For positive integers, unique factorisation makes $a\mid b$ equivalent to $v_p(a)\le v_p(b)$ for every prime $p$. ^[2] Since $v_p\bigl((i!)^{K}\bigr)=K\,v_p(i!)$, we get

$$(i!)^{K}\mid n!\iff v_p(n!)\ge K\,v_p(i!)$$

for every prime $p\le i$; primes above $i$ impose nothing, as they do not divide $i!$. For $i=10$ the factorisation $10!=3628800=2^{8}\cdot3^{4}\cdot5^{2}\cdot7$ gives four requirements: $v_2(n!)\ge8K$, $v_3(n!)\ge4K$, $v_5(n!)\ge2K$ and $v_7(n!)\ge K$.

2. Legendre's formula

Among $1,2,\dots,n$ exactly $\lfloor n/p\rfloor$ numbers are multiples of $p$, exactly $\lfloor n/p^{2}\rfloor$ are multiples of $p^{2}$, and so on; each multiple of $p^{j}$ contributes one factor $p$ beyond the $j-1$ already counted, so

$$v_p(n!)=\sum_{j\ge1}\Bigl\lfloor\frac{n}{p^{j}}\Bigr\rfloor,$$

a sum that stops as soon as $p^{j}\gt n$. ^[3] For $n=10$ and $p=2$ the terms are $5,2,1$, giving $v_2(10!)=8$; for $p=3$ they are $3,1$, giving $4$; and $v_5(10!)=2$, $v_7(10!)=1$, exactly the exponents above.

The same reference states the closed form that explains the shape of the whole problem. If $s_p(n)$ is the sum of the base-$p$ digits of $n$, then

$$v_p(n!)=\frac{n-s_p(n)}{p-1}.$$

3. The answer is a maximum over primes

For a fixed prime $p$, $v_p(n!)$ never decreases as $n$ grows, so the set of $n$ satisfying one requirement is an interval $[n_p,\infty)$, where

$$n_p=\min\{\,n:v_p(n!)\ge t_p\,\}$$

and $t_p=K\,v_p(i!)$ is its target. An intersection of such intervals is again an interval, so all requirements hold exactly when $n\ge\max_p n_p$, that is

$$N(i)=\max_{p\le i}n_p.$$

Which prime decides is not obvious, and more than one can attain the maximum; in Section 6 the winner is neither the smallest nor the largest prime.

4. A worked example small enough to check by hand

Replace $K$ by $5$ and take $i=10$. The targets are $t_2=40$, $t_3=20$, $t_5=10$ and $t_7=5$. Legendre's formula gives the four minima directly:

$v_2(44!)=22+11+5+2+1$${}=41\ge40$, while $v_2(43!)=39$, so $n_2=44$;
$v_3(45!)=15+5+1=21\ge20$, while $v_3(44!)=14+4+1=19$, so $n_3=45$;
$v_5(45!)=9+1=10\ge10$, while $v_5(44!)=8+1=9$, so $n_5=45$;
$v_7(35!)=5\ge5$, while $v_7(34!)=4$, so $n_7=35$.

Hence $N(10)=45$ for the exponent $5$, and multiplying out confirms it: $(10!)^{5}$ divides $45!$ but not $44!$. Note the tie between $3$ and $5$. The same exponent gives $N(9)=45$ and $N(11)=55$, again confirmed by dividing actual factorials.

5. How large $n_p$ is, and which prime wins

The digit form of Legendre's formula turns the requirement $v_p(n!)\ge t_p$ into

$$n-s_p(n)\ge K\bigl(i-s_p(i)\bigr),$$

because $t_p(p-1)=K\,v_p(i!)(p-1)$${}=K(i-s_p(i))$. Since $s_p(n)\ge1$, the requirement forces $n_p\gt K(i-s_p(i))$, and Section 8 bounds the excess by $(p-1)D\le51(p-1)$, which stays below $K$ for every $p\le10^6$. Lowering $s_p(i)$ by one raises $K(i-s_p(i))$ by $K$, more than any excess, so a strictly smaller digit sum means a strictly larger minimum: the winner is always a prime with the smallest digit sum at $i$, and the excesses decide between those.

Two consequences follow. If $i$ is a prime power, one $p$ has $s_p(i)=1$ and the strongest requirement of all, close to $K(i-1)$: at the prime $i=999983$ that prime itself wins. For $i/2\lt p\le i$ the digits of $i$ are $1$ and $i-p$, so $s_p(i)$ grows with the gap and such a prime matters only when it lies close to $i$. At $i=10^6$ the smallest digit sums are $s_2=7$ and $s_5=8$, so $2$ wins and the largest term is $N(10^6)=1234559248024796$.

6. The same example with the real exponent

Now keep $K=1234567890$ and follow $i=10,11,12$. At $i=10$ the minima are $n_2=9876543136$, $n_3=9876543138$, $n_5=9876543150$ and $n_7=7407407378$, so

$$N(10)=9876543150,$$

and the winner is $p=5$, not $p=2$: since $10=(1010)_2=(101)_3=(20)_5$ has digit sum $2$ in all three bases, the primes $2,3,5$ share the bound $K(10-2)=9876543120$, and the excess is largest for $5$.

At $i=11$ only the prime $11$ changes, with $v_{11}(11!)=1$ and target $K$. Its minimum is $n_{11}=12345678950$, which Legendre's formula confirms exactly: the terms $\lfloor n/11^{j}\rfloor$ are

$$1122334450,\ 102030404,\ 9275491,\ 843226,$$

$$76656,\ 6968,\ 633,\ 57,\ 5,$$

and they add up to $1234567890=K$, while $v_{11}(12345678949!)=K-1$. So $N(11)=12345678950$, well above the runner-up $9876543150$.

At $i=12=2^{2}\cdot3$ only $2$ and $3$ change: their minima rise to $12345678920$ and $12345678918$, both just below the leader, so $N(12)=12345678950$ again — the maximum can stand still for several steps. The sum $S(12)=34567901050$ is small enough to reproduce by hand.

7. Monotonicity and the incremental update

Two monotonicity facts make the computation cheap. First, $v_p(i!)=v_p((i-1)!)+v_p(i)$, so every target $t_p$ is non-decreasing in $i$ and changes only for the primes dividing $i$, growing there by $K\,v_p(i)$. Secondly, raising a target can only raise the minimum $n_p$ that answers it. Hence every $n_p$ is non-decreasing, and so is their maximum:

$$N(i)=\max\bigl(N(i-1),\ \max_{p\mid i}n_p\bigr).$$

A single running value therefore suffices: at step $i$ we recompute $n_p$ only for the prime divisors of $i$ and lift the maximum if one of them exceeds it.

This is what makes the problem tractable. With $\pi(i)$ the number of primes at most $i$, rescanning all of them for every $i$ would need $\sum_{i\le10^6}\pi(i)=40\,947\,676\,475$ inversions; the incremental version needs one per prime divisor, which is

$$\begin{aligned}\sum_{i\le10^{6}}\omega(i)&=\sum_{p\le10^{6}}\Bigl\lfloor\frac{10^{6}}{p}\Bigr\rfloor\\&=2\,853\,708,\end{aligned}$$

Here $\omega(i)$ counts the distinct prime divisors of $i$; both sides count the pairs $(i,p)$ with $p\mid i$, of which a fixed $p$ has $\lfloor10^6/p\rfloor$. ^[4] That saves a factor of about $14\,000$.

8. Inverting Legendre's formula

One primitive is still missing: given $p$ and a target $t$, find the least $n$ with $v_p(n!)\ge t$. Both ends of a search interval come from the formulas above; for the lower end, the geometric series gives

$$v_p(n!)\le\sum_{j\ge1}\frac{n}{p^{j}}=\frac{n}{p-1},$$

so $v_p(n!)\ge t$ forces $n\ge A$, where $A=t(p-1)$. For the upper end we use the targets of this problem, $t=K\,v_p(i!)$: then $n=K\,i$ passes, because $\lfloor Ki/p^{j}\rfloor\ge K\lfloor i/p^{j}\rfloor$ term by term; hence $n_p\le K\,i$, and any candidate has at most $D=1+\lfloor\log_p(Ki)\rfloor$ base-$p$ digits. Minimality means $v_p((n_p-1)!)\lt t$, in digit form $n_p-1-s_p(n_p-1)\lt A$, so

$$n_p\le A+s_p(n_p-1)\le A+(p-1)D.$$

The window above $A$ is therefore at most $(p-1)D$ wide, and $D$ is largest for $p=2$, where $D\le51$: the window is at most $51$ for $p=2$ and below $3\,000\,000$ for a prime near $10^6$, where $D=3$. The widest window met in the full run was $1\,984\,832$.

Inside the window $v_p(n!)$ is non-decreasing, so a binary search for the leftmost passing $n$ finds $n_p$ in about $\log_2$ of the window width evaluations ^[5] — at most $21$ even in the widest case, each the short division chain of Section 2. Finally $pt$ itself passes, since $v_p((pt)!)\ge\lfloor pt/p\rfloor=t$, so a failing probe stays below $pt$, and one step adds at most $t(p-1)$: every probe stays below $2pt\le4Ku$, because $t(p-1)\le K\,u$ and $p\le2(p-1)$.

9. Factorising every $i$ quickly

The update needs the prime divisors of each $i$, and trial division would dominate the cost. One sieve stores in $\operatorname{spf}[j]$ the smallest prime factor of $j$: start with $\operatorname{spf}[j]=j$ and, for each $j$ found prime, overwrite the still untouched entries of its multiples. Dividing $i$ by $\operatorname{spf}[i]$ repeatedly then factorises it in as many divisions as it has prime factors. ^[6]

10. The algorithm, and how large the sum is

Collecting the pieces:

build the smallest-prime-factor sieve up to $u$;
keep arrays $t_p$ and $n_p$, both starting at $0$, and a running maximum;
for $i=2,\dots,u$: factorise $i$; for each prime $p\mid i$ with exponent $e$ raise $t_p$ by $Ke$, invert Legendre's formula, and lift the running maximum if needed;
add the running maximum to the sum once $i\ge10$.

One caution about sizes. Each term is at most $K\,i$ by Section 8, so

$$S(10^{6})\lt K\sum_{i\le10^{6}}i\lt2^{70};$$

and since every $i\le10^6$ has at most $20$ binary digits, the prime $2$ alone gives $N(i)\gt K(i-20)$, whence

$$S(10^{6})\gt K\sum_{i=10}^{10^{6}}(i-20)\gt2^{69}.$$

The unreduced sum therefore needs exactly $70$ bits and leaves both $64$-bit ranges, although every single term fits comfortably; the reduction modulo $10^{18}$ happens only at the end. Inside an inversion the values are far smaller: targets and bounds stay below $K\,u\approx1.24\times10^{15}$ and every probe below $4Ku\lt2^{63}$. Three values can be checked before the full run: $S(12)=34567901050$, $S(20)=181481480341$ and the official $S(1000)=614538266565663$.

How the Code Works

valuation_factorial(n, p) is Section 2: it divides its local copy of $n$ by $p$ repeatedly and adds the quotients, costing $\lfloor\log_p n\rfloor+1$ divisions and never forming a power of $p$ that could overflow.

minimal_n_for_target(p, target, previous) is Section 8. It sets lower to the bound $A$, written target * (p - 1) and raised to previous if that is larger, then walks a variable hi upwards: while the valuation at hi falls short by $\delta$, it adds $\delta(p-1)$. By the digit form this lands exactly on $A+s_p(\mathrm{hi})$, so the gap closes only if the digit sum does not grow and the test must be repeated; each step adds at least $p-1$, so the loop ends. A binary search over [lower, hi] then returns the leftmost passing value.

The previous argument is a valid bound, since $n_p$ never decreases, but it never binds: between two updates of one prime the target grows by at least $K$, lifting $A$ by at least $K(p-1)$, more than the window $(p-1)D$. In all $2\,853\,708$ calls the arithmetic bound was indeed the larger.

build_spf(limit) fills the array with $\operatorname{spf}[j]=j$ and then, for each $i$ whose entry is still $i$, marks the multiples of $i$ from $i^{2}$ on where the entry is untouched: the sieve of Eratosthenes with the payload of Section 9.

solve_sum(upper) holds the arrays required (the targets $t_p$) and needed_n (the minima $n_p$), both indexed by the prime, together with current_max and the accumulator. Its inner while (x > 1) loop peels off $p=\operatorname{spf}[x]$, counts the exponent, raises required[p] by K * exp, calls the inversion and lifts current_max; the guard if (i >= 10) implements the lower limit of the sum. Indexing by the prime leaves the composite positions unused — $16$ MB in C++ at $u=10^{6}$ — in exchange for needing no prime list.

run_checkpoints() compares solve_sum(1000) with $614538266565663$; main runs that check unless --skip-checkpoints is given, rejects any other argument, and reduces solve_sum(1000000) modulo $10^{18}$. ^[9]

The three programs share this update but differ in naming and arithmetic. The names above are the C++ ones; Java has valuationFactorial, minimalNForTarget, buildSpf and solve, while Python nests its two helpers inside solve and sieves there. C++ keeps sieve entries, indices, primes and exponents in int, targets and minima in uint64_t, and the accumulator in unsigned __int128, a compiler extension it also uses for target * factor and the hi update. ^[7] Java has an int[] sieve, long targets and a BigInteger accumulator, whose sum.mod(outMod).longValue() is safe because the residue lies in $[0,10^{18})$. ^[8] Python uses arbitrary-precision integers. Only C++ has the $S(1000)$ checkpoint and the flag that skips it.

Complexity and Verification

The sieve costs $O(u\log\log u)$ time and, in C++ and Java, one $32$-bit word per entry. The main loop performs $\sum_{i\le u}\omega(i)$ inversions, $2\,853\,708$ at $u=10^{6}$, growing like $u(\ln\ln u+B)$ with $B\approx0.261497$, an estimate of $2\,887\,289$ within $1.2\%$ of the exact count. ^[4] Section 8 bounds one inversion: every probe stays below $4Ku$, hence has at most $D_p=1+\lfloor\log_p(4Ku)\rfloor$ base-$p$ digits, so one valuation evaluation costs at most $D_p$ divisions; the escalation probes lie in $[A,A+(p-1)D_p]$ and rise by at least $p-1$ each, so it needs $O(D_p)$ evaluations and the binary search $O(\log(pD_p))$. The total is therefore

$$O\Bigl(u\log\log u+\sum_{i\le u}\sum_{p\mid i}D_p\bigl(D_p+\log(pD_p)\bigr)\Bigr),$$

the first term covering the sieve and the factorisations. In plain numbers, the main pass performed $30\,360\,934$ valuation evaluations, $23\,677\,175$ binary-search steps — $8.3$ per inversion — and $448\,874\,683$ divisions in the valuation helper, with $3\,626\,619$ more in the factorisations; no escalation needed more than four steps. C++ and Java hold three arrays of $u+1$ entries of $4$, $8$ and $8$ bytes, about $20$ MB, while Python stores integer objects; the measured peak was $21.9$ MB against $77$ MB. On one machine, with Apple clang at -O2, OpenJDK 17 and CPython 3.11, C++ needs $0.85$ s including the checkpoint, Java $0.75$ s and Python $19$ s.

The following checks were run for this article. The digit-sum identity of Section 2 was compared with the sum form for seven primes and all $n\lt3000$. The formula $N(i)=\max_p n_p$ was checked against brute force on actual factorials for the pairs

$$\begin{aligned}(i,K)=\ &(4,1),(4,2),(5,3),(6,2),\\&(10,1),(10,2),(10,5),(12,3),(7,4),\end{aligned}$$

always agreeing. An independent implementation that rescans all primes $p\le i$ for every $i$, with no incremental update, reproduced $S(1000)=614538266565663$ and the values $N(10),N(11),N(12)$ above; the same rescan at $i=10^{6}$ over all $78\,498$ primes returned $N(10^{6})=1234559248024796$, the value the incremental loop carries. An instrumented pass then checked all $2\,853\,708$ inversions: every returned minimum passes its valuation test while its predecessor fails, no probe left the bounds of Section 8, and the maximum never decreased. The three published programs and that fourth implementation print the same final result.

Footnotes and References

The references below supply the background facts used above. The maximum over primes, the digit-sum description of the winner, the bounds framing the search and the operation counts are derived here.

Project Euler 320 — Factorials Divisible by a Huge Integer. The official statement defines N(i) as the smallest n with n! divisible by (i!)^1234567890, sets S(u) as the sum of N(i) for 10 ≤ i ≤ u, gives S(1000) = 614538266565663 and asks for S(1000000) modulo 10^18.
Wikipedia — Fundamental theorem of arithmetic. Every integer greater than 1 is a product of primes, uniquely up to order, and its canonical representation as a product of prime powers with distinct primes and positive exponents is unique. The exponent comparison of Section 1 follows from that uniqueness.
Wikipedia — Legendre's formula. It states both forms used here: the exponent of p in n! equals the sum of ⌊n/p^i⌋ over i ≥ 1, whose terms vanish once p^i exceeds n, and equals (n − s_p(n))/(p − 1) for the base-p digit sum s_p(n).
OEIS A013939 — partial sums of the number of distinct prime divisors. The sequence equals both the sum of ω(k) for k ≤ n and the sum of ⌊n/p⌋ over primes p ≤ n, the identity of Section 7; its term for n = 10 is 11 = 5 + 3 + 2 + 1. It also gives the asymptotics n(log log n + B) + o(n) with the Mertens constant B = 0.261497…
Wikipedia — Binary search. The article gives the leftmost-element procedure that narrows a lower and an upper bound until they meet and returns the lower one — the loop the code uses — with ⌊log₂(n) + 1⌋ iterations in the worst case.
cp-algorithms — Linear sieve. The page computes the minimum prime factor of every number in [2, n] and notes that this array factorises any such number in time proportional to the size of its factorisation. The solutions keep that payload but mark multiples like Eratosthenes, so their sieve costs O(n log log n).
GCC manual — 128-bit integers. The manual documents __int128 and unsigned __int128 as an extension on targets with a wide enough integer mode; the accumulator uses the unsigned form for the 70-bit sum.
Java API — BigInteger. The class provides immutable arbitrary-precision integers; mod always returns a result between 0 and the modulus minus one, which is why the final longValue conversion cannot overflow here.
C++ — Euler320.cpp. The linked immutable revision contains valuation_factorial, minimal_n_for_target, build_spf, solve_sum and run_checkpoints, and accepts one flag, --skip-checkpoints. It is the implementation record for this article.

Problem 320 source code

C++

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

namespace {

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

constexpr u64 K = 1234567890ULL;
constexpr u64 OUTPUT_MOD = 1000000000000000000ULL;

u64 valuation_factorial(u64 n, const int p) {
    u64 sum = 0;
    const u64 prime = static_cast<u64>(p);
    while (n > 0) {
        n /= prime;
        sum += n;
    }
    return sum;
}

u64 minimal_n_for_target(const int p, const u64 target, const u64 previous) {
    if (target == 0ULL) {
        return 0ULL;
    }

    const u64 factor = static_cast<u64>(p - 1);
    const u64 lower = std::max(previous, static_cast<u64>(static_cast<u128>(target) * factor));

    u64 hi = lower;
    while (true) {
        const u64 have = valuation_factorial(hi, p);
        if (have >= target) {
            break;
        }
        hi = static_cast<u64>(static_cast<u128>(hi) + static_cast<u128>(target - have) * factor);
    }

    u64 lo = lower;
    while (lo < hi) {
        const u64 mid = lo + (hi - lo) / 2ULL;
        if (valuation_factorial(mid, p) >= target) {
            hi = mid;
        } else {
            lo = mid + 1ULL;
        }
    }
    return lo;
}

std::vector<int> build_spf(const int limit) {
    std::vector<int> spf(static_cast<std::size_t>(limit + 1));
    std::iota(spf.begin(), spf.end(), 0);

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

u128 solve_sum(const int upper) {
    const std::vector<int> spf = build_spf(upper);

    std::vector<u64> required(static_cast<std::size_t>(upper + 1), 0ULL);
    std::vector<u64> needed_n(static_cast<std::size_t>(upper + 1), 0ULL);

    u64 current_max = 0;
    u128 sum = 0;

    for (int i = 2; i <= upper; ++i) {
        int x = i;
        while (x > 1) {
            const int p = spf[static_cast<std::size_t>(x)];
            int exp = 0;
            while (x % p == 0) {
                x /= p;
                ++exp;
            }

            required[static_cast<std::size_t>(p)] += K * static_cast<u64>(exp);
            const u64 updated = minimal_n_for_target(p, required[static_cast<std::size_t>(p)],
                                                     needed_n[static_cast<std::size_t>(p)]);
            needed_n[static_cast<std::size_t>(p)] = updated;
            if (updated > current_max) {
                current_max = updated;
            }
        }

        if (i >= 10) {
            sum += static_cast<u128>(current_max);
        }
    }

    return sum;
}

bool run_checkpoints() {
    const u128 sample = solve_sum(1000);
    if (sample != static_cast<u128>(614538266565663ULL)) {
        std::cerr << "Checkpoint failed: S(1000) 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 u128 total = solve_sum(1000000);
    const u64 answer = static_cast<u64>(total % static_cast<u128>(OUTPUT_MOD));
    std::cout << answer << '\n';
    return 0;
}

Python

def solve():
    K = 1234567890
    OUTPUT_MOD = 10**18
    upper = 1000000

    def valuation_factorial(n, p):
        s = 0
        while n > 0:
            n //= p
            s += n
        return s

    def minimal_n_for_target(p, target, previous):
        if target == 0:
            return 0
        factor = p - 1
        lower = max(previous, target * factor)
        hi = lower
        while True:
            have = valuation_factorial(hi, p)
            if have >= target:
                break
            hi += (target - have) * factor
        lo = lower
        while lo < hi:
            mid = lo + (hi - lo) // 2
            if valuation_factorial(mid, p) >= target:
                hi = mid
            else:
                lo = mid + 1
        return lo

    # Build SPF
    spf = list(range(upper + 1))
    for i in range(2, upper + 1):
        if i * i > upper:
            break
        if spf[i] == i:
            for j in range(i*i, upper + 1, i):
                if spf[j] == j:
                    spf[j] = i

    required = [0] * (upper + 1)
    needed_n = [0] * (upper + 1)
    current_max = 0
    total = 0

    for i in range(2, upper + 1):
        x = i
        while x > 1:
            p = spf[x]
            exp = 0
            while x % p == 0:
                x //= p
                exp += 1
            required[p] += K * exp
            updated = minimal_n_for_target(p, required[p], needed_n[p])
            needed_n[p] = updated
            if updated > current_max:
                current_max = updated
        if i >= 10:
            total += current_max

    return str(total % OUTPUT_MOD)

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

Java

import java.math.BigInteger;

public class Euler320 {
    static final long K = 1234567890L;
    static final long OUTPUT_MOD = 1000000000000000000L;

    static long valuationFactorial(long n, int p) {
        long sum = 0;
        long prime = p;
        while (n > 0) {
            n /= prime;
            sum += n;
        }
        return sum;
    }

    static long minimalNForTarget(int p, long target, long previous) {
        if (target == 0)
            return 0;

        long factor = p - 1;
        long lower = Math.max(previous, target * factor);

        long hi = lower;
        while (true) {
            long have = valuationFactorial(hi, p);
            if (have >= target) {
                break;
            }
            hi += (target - have) * factor;
        }

        long lo = lower;
        while (lo < hi) {
            long mid = lo + (hi - lo) / 2;
            if (valuationFactorial(mid, p) >= target) {
                hi = mid;
            } else {
                lo = mid + 1;
            }
        }
        return lo;
    }

    static int[] buildSpf(int limit) {
        int[] spf = new int[limit + 1];
        for (int i = 0; i <= limit; i++)
            spf[i] = i;
        for (int i = 2; i * i <= limit; ++i) {
            if (spf[i] == i) {
                for (int j = i * i; j <= limit; j += i) {
                    if (spf[j] == j)
                        spf[j] = i;
                }
            }
        }
        return spf;
    }

    public static String solve() {
        int upper = 1000000;
        int[] spf = buildSpf(upper);

        long[] required = new long[upper + 1];
        long[] neededN = new long[upper + 1];

        long currentMax = 0;
        BigInteger sum = BigInteger.ZERO;

        for (int i = 2; i <= upper; ++i) {
            int x = i;
            while (x > 1) {
                int p = spf[x];
                int exp = 0;
                while (x % p == 0) {
                    x /= p;
                    ++exp;
                }
                required[p] += K * exp;
                long updated = minimalNForTarget(p, required[p], neededN[p]);
                neededN[p] = updated;
                if (updated > currentMax) {
                    currentMax = updated;
                }
            }
            if (i >= 10) {
                sum = sum.add(BigInteger.valueOf(currentMax));
            }
        }

        BigInteger outMod = BigInteger.valueOf(1000000000000000000L);
        return String.valueOf(sum.mod(outMod).longValue());
    }

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