Fast Modular Multiplication

A modular multiplication requires a remainder operation, which is a slow operation if the modulus is a general integer. For example, contemporary processors can multiply two 64-bit integers, producing a 128-bit result, with a latency of 3 or 4 clock cycles. But, dividing a 128-bit integer by a 64-bit integer, producing a 64-bit quotient and a 64-bit remainder, is considerably slower (tens of clock cycles).

If the modulus is a power of two, say $2^n$,the remainder operation is very fast; the remainder is just the last n bits of the number being remaindered. In 1985, Peter Montgomery came up with a beautiful way to explore this to efficiently perform general remaindering operations without performing an expensive division.

The Greatest Common Divisor

Let $p_k$ be the $k$-th prime number, so that $p_1=2, p_2=3, p_3=5$, and so on,

Each positive integer can be factored into prime factors in a unique way (this is the fundamental theorem of arithmetic).

Let $a = \Pi^\infty_{k=1}p^{a_k}_k$, where $a_k$ is the number of times $p_k$ divides $a$. Since $a$ is a finite number, almost all of the $a_k$ values will be zero.

Likewise of $b$, let $b = \Pi^\infty_{k=1}p^{b_k}_k$.

Then $gcd(a,b) =\Pi^\infty_{k=1}p^{min(a_k,b_k)}_k$ and $lcm(a,b) =\Pi^\infty_{k=1}p^{max(a_k,b_k)}_k$

If $gcd(a,b) = 1$ then $a$ and $b$ are said to be relatively prime (or coprime).

The greatest common division can be generalized to a polynomial with integer coefficients!

Algorithm

Assume that $a \ge 0$ and that $b \ge 0$. Then:

• $gcd(a,b) = gcd(b,a)$, and so $gcd(a,b) = gcd(max(a,b), min(a,b))$. Thus, by exchanging $a$ with $b$ if necessary, we may assume that $a \ge b$.

• as any positive integer divides 0 we have $gcd(a,0) = a$ for $a \gt 0$. The mathematicians say that $gcd(0,0) = 0$, and so we can say that $gcd(a,0) = a$ as long as $a \ge 0$.

• If $a \ge b$ the $gcd(a,b) = gcd(a-b,b)$. We can keep subtracting $b$ from (the updated) $a$ until it becomes smaller than $b$, and so $gcd(a,b) = gcd(a mod b, b) = gcd(b, a mod b)$.

These observations give rise to the following so-called Euclid's algorithm (coded in C, but it can easily be translated to another programming language):

long gcd(long a, long b)
{
    while(b != 0) { long c = a % b; a = b; b = c; } return a;
}

The GNU MP library has a function, mpz_gcd, for this; pari-gp does this with the gcd function.