作用:
有时候我们想快速的知道一个数是不是素数,而这个数又特别的大导致 O( √n) 的算法不能通过,这时候我们可以对其进行 Miller-Rabin 素数测试,可以大概率测出其是否为素数。
两个定理:

费马小定理:如果p是素数,a是小于p的正整数,那么a^(p-1) mod p = 1。
证明:首先我们证明这样一个结论:如果p是一个素数的话,那么对任意一个小于p的正整数a,a, 2a, 3a, …, (p-1)a除以p的余数正好是一个1到p-1的排列。例如,5是素数,3, 6, 9, 12除以5的余数分别为3, 1, 4, 2,正好就是1到4这四个数。
反证法,假如结论不成立的话,那么就是说有两个小于p的正整数m和n使得na和ma除以p的余数相同。不妨假设n>m,则p可以整除a(n-m)。但p是素数,那么a和n-m中至少有一个含有因子p。这显然是不可能的,因为a和n-m都比p小。
用同余式表述,我们证明了:
(p-1)! ≡ a * 2a * 3a * … * (p-1)a (mod p)
也即:
(p-1)! ≡ (p-1)! * a^(p-1) (mod p)
两边同时除以(p-1)!,就得到了我们的最终结论:
1 ≡ a^(p-1) (mod p)

二次探测定理:如果p是素数,x是小于p的正整数,且x^2 mod p = 1,那么要么x=1,要么x=p-1。
这是显然的,因为x^2 mod p = 1相当于p能整除x^2-1,也即p能整除(x+1)(x-1)。由于p是素数,那么只可能是x-1能被p整除(此时x=1)或x+1能被p整除(此时x=p-1)。

我们下面来演示一下上面的定理如何应用在Fermat素性测试上。341可以通过以2为底的Fermat测试,因为2^340 mod 341=1。如果341真是素数的话,那么2^170 mod 341只可能是1或340;当算得2^170 mod 341确实等于1时,我们可以继续查看2^85除以341的结果。我们发现,2^85 mod 341=32,这一结果摘掉了341头上的素数皇冠,面具后面真实的嘴脸显现了出来。

这就是Miller-Rabin素性测试的方法。不断地提取指数n-1中的因子2,把n-1表示成d2^r(其中d是一个奇数)。那么我们需要计算的东西就变成了a的d2^r次方除以n的余数。于是,a^(d * 2^(r-1))要么等于1,要么等于n-1。如果a^(d * 2^(r-1))等于1,定理继续适用于a^(d * 2^(r-2)),这样不断开方开下去,直到对于某个i满足a^(d * 2^i) mod n = n-1或者最后指数中的2用完了得到的a^d mod n=1或n-1。这样,Fermat小定理加强为如下形式:
尽可能提取因子2,把n-1表示成d2^r,如果n是一个素数,那么或者a^d mod n=1,或者存在某个i使得a^(d2^i) mod n=n-1 ( 0<=i<r ) (注意i可以等于0,这就把a^d mod n=n-1的情况统一到后面去了)

Miller-Rabin素性测试同样是不确定算法,我们把可以通过以a为底的Miller-Rabin测试的合数称作以a为底的强伪素数(strong pseudoprime)。第一个以2为底的强伪素数为2047。第一个以2和3为底的强伪素数则大到1 373 653。

Miller-Rabin算法的代码也非常简单:计算d和r的值(可以用位运算加速),然后快速幂计算a^d mod n的值,最后把它平方r次。

算法流程
(1)对于偶数和 0,1,2 可以直接判断。

(2)设要测试的数为 ,我们取一个较小的质数 ,设 ,满足 (其中 是奇数)。

(3)我们先算出 ,然后不断地平方并且进行二次探测(进行 次)。

(4)最后我们根据费马小定律,如果最后 ,则说明 为合数。

(5)多次取不同的 进行 素数测试,这样可以使正确性更高

【备注】
(1)我们可以多选择几个 ,如果全部通过,那么 大概率是质数。

(2) 素数测试中,“大概率”意味着概率非常大,基本上可以放心使用。

(3)当 取遍小等于 30 的所有素数时,可以证明 范围内的数不会出错。

(4)另外,如果是求一个long long类型的平方,可能会爆掉,因此有时我们要用“快速积”,不能直接乘。

Miller-Rabin素性测试模板:

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll prime[10] = {2,3,5,7,11,13,17,19,23,29};
ll Quick_Multiply(ll a, ll b, ll c)
{
    ll ans = 0, res = a;
    while(b)
    {
        if(b & 1)
            ans = (ans + res) % c;
        res = (res + res) % c;
        b >>= 1;
    }
    return ans;
}
ll Quick_Power(ll a, ll b, ll m)
{
    ll ans = 1;
    while(b)
    {
        if(b & 1)
        {
            ans = ans * a % m;
        }
        a = a * a % m;
        b >>= 1;
    }
    return ans;
}
bool Miller_Rabin(ll x)    //判断素数
{
    ll i, j, k;
    ll s = 0, t = x-1;
    if(x == 2) return true; //2是素数
    if(x < 2 || !(x & 1)) return false;  //如果x是偶数或者是0,1,那它不是素数
    while(!(t & 1)) //将x-1分解成(2^s)*t的样子
    {
        s++;
        t >>= 1;
    }
    for(i = 0; i < 10 && prime[i] < x; i++) //随便选一个素数进行测试
    {
        ll a = prime[i];
        ll b = Quick_Power(a, t, x);   //先算出a^t
        for(j = 1; j <= s; j++) //然后进行s次平方
        {
            k = Quick_Multiply(b,b,x);  //求b的平方
            if(k == 1 && b != 1 && b != x-1)
                return false;
            b = k;
        }
        if(b != 1) return false;    //用费马小定理判断
    }
    return true;    //如果进行多次测试都是对的,那么x就很有可能是素数
}

int main()
{
    ll x;
    scanf("%lld",&x);
    if(Miller_Rabin(x)) printf("Yes");
    else printf("No");
    return 0;
}

参考博客:
http://www.matrix67.com/blog/archives/234
https://blog.csdn.net/forever_dreams/article/details/82314237