$2520$ is the smallest number that can be divided by each of the numbers from $1$ to $10$ without any remainder.
What is the smallest positive number that is evenly divisible by all of the numbers from $1$ to $20$?
Computing the LCM iteratively
We want to basically find the lowest common multiple (LCM) of 1, 2, ..., 20. Since $\operatorname{lcm}(a, b) = \operatorname{lcm}(b, a)$ and $\operatorname{lcm}(a, b, c) = \operatorname{lcm}(a, \operatorname{lcm}(b, c))$ we can iteratively compute the LCM using
function smallest_multiple(n)
result = 1
for i in 2:n
result = lcm(result, i)
end
return result
end
where the Julia lcm function computes the LCM using
and where $\operatorname{gcd}(a, b)$ is the greatest common divisor (GCD) of $a$ and $b$. Julia's gcd function computes the GCD using the binary GCD algorithm.
Computing smallest_multiple(20) is then done very quickly in just 249.069 ns.
We don't even need all of 1, 2, ..., 20. Any $k \leq 10$ can be doubled until it lands between 11 and 20, so $\operatorname{lcm}(1, 2, \dots, 20) = \operatorname{lcm}(11, 12, \dots, 20)$. We can't start any later than 11, since 11 has no other multiple up to 20. For general $n$ the latest starting point is OEIS A188666.
How big the answer gets
The answers to the general problem, $\operatorname{lcm}(1, 2, \dots, n) = 1, 2, 6, 12, 60, 60, 420, \dots$ are OEIS A003418, and they grow very quickly with smallest_multiple(42) being the largest that does not overflow Int64, returning 219060189739591200 in 975.391 ns.
smallest_multiple(Int128(88)) is the largest that does not overflow Int128, returning 8076030954443701744994070304101969600 in 7.824 μs. That's 8.02× longer than $n = 42$ for only about twice as many lcm calls, since 128-bit arithmetic is much slower than 64-bit arithmetic.
Why 42 and 88? The logarithm of $\operatorname{lcm}(1, 2, \dots, n)$ is the second Chebyshev function $\psi(n)$, and the prime number theorem is equivalent to $\psi(n) \sim n$, so $\operatorname{lcm}(1, 2, \dots, n) \approx e^n$. Int64 holds numbers up to $2^{63} \approx e^{43.7}$ and Int128 up to $2^{127} \approx e^{88.0}$.
We can keep going past 88 using BigInt. Going all the way to smallest_multiple(BigInt(100000)) returns a 43452-digit number in 582.811 ms.
Primorials
The primorial $n\#$ is the product of all the primes up to $n$,
so $20\# = 2 \cdot 3 \cdot 5 \cdot 7 \cdot 11 \cdot 13 \cdot 17 \cdot 19 = 9699690$. The primorials $n\# = 1, 2, 6, 6, 30, 30, 210, \dots$ are OEIS A034386, and the products of the first $n$ primes $2, 6, 30, 210, 2310, \dots$ are OEIS A002110. To compute $n\#$ we can sieve the primes up to $n$ and multiply them together.
Primorials are useful here because $\operatorname{lcm}(1, 2, \dots, n)$ is built from them. The exponent of a prime $p$ in $\operatorname{lcm}(1, 2, \dots, n)$ is the largest $k$ with $p^k \leq n$, which is the number of $k \geq 1$ with $p \leq n^{1/k}$. So we can collect the primes in layers. The first layer has every prime up to $n$, the second every prime up to $\sqrt{n}$, the third every prime up to $\sqrt[3]{n}$, and so on:
The GMP library, which Julia's BigInt is built on, already has a fast primorial function, mpz_primorial_ui. Its implementation sieves the primes up to $n$ and packs as many as fit into each 64-bit word. It then multiplies the words together with a binary splitting product. Computing $\operatorname{lcm}(1, 2, \dots, n)$ is then just a handful of calls to it and a few multiplications.
Computing the LCM with primorials
The layered product works, but we can simplify it. Every layer after the first only contains primes $p \leq \sqrt{n}$, and together those layers give each such prime its missing $p^{e-1}$, where $p^e$ is the largest power of $p$ up to $n$. So
where the second factor is OEIS A045948. For $n = 20$ it's $2^3 \cdot 3 = 24$, the extra factors that $16 = 2^4$ and $9 = 3^2$ need, and $20\# \cdot 24 = 9699690 \cdot 24 = 232792560$.
We can get $n\#$ straight from GMP's mpz_primorial_ui (Julia's BigInt is GMP's mpz_t, so we can pass it in directly)
function primorial(n::Integer)
n < 0 && throw(DomainError(n, "primorial is only defined for non-negative n"))
result = BigInt()
ccall((:__gmpz_primorial_ui, :libgmp), Cvoid, (Ref{BigInt}, Culong), result, Culong(n))
return result
end
and $p^{e-1}$ is the largest power of $p$ that is at most $n \div p$, which we can find with integer division so nothing overflows:
function smallest_multiple_primorial(n)
extra = big(1)
for p in sieve_of_eratosthenes(isqrt(n))
q = 1
while q <= n ÷ p ÷ p
q *= p
end
extra *= q
end
return primorial(n) * extra
end
This computes smallest_multiple_primorial(100000) in 401.156 μs, 1,450× faster than the 582.811 ms of the iterative loop. It also scales much better. The loop multiplies its whole growing result at every step, so its cost grows roughly like $n^2$, but going up to $n = 10^6$ takes 9.524 ms, only 23.7× longer for a 10× larger $n$. Going all the way to $n = 10^8$, a 43,428,686-digit number, takes 2.906 s, 305× longer than $n = 10^6$ for a 100× larger $n$.