Miller-Rabin随机法判定大质数

（转自：http://blog.csdn.net/z690933166/article/details/9860937）

在测试一些非常大的数时，如果用根号n的方法判断也会耗掉大量时间，实际上我们有更快的方法——Miller-Rabin随机法

我们知道，费马小定理：假如p是质数，且(a,p)=1，那么a^(p-1)≡1(mod p)。即假如p是质数，且a,p互质，那么a的(p-1)次方除以p的余数恒等于1。

在大多数情况下费马小定理的逆定理都是成立的，那么我们就得到了一个定理的直接应用，对于待验证的数p，我们不断取a∈［1，p-1]且a∈Z，验证a^(p-1) mod p是否等于1，不是则p果断不是素数，共取s次。其中a^(p-1) mod p可以通过把p-1写成二进制，由(a*b)mod c=(a mod c)*b mod c，可以在t=log(p-1)的时间内计算出解，如考虑整数相乘的复杂度，则一次计算的总复杂度为log³(p-1)。这个方法叫快速幂取模。

我们还有这样一个定理：对于0<x<p，x^2 mod p =1 => x=1或p-1

我们令p-1=(2^t)*u，即p-1为u二进制表示后面跟t个0。我们先计算出x[0]=a^u mod p ，再平方t次并在每一次模p，每一次的结果记为x[i]，最后也可以计算出a^(p-1) mod p。若发现x[i]=1而x[i-1]不等于1也不等于p-1，则发现p必然不是素数。

通过这个两个定理的应用，出错率大大降低，可以说基本不出错。

#include<cstdio>
#include<cstring>
#include<ctime>
#include<cstdlib>
long long n;
inline long long mul(long long a,long long b,long long p)//快速求a*b%p的答案 转换为二进制，一位一位的进行乘和取余运算。例如: b = 1011101那么a * b mod n = (a * 1000000 mod n + a * 10000 mod n + a * 1000 mod n + a * 100 mod n + a * 1 mod n) mod n    
{
	long long ans=0;
	while(b)
	 {
	 	if (b&1) ans=(ans+a)%p;
	 	a=(a+a)%p;
	 	b>>=1;
	 }
	return ans;
}
inline long long poww(long long x,long long num,long long p)//快速幂 
{
  long long sum=1,k=x%p;
  while (num)
   {
   	if (num&1) sum=mul(sum,k,p);
   	num=num>>1;
   	k=mul(k,k,p);
	   }	
  return sum;
}
inline bool miller_rabin(long long n)//miller_rabin算法，随机法判断大质数 
{
	if (n==2) return true;
	if (n<2||!(n&1)) return false;
	long long t=0,a,x,y,u=n-1;
	while (!(u&1)) t++,u>>=1;//若u为偶数，则右移一位，用t记录位数 
	for (int i=0;i<=20;++i)
	 {
	 	a=rand()*rand()%(n-1)+1;
	 	x=poww(a,u,n);
	 	for (int j=0;j<t;++j)//把移掉的位补回来 
	 	 {
	 	 	y=mul(x,x,n); 
	 	 	if (y==1&&x!=1&&x!=n-1) return false;//二次探测 
	 	 	x=y;
		  }
		if (x!=1) return false;//费马小定理 
	 }
	return true;
 }
int main()
{
	scanf("%I64d",&n);
	printf("%d",miller_rabin(n));
	return 0;
}