hdu1695 GCD 学习莫比乌斯反演

定理：f(x)为算术函数，F(x)为f(x)的和函数，那么有结论(这里d是n的倍数)。

还有另一种形式(这里n是d的倍数)。

其中就是莫比乌斯函数，它具有以下定义和性质：

(1)若d=1，那么μ(d)=1
(2)若d=p1p2...pk(p1...pk均为互异质数)，那么μ(d)=(−1)^k
(3)其他情况下，μ(d)=0

用处：当和函数F(x)易得时，求f(x)就可以通过莫比乌斯反演。

模板：

void Init(){//线性筛O(n) 
	mem(mu,0),mem(prime,0),mem(vis,0);
	mu[1]=1,cnt=0;
	for(int i=2;i<N;i++){
		if(!vis[i]){
			prime[cnt++]=i;
			mu[i]=-1;
		}
		for(int j=0;j<cnt&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j]) mu[i*prime[j]]=-mu[i];
			else{
				mu[i*prime[j]]=0;
				break;
			}
		}
	}
}

例题链接： http://acm.hdu.edu.cn/showproblem.php?pid=1695
题意：求满足a<=x<=b，c<=y<=d，gcd(x,y)==k的数对(x,y)的对数，(2,3)，(3,2)算一对。

题解：求gcd(x,y)==k，可以转化为gcd(x/k,y/k)==1。

设f(d)表示gcd(x,y)=d的对数，那么我们要求的就是f(1)。

设F(n)表示gcd(x,y)=(n的倍数)的对数。那么其实就满足了定理中的第二种形式，这样构造出来之后，有没有发现F(n)很好求？F(n)=(b/n)*(d/n)，(b=b/k，d=d/k)。

那么剩下的就可以根据结论：(其中d最大为min(b,d))。f(1)=mu(1)*F(1)+mu(2)*F(2)+...+mu(min(b,d))*F(min(b,d))。至于去重，因为重复的部分是在1<=x,y<=min(b,d)，减去这一部分即可，具体代码理解。

代码：

#include<set>
#include<map>
#include<stack>
#include<queue>
#include<vector>
#include<string>
#include<bitset>

#include<algorithm>
#include<cstring>
#include<cstdio>
#include<cmath>

#include<iomanip>
#include<iostream>

#define debug cout<<"aaa"<<endl
#define mem(a,b) memset(a,b,sizeof(a))
#define LL long long
#define lson l,mid,root<<1
#define rson mid+1,r,root<<1|1
#define MIN_INT (-2147483647-1)
#define MAX_INT 2147483647
#define MAX_LL 9223372036854775807i64
#define MIN_LL (-9223372036854775807i64-1)
using namespace std;

const int N = 100000 + 5;
const int mod = 1000000000 + 7;
int mu[N],prime[N],cnt,t,cas=0,a,b,c,d,k;
LL ans1,ans2;
bool vis[N];

void Init(){//线性筛O(n) 
	mem(mu,0),mem(prime,0),mem(vis,0);
	mu[1]=1,cnt=0;
	for(int i=2;i<N;i++){
		if(!vis[i]){
			prime[cnt++]=i;
			mu[i]=-1;
		}
		for(int j=0;j<cnt&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j]) mu[i*prime[j]]=-mu[i];
			else{
				mu[i*prime[j]]=0;
				break;
			}
		}
	}
}

int main(){
	Init();
	scanf("%d",&t);
	while(t--){
		scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
		if(k==0){
			printf("Case %d: 0\n",++cas);
			continue;
		}
		b/=k,d/=k;
		ans1=0,ans2=0;
		for(int i=1;i<=min(b,d);i++){
			ans1+=(LL)mu[i]*(b/i)*(d/i);
		}
		for(int i=1;i<=min(b,d);i++){
			ans2+=(LL)mu[i]*(min(b,d)/i)*(min(b,d)/i);
		}
		printf("Case %d: %lld\n",++cas,ans1-ans2/2);
	}
	return 0;
}