MbedTLS中的Montgomery算法实现解析（一）

MbedTLS中的Montgomery算法实现解析（一）

蒙哥马利算法模乘介绍

蒙哥马利模乘算法主要是为了进行大数运算a×b mod n，并且针对CPU进行了并行优化，来快速计算模乘。

在介绍蒙哥马利模乘之前，先让我们来了解蒙哥马利约减（Montgomery reduction）

蒙哥马利约减（Montgomery reduction）

首先，让我们来考虑一个简单的问题：

a mod n

如果a是一个2048位的整数，n是一个1024位的整数，那我们要如何来计算结果呢？

如果我们直接采用相除的方式，不论在时间还是在空间上都会产生非常大的开销，这是非常不划算的，因此我们需要采用一些手段来加速我们的运算。

在这里，我们考虑这样一种情况：

a + sn ≡ a mod n

这条等式是显然成立的，因为我们加上的是一个n的倍数，这不会对模除结果产生任何影响。

那么这里我们如果能够控制s的值，使a + sn恰好在低的1024位上全部为0，就能很好的计算a mod n，因为如果出现了这种情况就意味着两件事情：

1. 我们可以通过移位的手段来将a这个2048位的数字向右移动1024位，这就使得计算重新回到了2¹⁰²⁴的完全剩余系下。
2. CPU执行移位速度要比进行除法快得多，这能大大加速计算速度。

下面，我们就需要来求出这样一个s，令a + sn恰好在低的1024位上全部为0，也就是说

a+sn=k₁ × 2¹⁰²⁴

这样的s其实并不难求，首先我们将a拆成两个部分a₁和a₀，高的1024位作为a₁，低的1024位作为a₀，那么显然

a = a₁ × 2¹⁰²⁴ + a₀

所以

a+sn=a₁ × 2¹⁰²⁴ + a₀ + sn = k₁ × 2¹⁰²⁴

a₀ + sn = (k₁ - a₁) × 2¹⁰²⁴

我们记(k₁ - a₁)为k₂，那么显然：

a₀ + sn ≡ 0 (mod 2¹⁰²⁴)

a₀ ≡ - sn (mod 2¹⁰²⁴)

a₀ × n^-1 ≡ - sn × n^-1 (mod 2¹⁰²⁴)

[注]：这里的n^-1是n模2¹⁰²⁴的逆

a₀ × n^-1 ≡ -s (mod 2¹⁰²⁴)

s ≡ - a₀ × n^-1 (mod 2¹⁰²⁴)

在这里我们已经计算出了s，这个s能够使a+sn恰好在低的1024位上全部为0，下面我们将a+sn向右移动1024位，这相当于除以2¹⁰²⁴，得到的结果我们记为t。下面为了方便我们记2¹⁰²⁴为R，这里R只要满足2^m同时大于n即可，m并不一定是1024。

在这里我们先说结论t ≡ a × R^-1 (mod n)，这条算式就是蒙哥马利约减，下面我们来证明一下：

t = (a + sn) / R

t × R = a + sn

所以

tR ≡ a ( mod n )

t ≡ a × R^-1 (mod n)

蒙哥马利表示（Montgomery representation）

现在读者可能会问，为什么我们要大费周章的计算t？，主要的原因是我们在计算t的过程中，只会用到一次加法a+sn，一次乘法s×n，以及一次移位运算，这就计算出了t ≡ a × R^-1 (mod n)，这比进行除法要快得多得多。
至于n^-1，我们可以在预处理中计算出来，依然不会造成太大的计算量。

下面我们要介绍蒙哥马利表示，蒙哥马利表示也非常简单：

A = a × R mod n
B = b × R mod n

蒙哥马利模乘（Montgomery multiplication）

我们现在比如要计算一个X= a × b (mod n)，那么我们如何利用表示和约减来计算模乘的结果呢？

首先：

A = a × R (mod n)
B = b × R (mod n)

然后将A × B的结果计算一次约减：

Y = A × B × R^-1 (mod n)
Y = a × R × b × R × R^-1 (mod n)
Y = a × R × b (mod n)

最后再将Y进行一次约减：

X = Y × R^-1 (mod n)
X = a × R × b × R^-1 (mod n)
X= a × b (mod n)

实际上，除此之外，蒙哥马利还定义了加法运算X= a + b (mod n)

首先直接相加：

Y = A + B (mod n)
Y = a × R + b × R (mod n)
Y = (a + b) × R (mod n)

然后将Y进行一次约减：

X = Y × R^-1 (mod n)
X= a + b (mod n)

也就是说，乘法要经历两次约减，加法只用经历一次约减。