内蒙古建设工程造价管理网站,石青网站推广软件,织金网站建设,wordpress怎么制作响应式一、引言
在平时做题目或者进行运算的时候#xff0c;素数的出现次数总是十分频繁。这里总结了一些常见的判定素数和计算某个范围内素数个数的一些算法。部分代码来源于 kuangbin 的模板#xff0c;嗯毕竟都是跟着这个学的... 二、朴素判断素数算法
就判断素数而言#xf…一、引言
在平时做题目或者进行运算的时候素数的出现次数总是十分频繁。这里总结了一些常见的判定素数和计算某个范围内素数个数的一些算法。部分代码来源于 kuangbin 的模板嗯毕竟都是跟着这个学的... 二、朴素判断素数算法
就判断素数而言事实上是非常简单的了。根据定义判断一个整数n是否是素数只需要去判断在整数区间[2, n-1]之内是否具有某个数m使得n % m 0。代码可以这么写 int isPrime(int n) {int i;for (i 2; i n; i) {if (n % i 0) return 0;}return 1;
}
事实上这个算法是O(n)的感觉是很快了但是依旧无法满足需求。所以有一个算法是O(sqrt(n))的算法。代码可以这么写 int isPrime(long long n) {long long i;for (i 2; i * i n; i) {if (n % i 0) return 0;}return 1;
}
原理很巧妙也仅仅是把代码的i n变成了i * i n而已但是优化却是极其高的。可能你会注意到在上一份代码里面我定义的n为int类型而后面一份代码我却定义成了long long这是因为在1s内上一份代码能判断出来的数量级为1e8而后面一份代码判断出来的数量级却几乎可以达到1e16。而原因仅仅是因为a * b c的话一旦a是c的约数那么b也是如果a不是那么b也不是。 这个方法也可以称作试除法。 三、Miller_Rabin素性测试
尽管上面的O(sqrt(n))的算法已经非常优秀了但是面对更高数量级的“大数”却会显得力不从心。而这个时候Miller_Rabin的优越性就显示出来了。
Miller_Rabin的理论基础来源于费马小定理。值得一提的是费马小定理是数论四大定理之一。
费马小定理n是一个奇素数a是任何整数(1≤ a≤n-1) 则 a^(n-1)≡1(mod n)
要测试n是否是一个素数首先将n-1分解为(2^s) * d在每次测试开始时先随机选择一个介于[1, N - 1]的整数a之后如果对于所有的r∈[0, s - 1]若a^d mod N ≠ 1且a^((2^r) * d) mod N ≠ -1那么n就是一个合数否则n有3/4的几率是素数。
这也是为什么说这个算法只是素性测试了他不能完全保证结果正确但是当测试次数大约为20的时候基本可以确保正确率达到97%以上。
它的代码可以这么写 const int S 20;
LL mod_mul(LL a, LL b, LL n) {LL res 0;while (b) {if (b 1) res (res a) % n;a (a a) % n;b 1;}return res;
}
LL mod_exp(LL a, LL b, LL n) {LL res 1;while (b) {if (b 1) res mod_mul(res, a, n);a mod_mul(a, a, n);b 1;}return res;
}
bool miller_rabin(LL n) {if (n 2 || n 3 || n 5 || n 7 || n 11) return true;if (n 1 || !(n % 2) || !(n % 3) || !(n % 5) || !(n % 7) || !(n % 11)) return false;LL x, pre, u n - 1, k 0;while (!(u 1)) {k;u 1;}srand((LL)time(NULL));for (int i 0; i S; i) { //进行S次测试S越大结果越准确x rand() % (n - 2) 2; //在[2, n)中取随机数if (x % n 0) continue;x mod_exp(x, u, n); //计算x^u % npre x;for (int j 0; j k; j) {x mod_mul(x, x, n);if (x 1 pre ! 1 pre ! n - 1)return false;pre x;}if (x ! 1) return false;}return true;
} 四、筛选法
上面介绍的一些素数判断的算法在某些问题是基本可以适用的了。但是对于另外一类问题却十分尴尬。比如问你整数区间[1, n]内的素数个数是多少。这个时候如果一个个枚举一个个判断对于朴素方法来说最优也是O(nsqrt(n))即使是Miller_Rabin算法也无法在O(n)的时间内得到结果。于是就有了埃氏筛选法埃拉托斯特尼筛法。
对于筛选整数n以内的素数算法是这么描述的先把素数2的倍数全部删除剩下的数第一个为3再把素数3的倍数全部删除剩下的第一个数为5再把素数5的倍数全部删除······直到把n以内最后一个素数的倍数删除完毕得到n以内的所有素数。代码可以这么写 const int maxn 1e7 5;
int pri[maxn];
void getPrime(int n) {for (int i 0; i n; i) pri[i] i;pri[1] 0;for (int i 2; i n; i) {if (!pri[i]) continue;pri[pri[0]] i;for (int j 2; i * j n j n; j) {pri[i * j] 0;}}
}
而事实上观察可以发现上面的这种写法有很多次重复计算这样显然无法做到线性筛选而另外一种写法却可以得到线性筛选达到时间复杂度为O(n)代码可以这么写 const int MAXN 1e7 5;
void getPrime(){memset(prime, 0, sizeof(prime));for (int i 2;i MAXN;i) {if (!prime[i])prime[prime[0]] i;for (int j 1;j prime[0] prime[j] MAXN / i;j) {prime[prime[j] * i] 1;if (i%prime[j] 0) break;}}
}
来自kuangbin的模板。 五、容斥原理
从上面的代码可以发现显然这种筛法只能应付达到1e7这种数量级的运算即使是线性的筛选法也无法满足因为在ACM竞赛中1e8的内存是极有可能获得Memery Limit Exceed的。
于是可以考虑容斥原理。
以AHUOJ 557为例1e8的情况是筛选法完全无法满足的但是还是考虑a * b c的情况1e8只需要考虑10000以内的素数p[10000]然后每次先减去n / p[i]再加上n / (p[i] * p[j])再减去n / (p[i] * p[j] * p[k])以此类推...于是就可以得到正确结果了。
代码如下 #include cmath
#include cstdio
using namespace std;const int maxn 10005;
int sqrn, n, ans 0;
bool vis[maxn];
int pri[1500] {0};
void init(){vis[1] true;int k 0;for(int i 2; i maxn; i){if(!vis[i]) pri[k] i;for(int j 0; j k pri[j] * i maxn; j){vis[pri[j] * i] true;if(i % pri[j] 0) break;}}
}
void dfs(int num, int res, int index){for(int i index; pri[i] sqrn; i){if(1LL * res * pri[i] n){return;}dfs(num 1, res * pri[i], i1);if(num % 2 1){ans - n / (res * pri[i]);}else{ans n / (res * pri[i]);}if(num 1) ans;}
}
int main(){init();while(~scanf(%d,n) n){ans n;sqrn sqrt((double)n);dfs(1,1,0);printf(%d\n,ans-1);}return 0;
}
六、Meissel-Lehmer算法 最后介绍的这个算法可以说是黑科技级别的它能够在3s内求出1e11之内的素数个数。据说还有在300ms内求出1e11的个数的。可以参考wiki里面原理。然后给出来自Codeforces 665F题目里面的代码。 #define MAXN 100 // pre-calc max n for phi(m, n)
#define MAXM 10010 // pre-calc max m for phi(m, n)
#define MAXP 40000 // max primes counter
#define MAX 400010 // max prime
#define setbit(ar, i) (((ar[(i) 6]) | (1 (((i) 1) 31))))
#define chkbit(ar, i) (((ar[(i) 6]) (1 (((i) 1) 31))))
#define isprime(x) (( (x) ((x)1) (!chkbit(ar, (x)))) || ((x) 2))namespace pcf {long long dp[MAXN][MAXM];unsigned int ar[(MAX 6) 5] { 0 };int len 0, primes[MAXP], counter[MAX];void Sieve() {setbit(ar, 0), setbit(ar, 1);for (int i 3; (i * i) MAX; i, i) {if (!chkbit(ar, i)) {int k i 1;for (int j (i * i); j MAX; j k) setbit(ar, j);}}for (int i 1; i MAX; i) {counter[i] counter[i - 1];if (isprime(i)) primes[len] i, counter[i];}}void init() {Sieve();for (int n 0; n MAXN; n) {for (int m 0; m MAXM; m) {if (!n) dp[n][m] m;else dp[n][m] dp[n - 1][m] - dp[n - 1][m / primes[n - 1]];}}}long long phi(long long m, int n) {if (n 0) return m;if (primes[n - 1] m) return 1;if (m MAXM n MAXN) return dp[n][m];return phi(m, n - 1) - phi(m / primes[n - 1], n - 1);}long long Lehmer(long long m) {if (m MAX) return counter[m];long long w, res 0;int i, a, s, c, x, y;s sqrt(0.9 m), y c cbrt(0.9 m);a counter[y], res phi(m, a) a - 1;for (i a; primes[i] s; i) res res - Lehmer(m / primes[i]) Lehmer(primes[i]) - 1;return res;}
}
七、后记 写这么多...发现我会的也就这么点了....找了不少的模板...嘛...不过把这些知识梳理一遍思路也变得清晰许多。以后继续加油~~~
