FZU 1649 Prime number or not (MillerRabin素数测试)

题目链接：http://acm.fzu.edu.cn/problem.php?pid=1649 题目大意：很直接，判断一个数n(2<=n<=10^18)是不是素数.   当n达到long long的范围或者更大时，那么先筛好素数或者枚举1~sqrt(n)判断都行不通了，这时便要使用著名的素数测试算法---Miller_Rabin素数测试. 关于算法的学习Matrix67有一篇博客（“数论部分第一节：素数与素数测试”）讲的不错. 而且Miller-Rabin测试还要用到快速幂模计算a ^ b % c的知识，Here有介绍. 但是在Miller-Rabin测试中，判断的素数一般会达到64位，所以需要优化的（速度、不会long long溢出）快速幂模算法.也就是FZU 1752. 当模数也超过long long时，t * t乘法就会有溢出的危险，所以我们就要类似快速幂模设计一个基于加法的快速乘法模求a * b % m. 而且FZU 1752这道题的数据之强已经说明了非递归形式的速度要由于递归形式.   综合起来就是我们最后的 Miller-Rabin素数测试 or 计算a ^ b % c （a,b,c <= 10^18）的模版：  

#include 
using namespace std;

//return a * b % m
unsigned long long quick_add_mod(unsigned long long a, unsigned long long b, unsigned long long m){
    //为了防止long long型a * b溢出，有时需要把乘法变加法
    //且因为暴力加法会超时要使用二分快速乘法模（模仿二分快速幂模……）
    unsigned long long res = 0, tmp = a % m;
    while(b){
        if (b & 1)
        {
            res = res + tmp;
            res = (res >= m ? res - m : res);
        }
        b >>= 1;
        tmp <<= 1;
        tmp = (tmp >= m ? tmp - m : tmp);
    }
    return res;
}

//return a ^ b % m
long long exp_mod(long long a, long long b, long long m){
    long long res = 1 % m, tmp = a % m;
    while(b){
        if (b & 1){
            //如果m在int范围内直接用下一式乘就可以，否则需要用下二式把乘法化加法，用快速乘法模
            //res = (res * t) % m;
            res = quick_add_mod(res, tmp, m);
        }
        //同上
        //t = t * t % m;
        tmp = quick_add_mod(tmp, tmp, m);

        b >>= 1;
    }
    return res;
}

//Miller_Rabin素数测试, 素数return true.
bool Miller_Rabin(long long n){
    int a[5] = {2, 3, 7, 61, 24251};
    //一般Miller_Rabin素数测试是随机选择100个a,这样的错误率为0.25^100
    //但在OI&&ACM中，可以使用上面一组a，在这组底数下，10^16内唯一的强伪素数为46,856,248,255,981

    if (n == 2)
        return true;
    if (n == 1 || (n & 1) == 0)
        return false;

    long long b = n - 1;
    for (int i = 0; i < 5; i ++){
        if (a[i] >= n)
            break;
        while((b & 1) == 0)    b >>= 1;
        long long t = exp_mod(a[i], b, n);
        while(b != n - 1 && t != 1 && t != n-1){
            t = quick_add_mod(t, t, n);
            b <<= 1;
        }
        if (t == n - 1 || b & 1)
            continue;
        else
            return false;
    }
    return true;
}

int main(){
    long long n;
    while(cin >> n){
        if (Miller_Rabin(n))
            cout << "It is a prime number.\n";
        else
            cout << "It is not a prime number.\n";
    }
    return 0;
}