All Euler problems
Project Euler

Divisibility of Factorials

The smallest number m such that 10 divides m! is m = 5. The smallest number m such that 25 divides m! is m = 10. Let s(n) be the smallest number m such that n divides m!. So s(10) = 5 and s(25) = 1...

Source sync Apr 19, 2026
Problem #0549
Level Level 08
Solved By 3,203
Languages C++, Python
Answer 476001479068717
Length 309 words
number_theorysearchmodular_arithmetic

Problem Statement

This archive keeps the full statement, math, and original media on the page.

The smallest number \(m\) such that \(10\) divides \(m!\) is \(m=5\).

The smallest number \(m\) such that \(25\) divides \(m!\) is \(m=10\).

Let \(s(n)\) be the smallest number \(m\) such that \(n\) divides \(m!\).

So \(s(10)=5\) and \(s(25)=10\).

Let \(S(n)\) be \(\sum s(i)\) for \(2 \le i \le n\).

You are given \(S(100)=2012\).

Find \(S(10^8)\).

Problem 549: Divisibility of Factorials

Mathematical Analysis

Key Property: Multiplicativity over Prime Powers

If n=p1e1p2e2prern = p_1^{e_1} \cdot p_2^{e_2} \cdots p_r^{e_r}, then:

s(n)=max(s(p1e1),s(p2e2),,s(prer))s(n) = \max(s(p_1^{e_1}), s(p_2^{e_2}), \ldots, s(p_r^{e_r}))

Proof: The factorizations in m!m! for different primes are independent. For nm!n | m!, we need pieim!p_i^{e_i} | m! for every ii. Each condition is satisfied when ms(piei)m \geq s(p_i^{e_i}). Thus the minimum mm satisfying all conditions simultaneously is the maximum.

Computing s(pe)s(p^e) via Legendre’s Formula

By Legendre’s formula, the exponent of prime pp in m!m! is:

νp(m!)=j=1mpj\nu_p(m!) = \sum_{j=1}^{\infty} \left\lfloor \frac{m}{p^j} \right\rfloor

We need the smallest mm such that νp(m!)e\nu_p(m!) \geq e. Since mm must be a multiple of pp (otherwise νp(m!)=νp((m1)!)\nu_p(m!) = \nu_p((m-1)!)), we only check multiples of pp.

For a prime pp, s(p)=ps(p) = p. For s(pe)s(p^e), we search multiples of pp: p,2p,3p,p, 2p, 3p, \ldots until the factorial exponent reaches ee.

Editorial

For a prime power p^e, s(p^e) is the smallest m such that the exponent of p in m! >= e. By Legendre’s formula: v_p(m!) = sum(floor(m/p^k) for k=1,2,…). We use a sieve to find smallest prime factors, then factorize each n. We sieve primes** up to N=108N = 10^8 using a modified sieve. We then iterate over each prime pp: iterate through powers peNp^e \leq N and for each, compute s(pe)s(p^e) by finding the smallest multiple of pp whose factorial contains at least ee factors of pp. Finally, sieve approach: For each nn, compute s(n)=maxs(n) = \max of ss values over its prime power factorization. This is done during sieving: for each prime pp and each power pep^e, update all multiples.

Pseudocode

Sieve primes** up to $N = 10^8$ using a modified sieve
Initialize** an array `s[2..N]` where `s[i] = 0`
For each prime $p$**: iterate through powers $p^e \leq N$ and for each, compute $s(p^e)$ by finding the smallest multiple of $p$ whose factorial contains at least $e$ factors of $p$
Sieve approach**: For each $n$, compute $s(n) = \max$ of $s$ values over its prime power factorization. This is done during sieving: for each prime $p$ and each power $p^e$, update all multiples
Sum** all $s(i)$ for $2 \leq i \leq N$

Correctness

Theorem. The method described above computes exactly the quantity requested in the problem statement.

Proof. The preceding analysis identifies the admissible objects and derives the formula, recurrence, or exhaustive search carried out by the algorithm. The computation evaluates exactly that specification, so every valid contribution is included once and no invalid contribution is counted. Therefore the returned value is the required answer. \square

Complexity Analysis

  • Time: O(NloglogN)O(N \log \log N) for the sieve, plus work proportional to the number of prime powers.
  • Space: O(N)O(N) for the sieve and the ss array.

Answer

476001479068717\boxed{476001479068717}

Code

Each problem page includes the exact C++ and Python source files from the local archive.

C++ project_euler/problem_549/solution.cpp
#include <bits/stdc++.h>
using namespace std;

/*
 * Problem 549: Divisibility of Factorials
 * Find S(10^8) where s(n) = smallest m such that n | m!
 * Key: s(n) = max over prime power factors p^e of s(p^e)
 * Use sieve to compute s(n) for all n up to N.
 */

const int N = 100000000; // 10^8

// s[i] = smallest m such that i | m!
// We store as short values relative approach: use int array
static int s[N + 1];

// Compute the p-adic valuation of m!
// i.e., the exponent of p in m!
long long factorial_val(int m, int p) {
    long long v = 0;
    long long pk = p;
    while (pk <= m) {
        v += m / pk;
        if (pk > m / p) break; // prevent overflow
        pk *= p;
    }
    return v;
}

// Find smallest m such that p^e | m!
// m must be a multiple of p
int find_s_prime_power(int p, int e) {
    // Binary search: find smallest m (multiple of p) with factorial_val(m, p) >= e
    // Lower bound: e (since val(e!) >= e/p + ... but we need at least e)
    // Actually s(p^e) <= p * e (since val((pe)!) >= e + e/p + ... >= e)
    // But we can do binary search on multiples of p

    int lo = p, hi = (long long)p * e < N ? p * e : N;
    // Ensure hi is multiple of p
    hi = (hi / p) * p;
    if (hi < p) hi = p;

    while (lo < hi) {
        int mid = lo + (hi - lo) / 2;
        mid = (mid / p) * p;
        if (mid < lo) mid = lo;
        if (factorial_val(mid, p) >= e) {
            hi = mid;
        } else {
            mid += p;
            if (mid > hi) { lo = hi; break; }
            lo = mid;
        }
    }
    return lo;
}

int main() {
    ios_base::sync_with_stdio(false);

    // Initialize s array
    // s[0] = s[1] = 0 (not used)
    // For primes p, s[p] = p
    // For prime powers p^e, s[p^e] = find_s_prime_power(p, e)
    // For composite n, s[n] = max over prime power factors

    // Step 1: Sieve to find smallest prime factor
    vector<int> spf(N + 1, 0);
    for (int i = 2; i <= N; i++) {
        if (spf[i] == 0) { // i is prime
            for (int j = i; j <= N; j += i) {
                if (spf[j] == 0) spf[j] = i;
            }
        }
    }

    // Step 2: For each prime, precompute s(p^e) for all powers p^e <= N
    // We'll store these in a map for quick lookup, but it's faster to
    // just compute s(n) by factoring n using spf

    // Step 3: Compute s(n) for each n from 2 to N
    long long total = 0;

    for (int n = 2; n <= N; n++) {
        int tmp = n;
        int result = 0;
        while (tmp > 1) {
            int p = spf[tmp];
            int e = 0;
            while (tmp % p == 0) {
                tmp /= p;
                e++;
            }
            // compute s(p^e)
            int sp = find_s_prime_power(p, e);
            result = max(result, sp);
        }
        total += result;
    }

    cout << total << endl;

    return 0;
}