hihoCoder 92week Miller-Rabin质数测试

http://hihocoder.com/contest/hiho92/problem/1

这种质数算法是基于费马小定理的一个扩展。

首先我们要知道什么是费马小定理：

费马小定理：对于质数p和任意整数a，有a^p ≡ a(mod p)(同余)。反之，若满足a^p ≡ a(mod p)，p也有很大概率为质数。
将两边同时约去一个a，则有a^(p-1) ≡ 1(mod p)

也即是说：假设我们要测试n是否为质数。我们可以随机选取一个数a，然后计算a^(n-1) mod n，如果结果不为1，我们可以100%断定n不是质数。

否则我们再随机选取一个新的数a进行测试。如此反复多次，如果每次结果都是1，我们就假定n是质数。

该测试被称为Fermat测试。需要注意的是：Fermat测试不一定是准确的，有可能出现把合数误判为质数的情况。

Miller和Rabin在Fermat测试上，建立了Miller-Rabin质数测试算法。

与Fermat测试相比，增加了一个二次探测定理：

如果p是奇素数，则 x^2 ≡ 1(mod p)的解为 x ≡ 1 或 x ≡ p - 1(mod p)

如果a^(n-1) ≡ 1 (mod n)成立，Miller-Rabin算法不是立即找另一个a进行测试，而是看n-1是不是偶数。如果n-1是偶数，另u=(n-1)/2，并检查是否满足二次探测定理即a^u ≡ 1 或 a^u ≡ n - 1(mod n)。

举个 Matrix67 Blog 上的例子，假设n=341，我们选取的a=2。则第一次测试时，2^340 mod 341=1。由于340是偶数，因此我们检查2^170，得到2^170 mod 341=1，满足二次探测定理。同时由于170还是偶数，因此我们进一步检查2^85 mod 341=32。此时不满足二次探测定理，因此可以判定341不为质数。

将这两条定理合起来，也就是最常见的Miller-Rabin测试。

<span style="font-size:18px;">#include <stdio.h> 
#include <stdlib.h>

//求(x*y)%n
long long mod_pro(long long x,long long y,long long n) 
{ 
    long long ret=0,tmp=x%n; 
    while(y) 
    { 
        if(y&0x1)if((ret+=tmp)>n)ret-=n; 
        if((tmp<<=1)>n)tmp-=n; 
        y>>=1; 
    } 
    return ret; 
} 

//快速幂，求(a^b) % c
long long mod(long long a,long long b,long long c) 
{ 
    long long ret=1; 
    while(b) 
    { 
        if(b&0x1)ret=mod_pro(ret,a,c); 
        a=mod_pro(a,a,c); 
        b>>=1; 
    } 
    return ret; 
} 

long long ran() 
{ 
    long long ret=rand(); 
    return ret*rand(); 
} 

int is_prime(long long n,int t) 
{ 
    if(n<2)return 0; 
    if(n==2)return 1; 
    if(!(n&0x1))return 0; 
    long long k=0,m,a,i; 
    for(m=n-1;!(m&1);m>>=1,k++); 
    while(t--) 
    { 
        a=mod(ran()%(n-2)+2,m,n); 
        if(a!=1) 
        { 
            for(i=0;i<k&&a!=n-1;i++) 
                a=mod_pro(a,a,n); 
            if(i>=k)return 0; 
        } 
    } 
    return 1; 
} 
int main()
{
	int amount ;
	long long i;
	scanf("%d",&amount);
	while(amount--)
	{
		scanf("%lld",&i);
		if(is_prime(i,1)) printf("Yes\n");
		else printf("No\n");
	}
	return 0;
}</span>