第31章 数论算法：元素的幂以及模取幂

元素的幂：

在这里介绍求元素的幂X^n（X和n均为非负整数）的两种算法，一个是递归求解，另外一个是用反复平方法

递归求解：

当我们计算X^n时，当n=0，结果为1，这可以看做递归的基准情况，当n为偶数时，X^n=(X^2)^(N/2)，当n是基数时，X^n=(X^2)^((n-1)/2)*X。代码如下：

bool isEven(unsigned int n)
{
if(n%2==0) return true;
else return false;
}
unsigned long pow(unsigned long x,unsigned int n)
{
if(n==0) return 1;//基准情况
if(isEven(n))   return pow(x*x,n/2);
else return pow(x*x,n/2)*x;
}

不难得出，递归求解元素幂的时间复杂度为O（log n).

反复平方法：

我们是要求X^n，假设n可以用二进制表述出来，即n(k-1)…n(0)（n(i)=0或1），则n=2^(k-1)* n(k-1)+…+2^0* n(0);

这样的话，X^n=X^(2^0* n(0))* X^(2^1* n(1))…. *X^(2^k n(k))

我们迭代k次，即可得出X^n。不难得出k=O(log n)，所以反复平方法的时间复杂度为O（log n）。代码如下：

unsigned long pow(unsigned long x,unsigned int n)
{
unsigned long ret=1;
unsigned long tmp=x;

while(n!=0){
if(n%2!=0)
ret*=tmp;
tmp=tmp*tmp;
n/=2;
}

return ret;
}

模取幂：

模取幂是计算(a^b)%c这类问题（a，b 为非负整数，c为正整数），在这里介绍三种计算模取幂的三种算法。

第一种算法：

模运算具有这样的性质：(m*n)%d= ((m%d)*n)%d

所以a^b%c可写成((((((1*a)%c)*a)%c)*a)%c…*a)%c（有b个a）。因此依据这个等式，可以写成如下代码：

unsigned long  mod_1(unsigned long a, unsigned long b, unsigned long c)
{
unsigned long ret=1;

while(b--)
ret=(ret*a)%c;

return ret;
}

第二种算法：

模运算除了上面说的具有(m*n)%d= ((m%d)*n)%d这样的性质外，还具有这样的性质：当N=N1∗N2∗N3∗...∗Nn,(N mod m)=[1∗(N1 mod m)∗(N2 mod m)∗(N3 mod m)∗...(Nn mod m)] mod m

把上面两个性质结合在一起的话，N mod m可以写成如下形式：

(Nmodm)=((((((1)∗(N1 mod m)) mod m∗(N2 mod m))) mod m∗(N3 mod m)) mod m...∗(Nn mod m)) mod m

并且

if Ni=Ni−1∗Ni−1,(Ni mod m)={[Ni−1 mod m]∗[Ni−1 mod m]} mod m

当我们计算a^b%c时，b可以写成如下的形式：

b=p(0)*2^0+ p(1)*2^1+…p(i)*2^i+…+p(n-1)*2^n-1 (p(i)=0 or 1)

则a^b可以写成如下形式：

a^b=a^(p(0) * 2^0))* ..* a^(p(i) * 2^i) * … * a^(p(n-1)*2^n-1)。因此求解a^b%c的代码可以写成如下形式：

unsigned long mod_2(unsigned long a, unsigned long b, unsigned long c)
{
unsigned long ret=1;
unsigned long tmp=a%c;

while(b!=0){
if(b%2!=0)
ret=(ret*tmp)%c;
b/=2;
tmp=(tmp*tmp)%c;
}

return ret;
}

第三种算法：

第三种算法也就是第二种算法的递归版本，代码如下：

unsigned long mod_2Recursive(unsigned long a,unsigned long b,unsigned long c)
{
if(b==1)
return a%c;

unsigned long tmp=mod_2Recursive(a,b/2,c);
unsigned long ret=(tmp*tmp)%c;
if(b%2==1)
ret=(ret*a)%c;

return ret;
}