额，我们今天来讲一讲Miller-Rabin素性测试算法。
读者：怎么又是随机算法！！！(⊙o⊙)…
【好了，言归正传】
【费马小定理】
费马小定理只是个必要条件，符合费马小定理而非素数的数叫做Carmichael
Carmichael数是非常少的。
在1~100000000范围内的整数中，只有255个Carmichael数。
为此又有二次探测定理，以确保该数为素数。
这就构成了Miller-Rabin的基本原理 ╰(￣▽￣)╭
【二次探测定理】
二次探测定理 如果p是一个素数，0<x<p,则方程x^2≡1(mod p)的解为x=1,p-1
【Miller-Rabin】
让我们来模拟一下：
第一步：先计算出m、j，使得n-1=m*2^j，其中m是正奇数，j是非负整数
第二步：随机取一个b，2<=b
第三步：计算v=b^m mod n
第四步：如果v==1，通过测试，返回
第五步：令i=1
第六步：如果v=n-1，通过测试，返回
第七步：如果i==j，非素数，结束
第八步：v=v^2 mod n，i=i+1
第九步：循环到S（S根据数据进行改变）
【算法分析】
既然是随机算法，那么我们就来计算一下这种算法的正确率。
Miller-Rabin的正确率为75%。
怎么这么低(⊙o⊙)…，你可以多调用几次就好了！这样可以使正确概率提高为1-(1/4)^S
【代码】
#include<iostream>
#include<cstdio>
#include<cstring>
#include<windows.h>
#include<time.h>
using namespace std;
#define ll long long
const int S=5; //可以根据数据进行更改！！！
ll mod_mul(ll a, ll b, ll n)
{
ll res=0;
while(b)
{
if(b&1) res=(res + a)%n;
a=(a+a)%n;
b>>=1;
}
return res;
} ll mod_exp(ll a, ll b, ll n)
{
ll res=1;
while(b)
{
if(b&1) res=mod_mul(res,a,n);
a=mod_mul(a,a,n);
b>>=1;
}
return res;
}
bool miller_rabin(ll n)
{
if(n==2 || n==3 || n==5 || n==7 || n==11) return true;
if(n == 1 || !(n%2) || !(n%3) || !(n%5) || !(n%7) || !(n%11)) return false;
ll x,pre,u;
int i,j,k=0;
u=n-1;
//要求x^u % n
while(!(u&1))
{ //如果u为偶数则u右移，用k记录移位数
k++;
u >>= 1;
}
srand((ll)time(0));
for(i=0;i<S;++i)
{ //进行S次测试
x=rand()%(n-2) + 2;
//在[2, n)中取随机数
if((x%n)==0) continue;
x=mod_exp(x,u,n);
//先计算(x^u) % n，
pre=x;
for(j=0;j<k;++j)
{ //把移位减掉的量补上，并在这地方加上二次探测
x=mod_mul(x,x,n);
if(x==1 && pre!=1 && pre!=n-1) return false;
//二次探测定理，这里如果x = 1则pre 必须等于 1，或则 n-1否则可以判断不是素数
pre=x;
}
if(x!=1) return false;
//费马小定理
}
return true;
}
int main()
{
long long n;
scanf("%lld",&n);
cout<<miller_rabin(n);
system("pause");
return 0;
}
P.S. 感谢吉林大学的模板，当然这已经改过了，进行了优化。
【你以为这就完了？】
下面我们来研究一下Miller-Rabin的更高效算法
【代码】
#include<windows.h>
#include<iostream>
#include<time.h>
#include<cstdio>
#include<cmath>
using namespace std;
bool witness(long long a,long long n)
{
long long t,d,x;
d=1;
int i=ceil(log(n-1.0)/log(2.0)) - 1;
for(;i>=0;i--)
{
x=d; d=(d*d)%n;
if(d==1 && x!=1 && x!=n-1) return true;
if( ((n-1) & (1<<i)) > 0)
d=(d*a)%n;
}
return d==1? false : true;
}
bool miller_rabin(long long n)
{
int s[]={2,7,61};
if(n==2) return true;
if(n==1 || ((n&1)==0)) return false;
for(int i=0;i<3;i++)//注意这里的i<3不可以改，因为每次S都只为2,7,61！
if(witness(s[i], n)) return false;
return true;
}
int main()
{
long long n;
scanf("%lld",&n);
cout<<miller_rabin(n);
//system("pause");
return 0;
}
我们发现，它只调用witness函数3次，每次都a为2,7,61。时间很快啊。
但正确率……
:-D //迷之微笑