Miller-rabbin算法

这个算法用于快速判断一个数是否是质数。
那怎么判呢？
首先，既然要快速，我们就不能用朴素的$O(\sqrt n)$的试除法。我们会联想到质数的一些性质：

$$a^{p-1} \equiv 1~(mod~p)(0 < a < p)$$
也就是费马小定理了，我们可以由此想到一个简单的方法：随机一些a去判断$a^{p-1}$是否等于1。
然后按照套路的剧情发展，这个方法有毛病。
毛病在哪里呢？
你可以用这个方法试一下$561 = 11*51$，看看能不能被判断出来。
这就是毛病了：存在像561一样的强伪素数，无论a选多少，得到的结果都是1。
怎么改进呢？
看另一个简单的定理：
如果对于一个偶数x有下列同余式（p为质数）
$$a^x \equiv 1 ~(mod~p)~(0 < x < p)$$
那么
$$a^{\frac{x}{2}} \equiv \pm1~(mod~p)$$
为什么呢？
非常显然，由第一个式子可以得到
$$(a^{\frac{x}{2}} + 1)(a^{\frac{x}{2}} - 1) \equiv 0~(mod~p)\equiv p~(mod~p)$$
因为p是一个质数，所以$(a^{\frac{x}{2}} + 1)$和$(a^{\frac{x}{2}} - 1)$中必然有一项整除p，但是合数的话就不一定了。根据这个性质，我们就有了Miller-rabbin算法。
代码里有一些前人的经验，记不下来的话随机就好了。

Code

#include <iostream>
#include <cstring>
#include <cstdio>

typedef long long ll;

const int base[5] = {2, 3, 7, 61, 24251};

using namespace std;

inline ll read()
{
    ll X = 0; char ch = getchar();
    while (ch < '0' || ch > '9') ch = getchar();
    while (ch >= '0' && ch <= '9') X = X * 10 + ch - '0', ch = getchar();
    return X;
}

inline ll ksm(ll a, ll n, ll mod)
{
    ll b = 1;
    while (n) {
        if (n & 1) b = (__int128)a * b % mod;
        a = (__int128)a * a % mod;
        n >>= 1;
    }
    return b;
}

inline int miller_rabbin(ll n)
{
    if (n == 2 || n == 3 || n == 7 || n == 61 || n == 24251) return 1;
    if (n == 46856248255981ll || n < 2 || !(n & 1)) return 0;
    ll d = n - 1;
    while (!(d & 1)) d >>= 1;
    for (int i = 0; i < 5; i++) {
        ll t = ksm(base[i], d, n), k = d;
        for(; k != n - 1 && t != n - 1 && t != 1; k <<= 1)
            t = (__int128)t * t % n;
        if (t != n - 1 && (k & 1) != 1) return 0;
    }
    return 1;
}

int main(void)
{
    int t = read();
    while (t--) {
        ll n = read();
        printf(miller_rabbin(n) ? "Yes\n" : "No\n");
    }

    return 0;
}