Primorial in algorithms

In this post we present several practical algorithms based on primorials.

# Definition

The definition of the primorial, denoted as `n#`, is defined as the product of primes `p` less than or equal to `n`:

`n# = \prod_{p <= n} p`

Additionally, `p_n#` is defined as the product of the first `n` primes (where `p_n` is the `n`-th prime):

`p_n# = \prod_{k=1}^n p_k`

# Overview

Primorials can be used effectively in the following algorithms:

...and possibly in many others.

# Smoothness test

A positive integer `n` is said to be `B`-smooth if none of its prime factors exceed `B`.

Smooth numbers play an important role in modern general-purpose integer factorization algorithms, which use smooth numbers for creating a congruence of squares, which leads to the factorization of `n`.

Therefore checking if a large number `n` is `B`-smooth, without knowing the prime factorization of `n`, is a very important thing to do efficiently.

The clever idea is to use primorials, by computing `k = B#` and checking the greatest common divisor of `n` and `k`.

Let `n` be the number to be tested.
Let `k` be the primorial of `B`.
Compute the greatest common divisor: `g = gcd(n, k)`.
If `g` is greater than `1`, then `n = r g^e`, for some `e >= 1`:

if `r = 1`, then `n` is `B`-smooth.
otherwise, set `n = r` and go to step `3`.

If this step is reached, then `n` is not `B`-smooth.

This algorithm is efficient when `B` is fixed, as we can compute `k = B#` only once and reuse the same value of `k` for testing various values of `n` for `B`-smoothness.

# Squarefree pretest

A positive integer `n` is said to be squarefree if all of its prime factors have multiplicity `1`. For pre-checking if a number is squarefree, we propose the following algorithm, where `B` is a small bound, at most `B = floor(n^(1/3))`.

Let `n` be the number to be tested.
Let `k` be the primorial of `B`.
Return false if `n` is a perfect power > 1.
Compute the greatest common divisor: `g = gcd(n, k)`.
If `g` is greater than `1`, then `n = r g`:

if `r = 1`, then `n` is squarefree.
if `r` is a perfect power > 1, then `n` is not squarefree.
if `gcd(r, k) > 1`, then `n` is not squarefree.

If this step is reached, then `n` is probably squarefree.

Numbers of special form, such as `a^k ± 1`, very rarely have prime factors greater than `10^6` with multiplicity greater than 1. Therefore, by setting `B = 10^6` or `B = 10^7`, the above algorithm very efficiently identifies arbitrary large numbers as potential squarefree with high probability.

# Perfect power pretest

A composite positive integer `n` is said to be a perfect power if it can be represented as `n = a^b` for some `a >= 1` and `b >= 2`.

If `n` has the prime factorization:

`n = \prod_{k=1} p_k^(e_k)`

then `gcd(e_k)` represents the maximal possible value for `b` in `n = a^b`.

Based on this observation, we try to find small prime factors `p_k` of `n`, with multiplicity `e_k`.

If `gcd(e_k) = 1`, then `n` is not a perfect power.
On the other hand, if `gcd(e_k) = b`, for some `b >= 2` and `n = floor(n^(1/b))^b`, then `n` is a perfect power.

The most common algorithm for checking if a given number `n` is a perfect power or not, starts by iterating from `k = floor(log_2(n))` down to `k = 2`, at each step checking if `n = floor(n^(1/k))^k`. When such a value of `k` is found, then `n` is a `k`-th power. If no such value of `k` exists, then `n` is not a perfect power.

The following algorithm can be used as a fast pretest (with a small value of `B`) before checking if a number is really a perfect power:

Let `n` be the number to be tested.
Let `k` be the primorial of `B`.
Compute the greatest common divisor: `g = gcd(n,k)`.
If `g` is greater than `1`, then `n = r g^e`, for some `e >= 1`:

if `r = 1` and `e = 1`, then `n` is not a perfect power.
if `r = 1` and `e > 1`, then `n` is a perfect power.
Factorize `g` into primes `p_k` and compute multiplicity `e_k` in `n`:

if `gcd(e_k) = 1`, then `n` is not a perfect power.
if `gcd(e_k) = b`, with `b >= 2` and `n = floor(n^(1/b))^b` , then `n` is a perfect power.
if `n` divided by `\prod_{k=1} p_k^(e_k)` is prime, then `n` is not a perfect power. (optional)

If `n` is prime, then `n` is a not a perfect power. (optional)
If this step is reached, then `n` is probably a perfect power.

# Trial division

Sometimes is useful to efficiently find all the small prime factors `p <= B` of a given number `n`. The following algorithm does this efficiently for arbitrary large inputs, assuming that `B` is fixed and the primorial of `B` is pre-computed.

Let `n` be the number to be factorized.
Let `k` be the primorial of `B`.
Compute the greatest common divisor: `g = gcd(n,k)`.
If `g` is greater than `1`:

Factorize `g` into primes `p_k`, where `p_k <= B`, and compute multiplicity `e_k` in `n`.
This gives us `n = m \prod_{k=1} p_k^(e_k)`, for some `m >= 1`.

Otherwise, `n` is not divisible by any prime `p <= B`.

In practice, computing `gcd(n, B#)`, with `B#` pre-computed, is much more efficient than checking `n` for divisibility by `p` for every prime `p <= B`.

The worst case of this approach, is when `n` is composed entirely of small primes `p <= B`, as we have to compute `g = gcd(n, B#) = n`, followed by the factorization of `g` into primes, by checking `g` for divisibility by `p` for every prime `p <= B`, but this is unlikely to occur often for large `n`.

On the other hand, the best case of this algorithm, is when `n` has no small prime factors `p <= B`, as only one `gcd` is required to determine that.

On average, `g` is usually many times smaller than `n` and can be factorized quickly into primes, therefore this approach is considered of practical value.

# Primorial inflation/deflation

Another surprising application of primorials is in compressing a very large integer with many consecutive small prime factors into a relatively small integer.

The sequence of numbers that can be primoriarly deflated into an integer, is given by A025487, which are numbers `n` of the form `n = p_1^(e_1) * p_2^(e_2) * ... * p_k^(e_k)` with `e_(k-1) >= e_k` for all `k >= 2`. Among other terms, the sequence includes the factorial and the highly composite numbers.

Let `I(n)` be the primorial inflation of `n`, defined as:

`I(\prod_{k=1} p_k^(e_k)) = \prod_{k=1} (p_k#)^(e_k)`

For example:

`n = I(52900585920) = 2#^6 * 3#^3 * 5# * 7# * 13# * 61# * 1103#`

which is a highly-composite number with `503` digits, such that `n-1` and `n+1` are both prime numbers (cf. A321995).

Similarly, let `D(n)` be the primorial deflation of `n`, defined as:

`D(n) = ` the unique integer `m` such that `I(m) = n`.

For example:

`D(n!)` is given by A307035(n).
`D(H(n))` is given by A329902(n), where `H(n)` is the `n`-th highly-composite number.

A simple algorithm for finding `D(n)` would be to repeatedly divide `n` by the largest primorial `(p_k#)^(e_k)` that goes into `n`.
Then the product `\prod_{k=1} (p_k)^(e_k)` is the primorial deflation of `n`.

However, a more efficient algorithm would be to factorize `n` into primes `n = \prod_{k=1} (p_k)^(e_k)`, which allows us to compute the primorial deflation `D(n)` as:

`D(n) = \prod_{k=1} (p_k / f(p_k))^(e_k)`

where `f(p)` is the largest prime less than `p`, with `f(2) = 1`.

The above formula is a generalization of the concept of primorial deflation, due to Rémy Sigrist, which allows any positive integer to be primoriarly deflated. However, for numbers that are not terms of A025487, the value `D(n)` would be a fraction `D(n) = x/y`, where the inverse operation is `n = \frac{I(x)}{I(y)}` (cf. A319626, A329900).

Example:

`D(2^8 * 3^4 * 5^2 * 7) = 2^8 (3/2)^4 (5/3)^2 (7/5) = 5040`
`D(2 * 3^4 * 5 * 13) = 2 (3/2)^4 (5/3) (13/11) = 1755/88`

Additionally, the following two prime numbers:

`a = lcm(1..36497) - 1`
`b = lcm(1..55291) + 1`

can be deflated into:

`a = I(169905470978880) - 1`
`b = I(374772076136640) + 1`

and inflated as:

`a = 2#^6 * 3#^3 * 5# * 7# * 13# * 31# * 191# * 36497# - 1`
`b = 2#^6 * 3#^3 * 5# * 7# * 13# * 37# * 233# * 55291# + 1`

The terms "primorial inflation" and "primorial deflation" were first coined by Matthew Vandermast in A181815.

# Implementations

This section includes source-code implementation of the presented algorithms.

:: Perl

:: Sidef

Mathematics and computer science