素数是一个经常的涉及到得内容，所以有必要整理出有关解决素数相关问题的算法
学习资料：Eratosthenes筛法和欧拉筛法对比 一般筛法求素数+快速线性筛法求素数 数学技巧之素数筛选
素数与素性测试 〖数学算法〗素性测试 请看Miller-Rabin算法! Miller-Rabin素数测试 Pollard-rho大整数分解
1. 素数筛法
代码1(埃氏筛法)：
/*
埃氏筛法：返回n以内素数的个数， 复杂度O (nloglogn)
*/
int seive(int n) { //Eratosthenes (埃氏筛法)
int p = 0;
memset (is_prime, true, sizeof (is_prime));
is_prime[0] = is_prime[1] = false;
for (int i=2; i<=n; ++i) {
if (is_prime[i]) {
prime[++p] = i;
for (int j=2*i; j<=n; j+=i) is_prime[j] = false;
}
}
return p;
}
代码2(欧拉筛法)：
/*
欧拉筛法：返回n以内素数的个数， 复杂度O (n)
*/
int seive(int n) { //Euler (欧拉筛法)
int p = 0;
memset (is_prime, true, sizeof (is_prime));
is_prime[0] = is_prime[1] = false;
for (int i=2; i<=n; ++i) {
if (is_prime[i]) prime[++p] = i;
for (int j=1; j<=p && i*prime[j]<=n; ++j) {
is_prime[i*prime[j]] = false;
if (i % prime[j] == 0) break;
}
}
return p;
}
附上比较两种筛法的测试结果(新标签查看原图)
2. 素性测试
代码1：
/*
素性测试，输入正数，复杂度O (sqrt(n))
*/
bool is_prime(int n) {
for (int i=2; i*i<=n; ++i) {
if (n % i == 0) return false;
}
return n != 1; //1例外
}
代码2：
/*
素性测试，在小范围(1e5)内判素数个数以及单个数判素数有奇效，不适用于大范围判素数
*/
bool is_prime(int x) {
if (x == 2 || x == 3) return true;
if (x % 6 != 1 && x % 6 != 5) return false;
for (int i=5; i*i<=x; i+=6) {
if (x % i == 0 || x % (i + 2) == 0) return false;
}
return true;
}
代码3(Miller_Rabin 随机算法)：
/*
素性测试，Miller_Rabin 随机算法
可以判断< 2^63的数
素数：true，合数：false
*/
const int S = 20; //随机算法判定次数，S越大，判错概率越小
//非递归写法，递归写法可能会爆栈
ll GCD(ll a, ll b) {
if (a == 0) return 1;
if (a < 0) a = -a;
while (b) {
ll c = a % b;
a = b; b = c;
}
return a;
}
//计算 (a * b) % p，a，b是long long数，直接相乘可能会溢出
ll multi_mod(ll a, ll b, ll p) {
ll ret = 0;
a %= p; b %= p;
while (b) {
if (b & 1) {
ret += a;
if (ret >= p) ret -= p;
}
a <<= 1;
if (a >= p) a -= p;
b >>= 1;
}
return ret;
}
//计算 a ^ x % p
ll pow_mod(ll a, ll x, ll p) {
ll ret = 1;
a %= p;
while (x) {
if (x & 1) ret = multi_mod (ret, a, p);
a = multi_mod (a, a, p);
x >>= 1;
}
return ret;
}
/*
以a为基,n-1=x*2^t，a^(n-1) = 1(mod n) 验证n是不是合数
一定是合数返回true, 不一定返回false
*/
bool check(ll a, ll n, ll x, int t) {
ll ret = pow_mod (a, x, n);
ll last = ret;
for (int i=1; i<=t; ++i) {
ret = multi_mod (ret, ret, n);
if (ret == 1 && last != 1 && last != n - 1) return true; //合数
last = ret;
}
if (ret != 1) return true;
return false;
}
bool Miller_Rabin(ll n) {
if (n == 2) return true;
if (n < 2 || ! (n & 1)) return false; //偶数或1
ll x = n - 1; int t = 0;
while (! (x & 1)) {
x >>= 1; t++;
}
for (int i=1; i<=S; ++i) {
ll a = rand () % (n - 1) + 1; //需要cstdlib头文件
if (check (a, n, x, t)) return false; //合数
}
return true;
}
3. 整数分解
代码1：
/*
整数分解，试除法进行质因子分解，复杂度O(sqrt(n))
*/
map<int, int> factorize(int n) {
map<int, int> res;
for (int i=2; i*i<=n; ++i) {
while (n % i == 0) {
++res[i]; n /= i;
}
}
if (n != 1) res[n] = 1;
return res;
}
代码2：
/*
整数分解，先打个素数表优化试除法
*/
vector<int> factorize(int n) {
int p = seive (100000);
vector<int> ret;
for (int i=1; i<=p && prime[i]<=n/prime[i]; ++i) {
while (n % prime[i] == 0) {
n /= prime[i];
ret.push_back (prime[i]);
}
}
if (n != 1) ret.push_back (n);
return ret;
}
代码3(Pollard_rho 随机算法)：
/*
大整数分解，Pollard_rho 随机算法
factorize ()保存质因数在vector
*/
ll Pollard_rho(ll x, ll c) {
ll i = 1, k = 2;
ll a = rand () % x;
ll b = a;
while (1) {
i++;
a = (multi_mod (a, a, x) + c) % x;
ll d = GCD (b - a, x);
if (d != 1 && d != x) return d;
if (b == a) return x;
if (i == k) b = a, k += k;
}
}
void factorize(ll n, vector<ll> &ret) { //ret保存质因数，无序
if (Miller_Rabin (n)) { //素数
ret.push_back (n); return ;
}
ll p = n;
while (p >= n) p = Pollard_rho (p, rand () % (n - 1) + 1);
factorize (p, ret);
factorize (n / p, ret);
}
int main(void) {
srand (time (NULL)); //需添加ctime头文件
int T; scanf ("%d", &T);
while (T--) {
ll n; scanf ("%I64d", &n);
if (Miller_Rabin (n)) puts ("Prime");
else {
if (n <= 1) { //1的情况特判，否则RE
puts ("-1"); continue;
}
vector<ll> ans;
factorize (n, ans);
sort (ans.begin (), ans.end ());
for (int i=0; i<ans.size (); ++i) {
printf ("%I64d%c", ans[i], (i == ans.size ()-1) ? '
' : ' ');
}
}
}
return 0;
}
其他：
代码：
/*
约数枚举，复杂度O (sqrt(n))
*/
vector<int> divisor(int n) {
vector<int> ret;
for (int i=1; i*i<=n; ++i) {
if (n % i == 0) {
ret.push_back (i);
if (n / i != i) ret.push_back (n / i);
}
}
return res;
}