Miller-Rabin

算法简介

Miller-Rabin 算法可以较快地判断一个数是不是质数。
但是，这个算法是概率性算法，即有一定的概率会误判。
所以，一般会使用较多的数检验，从而降低误判的概率。

前置定理

费马小定理：(a^{p-1} equiv 1 (mod p))
注意：费马小定理的逆命题并不一定成立。
二次探测定理：若 (p) 为奇素数且 (x^2 equiv 1 (mod p)) ， 则 (x equiv 1 (mod p)) 或 (x equiv p-1 (mod p))

前置算法

快速乘，快速幂

判断过程

设当前在判断的数是 (x) 。
先判掉 (x) 为偶数的情况。
此时， (x) 为奇数，则 (x-1) 为偶数，
所以 (x-1) 必然能表示成 (2^q*p) （(q) ， (p) 为正整数） 。
先在 (1) ~ (n-1) 中随机一个 (a) ， 若 (a^{x-1} mod x eq 1) ，根据费马小定理， (x) 为合数。
因为费马小定理的逆定理并不一定成立，所以现在还并不能判定 (x) 就是质数。
接下来，枚举一个 (y) ， (0 leq y leq q-1) 。
若 (a^{2^y*p} mod x=1) 且 (a^{2^{y-1} *p} mod x eq 1) 且 (a^{2^{y-1} *p} mod x eq x-1) ，根据二次探测定理， (x) 为合数。
若到最后都没有判断 (x) 为合数，则认为 (x) 为质数。
多用几个 (a) 判断以提高正确率。

code

#include<iostream>
#include<cstdlib>
#include<cstdio>
#include<ctime>
#define int long long
#define LD long double
#define ull unsigned long long
using namespace std;
int MUL(int a,int b,int p) //a*b%p
{
	int x=(LD)a/p*b;
	return ((ull)a*b-(ull)x*p+p)%p;
}
int POW(int a,int b,int p) //a^b%p
{
	if(!b) return 1;
	if(b==1) return a;
	int sum=POW(a,b/2,p);
	if(b%2) return MUL(MUL(sum,sum,p),a,p);
	return MUL(sum,sum,p);
}
bool check(int x)
{
	if(x==0||x==1) return false;
	if(x==2) return true;
	if(x%2==0) return false;
	int p=x-1,q=0;
	while(p%2==0) q++,p/=2;
	for(int i=1;i<=10;i++)
	{
		int a=rand()%(x-1)+1;
		if(POW(a,x-1,x)!=1) return false;
		int lst=1;
		for(int j=0;j<q;j++)
		{
			int t=POW(a,(1ll<<j)*p,x);
			if(t==1&&lst!=1&&lst!=x-1) return false;
			lst=t;
		}
		if(lst!=1&&lst!=x-1) return false;
	}
	return true;
}
int x;
signed main()
{
	srand(time(0));
	scanf("%lld",&x);
	puts(check(x)?"YES":"NO");
	return 0;
}

一些提醒

此算法的时间复杂度为 (O(k log^2 n)) （其中， (k) 为用于判断的质数个数），在数据较大时非常优秀。
但是，在数据比较小的时候此算法的效率甚至不如 (O(sqrt{n})) 的暴力。所以，在数据小时要慎用。