Miller Rabin

Miller-Rabin 算法是一种高效的素性测试算法，对于正整数 \(n\)，可以以 \(O(k \log n)\) 的时间复杂度对其进行素性测试，其中 \(k\) 表示测试的次数。当选取前 12 个素数为底数（\(k=12\)）时，该算法对于 long long 范围内的正整数可以保证正确性。

朴素的试除法时间复杂度是 \(O(\sqrt{n})\)，在预处理过素数之后可以做到 \(O(\frac{\sqrt{n}}{\ln n})\)。此处从略。

换一个角度理解试除法，我们发现对每一个试除的除数 \(i\)，\(i \nmid n\) 是 \(n\) 为素数的必要条件。我们实际上就是在验证一系列 \(n\) 为素数的必要条件是否成立。从 \(2\) 尝试到 \(\lfloor\sqrt{n}\rfloor\) 有足够多个，已经具有了充分性。但我们是不是可以寻找一些更高效的必要条件，使得每一次测试都能排除尽可能多的合数，使得我们只需要较少的测试次数就能得到 \(n\) 的素性？

Fermat 素性测试

费马小定理：对素数 \(p\)，若整数 \(a\) 不是 \(p\) 的倍数，则有 \(a^{p-1} \equiv 1 \pmod p\)。

我们会马上想到费马小定理。由费马小定理可知，对待测正整数 \(n\) 和 \([2,n-1]\) 的正整数 \(a\)，\(a^{n-1} \equiv 1 \pmod n\) 成立是 \(n\) 为素数的必要条件。此处有 \(a<n\)，自动满足 \(a\) 不是 \(n\) 的倍数。

那我们自然会选取几个 \([2,n-1]\) 之间的 \(a\)，检查是否有 \(a^{n-1} \equiv 1 \pmod n\) 成立。如果发现没有矛盾，就认为 \(n\) 很可能是素数，或者说 \(n\) 通过了这次 Fermat 测试。

但有这样一类数，\([2,n-1]\) 内的所有与 \(n\) 互素的底数 \(a\) 都满足 \(a^{n-1} \equiv 1 \pmod n\)，但 \(n\) 却是合数。一个例子是 \(561=3 \times 11 \times 17\)。这种数被称为 Carmichael 数。这也提醒我们费马小定理的逆定理不成立。

二次探测定理

二次探测定理：若 \(n\) 为奇素数，同余方程 \(x^2 \equiv 1 \pmod n\) 仅存在平凡解 \(x=1\) 和 \(x=n-1\)。

不难发现，这又可以作为一个 \(n\) 为素数的必要条件。那我们把费马小定理和二次探测定理这两个必要条件拼一下，就得到了正确性有保证的 Miller-Rabin 算法。

Miller-Rabin 算法

我们被给到一个正整数 \(n\)，现在需要测试其素性。考虑 \(a^{n-1}\)。需要检测 \(n\) 作为模数的性质（二次探测定理）和 \(n\) 作为指数和模数的性质（费马小定理）。

因为 \(n\) 是奇数，所以 \(n-1\) 是偶数。不断将 \(n-1\) 除以 \(2\)，直到得到奇数 \(u\)。记录除以 \(2\) 的次数，记为 \(t\)。可知 \(n-1=u \times 2^t\)。指数的每个 \(2\) 因子都相当于一次平方。

接下来选取 \(a\)。选前 12 个质数可以保证 long long 内的正确性。也可以随机一些 \([2,n-1]\) 的整数。

我们用快速幂计算 \(x_0=a^u \mod n\)，递推 \(x_{i+1}=x_i^2 \mod n\)，我们用 \(\left\{ x_i \right\}\) 验证二次探测定理。每次只需要检测是否有 \(x_i=1\) 或 \(x_i=n-1\)，有以下几种情况：

• 初始就有 \(x_0=1\) 或 \(x_0=n-1\)：\(i \ge 1\) 的每一个 \(x_i\) 都为 \(1\)。该底数通过。
• 中途出现 \(x_i=n-1\)：下一步就会出现 \(x_{i+1}=1\)。该底数通过。
• 中途出现 \(x_i=1\)：上一步未跳出，肯定不是 \(x_{i-1}=n-1\)。不通过。

注意递推只能进行 \(t-1\) 轮，这是因为第 \(t\) 轮的 \(x_t=a^{n-1} \mod n\)。这时我们需要验证的是费马小定理，也就是 \(x_t \equiv 1 \pmod n\)，这需要 \(x_{t-1}=n-1\)。如果在第 \(t-1\) 轮仍没有跳出，这说明费马小定理不成立，也不通过。

只有每一个底数都通过，才能判断 \(n\) 通过了 Miller-Rabin 素性检验。

实现

#include<stdio.h>
#include<stdbool.h>
#define ll long long
#define lll __int128
ll power(ll x,ll y,ll mod)
{
	ll ret=1;
	while(y)
	{
		if(y%2==1) ret=(lll)ret*x%mod;
		x=(lll)x*x%mod; y/=2;
	}
	return ret;
}
const int p[12]={2,3,5,7,11,13,17,19,23,29,31,37};
bool miller(ll n)
{
	if(n<3||n%2==0) return n==2;
	if(n%3==0) return n==3;
	ll u=n-1,t=0;
	while(u%2==0) u/=2, t++;
	for(int i=0;i<12;i++)
	{
		if(p[i]>=n) break;
		ll x=power(p[i],u,n);
		if(x==1||x==n-1) continue;
		bool flag=0;
		for(int j=1;j<=t-1;j++)
		{
			x=(lll)x*x%n;
			if(x==1) return 0;
			if(x==n-1) {flag=1; break;}
		}
		if(!flag) return 0;
	}
	return 1;
}

int main()
{
	int T; scanf("%d",&T);
	while(T--)
	{
		ll x; scanf("%lld",&x);
		if(miller(x)) printf("YES\n");
		else printf("NO\n");
	}
	return 0;
}