有关素数判断的一些算法（总结&&对比）
素性测试是数论题中比较常用的一个技巧。它可以很基础，也可以很高级（哲学）。这次主要要介绍一下有关素数判断的奇技淫巧
素数的判断主要分为两种：范围筛选型&&单个判断型
我们先从范围筛选型这种常用的开始讲起，这里采用模板题Luogu P3383 【模板】线性筛素数来进行测试
1.埃氏筛
这是最常用的筛法了，思路也很简单：任何一个素数的倍数都是合数
然后我们O(n)扫一遍，同时筛去素数的倍数
但是有一些数如6，会被2和3都筛去一次，就造成了效率上的浪费，所以复杂度经证明为**O(n log log n)
CODE
#include<cstdio>
using namespace std;
const int N=10000005;
bool vis[N];
int n,m,x;
inline char tc(void)
{
static char fl[100000],*A=fl,*B=fl;
return A==B&&(B=(A=fl)+fread(fl,1,100000,stdin),A==B)?EOF:*A++;
}
inline void read(int &x)
{
x=0; char ch=tc();
while (ch<'0'||ch>'9') ch=tc();
while (ch>='0'&&ch<='9') x=x*10+ch-'0',ch=tc();
}
inline void get_prime(int m)
{
register int i,j;
for (vis[1]=1,i=2;i<=m;++i)
if (!vis[i]) for (j=i<<1;j<=m;j+=i) vis[j]=1;
}
int main()
{
//freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
read(n); read(m); get_prime(n);
while (m--)
{
read(x);
puts(vis[x]?"No":"Yes");
}
return 0;
}
2.线性筛（欧拉筛）
这其实是对上者的优化，我们意识到一个数应该只有它的最小质因数删去，所以我们可以一边筛数的同时一边记录素数，这就是真正的O(n)复杂度
CODE
#include<cstdio>
using namespace std;
const int N=10000005;
int prime[N],n,m,x,cnt;
bool vis[N];
inline char tc(void)
{
static char fl[100000],*A=fl,*B=fl;
return A==B&&(B=(A=fl)+fread(fl,1,100000,stdin),A==B)?EOF:*A++;
}
inline void read(int &x)
{
x=0; char ch=tc();
while (ch<'0'||ch>'9') ch=tc();
while (ch>='0'&&ch<='9') x=x*10+ch-'0',ch=tc();
}
inline void Euler(int n)
{
register int i,j;
for (vis[1]=1,i=2;i<=n;++i)
{
if (!vis[i]) prime[++cnt]=i;
for (j=1;j<=cnt&&i*prime[j]<=n;++j)
{
vis[i*prime[j]]=1;
if (!(i%prime[j])) break;
}
}
}
int main()
{
//freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
read(n); read(m);
Euler(n);
while (m--)
{
read(x);
puts(vis[x]?"No":"Yes");
}
return 0;
}
注意上面的那句话：
if (!(i%prime[j])) break;
这保证了线性筛的效率，不会产生重复，因为当i%prime[j]==0时这个数就是让后面的数删去。
3.基础素性测试
这是最基本的素数判定法了吧。从2到sqrt(x)枚举是否有数能够整除x
证明的话很简单，因为如果这个数是素数，那么它的因数必定为1和x，若其因数大于sqrt(x)，那么平方后就大于x，这显然不可能。
所以我们O(sqrt(x))判断一次
CODE
#include<cstdio>
#include<cmath>
using namespace std;
int n,m,x;
inline char tc(void)
{
static char fl[100000],*A=fl,*B=fl;
return A==B&&(B=(A=fl)+fread(fl,1,100000,stdin),A==B)?EOF:*A++;
}
inline void read(int &x)
{
x=0; char ch=tc();
while (ch<'0'||ch>'9') ch=tc();
while (ch>='0'&&ch<='9') x=x*10+ch-'0',ch=tc();
}
inline bool check(int x)
{
if (!(x^1)) return 0;
register int i; int bound=(int)sqrt(x);
for (i=2;i<=bound;++i)
if (!(x%i)) return 0;
return 1;
}
int main()
{
//freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
read(n); read(m);
while (m--)
{
read(x);
puts(check(x)?"Yes":"No");
}
return 0;
}
4.对于算法3的优化
首先我们看一个结论：
大于等于5的质数一定和6的倍数相邻。
证明等参考：dalao's blog
然后同3，我们只不过每次快进6个单位，然后常数就得到了难以言喻都优化（一跃成为此题最快的算法）
CODE
#include<cstdio>
#include<cmath>
using namespace std;
int n,m,x;
inline char tc(void)
{
static char fl[100000],*A=fl,*B=fl;
return A==B&&(B=(A=fl)+fread(fl,1,100000,stdin),A==B)?EOF:*A++;
}
inline void read(int &x)
{
x=0; char ch=tc();
while (ch<'0'||ch>'9') ch=tc();
while (ch>='0'&&ch<='9') x=x*10+ch-'0',ch=tc();
}
inline bool check(int x)
{
if (!(x^1)) return 0;
if (!(x^2)||!(x^3)) return 1;
if ((x%6)^1&&(x%6)^5) return 0;
register int i; int bound=(int)sqrt(x);
for (i=5;i<=bound;i+=6)
if (!(x%i)||!(x%(i+2))) return 0;
return 1;
}
int main()
{
//freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
read(n); read(m);
while (m--)
{
read(x);
puts(check(x)?"Yes":"No");
}
return 0;
}
5.Miller-Rabin算法
这是历史上判断素数最快的方法了吧（但在此题中被算法4吊打了）
首先，这个算法基于费马小定理和二次探测定理：
二次探测定理：如果p是奇素数，则 x2≡1(modp)的解为x = 1或x = p - 1(mod p)
所以我们可以把x变成r*2^t的形式，其中r是一个奇数
然后我们结合两种算法&&快速幂就可以稳定O(log x)进行单次判断了
但是这个算法是一个非完美算法，它每一次都25%的概率是错的，所以我们可以多选择几个数多弄几次
但是偶然在网上看到一段话：
对于大数的素性判断，目前Miller-Rabin算法应用最广泛。一般底数仍然是随机选取，但当待测数不太大时，选择测试底数就有一些技巧了。比如，如果被测数小于4 759 123 141，那么只需要测试三个底数2, 7和61就足够了。当然，你测试的越多，正确的范围肯定也越大。如果你每次都用前7个素数(2, 3, 5, 7, 11, 13和17)进行测试，所有不超过341 550 071 728 320的数都是正确的。如果选用2, 3, 7, 61和24251作为底数，那么10^16内唯一的强伪素数为46 856 248 255 981。这样的一些结论使得Miller-Rabin算法在OI中非常实用。通常认为，Miller-Rabin素性测试的正确率可以令人接受，随机选取k个底数进行测试算法的失误率大概为4^(-k)。
所以对于这一题n=10000000的范围就只需要选择2，7，61即可
CODE
#include<cstdio>
using namespace std;
typedef long long LL;
const int prime[3]={2,7,61};
int n,m,x;
inline char tc(void)
{
static char fl[100000],*A=fl,*B=fl;
return A==B&&(B=(A=fl)+fread(fl,1,100000,stdin),A==B)?EOF:*A++;
}
inline void read(int &x)
{
x=0; char ch=tc();
while (ch<'0'||ch>'9') ch=tc();
while (ch>='0'&&ch<='9') x=x*10+ch-'0',ch=tc();
}
inline int quick_pow(int x,int p,int mod)
{
int tot=1;
while (p)
{
if (p&1) tot=((LL)tot*x)%mod;
x=((LL)x*x)%mod; p>>=1;
}
return tot;
}
inline bool Miller_Rabin(int x)
{
if (!(x^2)) return 1;
if (x<2||!(x&1)) return 0;
int t=0,u=x-1;
while (!(u&1)) ++t,u>>=1;
for (register int i=0;i<3;++i)
{
if (!(x^prime[i])) return 1;
if (!(x%prime[i])) return 0;
int lst=quick_pow(prime[i],u,x);
for (register int j=1;j<=t;++j)
{
int now=((LL)lst*lst)%x;
if (!(now^1)&&lst^1&&lst^(x-1)) return 0; lst=now;
}
if (lst^1) return 0;
}
return 1;
}
int main()
{
//freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
read(n); read(m);
while (m--)
{
read(x);
puts(Miller_Rabin(x)?"Yes":"No");
}
return 0;
}
最后给出5个算法的运行结果（无O2）
