BZOJ2301 莫比乌斯反演
题意:a<=x<=b,c<=y<=d,求满足gcd(x,y)=k的数对(x,y)的数量 ((x,y)和(y,x)不算同一个)
比hdu1695多加了个下界,还有顺序不一样的也算上了。
因为G(x,y)本来就是顺序不一样的算不同方案,所以这题的公式就是:
Ans=G(b/k,d/k)-G((a-1)/k,d/k)-G(b/k,(c-1)/k)+G((a-1)/k,(c-1)/k)
但是本题数据很大,直接计算会TLE,
有一个优化:http://www.cnblogs.com/zhsl/p/3269288.html
//calc G(a,b) LL _G(int a,int b) //朴素算法 { int tx=min(a,b),ty=max(a,b); LL ans = 0; for(int i = 1; i <= tx; i++) ans += (LL)mu[i]*(tx/i)*(ty/i); return ans; } LL G(int n,int m) //加分块优化 { LL ans = 0; if(n > m) swap(n,m); for(int i = 1, la = 0; i <= n; i = la+1) { la = min(n/(n/i),m/(m/i)); ans += (LL)(msum[la] - msum[i-1])*(n/i)*(m/i); //事先预处理:msum[n]=SUM(mu[1..n]) } return ans; }
不知为什么自己测AC了结果bzoj上一直RuntimeError......要完蛋了orz,改成scanf和printf就过了
1 #include <iostream> 2 #include <cstring> 3 #include <cmath> 4 #include <cstdio> 5 using namespace std; 6 #define LL long long 7 #define MMX 50100 8 LL mu[MMX],msum[MMX]; 9 bool check[MMX]; 10 int prime[MMX]; 11 int a,b,c,d,k,T; 12 13 void Moblus() 14 { 15 memset(check,false,sizeof(check)); 16 mu[1] = 1; 17 int tot = 0; 18 for(int i = 2; i <= MMX; i++) 19 { 20 if( !check[i] ) 21 { 22 prime[tot++] = i; 23 mu[i] = -1; 24 } 25 for(int j = 0; j < tot; j++) 26 { 27 if(i * prime[j] > MMX) break; 28 check[i * prime[j]] = true; 29 if( i % prime[j] == 0) 30 { 31 mu[i * prime[j]] = 0; 32 break; 33 } 34 else 35 { 36 mu[i * prime[j]] = -mu[i]; 37 } 38 } 39 } 40 msum[1]=mu[1]; 41 for (int i=2;i<=MMX;i++) 42 msum[i]=msum[i-1]+mu[i]; 43 } 44 45 //calc G(a,b) 46 LL _G(int a,int b) //朴素算法 47 { 48 int tx=min(a,b),ty=max(a,b); 49 LL ans = 0; 50 for(int i = 1; i <= tx; i++) 51 ans += (LL)mu[i]*(tx/i)*(ty/i); 52 return ans; 53 } 54 55 LL G(int n,int m) //加分块优化 56 { 57 LL ans = 0; 58 if(n > m) swap(n,m); 59 for(int i = 1, la = 0; i <= n; i = la+1) 60 { 61 la = min(n/(n/i),m/(m/i)); 62 ans += (LL)(msum[la] - msum[i-1])*(n/i)*(m/i); //事先预处理:msum[n]=SUM(mu[1..n]) 63 } 64 return ans; 65 } 66 67 int main() 68 { 69 //freopen("in.txt","r",stdin); 70 //freopen("out.txt","w",stdout); 71 72 cin>>T; 73 Moblus(); 74 while (T--) 75 { 76 scanf("%d%d%d%d%d",&a,&b,&c,&d,&k); 77 //cout<<G(b/k,d/k)<<" "<<G((a-1)/k,d/k)<<" "<<G(b/k,(c-1)/k)<<" "<<G((a-1)/k,(c-1)/k)<<endl; 78 LL res=G(b/k,d/k)-G((a-1)/k,d/k)-G(b/k,(c-1)/k)+G((a-1)/k,(c-1)/k); 79 printf("%lld\n",res); 80 } 81 return 0; 82 }
实测:
朴素算法:
1005 | root | 1002 | Accepted | 380K | 47400MS | G++ | 1.6K | 2014-11-15 15:37:03 |
优化算法:
1014 | root | 1002 | Accepted | 952K | 6011MS | G++ | 1.81K | 2014-11-15 17:17:59 |
posted on 2014-11-15 15:57 Pentium.Labs 阅读(202) 评论(0) 编辑 收藏 举报