转载自：http://www.dxmtb.com/blog/miller-rabbin/
普通的素数测试我们有O(√ n)的试除算法。事实上，我们有O(slog³n)的算法。
定理一：假如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果断不是素数。
可以证明，使用以上两个定理以后，检验s次出错的概率至多为2^(-s)，所以这个算法是很可靠的。
需要注意的是，为了防止溢出（特别大的数据），a*b mod c 也应用类似快速幂取模的方法计算。当然，数据不是很大就可以免了
codevs 1702 素数判定 2
一个数，他是素数么？
设他为P满足(P<=2^63-1)
P
Yes|No
2
Yes
算法导论——数论那一节
注意Carmichael Number
分类标签 Tags 点此展开
/*数据确实没有错，我的错误点是，在两个long long相乘的时候，直接进行了乘法运算，很有可能溢出，所以错了后面的点，long long数据要再写一个相乘%的子函数*/
#include<iostream>
using namespace std;
#include<cstdio>
#define S 10
#include<cstdlib>
#include<ctime>
#define ll long long
ll cas, maxz;
ll read()
{
ll ans=;char c;
c=getchar();
while(c<''||c>'') c=getchar();
while(c>=''&&c<='')
{
ans=ans*+c-'';
c=getchar();
}
return ans;
}
ll quick_mul_mod(ll a,ll b,ll c)//a*b%c
{
ll ret=;
a%=c;b%=c;
while(b)
{
if(b&)
{
ret+=a;/*相当于加了b次a*/
ret%=c;
b--;
}
a<<=;
a%=c;
b>>=;
}
return ret;
}
ll quick_mod(ll a,ll b,ll c)//ji suan a^b%c
{
ll ans=;
a%=c;
while(b)
{
if(b&)
{
b--;
ans=quick_mul_mod(ans,a,c);/*注意只要是long long的乘法，很有可能溢出的*/
}
b>>=;
a=quick_mul_mod(a,a,c);
}
return ans;
}
bool Miller_rabin(ll n)
{
if(n==) return true;
if(n<=||!(n&)) return false;
ll u=n-,t=;
while(!(u&))
{
u>>=;
t++;
}
for(int i=;i<S;++i)
{
ll x=rand()%(n-)+;
x=quick_mod(x,u,n);
for(int i=;i<=t;++i)
{
ll y=quick_mul_mod(x,x,n);
if(y==&&x!=&&x!=n-)
return false;
x=y;
}
if(x!=) return false;
}
return true;
} int main()
{
ll n=read();
if(Miller_rabin(n))
printf("Yes\n");
else printf("No\n");
return ;
}