Miller_Rabin算法

为了测试一个大整数是不是素数，我们不能够使用传统的测试是否有因子的方法，因为那样的时间复杂度至少也是O(n)O(n)O(n)，空间复杂度是O(n)O(n)O(n)（使用线性筛数法），时间复杂度还好说，空间复杂度在n比较大的时候完全不可以接受，这就需要我们使用Miller_RabinMiller\_RabinMiller_Rabin算法进行解决。

算法的核心使用了费马小定理和二次探测定理

费马小定理：若ppp为质数，则ap−1≡1(modp)a^{p-1}\equiv 1\pmod pap−1≡1(modp)

这是一个质数的必要条件，但是有一类特殊的数也满足这个等式并且是合数（卡米歇尔数），为了剔除这些特殊的数，我们还需要二次探测定理帮助判断

二次探测定理：如果ppp为素数，x2≡1(modp)x^2\equiv 1 \pmod px2≡1(modp)，则x=1或x=p−1x=1 或 x=p-1x=1或x=p−1

因此对于满足费马小定理的数ppp，即ap−122≡1(modp){a^{\frac{p-1}{2}}}^2\equiv 1\pmod pa2p−12≡1(modp)，我们判断ap−12a^{\frac{p-1}{2}}a2p−1是否为1或者p−1p-1p−1，如果不是说明ppp不是素数，如果是我们可以继续对ap−12a^{\frac{p-1}{2}}a2p−1进行探测

每一次检测都无法100%确定这是一个素数，但是判断8次几乎不可能出错。

这里的mult函数是乘积函数，因为接近long long 的数字相乘会爆long long,这里按位相乘。

在研究的很多人的实现后我实现了以下，觉得我这个复杂度低一些(因为自己写的缘故看起来很亲切)。实际过题时间也比之前看到的快。

typedef long long ll;
ll mult(ll x,ll y,ll p)
{long double d=1; d=d*x/p*y;return ((x*y-((ll)d)*p)%p+p)%p;
}ll quick_pow(ll a,ll b,ll p)
{a%=p; ll ret=1;while(b){if(b&1) ret=mult(ret,a,p);a=mult(a,a,p); b>>=1;}return ret;
}bool Miller_Rabin(ll n)
{if(n<2) return false;const int times=8;const ll prime[times]={2,3,5,7,11,13,17,61};for(int i=0;i<times;i++)if(n==prime[i]) return true;else if(n%prime[i]==0) return false;ll x=n-1; while(~x&1) x>>=1;for(int i=0;i<times;i++){ll a=prime[i]; ll now=quick_pow(a,x,n); ll last;if(x==n-1){if(now!=1) return false;}else{while(x!=n-1){x<<=1; last=now; now=mult(now,now,n);if(now==1){if(last!=1 && last!=n-1) return false;break;}}}}return true;
}

之前看到的别人的板子

ll mult(ll a, ll b, ll p) {a %= p;b %= p;ll res = 0, tmp = a;while(b) {if (b & 1) {res += tmp;if (res > p) res -= p;}tmp <<= 1;if (tmp > p) tmp -= p;b >>= 1;}return res;
}ll quick_pow(ll a, ll b, ll p) {ll res = 1;ll tmp = a % p;while(b) {if (b & 1) res = mult(res, tmp, p);tmp = mult(tmp, tmp, p);b >>= 1;}return res;
}bool check(ll a, ll n, ll x, ll t) {ll res = quick_pow(a, x, n);ll last = res;for (ll i = 1; i <= t; i++) {res = mult(res, res, n);if (res == 1 && last != 1 && last != n-1) return true;last = res;}return res != 1;
}bool Miller_Rabin(ll n) {if (n < 2) return false;if (n == 2) return true;if ((n & 1) == 0) return false;ll x = n-1;ll t = 0;while ((x & 1) == 0) {x >>= 1;t++;}srand(time(NULL));const ll tims = 8;for (ll i = 0; i < tims; i++) {ll a = rand() % (n-1) + 1;if (check(a, n, x, t)) return false;}return true;
}