Currently there may be errors shown on top of a page, because of a missing Wiki update (PHP version and extension DPL3). |
Topics | Help • Register • News • History • How to • Sequences statistics • Template prototypes |
Montgomery multiplication
Montgomery multiplication is an algorithm used to perform multiple precision modular multiplication quickly by replacing division (which is a slow operation) by multiplications.
In this article the modular multiplication to be done will be [math]\displaystyle{ ab\,\equiv \,c\,\pmod{m} }[/math]. We will also assume that [math]\displaystyle{ m\,\lt \,2^n }[/math].
Before performing the actual multiplication, both factors must be converted to Montgomery representation. The algorithm will find the Montgomery representation of the product which will be converted back to the actual product. All these steps are slower than a naive modular multiplication, using a multiplication followed by a division, but when several multiplications and additions can be done without converting to and from Montgomery representation (for example: modular exponentiation, ECM, p-1, etc.) this method is really fast.
The Montgomery representation of a number [math]\displaystyle{ a }[/math] is the value
- [math]\displaystyle{ a'=2^n\,a\,\bmod{m} }[/math].
From this formula we can see that [math]\displaystyle{ p+q\,\equiv \, r }[/math] implies [math]\displaystyle{ p'+q'\,\equiv \, r' }[/math].
The Montgomery multiplication computes [math]\displaystyle{ c'\,=\,2^{-n}\,a'b'\bmod{m} }[/math]. Notice that [math]\displaystyle{ pq\,\equiv \, r }[/math] implies [math]\displaystyle{ p'q' \equiv \, r' }[/math]. The sections below show how to compute the Montgomery multipication.
In order to convert back from Montgomery representation to normal, just perform a Montgomery multiplication using the number 1 as the second factor.
Depending on the length of the factors, there are two ways the Montgomery multiplication can be applied.
Large factors
When the factors are very large and the Karatsuba multiplication or Fourier transform multiplication is faster than long multiplication, the division implied in the modular multiplication can be replaced by two multiplications.
First the negative of the modular inverse of [math]\displaystyle{ m \,\pmod{2^n} }[/math] is computed. Let's call this number [math]\displaystyle{ R }[/math]. Notice that this assumes that [math]\displaystyle{ m }[/math] is odd (which is almost always true).\n\nNow we perform the actual Montgomery multiplication [math]\displaystyle{ c' \,\equiv \,2^{-n}a'b' }[/math].
- [math]\displaystyle{ x = a'b' }[/math]
- [math]\displaystyle{ s = (x \,\bmod{2^n})R\,\bmod{2^n} }[/math]
- [math]\displaystyle{ t = (x + sm)/(2^n) }[/math]
- if [math]\displaystyle{ t \,\lt \, m }[/math] then [math]\displaystyle{ c' \,=\, t }[/math] else [math]\displaystyle{ c' \,=\, t - m }[/math].
Example
- Let m=5 , n=3 , a = b = 3.
- We get: a' = b' = 3 × 8 mod 5 = 4
- R = -5-1 mod 8 = 3
- x = a' b' = 4 × 4 = 16
- s = (16 mod 8) × 3 mod 8 = 0
- t = ((x + sm)/(2n) = (16 + 0 × 5) / 8 = 2
- c' = 2.
This is right since c = 9 implies that c' = 9 × 8 mod 5 = 2.
Small factors
When the factors are not large, so the long multiplication is faster, the division can be replaced by a little more than one multiplication.
In this method we will work with the multiprecision factors directly as two arrays of integers.
Let [math]\displaystyle{ a' = a'[0] \,+\, a'[1] * 2^k \,+\, a'[2] * 2^{2k} \,+\, ... \,+\, a'[s-1] * 2^{(s-1)k} }[/math].
and similarly for b' and c'. The number [math]\displaystyle{ k }[/math] is the number of bits of the integer, and [math]\displaystyle{ n = ks }[/math].
Let [math]\displaystyle{ R }[/math] be the modular inverse of [math]\displaystyle{ m\,\pmod{2^k} }[/math].
- [math]\displaystyle{ T = 0 }[/math]
- for i from 0 to s-1 perform the following:
- [math]\displaystyle{ T = T + a'[i] * b' * 2^{ki} }[/math]
- [math]\displaystyle{ h = T[i] * R \bmod{2^k} }[/math]
- [math]\displaystyle{ T = T + h * m * 2^{ki} }[/math]
- [math]\displaystyle{ T = T/2^n }[/math]
- if [math]\displaystyle{ T \,\lt \, m }[/math] then [math]\displaystyle{ c' \,=\, T }[/math] else [math]\displaystyle{ c' \,=\, T-m }[/math]