使用试除法得到质因数的复杂度是\(O(\sqrt{n})\),如果数字较大比如\(10^{18}\)就束手无策了。
在上一篇文章有写到Miller Rabin可以在\(O(klog^{2}{n})\)的复杂度下判断一个数字是不是素数,这次要介绍的Pollard Rho(普拉德柔)算法就是结合Miller Rabin的快速素性检验来进行的一个复杂度大约为\(O(n^{\frac{1}{4}})\)来分解一个数字的所有质因数的算法。为什么说是大约呢,因为真不好算哈哈哈。
介绍 / Introduction
同样不深究原理,会用好模板就行了。
如果想要了解原理的请看:Pollard_Rho - ZenyZ - 博客园 (cnblogs.com)
简单来说就是随机找到一个因子(不一定是质因子),利用生日悖论来提高找到因子的概率,用倍增来提高找因子的速度,用累乘gcd来优化常数,不断地分治进行,如果到了一个质数就停止递归并将其加入到因子列表中。
例题:P4718 【模板】Pollard-Rho算法 - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)
模板 / Templete
#include <bits/stdc++.h> #define int long long using namespace std; //四个工具函数 int qmi(int a, int b, int p) // 快速幂 { int res = 1; while (b)// 注意!中间结果可能溢出,需要使用__int128过渡 { if (b & 1)res = (__int128)res * a % p; a = (__int128)a * a % p, b >>= 1; } return res; } int getabs(int x){return x < 0 ? -x : x;}//绝对值 int f(int t,int c,int x){return ((__int128)t * t + c) % x;}//随机化函数 int gcd(int a,int b){return b == 0 ? a : gcd(b, a % b);}//求最大公因数 //判断x是否是质数 bool isprime(int x)//Miller_Rabin { //特判 if (x < 3)return x == 2; if (x % 2 == 0)return false; //用下面这几个数字判断不会出错 int A[] = {2, 325, 9375, 28178, 450775, 9780504, 1795265022}, d = x - 1, r = 0; while (d % 2 == 0)d >>= 1, ++r; for (auto a : A) { int v = qmi(a, d, x); if (v <= 1 || v == x - 1)continue; for (int i = 0; i < r; ++i) { v = (__int128)v * v % x; if (v == x - 1 && i != r - 1) { v = 1; break; } if (v == 1)return false; } if (v != 1)return false; } return true; } //找出x的一个因子 int Pollard_Rho(int x) { int s = 0, t = 0, c = 1LL * rand() % (x - 1) + 1;//[1, x - 1] for(int i = 1; true; i <<= 1, s = t) { int mul_sum = 1, g = 1; for(int j = 1;j <= i; ++ j) { t = f(t, c, x); mul_sum = (__int128)mul_sum * getabs(t - s) % x; if(j % 127 == 0) { g = gcd(mul_sum, x); if(g > 1)return g; } } g = gcd(mul_sum, x); if(g > 1)return g; } return -1; } //分治主体 void find_factor(int x, vector<int> &v) { if(isprime(x))return v.push_back(x), void(); int yz = x; //一定要确保yz是小于x的,否则重新找 while(yz >= x)yz = Pollard_Rho(x); find_factor(yz, v), find_factor(x / yz, v); } signed main() { int N;cin >> N; for(int i = 1;i <= N; ++ i) { int x;cin >> x; if(isprime(x))cout << "Prime" << '\n'; else { vector<int> v; find_factor(x, v); //找完之后所有的质数都在v这个数组中 //v的大小不会超过70,毕竟1e18也才2^63 sort(v.begin(), v.end()); cout << v.back() << '\n'; } } return 0; }
Comments NOTHING