洲阁筛

洲阁筛真是个不错的暴力啊。。

简单的写了个求1~n质数个数

-O2  N=1e11要2.5s

 1 #pragma GCC optimize(2)
 2 #include <bits/stdc++.h>
 3 using namespace std;
 4 #define LL long long
 5 int L0[1000005],L1[1000005],sn,m,P[1000005],x,ga[1000005],gb[1000005];
 6 LL n,a[1000005],b[1000005];
 7 int main(){
 8     n=100000000000;
 9 //    scanf("%lld",&n);
10     sn=sqrt(n);
11     for (int i=2;i<=sn;++i){
12         if (!P[i]) P[++m]=i;
13         for (int j=1;j<=m;++j){
14             x=P[j]*i;
15             if (x>sn) break;
16             P[x]=1;
17             if (i%P[j]==0) break;
18         }
19     }
20     for (int i=1;i<=sn;++i) a[i]=i,b[i]=n/i;
21     for (int i=1,j=1;i<=sn;++i){
22         while (j<=m&&P[j]*P[j]<=i) ++j;
23         L0[i]=j;
24     }
25     for (int i=sn,j=L0[sn];i;--i){
26         while (j<=m&&P[j]*P[j]<=b[i]) ++j;
27         L1[i]=j;
28     }
29     for (int i=1;i<=m;++i){
30         for (int j=1;j<=sn&&i<L1[j];++j){
31             LL y=j*P[i]; gb[j]=i;
32             if (y<=sn) b[j]-=b[y]+gb[y]+1-i;
33             else b[j]-=a[n/y]+ga[n/y]+1-i;
34         }
35         for (int j=sn;j&&i<L0[j];--j){
36             LL y=j/P[i]; ga[j]=i;
37             a[j]-=a[y]+ga[y]+1-i;
38         }
39     }
40     printf("%lld\n",b[1]+m-1);
41     return 0;
42 }
go ahead

 

外加一个模板题 http://www.spoj.com/problems/DIVCNT3/

代码:

带上求f的那一部分,N=1e11要12s多

不开O2可能是要死的吧。。

 1 #pragma GCC optimize(2)
 2 #include <bits/stdc++.h>
 3 using namespace std;
 4 #define LL long long
 5 int T,W[340005],L0[340005],L1[340005],L2[340005],sn,m,P[340005],x,ga[340005],gb[340005];
 6 LL n,a[340005],b[340005],F[340005];
 7 void hh(){
 8     scanf("%lld",&n);
 9     sn=sqrt(n); F[1]=1; m=0;
10     for (int i=2;i<=sn;++i){
11         if (!P[i]) P[++m]=i,F[i]=4,W[i]=1;
12         for (int j=1;j<=m;++j){
13             x=P[j]*i;
14             if (x>sn) break;
15             P[x]=1;
16             if (i%P[j]==0){
17                 W[x]=W[i]+1;
18                 F[x]=F[i]/(3*W[i]+1)*(3*W[x]+1);
19                 break;
20             }else{
21                 W[x]=1; F[x]=F[i]*4;
22             }
23         }
24     }
25     for (int i=1;i<=sn;++i)
26         a[i]=i,b[i]=n/i,ga[i]=gb[i]=0;
27     for (int i=1,j=1;i<=sn;++i){
28         while (j<=m&&P[j]*P[j]<=i) ++j;
29         L0[i]=j;
30     }
31     for (int i=sn,j=L0[sn];i;--i){
32         while (j<=m&&1ll*P[j]*P[j]<=b[i]) ++j;
33         L1[i]=j;
34     }
35     for (int i=1,j=1;i<=sn;++i){
36         while (j<=m&&P[j]<=i) ++j;
37         L2[i]=j;
38     }
39     for (int i=1;i<=m;++i){
40         for (int j=1;j<=sn&&i<L1[j];++j){
41             LL y=j*P[i]; gb[j]=i;
42             if (y<=sn) b[j]-=b[y]+gb[y]+1-i;
43             else b[j]-=a[n/y]+ga[n/y]+1-i;
44         }
45         for (int j=sn;j&&i<L0[j];--j){
46             int y=j/P[i]; ga[j]=i;
47             a[j]-=a[y]+ga[y]+1-i;
48         }
49     }
50     LL ans=0;
51     for (int i=1;i<=sn;++i)
52         ans+=F[i]*((b[i]-m+gb[i]-1)*4+1);
53     for (int i=1;i<=sn;++i)
54         a[i]=b[i]=1,ga[i]=gb[i]=m+1;
55     for (int i=m;i;--i){
56         for (int j=1;j<=sn&&i<L1[j];++j){
57             LL y=n/j;
58             if (gb[j]!=i+1) b[j]=(m-i)*4+1;
59             y/=P[i]; gb[j]=i;
60             for (int c=1;y;++c,y/=P[i])
61             if (y<=sn)
62                 if (ga[y]==i+1) b[j]+=a[y]*(3*c+1);
63                 else b[j]+=(max(0,L2[y]-i-1)*4+1)*(3*c+1);
64             else
65                 if (gb[n/y]==i+1) b[j]+=b[n/y]*(3*c+1);
66                 else b[j]+=((m-i)*4+1)*(3*c+1);
67         }
68         for (int j=sn;j&&i<L0[j];--j){
69             int y=j;
70             if (ga[j]!=i+1) a[j]=max(0,L2[y]-i-1)*4+1;
71             y/=P[i]; ga[j]=i;
72             for (int c=1;y;++c,y/=P[i])
73             if (ga[y]==i+1) a[j]+=a[y]*(3*c+1);
74             else a[j]+=(max(0,L2[y]-i-1)*4+1)*(3*c+1);
75         }
76     }
77     ans+=b[1]-(ga[sn]==1?a[sn]:m*4+1);
78     for (int i=1;i<=sn;++i) P[i]=0;
79     printf("%lld\n",ans);
80 }
81 int main(){
82     scanf("%d",&T);
83     while (T--) hh();
84     return 0;
85 }
鸣狐

 

posted @ 2017-12-30 12:10  cyz666  阅读(966)  评论(0编辑  收藏  举报