素数判定
素数不需要解释了,那么素数如何判定?
- 最简单的算法,暴力测试,就是最简单的,从2枚举到 sqrt(n) 就可以知道是不是素数了。
Fermat
小定理费马小定理:对于质数p和任意整数a,有 ap≡a(mod p) (同余)。反之,若满足 ap≡a(mod p) , p 也有很大概率为质数。将两边同时约去一个a,则有
a(p−1)≡1(mod p) Mfiller-Rabin素数判定,在Fermat基础上增加了二次判定:
如果p是奇素数,则 x2≡1(mod p) 的解为 x≡1 或 x≡p−1(mod p)
Fermat小定理
这个肯定是由费马提出来的,这个定理我们看到之后,其实发现是来判断一个数字是不是合数的,而不是素数。 为什么?
因为他说的是,所有素数p一定满足a^(p-1) % p = 1
.所以我们对一个数字如果判断发现不满足了,就一定是合数。那么满足一定是素数么? 不一定啊,因为有可能比较特殊的合数也是有可能有这个效果的。所以要多测试。
结论:
如果不符Fermat小定理,一定是合数, 如果符合Fermat小定理不一定是素数(但是是素数的可能性比较高)。
Miller-Rabin算法
其实二次判定是包含Fermat小定理的,我最开始没看懂这个x是什么得来的,其实这个x是:
所以
那么x自然也知道了: 这个x是一系列的数字,让在指数
2^u
哪里拿出一个
2
就编程
2^(u-1) * 2
,这样不就有平方了么,剩下的底数就是
x
了,并且在还可以继续不停的这样操作,到 u = 0为止。
注意Miller-Rabin算法同样不保证能够一定是素数,所以多次检测,会让错误率变低:
1/(4^S)
,这里S是检测次数。
举例说明
举个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不为质数。
可以发现这个例子里面:x = 2^85 , x= 2^170
(看到了没有出现x = 2^340
的情况,不然一平方操作就是2^680 mod 341
了)
代码实现
先上一下伪代码:参考自:hihocoder 1287
我自己做了一点更改,因为原版的伪代码,在while循环哪里,有一个比较误导的操作:
Miller-Rabin(n):
If (n <= 2) Then
If (n == 2) Then
Return True
End If
Return False
End If
If (n mod 2 == 0) Then
// n为非2的偶数,直接返回合数
Return False
End If
// 我们先找到的最小的a^u,再逐步扩大到a^(n-1)
u = n - 1; // u表示指数
while (u % 2 == 0)
u = u / 2
End While // 提取因子2
For i = 1 .. S // S为设定的测试次数
a = rand_Number(2, n - 1) // 随机获取一个2~n-1的数a
x = a^u % n
tu = u
While (tu < n)
// 依次次检查每一个相邻的 a^u, a^2u, a^4u, ... a^(2^k*u)是否满足二次探测定理
y = x^2 % n
If (y == 1 and x != 1 and x != n - 1) // 二次探测定理
// 若y = x^2 ≡ 1(mod n)
// 但是 x != 1 且 x != n-1
Return False
End If
x = y
tu = tu * 2
End While
If (x != 1) Then // Fermat测试
Return False
End If
End For
Return True
看到伪代码发现是不是挺简单的?But! 这里要时刻注意,不要整数溢出(也就是超出数字范围)。
用long long
类型也不能单纯解决办法,尤其在计算2^340这种操作的时候,复杂度也不低,所以需要用 我之前写过的一个内容:快速幂
AC代码(hihocoder 1287)
/*
ID: chenfus2
PROG: pprime
LANG: C++
*/
#include <iostream>
#include <set>
#include <vector>
#include <string>
#include <cstring>
#include <cstdio>
#include <cstdlib>
#include <time.h>
using namespace std;
#ifndef ll
#define ll long long
#endif
// rand 2 - e
ll rand_Number(ll e){
return rand()%e+1;;
}
// a*b % c
ll mult_mod(ll a,ll b,ll c){
a %= c;
b %= c;
ll res = 0;
while(b){
if(b & 1){
res += a;
res %= c;
}
a <<= 1;
if(a >= c)
a %= c;
b >>= 1;
}
return res;
}
//a^u % num
ll getRemainder(ll a,ll u, ll num){
ll cur = 1;
ll nxt = a;
while (u) {
if ((u&1) > 0) {
cur = mult_mod(cur, nxt, num);
}
nxt = mult_mod(nxt, nxt, num);
u = u>>1;
}
return cur%num;
}
bool checkPrime(ll num){
if (num == 2) {
return true;
}
if (num < 2 || num % 2 == 0) {
return false;
}
ll u = num-1;
int S = 20; //检测次数
while (u % 2 == 0) {
u /= 2;
}
for (int i = 0; i < S; i++) {
ll a = rand_Number(num-1);
ll x = getRemainder(a, u, num);
ll y = x;
ll tu = u;
while (tu < num) {
y = mult_mod(x, x, num);;
if (y == 1 && x != 1 && x != num-1) {
return false;
}
x = y;
tu *= 2;
}
if (x != 1) {
return false;
}
}
return true;
}
int main(){
ll m,n;
srand((unsigned)time(NULL));
scanf("%lld", &m);
for (int i = 0; i < m; i++) {
scanf("%lld", &n);
if (checkPrime(n)) {
cout<<"Yes"<<endl;
continue;
}else{
cout<<"No"<<endl;
continue;
}
}
return 0;
}
时间复杂度:S(logN)