学习笔记 - Miller_Rablin素性测试 | Lucky_Glass's Blog

以心为名，名为爱惜

0%

学习笔记 - Miller_Rablin素性测试

发表于 更新于 分类于 Valine：

被逼无奈……开始自学 Millar_Rablin……

虽然这个算法现在的初二都讲过，我也不知道为什么我们什么都没讲 QwQ

Part1. 引入

我们如何证明一个数 $n$ 是素数？

首先我们有一个 $O(\sqrt n)$ ，即暴力枚举可能是 $n$ 的因数的数 的算法。其次，如果我们要判断多个数，且数的范围不大（数组存得下），我们可以想到线性筛。

However，假如数 $n$ 在 long long 的范围，且要判断多个数是否是素数。怎样高效搞笑地解决这个问题？于是我们了解到—— Miller_Rablin！也就是素性测试。

Part2. 前置

Part2/1. 费马测试

费马小定理：如果 $p$ 是素数，对于任意 $a$，满足 $a^{p-1}\equiv1\pmod p$。

但是并不是说对于某个 $a$，满足 $a^{p-1}\equiv1\pmod p$，$p$ 就是素数（充分但不必要）。

尽管如此，这还是一个强有力的检测方法——只要我们多取几个 $a$，如果都成立，则 $p$ 是质数的可能性就很大（很大而已）。比如有一些合数（被称作Carmichael数），它可以通过所有的费马测试，于是就不能仅仅通过费马测试来测试素数。

显然，质数一定能通过费马测试，但是通过费马测试的可能是合数。

Part2/2. 二次探测

我们先证明这样一个结论：

如果 $p$ 为大于 $2$ 的素数，则对于同余方程 $x^2\equiv 1\pmod p$ ，仅有解 $x\equiv \pm1$。

证明如下：

由原方程变形可得 $x^2-1\equiv 0\pmod p$ ，则 $(x-1)(x+1)\equiv 0$；

因为 $p$ 是素数，所以只有以下三种情况：① $p|(x-1)$ ② $p|(x+1)$ ③ $p|(x-1)$ 且 $p|(x+1)$；
我们可以证明 ③ 是不成立的：

假设 ③ 成立，则令 $\begin{cases}x-1=ap\\x+1=bp\end{cases}$ （$a,b$ 均为整数）；

相减可得 $(b-a)p=2$；

因为 $p$ 是大于 $2$ 的素数，可得矛盾。则 ③ 不成立。

对于 ①，$x-1\equiv 0\pmod p$，$x\equiv1$；
对于 ②，$x+1\equiv 0\pmod p$，$x\equiv -1$；

于是我们成功地证明了这个命题。

现在再考虑它的逆否命题：

如果 $x^2\equiv1\pmod p$ 有一个解 $x’\bmod p\not\in\{1,-1\}$，则 $p$ 不是素数。

因为原命题与它的逆否命题等价，所以这个命题也是成立的。

Part3. 素性测试

我们先假设 $p$ 是一个大于 $2$ 的素数（$p\leq2$ 时特判一下）。

首先 $p$ 一定能通过费马测试，也就是说对于任意整数 $a^{p-1}\equiv1\pmod p$；
根据这一点，我们先随机化一些 $a$，如果有一个 $a$ 不满足 $a^{p-1}\equiv 1\pmod p$，那么 $p$ 就不是素数。

但是不能仅仅用费马测试，对于一个整数 $a$，还要通过二次探测——$x^2\equiv 1\pmod p$ 的解仅有 $x\equiv \pm1$。
我们考虑这样实现二次探测：

因为 $p-1$ 是偶数，所以可以令 $p-1=2^{s}d$ （$s\ge1$，$d$ 是奇数）。

首先计算 $a^d$，如果 $a^d\equiv\pm1$，则有 $a^{2d}\equiv1$，也就是说 $a^d$ 是原方程的一个解。又因为 $a^{2^id}\equiv\left(a^{d}\right)^{2^i}\equiv1$，所以 $a^{2^id}$ 都是原方程的解。所以原方程的解只包含解 $x\equiv \pm1$。通过测试。

然后枚举 $i$ 从 $1$ 到 $s-1$：
① 如果 $a^{2^id}\equiv-1$，则有 $\left(a^{2^id}\right)^2\equiv1$，也就是说 $a^{2^id}$ 本身就是原方程的一个解，那么剩下的解就是 $\left(a^{2^id}\right)^{2^j}\equiv1$ 。通过测试。

② （首先你要知道 当你枚举到当前这个 $a^{2^id}$ 时，之前不存在 $a^{2^id}\equiv \pm1$，不然就判断出是否通过测试了）如果 $a^{2^id}\equiv1$ ，则有 $\left(a^{2^{i-1}d}\right)^2\equiv1$ ，也就是说 $a^{2^{i-1}d}$ 是原方程的解。但是 $a^{2^{i-1}d}\not \equiv\pm1$，所以就存在一个解 $x\not \equiv\pm1$ 。也就不通过测试。

最后如果没有找到解，也不通过测试。

Part4. 源代码

例题是 hihocoder 的 “数论一·Miller-Rabin质数测试”

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
/*Lucky_Glass*/
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
using namespace std;

typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;

inline ll Times(ll Vx,ll Vy,ll Vm){ //这里两个 long long 作乘法会爆 long long……
	ll Vz=(ld)Vx/Vm*Vy;
	ll res=(ull)Vx*Vy-(ull)Vz*Vm;
	return (res+Vm)%Vm;
}
ll QPow(ll Vx,ll Vy,ll Vm){
	ll Rr=1;Vx%=Vm;
	while(Vy){
		if(Vy&1) Rr=Times(Rr,Vx,Vm);
		Vx=Times(Vx,Vx,Vm);
		Vy>>=1;
	}
	return Rr;
}
bool MillerRabin(ll Vd,ll n){
	ll num=2+1ll*rand()*rand()%(n-4); //2<=a<=n-1
	ll Vx=QPow(num,Vd,n);
	if(Vx==1 || Vx==n-1) return true;
	while(Vd!=n-1){
		Vx=Times(Vx,Vx,n);
		Vd*=2;
		if(Vx==1) return false;
		if(Vx==n-1) return true;
	}
	return false;
}
bool CheckPrime(ll n,int tims){
	if(n<=1 || n==4) return false;
	if(n<=3) return true;
	if(n%2==0) return false;
	ll Vd=n-1;
	while(Vd%2==0) Vd/=2;
	while(tims--)
		if(!MillerRabin(Vd,n))
			return false;
	return true;
}
int main(){
	bool fir=true;
	ll n,Ct;
	scanf("%lld",&Ct);
	while(Ct--){
		scanf("%lld",&n);
		if(fir) fir=false,srand(n);
		if(CheckPrime(n,20))
			printf("Yes\n");
		else
			printf("No\n");
	}
	return 0;
}

The End

Thanks for reading!

Email: lucky_glass@foxmail.com，欢迎提问~