Miller_Rabin 素数测试

费马定理的逆定理几乎可以用来判断一个数是否为素数，但是有一些数是判断不出来的，因此，Miller_Rabin测试方法对费马的测试过程做了改进，克服其存在的问题。

推理过程如下（摘自维基百科）：

摘自另一篇博文（手动滑稽）：

原理明白了，就直接上代码了（KuangBin大神的板子）：

代码思路是，

• Miller_Rabin()函数随机选取 s 个a，a用做“基底”
• check() 函数是用来判断x是否等于1,也就是判断a是否是n的凭证。
• Mul_mod()函数是 快速乘 ，求 a^t % n 之后的值是否为正负一，因为两个数直接乘的话会很大，可能会爆Long Long, 因此用快速乘边乘边mod。
• pow_mod()函数是 快速幂 ，在刚开始第一次判断 a^(n-1) % n 时会用到。 

/* *************************************************
* Miller_Rabin 算法进行素数测试
* 速度快可以判断一个 < 2^63 的数是不是素数
*
**************************************************/
const int S = 8; //随机算法判定次数一般 8 ∼ 10 就够了
// 计算 ret = (a*b)%c
//a,b,c < 2^63
long long mult_mod(long long a,long long b,long long c)
{
    a %= c;
    b %= c;
    long long ret = 0;
    long long tmp = a;
    while(b)
    {
        if(b & 1)
        {
            ret += tmp;
            if(ret > c)ret -= c;//直接取模慢很多
        }
        tmp <<= 1;
        if(tmp > c)tmp -= c;
        b >>= 1;
    }
    return ret;
}
// 计算 ret = (a^n)%mod
long long pow_mod(long long a,long long n,long long mod)
{
    long long ret = 1;
    long long temp = a%mod;
    while(n)
    {
        if(n & 1)ret = mult_mod(ret,temp,mod);
        temp = mult_mod(temp,temp,mod);
        n >>= 1;
    }
    return ret;
}
// 通过 a^(n − 1)=1(mod n)来判断 n 是不是素数
// n − 1 = x ∗ 2 t 中间使用二次判断
// 是合数返回 true, 不一定是合数返回 false
bool check(long long a,long long n,long long x,long long t)
{
    long long ret = pow_mod(a,x,n);
    long long last = ret;
    for(int i = 1; i <= t; i++)
    {
        ret = mult_mod(ret,ret,n);
        if(ret == 1 && last != 1 && last != n - 1)return true;//合数
        last = ret;
    }
    if(ret != 1)return true;
    else return false;
}
//**************************************************
// Miller_Rabin 算法
// 是素数返回 true,(可能是伪素数)
// 不是素数返回 false
//**************************************************
bool Miller_Rabin(long long n)
{
    if( n < 2)return false;
    if( n == 2)return true;
    if( (n&1) == 0)return false;//偶数
    long long x = n - 1;
    long long t = 0;
    while( (x&1)==0 )
    {
        x >>= 1;
        t++;
    }
    srand(time(NULL));/* *************** */
    for(int i = 0; i < S; i++)
    {
        long long a = rand()%(n - 1) + 1;
        if( check(a,n,x,t) )
            return false;
    }
    return true;
}