首页 > 代码库 > 蒙哥马利算法详解

蒙哥马利算法详解

这篇文章为大家梳理一下整个蒙哥马利算法的本质,蒙哥马利算法并不是一个独立的算法,而是三个相互独立又相互联系的算法集合,其中包括

  • 蒙哥马利乘模,是用来计算x?y (mod N)<script type="math/tex" id="MathJax-Element-1">x\cdot y\ (mod\ N)</script>
  • 蒙哥马利约减,是用来计算t?ρ?1 (mod N)<script type="math/tex" id="MathJax-Element-2">t\cdot \rho^{-1}\ (mod\ N)</script>
  • 蒙哥马利幂模,是用来计算xy (mod N)<script type="math/tex" id="MathJax-Element-3">x^y\ (mod\ N)</script>

其中蒙哥马利幂乘是RSA加密算法的核心部分。

基本概念

梳理几个概念,试想一个集合是整数模N之后得到的
ZN={0,1,2,?,N?1}<script type="math/tex" id="MathJax-Element-4">Z_N=\left\{0,1,2,\cdots,N-1\right\}</script>

注:N在base-b进制下有lN<script type="math/tex" id="MathJax-Element-5">l_N</script>位。 比如10进制和100进制,都属于base-10进制,因为100=102<script type="math/tex" id="MathJax-Element-6">100=10^2</script>,所以b=10。在10进制下,667的lN=3<script type="math/tex" id="MathJax-Element-7">l_N=3</script>

这样的集合叫做N的剩余类环,任何属于这个集合Z的x满足以下两个条件:
1. 正整数
2. 最大长度是lN<script type="math/tex" id="MathJax-Element-8">l_N</script>

这篇文章中讲到的蒙哥马利算法就是用来计算基于ZN<script type="math/tex" id="MathJax-Element-9">Z_N</script>集合上的运算,简单讲一下原因,因为RSA是基于大数运算的,通常是1024bit或2018bit,而我们的计算机不可能存储完整的大数,因为占空间太大,而且也没必要。因此,这种基于大数运算的加密体系在计算的时候都是基于ZN<script type="math/tex" id="MathJax-Element-10">Z_N</script>集合的,自然,蒙哥马利算法也是基于ZN<script type="math/tex" id="MathJax-Element-11">Z_N</script>。

在剩余类环上,有两种重要的运算,一类是简单运算,也就是加法和减法,另一类复杂运算,也就是乘法。我们比较熟悉的是自然数集上的运算,下面看下怎么从自然数集的运算演变成剩余类环上的运算。

对于加法运算,如果计算x±y (mod N)<script type="math/tex" id="MathJax-Element-12">x\pm y\ (mod\ N)</script> (0?x,y<N<script type="math/tex" id="MathJax-Element-13">0\leqslant x,yx±y<script type="math/tex" id="MathJax-Element-14">x\pm y</script>

0?x+y?2?(N?1)<script type="math/tex" id="MathJax-Element-15">\qquad 0\leqslant x+y\leqslant 2\cdot(N-1)</script>

?(N?1)?x?y?(N?1)<script type="math/tex" id="MathJax-Element-16">-(N-1)\leqslant x-y\leqslant (N-1)</script>

我们可以简单的通过加减N来实现从自然数到剩余类集的转换

另外一类是乘法操作,也就是x?y (mod N)<script type="math/tex" id="MathJax-Element-17">x\cdot y\ (mod\ N)</script>(0?x,y<N<script type="math/tex" id="MathJax-Element-18">0\leqslant x,y

0?x?y?(N?1)2<script type="math/tex" id="MathJax-Element-19">0\leqslant x\cdot y\leqslant (N-1)^2</script>

如果在自然数集下,令t=x?y<script type="math/tex" id="MathJax-Element-20">t=x\cdot y</script>,那么对于modN<script type="math/tex" id="MathJax-Element-21">\mod N</script>我们需要计算

t?N??tN?<script type="math/tex" id="MathJax-Element-22">t-(N\cdot \lfloor\frac{t}{N}\rfloor)</script>

加减操作很简单,具体的算这里就不细说了,我们用ZN?ADD<script type="math/tex" id="MathJax-Element-23">Z_N-ADD</script> 来代表剩余类环上的加法操作。既然我们可以做加法操作,那么我们就可以扩展到乘法操作,算法如下

技术分享

但是这并不是一个好的解决方案,因为通常来说,我们不会直接做w位乘w位的操作,这个后面会用蒙哥马利的乘法来代替解决。

对于取模操作,一般有以下几种方法

1,根据以下公式,来计算取模操作

t?N??tN?<script type="math/tex" id="MathJax-Element-24">t-(N\cdot \lfloor\frac{t}{N}\rfloor)</script>

这种解法有以下特征

  • 整个计算过程是基于标准的数字表示
  • 不需要预计算(也就是提前计算一些变量,以备使用)
  • 涉及到一个除法操作,非常费时和复杂

2,用Barrett reduction算法,这篇文章不细说,但是有以下特征

  • 基于标准的数字表示
  • 不需要预计算
  • 需要2?(lN+1)?(lN+1)<script type="math/tex" id="MathJax-Element-25">2 \cdot (l_N+1) \cdot (l_N+1)</script> 次数乘运算

3,用蒙哥马利约减,也就是下面要讲的算法,有以下特征

  • 不是基于标准的数字表示(后文中有提到,是基于蒙哥马利表示法)
  • 需要预计算
  • 需要2?(lN)?(lN)<script type="math/tex" id="MathJax-Element-26">2 \cdot (l_N) \cdot (l_N)</script> 次数乘运算

蒙哥马利预备知识

在将蒙哥马利算法之前,先看一下在自然数下的乘法公式

计算x?y<script type="math/tex" id="MathJax-Element-932">x\cdot y</script>,想象一下我们常用的计算乘法的方法,用乘数的每一位乘上被乘数,然后把得到的结果相加,总结成公式,可以写成如下的形式。

x?y=x?sumly?1i=0yi?bi<script type="math/tex" id="MathJax-Element-933">x\cdot y=x\cdot sum_{i=0}^{l_y-1}y_i \cdot b^i</script>

=sumly?1i=0yi?x?bi<script type="math/tex" id="MathJax-Element-934">\qquad=sum_{i=0}^{l_y-1}y_i \cdot x \cdot b^i</script>

尝试下面一个例子,10进制下(也就是b=10),y=456(也就是ln=3<script type="math/tex" id="MathJax-Element-935">l_n=3</script>),计算x?y<script type="math/tex" id="MathJax-Element-936">x\cdot y</script>,公式可演变如下:

x?y=(y0?x?100)+(y1?x?101)+(y2?x?102)<script type="math/tex" id="MathJax-Element-937">x\cdot y=(y_{0}\cdot x\cdot 10^{0})+(y_{1}\cdot x\cdot 10^{1})+(y_{2}\cdot x\cdot 10^{2})</script>
=(y0?x?0)+(y1?x?10)+(y2?x?100)<script type="math/tex" id="MathJax-Element-938">\qquad=(y_{0}\cdot x\cdot 0)+(y_{1}\cdot x\cdot 10)+(y_{2}\cdot x\cdot 100)</script>
=(y0?x)+10?(y1?x+10?(y2?x?+10?(0)))<script type="math/tex" id="MathJax-Element-939">\qquad=(y_{0}\cdot x)+10\cdot(y_{1}\cdot x+10\cdot(y_{2}\cdot x\cdot +10\cdot(0)))</script>

最后一次演变其实就是用霍纳法则(Horner’s rule)所讲的规则,关于霍纳法则,可以自行百度。

这个计算过程写成代码实现的算法应该是这样的:
技术分享

接下来我们来看下面这样的计算,计算(x?y)/1000<script type="math/tex" id="MathJax-Element-940">(x\cdot y)/1000</script>,由前面可以知道

x?y=(y_0?x)+10?(y_1?x+10?(y_2?x?+10?(0)))<script type="math/tex" id="MathJax-Element-941">x\cdot y=(y\_{0}\cdot x)+10\cdot(y\_{1}\cdot x+10\cdot(y\_{2}\cdot x\cdot +10\cdot(0)))</script>

由此可以知道:

x?y1000=(y0?x?100)+(y1?x?101)+(y2?x?102)1000<script type="math/tex" id="MathJax-Element-942">\frac{x\cdot y}{1000}=\frac{(y_{0}\cdot x\cdot 10^{0})+(y_{1}\cdot x\cdot 10^{1})+(y_{2}\cdot x\cdot 10^{2})}{1000}</script>

=(y0?x?0)+(y1?x?10)+(y2?x?100)1000<script type="math/tex" id="MathJax-Element-943">\qquad=\frac{(y_{0}\cdot x\cdot 0)+(y_{1}\cdot x\cdot 10)+(y_{2}\cdot x\cdot 100)}{1000}</script>

=(y0?x)1000+(y1?x)100+(y2?x)10<script type="math/tex" id="MathJax-Element-944">\qquad=\frac{(y_{0}\cdot x)}{1000}+\frac{(y_{1}\cdot x)}{100}+\frac{(y_{2}\cdot x)}{10}</script>

=(((((y0?x)/10)+y1?x)/10)+y2?x)/10<script type="math/tex" id="MathJax-Element-945">\qquad=(((((y_0\cdot x)/10)+y_1\cdot x)/10)+y_2\cdot x)/10</script>

这个计算过程写成代码实现的算法是这样的:
技术分享

接下来我们再来看在剩余类集合下的乘法操作 x?y/1000 (mod 667)<script type="math/tex" id="MathJax-Element-946">x\cdot y/1000\ (mod\ 667)</script>

我们知道剩余类集合Z667={0,1?666}<script type="math/tex" id="MathJax-Element-947">Z_{667}=\left\{0,1 \cdots 666\right\}</script>,是不存在小数的,而如果我们采用自然数集的计算方式的话,就会出现小数,比如前面的例子,除10就会有小数。

这个问题是这样的,我们知道u·6670(mod667)<script type="math/tex" id="MathJax-Element-948">u·667 \equiv 0 (mod 667)</script>(<script type="math/tex" id="MathJax-Element-949">\equiv</script>表示取模相等),所以我们可以选择一个合适的u,用u乘667再加上r,使得和是一个可以除10没有小数,这样在mod 667之后依然是正确的结果。至于u怎么算出来,这篇文章会在后面的章节说明。

这个过程之后x?y/1000 (mod 667)<script type="math/tex" id="MathJax-Element-950">x\cdot y/1000\ (mod\ 667)</script> 的计算算法可以写成如下的形式
技术分享

至此,你可能还不明白上面说这一堆演变的原因,其实很简单,原来是一个(x?y) (mod 667)<script type="math/tex" id="MathJax-Element-951">(x\cdot y)\ (mod\ 667)</script>的运算,这个运算中的模操作,正常情况下是要通过除法实现的,而除法是一个特别复杂的运算,要涉及到很多乘法,所以在大数运算时,我们要尽量避免除法的出现。而通过以上几个步骤,我们发现(x?y)/1000 (mod 667)<script type="math/tex" id="MathJax-Element-952">(x\cdot y)/1000\ (mod\ 667)</script>这个操作是不用除法的。等等,算法中明明有个除10的操作,你骗谁呢。不知道你有没有发现,除数其实是我们的进制数,除进制数在计算机中是怎么做呢,其实很简单,左移操作就ok了。所以这个计算方法是不涉及到除法操作的。

但是我们要计算的明明是(x1?y1) (mod 667)<script type="math/tex" id="MathJax-Element-953">(x_1\cdot y_1)\ (mod\ 667)</script>,怎么现在变成了(x2?y2)/1000 (mod 667)<script type="math/tex" id="MathJax-Element-954">(x_2\cdot y_2)/1000\ (mod\ 667)</script>,所以在下一步,我们要思考的是怎么样让(x1?y1) (mod 667)<script type="math/tex" id="MathJax-Element-955">(x_1\cdot y_1)\ (mod\ 667)</script>转变成(x2?y2)/1000 (mod 667)<script type="math/tex" id="MathJax-Element-956">(x_2\cdot y_2)/1000\ (mod\ 667)</script>这种形式。

考虑这样两个算法
- 第一个是输入x和y,计算x?y?ρ?1<script type="math/tex" id="MathJax-Element-957">x \cdot y \cdot \rho^{-1}</script>
- 第二个算法,输入一个t,计算t?ρ?1<script type="math/tex" id="MathJax-Element-958">t \cdot \rho^{-1}</script>。

x?y (mod 667)=((x?1000)?(y?1000)/1000)/1000 (mod 667)<script type="math/tex" id="MathJax-Element-959">x\cdot y\ (mod\ 667)=((x\cdot1000)\cdot(y\cdot1000)/1000)/1000\ (mod\ 667)</script>

是不是变成了我们需要的(x?y)/1000 (mod 667)<script type="math/tex" id="MathJax-Element-960">(x\cdot y)/1000\ (mod\ 667)</script>模式,而且这个转变过程是不是可以通过上面两个算法来实现,输入值如果是(x?1000)<script type="math/tex" id="MathJax-Element-961">(x\cdot1000)</script>和(y?1000)<script type="math/tex" id="MathJax-Element-962">(y\cdot1000)</script>,则通过第一个算法可以得到((x?1000)?(y?1000)/1000)<script type="math/tex" id="MathJax-Element-963">((x\cdot1000)\cdot(y\cdot1000)/1000)</script>,把结果作为第二个算法的输入值,则通过第二个算法可以得到((x?1000)?(y?1000)/1000)/1000<script type="math/tex" id="MathJax-Element-964">((x\cdot1000)\cdot(y\cdot1000)/1000)/1000</script>。

扯了一大顿,终于引出了今天文章的主角,前面讲到的两个算法,第一个就是蒙哥马利乘模,第二个就是蒙哥马利约减。下面我们来讲这两个算法的详解。

正如前面提到的蒙哥马利算法的三个特性之一是,不是基于普通的整数表示法,而是基于蒙哥马利表示法。什么是蒙哥马利表示法呢,其实也很简单,上面我们提到,要让(x1?y1) (mod 667)<script type="math/tex" id="MathJax-Element-965">(x_1\cdot y_1)\ (mod\ 667)</script>转变成(x2?y2)/1000 (mod 667)<script type="math/tex" id="MathJax-Element-966">(x_2\cdot y_2)/1000\ (mod\ 667)</script>这种形式,我们需要将输入参数变成(x?1000)<script type="math/tex" id="MathJax-Element-967">(x\cdot1000)</script>和(y?1000)<script type="math/tex" id="MathJax-Element-968">(y\cdot1000)</script>,而不是x和y本身,而(x?1000)<script type="math/tex" id="MathJax-Element-969">(x\cdot1000)</script>和(y?1000)<script type="math/tex" id="MathJax-Element-970">(y\cdot1000)</script> 其实就是蒙哥马利表示法。

所以我们先定义几个概念:

  • 蒙哥马利参数
    给定一个N,N在b进制(例如,二进制时,b=2)下共有l位,gcd(N,b)=1<script type="math/tex" id="MathJax-Element-971">gcd(N,b)=1</script>,先预计算以下几个值(这就是前面提到的特性之一,需要预计算):
    • ρ=bk<script type="math/tex" id="MathJax-Element-972">\rho = b^k</script> 指定一个最小的k,使得bk>N<script type="math/tex" id="MathJax-Element-973">b^k>N</script>
    • ω=?N?1(mod ρ)<script type="math/tex" id="MathJax-Element-974">\omega = -N^{-1} (mod\ \rho)</script>
      这两个参数是做什么用的呢,你对照前面的演变过程可以猜到ρ<script type="math/tex" id="MathJax-Element-975">\rho </script> 就是前面演变中的1000,而ω<script type="math/tex" id="MathJax-Element-976">\omega</script> 则是用来计算前面提到的u的。
  • 蒙哥马利表示法
    对于x,0?x?N?1<script type="math/tex" id="MathJax-Element-977">0\leqslant x\leqslant N-1</script>,x的蒙哥马利表示法表示为x=x?ρ (mod N)<script type="math/tex" id="MathJax-Element-978">x=x\cdot \rho\ (mod\ N)</script>

蒙哥马利约减

蒙哥马利约减的定义如下
给定一些整数t,蒙哥马利约减的计算结果是t?ρ?1 (mod N)<script type="math/tex" id="MathJax-Element-74">t\cdot \rho^{-1}\ (mod\ N)</script>

蒙哥马利约减的算法可表示为
技术分享

蒙哥马利约减可以算作是下面要说的蒙哥马利模乘当x=1<script type="math/tex" id="MathJax-Element-75">x=1</script>时的一种特殊形式,。同时它又是蒙哥马利乘模要用到的一部分,这在下一部分讲蒙哥马利乘模的时候有讲到。

蒙哥马利约减可以用来计算某个值得取模操作,比如我们要计算m(mod N)<script type="math/tex" id="MathJax-Element-76">m(mod\ N)</script>,只需要将m
的蒙哥马利表示法m?ρ<script type="math/tex" id="MathJax-Element-77">m\cdot \rho</script>作为参数,带入蒙哥马利约减,则计算结果就是m(mod N)<script type="math/tex" id="MathJax-Element-78">m(mod\ N)</script>。

蒙哥马利乘模

一个蒙哥马利乘模包括整数乘法和蒙哥马利约减,现在我们有蒙哥马利表示法:

x=x?ρ (mod N)<script type="math/tex" id="MathJax-Element-79"> \hat{x}=x\cdot\rho\ (mod\ N)</script>
y=y?ρ (mod N)<script type="math/tex" id="MathJax-Element-80"> \hat{y}=y\cdot\rho\ (mod\ N)</script>

它们相乘的结果是

t=x?y<script type="math/tex" id="MathJax-Element-81">t=\hat{x}\cdot\hat{y}</script>
 =(x?ρ)?(y?ρ)<script type="math/tex" id="MathJax-Element-82">\ =(x\cdot\rho)\cdot(y\cdot\rho)</script>
 =(x?y)?ρ2<script type="math/tex" id="MathJax-Element-83">\ =(x\cdot y)\cdot\rho^2</script>

最后,用一次蒙哥马利约减得到结果

t=(x?y)?ρ (mod N)<script type="math/tex" id="MathJax-Element-84"> \hat{t}=(x \cdot y) \cdot \rho\ (mod\ N)</script>

上面我们可以看出,给出的输入参数是x<script type="math/tex" id="MathJax-Element-85"> \hat{x}</script> 和y<script type="math/tex" id="MathJax-Element-86"> \hat{y}</script>, 得到的结果是(x?y)?ρ (mod N)<script type="math/tex" id="MathJax-Element-87">(x \cdot y) \cdot \rho\ (mod\ N)</script>,所以蒙哥马利乘法也可以写成如下形式,已知输入参数x和y,蒙哥马利乘法是计算(x?y)?ρ?1 (mod N)<script type="math/tex" id="MathJax-Element-88">(x \cdot y) \cdot \rho ^ {-1}\ (mod\ N)</script>

举个例子:
b=10,也就是说在10进制下,N=667
bk>N<script type="math/tex" id="MathJax-Element-89">b^k>N</script>的最小的k是3,所以ρ=bk=103=1000<script type="math/tex" id="MathJax-Element-90">\rho=b^k=10^3=1000</script>
ω=?N?1 (mod ρ)=?667?1 (mod ρ)=997<script type="math/tex" id="MathJax-Element-91">\omega=-N^{-1}\ (mod\ \rho)=-667^{-1}\ (mod\ \rho)=997</script>

因为x=421<script type="math/tex" id="MathJax-Element-92">x=421</script>,所以x=x?ρ(mod N)=421?1000(mod 667)=123<script type="math/tex" id="MathJax-Element-93">\hat{x}=x\cdot\rho(mod\ N)=421\cdot1000(mod\ 667)=123</script>
因为y=422<script type="math/tex" id="MathJax-Element-94">y=422</script>,所以y=y?ρ(mod N)=422?1000(mod 667)=456<script type="math/tex" id="MathJax-Element-95">\hat{y}=y\cdot\rho(mod\ N)=422\cdot1000(mod\ 667)=456</script>

所以计算x<script type="math/tex" id="MathJax-Element-96">\hat{x}</script>和y<script type="math/tex" id="MathJax-Element-97">\hat{y}</script>蒙哥马利乘结果是

x?y?ρ?1=(421?1000?422?1000)?1000?1 (mod 667)<script type="math/tex" id="MathJax-Element-98">\hat{x}\cdot\hat{y}\cdot\rho^{-1}=(421\cdot1000\cdot422\cdot1000)\cdot1000^{-1}\ (mod\ 667)</script>
(421?422)?1000 (mod 667)<script type="math/tex" id="MathJax-Element-99">\qquad\qquad(421\cdot422)\cdot1000\ (mod\ 667)</script>
(240)?1000 (mod 667)<script type="math/tex" id="MathJax-Element-100">\qquad\qquad(240)\cdot1000\ (mod\ 667)</script>
547 (mod 667)<script type="math/tex" id="MathJax-Element-101">\qquad\qquad547\ (mod\ 667)</script>

然后总结一下蒙哥马利约减和蒙哥马利乘法的伪代码实现,这个算法其实就是从蒙哥马利预备知识讲到的算法演变来的。
技术分享

上面的例子用这个算法可以描述为
技术分享

蒙哥马利算法是一套很完美的算法,为什么这么说呢,你看一开始已知x<script type="math/tex" id="MathJax-Element-102">x</script>,我们要求x=x?ρ<script type="math/tex" id="MathJax-Element-103">\hat{x}=x \cdot \rho</script>,这个过程可以通过蒙哥马利乘法本身来计算,输入参数x<script type="math/tex" id="MathJax-Element-104">x</script>和ρ2<script type="math/tex" id="MathJax-Element-105">\rho^2</script>,计算结果就是x=x?ρ<script type="math/tex" id="MathJax-Element-106">\hat{x}=x \cdot \rho</script>。然后在最后,我们知道x=x?ρ<script type="math/tex" id="MathJax-Element-107">\hat{x}=x \cdot \rho</script>,要求得x<script type="math/tex" id="MathJax-Element-108">x</script>的时候,同样可以通过蒙哥马利算法本身计算,输入参数x<script type="math/tex" id="MathJax-Element-109">\hat{x}</script>和1<script type="math/tex" id="MathJax-Element-110">1</script>,计算结果就是x<script type="math/tex" id="MathJax-Element-111">x</script>。有没有一种因就是果,果就是因的感觉,这就是为什么说蒙哥马利算法是一套很完美的算法。

蒙哥马利幂模

最后,才说到我们最开始提到的RSA的核心幂模运算,先来看一下普通幂运算的算法是怎么得出来的。

以下资料来自于百度百科快速模幂运算
针对快速模幂运算这一课题,西方现代数学家提出了大量的解决方案,通常都是先将幂模运算转化为乘模运算。
例如求D=C^15%N
由于:a*b % n = (a % n)*(b % n) % n
所以令:
C1 =C*C % N =C^2 % N
C2 =C1*C % N =C^3 % N
C3 =C2*C2 % N =C^6 % N
C4 =C3*C % N =C^7 % N
C5 =C4*C4 % N =C^14 % N
C6 =C5*C % N =C^15 % N
即:对于E=15的幂模运算可分解为6 个乘模运算,归纳分析以上方法可以发现:
对于任意指数E,都可采用以下算法计算D=C**E % N:
D=1
WHILE E>0
IF E%2=0
C=C*C % N
E=E/2
ELSE
D=D*C % N
E=E-1
RETURN D
继续分析会发现,要知道E 何时能整除 2,并不需要反复进行减一或除二的操作,只需验证E 的二进制各位是0 还是1 就可以了,从左至右或从右至左验证都可以,从左至右会更简洁,
设E=Sum[i=0 to n](E*2**i),0<=E<=1
则:
D=1
FOR i=n TO 0
D=D*D % N
IF E=1
D=D*C % N
RETURN D这样,模幂运算就转化成了一系列的模乘运算。

算法可以写成如下的形式
技术分享

如果我们现在用蒙哥马利样式稍作改变,就可以变成如下的形式:
技术分享

以上就是蒙哥马利算法的全部,通过蒙哥马利算法中的约减运算,我们将大数运算中的模运算变成了移位操作,极大地提高了大数模乘的效率。

但是在以上的算法,可以发现还有两个变量的计算方式不是很清楚,一个是ω<script type="math/tex" id="MathJax-Element-112">\omega</script>,前面说过ω=?N?1(modN)<script type="math/tex" id="MathJax-Element-113">\omega = -N^{-1} (mod N)</script> ,其实在算法中,我们看到,omega<script type="math/tex" id="MathJax-Element-114">omega</script>仅仅被用来做modb<script type="math/tex" id="MathJax-Element-115">\mod b</script>操作,所以事实上,我们只需要计算modb<script type="math/tex" id="MathJax-Element-116">\mod b</script>即可。

尽管N有可能是合数(因为两个素数的乘积不一定是素数),但通常N和ρ<script type="math/tex" id="MathJax-Element-117">\rho</script>(也就是N和b)是互质的,也就是说N?(b)=1(mob b)<script type="math/tex" id="MathJax-Element-118">N^{\phi(b)}=1(mob\ b)</script>(费马定理),N?(b)?1=N?1(mob b)<script type="math/tex" id="MathJax-Element-119">N^{\phi(b)-1}=N^{-1}(mob\ b)</script>
因为b=2ω<script type="math/tex" id="MathJax-Element-120">b=2^\omega</script>,所以?(b)=2(ω?1)<script type="math/tex" id="MathJax-Element-121">\phi(b)=2^{(\omega-1)}</script>,写成算法是这样的
技术分享

还有一个参数是ρ2<script type="math/tex" id="MathJax-Element-122">\rho^2</script>,还记得前面说过ρ<script type="math/tex" id="MathJax-Element-123">\rho</script>是怎么得出来吗,选定一个最小的k<script type="math/tex" id="MathJax-Element-124">k</script>,使得bk>N<script type="math/tex" id="MathJax-Element-125">b^k>N</script>,我们还知道N<script type="math/tex" id="MathJax-Element-126">N</script>在b<script type="math/tex" id="MathJax-Element-127">b</script>进制下是lN<script type="math/tex" id="MathJax-Element-128">l_N</script>位,所以当k=lN<script type="math/tex" id="MathJax-Element-129">k=l_N</script>的时候肯定是符合要求。

b=2ω<script type="math/tex" id="MathJax-Element-130">b=2^{\omega}</script> 所以ρ=bk=(2ω)k<script type="math/tex" id="MathJax-Element-131">\rho=b^k=({2^{\omega}})^k</script>

ρ2=(2w)k)2=22?k?ω=22?lN?ω<script type="math/tex" id="MathJax-Element-132">\rho^2={({2^w})^k)}^2=2^{2\cdot k\cdot \omega}=2^{2\cdot l_N\cdot \omega}</script>,算法如下

技术分享

至此整个蒙哥马利算法就全部说完了。通过这个算法,我们可以实现快速幂模。

<script type="text/javascript"> $(function () { $(‘pre.prettyprint code‘).each(function () { var lines = $(this).text().split(‘\n‘).length; var $numbering = $(‘
    ‘).addClass(‘pre-numbering‘).hide(); $(this).addClass(‘has-numbering‘).parent().append($numbering); for (i = 1; i <= lines; i++) { $numbering.append($(‘
  • ‘).text(i)); }; $numbering.fadeIn(1700); }); }); </script>

    蒙哥马利算法详解