【转】线性筛素数，线性筛欧拉函数，求前N个数的约数个数

转自http://blog.csdn.net/ji414341055/article/details/5771066

先来最基本的线性筛素数，以后的算法其实都是基于这个最基本的算法：

[cpp]  view plain copy

1. #include<stdio.h>  
2. #include<string.h>  
3. #define M 10000000  
4. int prime[M/3];  
5. bool flag[M];  
6. void get_prime()  
7. {  
8.     int i,j,k;  
9.         memset(flag,false,sizeof(flag));  
10.         k=0;  
11.         for(i=2;i<M;i++){  
12.             if(!flag[i])                              
13.             prime[k++]=i;  
14.             for(j=0;j<k&&i*prime[j]<M;j++){  
15.                     flag[i*prime[j]]=true;              
16.                     if(i%prime[j]==0)               
17.                         break;  
18.             }  
19.         }  
20. }  
21. int main()  
22. {}  

利用了每个合数必有一个最小素因子，每个合数仅被它的最小素因子筛去正好一次，所以是线性时间。
代码中体现在： if(i%prime[j]==0) break;
----------------------------------------------------------------------- 我是低调的分割线 ------------------------------------------------------------------------------------------
然后可以利用这种线性筛法求欧拉函数，需要用到以下几个性质：
//(1) 若(N%a==0 && (N/a)%a==0) 则有:E(N)=E(N/a)*a;
//(2) 若(N%a==0 && (N/a)%a!=0) 则有:E(N)=E(N/a)*(a-1);  
其中a是N的质因数。
关于欧拉函数还有以下性质：
(1) phi[p]=p-1;  (p为素数);
(2)若N=p^n(p为素数)，则 phi[N]=(p-1)*p^(n-1);
关于欧拉函数，Wiki有很详细的介绍。

[cpp]  view plain copy

1. #include<stdio.h>  
2. #include<string.h>  
3. #define M 10000000  
4. int prime[M/3],phi[M];  
5. bool flag[M];  
6. void get_prime()  
7. {  
8.     int i,j,k;  
9.     memset(flag,false,sizeof(flag));  
10.     k=0;  
11.     for(i=2;i<M;i++){  
12.         if(!flag[i]){                              
13.             prime[k++]=i;  
14.             phi[i]=i-1;  
15.         }  
16.         for(j=0;j<k&&i*prime[j]<M;j++){  
17.             flag[i*prime[j]]=true;              
18.             if(i%prime[j]==0){  
19.                 phi[i*prime[j]]=phi[i]*prime[j];  
20.                     break;  
21.             }  
22.             else  
23.                 phi[i*prime[j]]=phi[i]*(prime[j]-1);  
24.         }  
25.         }  
26. }  
27. int main()  
28. {}  

-----------------------------------------------------------------------我 是低调的分割线 -----------------------------------------------------------------------------------------
求约数个数略微复杂一点，但大体还是那个意思。
约数个数的性质，对于一个数N，N=p1^a1 + p2^a2 + ... + pn^an。其中p1 ,p2, p3... pn是N的质因数，a1 ,a2, a2,...an为相应的指数，则
                                                           div_num[N]=(p1+1)*(p2+1)*(p3+1)* ... *(pn+1);
结合这个算法的特点，在程序中如下运用：
  对于div_num：

(1)如果i|prime[j] 那么 div_num[i*prime[j]]=div_sum[i]/(e[i]+1)*(e[i]+2)                  //最小素因子次数加1
(2)否则 div_num[i*prime[j]]=div_num[i]*div_num[prime[j]]                                     //满足积性函数条件

  对于e：

(1)如果i|pr[j]  e[i*pr[j]]=e[i]+1; //最小素因子次数加1
(2)否则 e[i*pr[j]]=1;              //pr[j]为1次

 
 

  
  

   
   

    
    [cpp] 
    
    view plain
    
    copy
    
    

   
   
  
  
  
  
#include<stdio.h>  
#include<string.h>  
#define M 10000000  
int prime[M/3],e[M/3],div_num[M];           // e[i]表示第i个素数因子的个数  
bool flag[M];  
void get_prime()  
{  
    int i,j,k;  
    memset(flag,false,sizeof(flag));  
    k=0;  
    for(i=2;i<M;i++){  
        if(!flag[i]){                              
                    prime[k++]=i;  
                    e[i]=1;  
                div_num[i]=2;                       //素数的约数个数为2  
            }  
            for(j=0;j<k&&i*prime[j]<M;j++){  
                    flag[i*prime[j]]=true;              
                    if(i%prime[j]==0){  
                        div_num[i*prime[j]]=div_num[i]/(e[i]+1)*(e[i]+2);  
                        e[i*prime[j]]=e[i]+1;  
                        break;  
                    }  
                    else{  
                        div_num[i*prime[j]]=div_num[i]*div_num[prime[j]];  
                        e[i*prime[j]]=1;  
                    }  
            }  
        }  
}  
int main()  
{}