转:http://blog.chinaunix.net/uid-9162199-id-111828.html
蒙哥马利对自然是的快速求幂的原理在北邮人论坛上面找到了,通俗易懂。
"
快速幂运算的好处是减少了运算量,极大地提高了速度。例如A^65535(A的65535次幂),原始算法要做65535-1=65534次乘法,而快速幂运算只需要做(16-1)×2=30次乘法。
原理其实很好懂,假设要计算A^B,即底数是A,指数是B。把B写成二进制形式,拿4位来举例:B=b4b3b2b1(二进制)。
先用B=1111(二进制)来做解释。显然A^B = A × A^2 × A^4 × A^8。
又显然,
A^2 = A × A,
A^4 = A^2 × A^2,
A^8 = A^4 × A^4。
假如B的二进制位数更多,则依此类推。
上面这段看懂了吗?如果看懂的,就应该能够写出B=1111(二进制)情况下,A^B的快速幂运算程序。
然后再来看B=b4b3b2b1(二进制)的情况。这时候A^B = A1 × A2 × A3 × A4。
其中,
A1 = A(如果b1=1)或者1(如果b1=0),
A2 = A^2(如果b2=1)或者1(如果b2=0),
A3 = A^4(如果b3=1)或者1(如果b3=0),
A4 = A^8(如果b4=1)或者1(如果b4=0)。
假如B的二进制位数更多,则依此类推。
这一段看懂了吗?如果看懂的,我不用继续解释,能不能把上面B=1111(二进制)的程序改一改,加入对b4、b3、b2、b1等于0或者1的条件分支判断,让它变成A^B的快速幂运算程序?
"
根据上面的解释我实现了算法:
#include
/* a^b蒙哥马利快速幂运算算法 */
long long power(long long a, long long b)
{
if (b == 0)
return 1;
long long total = 1;
while (b != 0)
{
if ((b & 1) != 0)
{
total *= a;
}
a *= a;
b >>= 1;
}
return total;
}
int main()
{
long long n, t, k;
while (scanf("%lld%lld", &n, &t) != 2)
{
k = power(n, t);
printf("%lld\n", k);
}
return 0;
}
解释如下(没看懂,数学功底不够):
“
如果直接计算C=A*B mod N,必须做很费时的mod N运算。从我前面讲的快速幂运算算法可知,乘法运算要进行多次,所以mod N运算也需要多次,比较慢。
假如我们做个转换,A-->A',B-->B',然后计算C'=A'*B' mod N,这时mod N运算可以快得多(当然这个转换是需要技巧的)。等到所有乘法运算都做完之后,最后结果再做一次转换,C'-->C,就转回到原本需要的结果了。
这个跟快速乘法运算有点类似:先用快速傅立叶变换,A-->A',B-->B',然后乘法变加法,C'=A'+B',最后来个傅立叶逆变换,C'-->C,那么得到的C=A*B。加法当然比乘法快得多。
然后就是蒙哥马利算法之中,那个神奇的转换了:
A'=A*2^n mod N,
B'=B*2^n mod N,
显然,乘以2^n就是左移n位。至于mod N嘛……我相信n足够大,并且取某些特定值的时候,mod N可以变成不是除法的运算。例如,变成少量乘法、加减法。
由于N是奇数(通常还是大素数或者两个素数的乘积),所以2^n与N的最大公约数是1。2^n有乘法逆元(mod N域内),表示成2^-n(mod N),实际上逆元也是一个正整数,2^n乘这个整数会得1(mod N)。求逆元的算法不讨论了,欧几里得算法是最常规的一种。
显然,上述逆转换就是:
C=C'*2^n-1 mod N
我是现炒现卖的,原理大致是这些。具体的实现细节肯定还有很多问题,例如n值怎样选取,我就没有查到。貌似这里就衍生出不同的技巧性算法。
”
本站仅提供存储服务,所有内容均由用户发布,如发现有害或侵权内容,请
点击举报。