soj 1454. Pseudoprime numbers

描述：

费马小定理指出 for any prime number p and for any integer a > 1, a^p == a (mod p)。

对于合数p，如果存在1 < a < p满足上式，那么就称 p 为以 a 为底的伪素数。

给出p和a，问p是不是以a为底的伪素数。2 < p ≤ 1,000,000,000 and 1 < a < p

思路：

先用Miller-Rabin测试p是不是素数，再测试a^p%p是否等于a。a^p%p要用快速指数幂和二分相乘防止溢出。

代码：

#include <cstdio>
#include <ctime>
#include <cstdlib>
using namespace std;
#define Times 12
typedef long long LL;

LL random(LL n)
{
	return (LL)((double)rand()/RAND_MAX * n + 0.5);
}
LL multi(LL a, LL b, LL &m) // (a * b) % m
{
	LL ret = 0;
	while (b)
	{
		if (b & 1) ret = (ret + a) % m;
		b >>= 1;
		a = (a<<1) % m;
	}
	return ret;
}
LL pow(LL a, LL b, LL &m) // (a ^ b) % m
{
	LL ans = 1;
	a %= m;
	while (b)
	{
		if (b & 1) ans = multi(ans, a, m);
		b >>= 1;
		a = multi(a,a,m);
	}
	return ans;
}
bool Witness(LL a, LL n) // if n is base-a pseudoprime
{
	LL m = n-1;
	int j = 0;
	while (! (m & 1))
	{
		j ++;
		m >>= 1;
	}
	LL x = pow(a, m, n);
	if (x == 1 || x == n-1) return false;
	while (j --)
	{
		x = multi(x, x, n);
		if (x == n-1) return false;
	}
	return true;
}
bool Miller_Rabin(LL n) // if n is a prime
{
	if (n < 2) return false;
	if (n == 2) return true;
	if (! (n&1)) return false;
	for (int i = 1; i <= Times; ++ i)
	{
		LL a = random(n-2) + 1;
		if (Witness(a, n)) return false;
	}
	return true;
}
bool PseudoPrime(LL a, LL n)
{
	return pow(a, n, n) == a;
}
int main()
{
	LL n, a;
	while (scanf("%lld %lld", &n, &a))
	{
		if (n == 0 && a == 0) break;
		if (!Miller_Rabin(n) && PseudoPrime(a, n)) printf("yes\n");
		else printf("no\n");
	}
}