Fermat素性检验
费马小定理:
n
n
n是质数,
g
c
d
(
a
,
n
)
=
1
gcd\left(a,n\right) = 1
gcd(a,n)=1,则
a
n
−
1
≡
1
(
m
o
d
n
)
a^{n-1} \equiv 1\left(\mathop{mod}n\right)
an−1≡1(modn)
反过来不一定成立,即:如果
a
n
−
1
≢
1
(
m
o
d
n
)
a^{n-1} \not\equiv 1\left(\mathop{mod} n\right)
an−1≡1(modn),则
n
n
n也不一定是合数,例如
2
340
≡
1
(
m
o
d
341
)
2^{340} \equiv 1\left(\mathop{mod} 341\right)
2340≡1(mod341)
所以这里的思想就是,在
[
2
,
n
−
1
]
\left[2,n-1\right]
[2,n−1]中随机挑选
a
a
a,检测是否每一次都有
a
n
−
1
≡
1
(
m
o
d
n
)
a^{n-1} \equiv 1\left(\mathop{mod} n\right)
an−1≡1(modn)
如果不是,则
n
n
n是合数(因为此时
∀
a
∈
[
2
,
n
−
1
]
,
s
.
t
.
g
c
d
(
a
,
n
)
=
1
\forall a\in\left[2,n-1\right], s.t.\ gcd\left(a,n\right) =1
∀a∈[2,n−1],s.t. gcd(a,n)=1, 如果
n
n
n是质数,则由费马小定理,矛盾)
#include <cstdio>
#include <stdlib.h>
typedef long long LL;
LL quick_pow(LL a, LL b, LL p) {
LL res = 1;
while (b) {
if (b & 1)res = (res * a) % p;
a = (a * a) % p;
b >>= 1;
}
return res;
}
const int test_time = 500;
bool judge(int n) {
if (n < 3)return n == 2;
for (int i = 1; i <= test_time; ++i) {
int a = rand() % (n - 2) + 2;
if (quick_pow(a, n - 1, n) != 1)return false;
}
return true;
}
Carmichael number
n n n是一个正合数,且 ∀ b , g c d ( b , n ) = 1 \forall b, gcd\left(b,n\right) = 1 ∀b,gcd(b,n)=1,有 b n − 1 ≡ 1 ( m o d n ) b^{n-1} \equiv 1\left(\mathop{mod} n\right) bn−1≡1(modn),则称 n n n是一个卡迈克尔数/费马伪素数/Carmichael number
性质
性质1:Carmichael number必定是至少三个素数的乘积
Korselt’s criterion
n
n
n是一个正合数,当且仅当对于所有的素因数
p
p
p,
p
2
∤
n
p^2 \nmid n
p2∤n且
(
p
−
1
)
∣
(
n
−
1
)
(p-1)\mid (n-1)
(p−1)∣(n−1)
其实这也说明了Carmichael number必然是一个奇数
例题
#include<cstdio>
const int N = 65000;
bool visit[N + 5];
int primes[N / 3], cnt;
void init(int n) {
primes[0] = 2;
cnt = 1;
for (int i = 4; i <= n; i += 2) {
visit[i] = true;
}
for (int i = 3; 3 * i <= n; i += 2) {
if (!visit[i]) {
primes[cnt] = i;
++cnt;
}
for (int j = i * i; j <= n; j += i) {
visit[j] = true;
}
}
}
bool judge(int n) {
//n是一个奇合数
if ((n & 1) == 0 || !visit[n])return false;
// 遍历n的素因子
// n是一个奇数,所以最大的因子不超过n/3
for (int i = 1; primes[i] * 3 <= n; ++i) {
if (n % primes[i] != 0) {
continue;
}
//p^2 \not\mid n and (p-1)|(n-1)
if (n % (primes[i] * primes[i]) == 0 || (n - 1) % (primes[i] - 1) != 0) {
return false;
}
}
return true;
}
int main() {
init(N);
int n;
while (scanf("%d", &n), n) {
if (judge(n))printf("The number %d is a Carmichael number.\n", n);
else printf("%d is normal.\n", n);
}
return 0;
}
Miller-Rabin素性测试
引理
如果
p
p
p是质数,则
x
2
≡
1
(
m
o
d
p
)
x^2 \equiv 1\left(\mathop{mod} p\right)
x2≡1(modp)的解为
x
≡
1
(
m
o
d
p
)
x \equiv 1\left(\mathop{mod} p\right)
x≡1(modp)或
x
≡
−
1
(
m
o
d
p
)
x \equiv -1\left(\mathop{mod} p\right)
x≡−1(modp)
证明:
p
∣
(
x
2
−
1
)
⇒
p
∣
(
x
−
1
)
(
x
+
1
)
p \mid\left(x^2-1\right) \Rightarrow p \mid \left(x-1\right)\left(x+1\right)
p∣(x2−1)⇒p∣(x−1)(x+1)
又因为
p
p
p是质数,
p
∣
(
x
−
1
)
p\mid \left(x-1\right)
p∣(x−1)或
p
∣
(
x
+
1
)
p\mid \left(x+1\right)
p∣(x+1)
原理
结合引理和Fermat素性检验
n
n
n是奇质数
考虑
a
<
n
a<n
a<n,
n
−
1
n-1
n−1是一个偶数,设
n
−
1
=
2
s
u
n-1 = 2^s u
n−1=2su
如果
a
u
≡
1
(
m
o
d
n
)
a^{u} \equiv 1\left(\mathop{mod} n\right)
au≡1(modn),或者
∃
t
<
s
,
\exists t <s,
∃t<s,满足
a
2
t
u
≡
−
1
(
m
o
d
n
)
a^{2^{t} u}\equiv -1\left(\mathop{mod}n\right)
a2tu≡−1(modn) ,则
n
n
n可能是质数,否则就是合数
证明:
当
a
u
≡
±
1
(
m
o
d
n
)
a^u \equiv \pm1\left(\mathop{mod} n\right)
au≡±1(modn),则
a
n
−
1
=
(
a
u
)
2
s
≡
1
(
m
o
d
n
)
,
则
a^{n-1} =\left(a^{u}\right)^{2^s}\equiv1\left(\mathop{mod} n\right),则
an−1=(au)2s≡1(modn),则n$有可能是质数,返回
当
a
u
≢
±
1
(
m
o
d
p
)
a^{u} \not \equiv \pm 1\left(\mathop{mod} p\right)
au≡±1(modp),我们接着平方
如果
a
2
u
≡
1
(
m
o
d
n
)
a^{2u} \equiv 1\left(\mathop{mod} n\right)
a2u≡1(modn),因为
a
u
≢
±
1
(
m
o
d
p
)
a^{u} \not \equiv \pm 1\left(\mathop{mod} p\right)
au≡±1(modp),根据引理,
n
n
n是合数
如果
a
2
u
≡
−
1
(
m
o
d
n
)
a^{2u} \equiv -1 \left(\mathop{mod}n\right)
a2u≡−1(modn),如果
2
u
=
n
−
1
2u = n-1
2u=n−1,则根据费马小定理
n
n
n是合数,否则
n
n
n有可能是质数,返回
如果
a
2
u
≢
±
1
(
m
o
d
n
)
a^{2u} \not\equiv \pm 1 \left(\mathop{mod}n\right)
a2u≡±1(modn),我们就接着平方
同理对于
t
<
s
t<s
t<s
如果
a
2
t
u
≡
1
(
m
o
d
n
)
a^{2^t u} \equiv 1\left(\mathop{mod} n\right)
a2tu≡1(modn),因为
a
2
t
−
1
u
≢
±
1
(
m
o
d
p
)
a^{2^{t-1}u} \not \equiv \pm 1\left(\mathop{mod} p\right)
a2t−1u≡±1(modp),根据引理,
n
n
n是合数
如果
a
2
t
u
≡
−
1
(
m
o
d
n
)
a^{2^tu} \equiv -1 \left(\mathop{mod}n\right)
a2tu≡−1(modn),
n
n
n有可能是质数,返回
如果
a
2
t
u
≢
±
1
(
m
o
d
n
)
a^{2^t u} \not\equiv \pm 1 \left(\mathop{mod}n\right)
a2tu≡±1(modn),我们就接着平方
到最后 a n − 1 ≡ 1 ( m o d n ) a^{n-1}\equiv 1\left(\mathop{mod} n\right) an−1≡1(modn)则 n n n是合数
此外
针对
[
1
,
2
32
)
\left[1,2^{32}\right)
[1,232)范围内的数,只需要选择
{
2
,
7
,
61
}
\left\{2,7,61\right\}
{2,7,61}作为基底进行Miller-Rabin素性测试就可以确定素性
针对
[
1
,
2
64
)
\left[1,2^{64}\right)
[1,264)范围内的数,只需要选择
{
2
,
325
,
9375
,
28178
,
450775
,
9780504
,
1795265022
}
\left\{2, 325, 9375, 28178, 450775, 9780504, 1795265022\right\}
{2,325,9375,28178,450775,9780504,1795265022}作为基底进行Miller-Rabin素性测试就可以确定素性
针对
[
1
,
2
78
)
\left[1,2^{78}\right)
[1,278)范围内的数,只需要选择
{
2
,
3
,
5
,
7
,
11
,
13
,
17
,
19
,
23
,
29
,
31
,
37
}
\left\{2,3,5,7,11,13,17,19,23,29,31,37\right\}
{2,3,5,7,11,13,17,19,23,29,31,37}(前12个质数)作为基底进行Miller-Rabin素性测试就可以确定素性
注意:
如果用这种方法,
1.所有数都要取一遍,不能只选择小于
n
n
n的
2.
a
a
a换成
a
m
o
d
n
a\mathop{mod} n
amodn
3.如果
a
≡
0
(
m
o
d
n
)
a\equiv 0\left(\mathop{mod} n\right)
a≡0(modn)则直接通过该轮测试
另外,假设 广义 Riemann 猜想(generalized Riemann hypothesis, GRH)成立,则对数 n n n最多只需要测试 [ 2 , min { n − 2 , ⌊ 2 ln 2 n ⌋ } ] \left[2,\min\left\{n-2,\left\lfloor 2\ln^2n\right\rfloor\right\}\right] [2,min{n−2,⌊2ln2n⌋}]中的全部整数即可确定 n n n的素性
#include <cstdio>
#include <stdlib.h>
typedef long long LL;
LL quick_pow(LL a, LL b, LL p) {
LL res = 1;
while (b) {
if (b & 1)res = (res * a) % p;
a = (a * a) % p;
b >>= 1;
}
return res;
}
const int test_time = 500;
bool Rabin_Miller(int n) {
if (n < 3 || (n & 1) == 0)return n == 2;
int u = n - 1, t = 0;
// n - 1 = u * (2^t)
while ((u & 1) == 0) {
u >>= 1;
++t;
}
for (int i = 1; i <= test_time; ++i) {
int a = rand() % (n - 2) + 2, v = quick_pow(a, u, n);
if (v == 1)continue;
int s;
for (s = 0; s < t; ++s) {
//a^{u * 2^s}= -1 (mod n)
if (v == n - 1)break;
v = 1LL * v * v % n;
}
if (s == t)return false;
}
return true;
}
第二种: [ 1 , 2 32 ) \left[1,2^{32}\right) [1,232)
#include <cstdio>
#include <cstdlib>
#include <vector>
using namespace std;
typedef long long LL;
LL quick_pow(LL a, LL b, LL p) {
LL res = 1;
while (b) {
if (b & 1)res = (res * a) % p;
a = (a * a) % p;
b >>= 1;
}
return res;
}
vector<unsigned> bases1 = { 2,7,61 };
bool Rabin_Miller2(unsigned n) {
if (n < 3 || (n & 1) == 0)return n == 2;
unsigned u = n - 1, t = 0;
// n - 1 = u * (2^t)
while ((u & 1) == 0) {
u >>= 1;
++t;
}
for (unsigned a : bases1) {
a %= n;
if (a == 0)continue;
unsigned v = quick_pow(a, u, n);
if (v == 1)continue;
unsigned s;
for (s = 0; s < t; ++s) {
//a^{u * 2^s}= -1 (mod n)
if (v == n - 1)break;
v = 1ULL * v * v % n;
}
if (s == t)return false;
}
return true;
}
第三种:
[
1
,
2
64
)
\left[1,2^{64}\right)
[1,264)
__int128要gcc或者g++才有,如果你能确保不溢出,换成long long也可以
#include <cstdio>
using namespace std;
typedef long long LL;
LL quick_pow(LL a, LL b, LL p) {
LL res = 1;
while (b) {
if (b & 1)res = (res * a) % p;
a = (a * a) % p;
b >>= 1;
}
return res;
}
const LL bases[] = {2, 325, 9375, 28178, 450775, 9780504, 1795265022};
bool Rabin_Miller(LL n) {
if (n < 3 || (n & 1) == 0)return n == 2;
LL u = n - 1, t = 0;
// n - 1 = u * (2^t)
while ((u & 1) == 0) {
u >>= 1;
++t;
}
for (int i = 0; i < 7; ++i) {
LL a = bases[i] % n;
if (a == 0)continue;
LL v = quick_pow(a, u, n);
if (v == 1)continue;
LL s;
for (s = 0; s < t; ++s) {
//a^{u * 2^s}= -1 (mod n)
if (v == n - 1)break;
v = (__int128)v * v % n;
}
if (s == t)return false;
}
return true;
}
第四种:
[
1
,
2
78
)
\left[1,2^{78}\right)
[1,278)
__int128要gcc或者g++才有,如果你能确保不溢出,换成long long也可以
#include <cstdio>
using namespace std;
typedef long long LL;
LL quick_pow(LL a, LL b, LL p) {
LL res = 1;
while (b) {
if (b & 1)res = (res * a) % p;
a = (a * a) % p;
b >>= 1;
}
return res;
}
const LL bases[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
bool Rabin_Miller(LL n) {
if (n < 3 || (n & 1) == 0)return n == 2;
LL u = n - 1, t = 0;
// n - 1 = u * (2^t)
while ((u & 1) == 0) {
u >>= 1;
++t;
}
for (int i = 0; i < 12; ++i) {
if(bases[i] == n) return true;
LL a = bases[i] % n;
if (a == 0)return false;
LL v = quick_pow(a, u, n);
if (v == 1)continue;
LL s;
for (s = 0; s < t; ++s) {
//a^{u * 2^s}= -1 (mod n)
if (v == n - 1)break;
v = (__int128)v * v % n;
}
if (s == t)return false;
}
return true;
}
参考:
https://www.ams.org/journals/mcom/1990-55-191/S0025-5718-1990-1023756-8/
https://miller-rabin.appspot.com/#