BZOJ 1101 Zap(莫比乌斯反演)
http://www.lydsy.com/JudgeOnline/problem.php?id=1101
给定a,b,d,求有多少gcd(x,y)==d(1<=x<=a&&1<=y<=b)
思路:
Σgcd(x,y)==d (1<=x<=a,1<=y<=b)
=
Σgcd(x,y)==1 (1<=x<=a/d,1<=y<=b/d)
令G(i)=num(i|gcd(x,y))=n/i*m/i
g(i)=num(i=gcd(x,y))
g(i)=ΣG(d)*u(d/i) (i|d)
则答案就是g(1)
g(1)=ΣG(i)*u(i) (1<=i<=min(n,m))
=Σ(n/i)*(m/i)*u(i)
因此做出u(i)的前缀和,这样就可以一起处理n/i和m/i相同的i
1 #include<algorithm> 2 #include<cstdio> 3 #include<cmath> 4 #include<cstring> 5 #include<iostream> 6 int mul[200005],mark[200005],sum[200005],p[200005]; 7 int read(){ 8 char ch=getchar();int t=0,f=1; 9 while (ch<'0'||ch>'9'){if (ch=='-') f=-1;ch=getchar();} 10 while ('0'<=ch&&ch<='9'){t=t*10+ch-'0';ch=getchar();} 11 return t*f; 12 } 13 void init(){ 14 mul[1]=1; 15 for (int i=2;i<=50000;i++){ 16 if (!mark[i]){ 17 p[++p[0]]=i; 18 mul[i]=-1; 19 } 20 for (int j=1;j<=p[0]&&i*p[j]<=50000;j++){ 21 mark[p[j]*i]=1; 22 if (i%p[j]) mul[p[j]*i]=mul[i]*(-1); 23 else { 24 mul[p[j]*i]=0; 25 break; 26 } 27 } 28 } 29 sum[0]=0; 30 for (int i=1;i<=50000;i++) 31 sum[i]=sum[i-1]+mul[i]; 32 } 33 int cal(int a,int b){ 34 if (a>b) std::swap(a,b); 35 int res=0,n=a,m=b; 36 for (int i=1,j;i<=a;i=j+1){ 37 j=std::min(n/(n/i),m/(m/i)); 38 res+=(a/i)*(b/i)*(sum[j]-sum[i-1]); 39 } 40 return res; 41 } 42 int main(){ 43 int Q=read(); 44 init(); 45 while (Q--){ 46 int a=read(),b=read(),d=read(); 47 printf("%d\n",cal(a/d,b/d)); 48 } 49 }