如何快速求出区间范围[0,x]内的素数?
Brute Force 总是那么的让人觉得亲近。
枚举每一个数,判断该数是否为素数?
bool isPrime( int x ) { if(x < 2) return 0; if(x == 2) return 1; for ( int i=2; i*i<=x; ++i ) if(0 == x%i) return 0; return 1; }
// 用prime[0]存素数个数 for ( int i=0; i<=x; ++i ) if( isPrime(i) ) prime[++prime[0] ] = i;
暴力解决问题在很多时候都会存在效率低下这个问题。那么对于求素数这个问题如何提高效率呢?
其实很多数都是有公共因子的(公约数)。如果直接把这个公约数拿出来去踢出其倍数,被踢出的数都不会是素数的。例如:
2 可以踢出: 4 6 8 10 12 。。。2*x (x=2,3,4,5,6,7,8....)
3 可以踢出: 6 9 12 15 18。。。3*x(x=2,3,4,5,6,7,8....)
5 可以踢出: 10 15 20 25。。。5*x(x=2,3,4,5,6,7,8.....)
......... 依次做下去,剩余的数就是素数了。
int prime[3438]; bool isPrime[40007]; for ( int i=1; i<=x; ++i ) // 初始化所有的数都是素数 isPrime[i] = 1; for ( int i=2; i<=x; ++i ) if(isPrime[i] ) { // 该数是素数 for ( int j=i+i; j<x; j+=i ) //踢出不是素数的数 isPrime[j] = 0; } for ( int i=2; i<=x; ++i ) // 存取素数 if(isPrime[i] ) prime[++prime[0] ] = i;
这样看貌似效率是提高了不少额。可是仔细分析上面的踢出非素数的过程会发现,有些数被踢出了好几次,为什么会出现这个问题呢?
x 可以踢出: x*2 x*3 x*4 x*(x-1) x*x x*(x+1).....
不难发现其实x*2 x*3 x*4....x*(x-1)这些被更小的因子给踢出了。也就是说对于x应该踢出的数是从
x*x开始的。那么是不是这就已经优化到顶了呢? 没有哦。
开始考虑x*x是偶数,还是奇数。如果是偶数的话,x必然是偶数的。这中情况可以直接用素因子2直接把这些数踢出掉。
x*x是奇数,x必然是奇数了。 x*x+x是什么数呢? 是偶数。奇数+奇数=偶数。因为奇数%2=1,两个1相加就是偶数了。这样的话又会造成部分偶数被踢出多次。如何解决这个问题。呵呵,只需对x*x+x+x就好了,因为x*x+x是偶数 x*x+x+x,就是奇数了, 接下来的奇数就是x*x+2*(x+x)了。优化到此结束。
int prime[3438]; bool isPrime[40007]; void getPrime( int x ){ for ( int i=1; i<x; i+=2 ) isPrime[i] = 1, isPrime[i-1] = 0; // 首先踢出2的倍数 prime[prime[0]=1 ] = 2; // 存取第一个素数 for ( int i=3; ; i+=2 ) if(isPrime[i] ) { // 只考虑奇数了 int j = i*i, k = i+i; if(j >= x) break; while(j < x ) { isPrime[j] = 0; j += k; } }// 存素数 for ( int i=3; i<x; i += 2 ) if(isPrime[i]) prime[++prime[0] ] = i; }
对次求解素数已经完毕了。
接下来是将x分解成素因子的乘积
求x = p[1]^cnt[1]*p[2]^cnt[2]*p[3]^cnt[3]。。。。p[...]^cnt[...] p[]是素因子
int p[3438], cnt[3438]; void getPrimeDivisor( int x ) { p[0] = cnt[0] = 0; int t; for ( int i=1; prime[i]*prime[i]<=x && i<=prime[0]; ++i ) { t = 0; while( x%prime[i] == 0 ) { ++t; x /= prime[i]; } if( t ) p[++p[0] ] = prime[i], cnt[++cnt[0] ] = t; } if(x > 1) p[++p[0] ] = x, cnt[++cnt[0] ] = 1;
// 此时x必定是个素数。 因为在sqrt(x)内都没有找到其因子 };
此时x的因子个数=(cnt[1]+1)*(cnt[2]+1)*(cnt[3]+1)*(cnt[4]+1)....
额,只知道因子个数而已。不是很好玩ing。能不能求出x的所有约数呢? 恭喜你,这是可以滴。
利用乘法原理对约数集合乘以一个约数得到新的集合,
int divisor[1500]; void getDivisor( int x ) { // 或得x的所有因子 getPrimeDivisor(x); // 对x进行素因子分解 divisor[0] = 1; // 存取因子个数 divisor[1] = 1; for ( int i=1; i<=p[0]; ++i ) { int nowNum = divisor[0]; int base = 1; for ( int j=1; j<=cnt[i]; ++j ) { base *= p[i]; for ( int k=1; k<=divisor[0]; ++k ) divisor[++nowNum ] = divisor[k]*base; } divisor[0] = nowNum; } }
下面通过一道题来测试代码的正确性:
hdoj 4432 http://acm.hdu.edu.cn/showproblem.php?pid=4432
题意是求n的所有约数,并把每个约数转化为m进制,把每个进制为的平方加起来。再把最后的结果转化为m进制。简单吧。
#include <iostream> #include <cstdio> #include <cstring> #include <algorithm> using namespace std; int prime[3438]; bool isPrime[40007]; void getPrime( int x ){ for ( int i=1; i<x; i+=2 ) isPrime[i] = 1, isPrime[i-1] = 0; prime[prime[0]=1 ] = 2; for ( int i=3; ; i+=2 ) if(isPrime[i] ) { int j = i*i, k = i+i; if(j >= x) break; while(j < x ) { isPrime[j] = 0; j += k; } } for ( int i=3; i<x; i+=2 ) if(isPrime[i]) prime[++prime[0] ] = i; } int p[3438], cnt[3438]; void getPrimeDivisor( int x ) { p[0] = cnt[0] = 0; int t; for ( int i=1; prime[i]*prime[i]<=x && i<=prime[0]; ++i ) { t = 0; while( x%prime[i] == 0 ) { ++t; x /= prime[i]; } if(t ) p[++p[0] ] = prime[i], cnt[++cnt[0] ] = t; } if(x > 1) p[++p[0] ] = x, cnt[++cnt[0] ] = 1; }; int divisor[1500]; void getDivisor( int x ) { getPrimeDivisor(x); divisor[0] = 1; divisor[1] = 1; for ( int i=1; i<=p[0]; ++i ) { int nowNum = divisor[0]; int base = 1; for ( int j=1; j<=cnt[i]; ++j ) { base *= p[i]; for ( int k=1; k<=divisor[0]; ++k ) divisor[++nowNum ] = divisor[k]*base; } divisor[0] = nowNum; } } int get( int x, int bit ){ int ret = 0; while(x ) { ret += (x%bit)*(x%bit); x /= bit; } return ret; }; void print( int x, int bit ) { int a[100]={0}; while( x ) { a[++a[0] ] = x%bit; x /= bit; } for ( int i=a[0]; i>0; --i) if(a[i] > 9) printf("%c", a[i]-10+'A'); else printf("%d", a[i]); puts(""); } int main(){ getPrime(32000); int n, m; while( ~scanf("%d%d", &n, &m) ) { getDivisor(n); int ans = 0; for ( int i=1; i<=divisor[0]; ++i ) { ans += get(divisor[i], m); } print(ans, m); } }
关于素数的运用还有很多的,比如求数x的所有约数的和。初等数论还是挺有意思的,可惜没时间整了,悲剧啊。