算法原理
模重复平方算法是用来快速计算 b n m o d m b^n \mod m bnmodm的一个算法。
考虑直接计算 b n m o d m b^n \mod m bnmodm,需要 n − 1 n-1 n−1次乘法,也就是递归计算 b n ≡ ( b n − 1 m o d m ) ⋅ b m o d m . b^n \equiv (b^{n-1} \mod m) \cdot b \mod m. bn≡(bn−1modm)⋅bmodm.不过,当 n n n很大的时候,计算会非常耗时。
现在,考虑将 n n n的二进制表示 n = n 0 + n 1 2 + n 2 2 2 + ⋯ + n k − 1 2 k − 1 n=n_0+n_12+n_22^2+\cdots +n_{k-1}2^{k-1} n=n0+n12+n222+⋯+nk−12k−1,其中 n 0 , n 1 , ⋯ , n k − 1 n_0,n_1,\cdots,n_{k-1} n0,n1,⋯,nk−1为0或者1.
则:
b n ≡ b n 0 + n 1 2 + n 2 2 2 + ⋯ + n k − 1 2 k − 1 m o d m ≡ b n 0 ⋅ ( b 2 ) n 1 ⋅ ( ( b 2 ) 2 ) n 2 ⋯ ( b 2 k − 1 ) n k − 1 m o d m \begin{aligned} b^n& \equiv b^{n_0+n_12+n_22^2+\cdots +n_{k-1}2^{k-1}} \mod m\\ &\equiv b^{n_0} \cdot (b^2)^{n_1}\cdot ((b^2)^2)^{n_2} \cdots (b^{2^{k-1}})^{n_{k-1}} \mod m \end{aligned} bn≡bn0+n12+n222+⋯+nk−12k−1modm≡bn0⋅(b2)n1⋅((b2)2)n2⋯(b2k−1)nk−1modm
通过观察,我们可以知道,我们需要 b b b的偶次幂。设 b 0 = b b_0=b b0=b的话,下一个偶次幂 b 1 = b 0 2 m o d m b_1=b_0^2 \mod m b1=b02modm,以此类推, b 2 = b 1 2 m o d m , ⋯ , b k − 1 = b k − 2 2 m o d m b_2=b_1^2 \mod m, \cdots ,b_{k-1}=b_{k-2}^2 \mod m b2=b12modm,⋯,bk−1=bk−22modm,然后代入上面的公式,有:
b n ≡ b 0 n 0 ⋅ b 1 n 1 ⋅ b 2 n 2 ⋯ b k − 1 n k − 1 m o d m . b^n \equiv b_0^{n_0} \cdot b_1^{n_1} \cdot b_2^{n_2} \cdots b_{k-1}^{n_{k-1}} \mod m. bn≡b0n0⋅b1n1⋅b2n2⋯bk−1nk−1modm.
又由于 n 0 , n 1 , ⋯ , n k − 1 n_0,n_1,\cdots,n_{k-1} n0,n1,⋯,nk−1为0或者1,所以,我们只需要通过不断平方计算出 b i b_i bi,令初始结果为1,然后,如果 n i n_i ni为1则乘上 b i m o d m b_i \mod m bimodm,否则,不管,继续计算 b i + 1 b_{i+1} bi+1.
可以看出,通过上面的计算,我们只需要最多 2 [ l o g 2 n ] 2[log_2n] 2[log2n]次乘法,和一次 n n n的分解。
c++实现代码
#include<iostream>
using namespace std;
int msa(int b,int n,int m){
int b1=b;
int res=1;
while(n>0){
if(n%2==1){
res=res*b1%m;
}
b1=b1*b1%m;
n=n/2;
}
return res;
}
int main(){
int b,n,m;
cout<<"Input b:"<<endl;
cin>>b;
cout<<"Input n:"<<endl;
cin>>n;
cout<<"Input m:"<<endl;
cin>>m;
cout<<"The result of b^n mod m is:"<<endl;
cout<<msa(b,n,m)<<endl;
}
gmp函数
在gmp库中,实现了大整数的幂运算,其函数原型如下:
void mpz_powm (mpz_t rop, const mpz_t base, const mpz_t exp, const mpz_t mod)
rop是返回的结果。