一、思路
对于输入的一个数字n和a，我们用快速幂判断 n ^ a % n 是否等于n，如果不等于直接输出no，等于的话，再判断n是否为素数，如果n为素数，输出no，否则输出yes。
判断素数的话，对于30000以下的，我采用埃氏筛法打表，对于3000以上的，我用米勒拉宾判断。
米勒拉宾算法，就是求出一组r和d，使得 d * 2 * r = n-1，产生小于n大于0的随机数m，判断m ^ (n - 1) == 1 (mod n)是否成立，成立的话，再定义i从0到r-1循环，判断每一个 (m ^ d)^(2^i)是否为1或者p-1，如果是的话，直接判断过，如果都不是，那么n不是素数。
二、代码
#include <iostream>
#include <algorithm>
using namespace std;
typedef unsigned long long ll;
ll r, d;
int primeLimit = 30000;
bool isPrime[30007];
void sieve()
{
for (int i = 0; i <= primeLimit; i++)
{
isPrime[i] = true;
}
isPrime[0] = false;
isPrime[1] = false;
for (int i = 1; i * i <= primeLimit; i++)
{
if (!isPrime[i])
{
continue;
}
for (int j = 2 * i; j <= primeLimit; j += i)
{
isPrime[j] = false;
}
}
}
void getRd(ll p)
{
d = p - 1;
r = 0;
while (!(d & 1))
{
d = d >> 1;
r++;
}
}
ll getLongRand(ll p)
{
int x = 0;
while (x == 0)
{
x = rand();
if (x < 0)
{
x = x * (-1);
}
if (x >= p)
{
x = x % p;
}
}
return (ll)x;
}
ll mulMod(ll a, ll b, ll mod)
{
ll res = 0;
while (b)
{
if (b & 1)
{
res = (res + a) % mod;
}
a = (a << 1) % mod;
b = b >> 1;
}
return res;
}
ll powMod(ll a, ll b, ll mod)
{
ll res = 1;
while (b)
{
if (b & 1)
{
res = mulMod(res, a, mod);
}
a = mulMod(a, a, mod);
b = b >> 1;
}
return res;
}
bool millerRabin(ll p)
{
ll a = getLongRand(p);
if (powMod(a, p - 1, p) != 1)
{
return false;
}
ll k = powMod(a, d, p);
for (int i = 0; i < r; i++)
{
ll parameter = powMod(k, powMod(2, i, p), p);
if (parameter == 1 || parameter == p - 1)
{
return true;
}
}
return false;
}
bool multipleMillerRabin(ll p)
{
getRd(p);
for (int i = 0; i < 20; i++)
{
if (!millerRabin(p))
{
return false;
}
}
return true;
}
bool judgePrime(ll p)
{
if (p <= primeLimit)
{
return isPrime[(int)p];
}
else if (!(p & 1))
{
return false;
}
else
{
return multipleMillerRabin(p);
}
}
bool fermet(ll a, ll p)
{
return powMod(a, p, p) == (a % p);
}
int main()
{
sieve();
ll a, p;
while (true)
{
scanf("%lld%lld", &p, &a);
if (a == 0 && p == 0)
{
break;
}
if (fermet(a, p) && !judgePrime(p))
{
printf("yes\n");
}
else
{
printf("no\n");
}
}
return 0;
}三、我对米勒拉宾算法个人小小理解
// 这里引入费马小定理， 当 p 时一个素数时，对任意的整数a，一定有 pow ( a , p ) % p == a % p
// 当 1 <= a <= p - 1 我们对等式两边同时除以 a，得到 pow ( a , p - 1 ) % p == 1
// 所以只要对于大整数 p 和随机数 a ，如果不满足这个规则，那么p一定不是素数
// 为了提高算法的准确性，我们把 pow ( a , p - 1 ) % p == 1 进行推导，这里用 a ^ (p-1) 代表a的 p - 1 次方
// a ^ ( p - 1 ) % p == 1
// 定义 r d，使得 ( 2 ^ r ) * d = p - 1
// a ^ ( ( 2 ^ r ) * d ) % p == 1
// ( a ^ ( ( 2 ^ r ) * d ) - 1 ) % p == 0
// 这里可以使用平方差公式进行推导
// (a ^ ( (2 ^ (r - 1) ) * d) - 1) * (a ^ ( (2 ^ (r - 1) ) * d) + 1 ) % p == 0
// (a ^ ( (2 ^ (r - 2) ) * d) - 1) * (a ^ ( (2 ^ (r - 2) ) * d) + 1 ) * (a ^ ( (2 ^ (r - 1) ) * d) + 1) % p == 0
// ....
// (a ^ ( (2 ^ (r - r) ) * d) - 1) * (a ^ ( (2 ^ (r - r) ) * d) + 1) * ... (a ^ ( (2 ^ (r - 1) ) * d) + 1) % p == 0
// 同时p是素数，如果这些项相乘 % p == 0，那么其中的某一项一定等于0，那么就推出了另外一个判断依据
// 在 0 <= i < r 必定存在 a ^ ( (2 ^ i) * d) 等于 p 或等于 1，否则 p 一定不是素数
// a ^ ( (2 ^ i) * d) = (a ^ d) ^ (2 ^ i)
四、运行情况
