bzoj2301:[HAOI2011]Problem b(容斥+莫比乌斯反演+分块)

www.cnblogs.com/shaokele/

---

2301: [HAOI2011]Problem b##

　　Time Limit: 50 Sec Memory Limit: 256 MB

Description###

　　对于给出的n个询问，每次求有多少个数对(x,y)，满足a≤x≤b，c≤y≤d，且gcd(x,y) = k，gcd(x,y)函数为x和y的最大公约数。

Input###

　　第一行一个整数n，接下来n行每行五个整数，分别表示a、b、c、d、k

Output###

　　共n行，每行一个整数表示满足要求的数对(x,y)的个数

Sample Input###

　　2
　　2 5 1 5 1
　　1 5 1 5 2

Sample Output###

　　14
　　3

HINT###

　　100%的数据满足：1≤n≤50000，1≤a≤b≤50000，1≤c≤d≤50000，1≤k≤50000

题目地址 bzoj 2301: [HAOI2011]Problem b
题目已经很简洁明了了，不再复述：）懒一下www.cnblogs.com/shaokele/

---

题解：####

　　莫比乌斯反演最常用的几个等式：
　　　　 　　(gcd(a,b)=sum_{d|gcd(a,b)}^{}μ(d))
　　
　　　　 　　([gcd(a,b))==(1]=sum_{d|a,d|b}^{}μ(d))

　　首先容斥
　　(f(n,m)=sum_{i=1}^{n}sum_{j=1}^m[gcd(i,j)=k])
　　(∴ans=f(b,d)-f(a-1,d)-f(b,c-1)+f(a-1,c-1))
　　
　　化简:
　　(f(n,m)⇒sum_{i=1}^{⌊frac{n}{k}⌋}sum_{j=1}^{⌊frac{m}{k}⌋}[gcd(i,j)=1])
　　　　　　(⇒sum_{i=1}^{⌊frac{n}{k}⌋}sum_{j=1}^{⌊frac{m}{k}⌋}sum_{d|(i,j)}μ(d))
　　　　　　(⇒sum_dμ(d)sum_{d|i}sum_{d|j}1)
　　　　　　(⇒sum_dμ(d)⌊frac{n}{kd}⌋⌊frac{m}{kd}⌋)
　　由于(⌊frac{n}{kd}⌋⌊frac{m}{kd}⌋)是可以分块处理的
　　所以每个询问都只要(O(sqrt{n}))就能得到答案

---

AC代码

#include <cstdio>
#include <algorithm>
using namespace std;
const int N=50005;
int Q,a,b,c,d,k,tot,ans;
int mu[N],p[N],sum[N];
bool vis[N];
inline int read(){
	int x=0,f=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();}
	return x*f;
}
void Get_mu(){
	mu[1]=1;
	for(int i=2;i<=50000;i++){
		if(!vis[i]){
			p[++tot]=i;
			mu[i]=-1;
		}
		for(int j=1;j<=tot;j++)
			if(i*p[j]<=50000){
				vis[i*p[j]]=1;
				if(i%p[j]==0){
					mu[i*p[j]]=0;
					break;
				}
				mu[i*p[j]]=-mu[i];
			}else break;
	}
}
void init(){
	Get_mu();
	for(int i=1;i<=50000;i++)
		sum[i]=sum[i-1]+mu[i];
}
int solve(int n,int m){
	n/=k;m/=k;
	if(n>m)swap(n,m);
	int res=0,last;
	for(int i=1;i<=n;i=last+1){
		last=min(n/(n/i),m/(m/i));
		res+=(n/i)*(m/i)*(sum[last]-sum[i-1]);
	}
	return res;
}
int main(){
	Q=read();
	init();
	while(Q--){
		a=read();b=read();c=read();d=read();k=read();
		ans=solve(b,d)-solve(a-1,d)-solve(b,c-1)+solve(a-1,c-1);
		printf("%d
",ans);
	}
	return 0;
}