BZOJ2301 容斥原理,莫比乌斯反演

来源：互联网发布：ubuntu自带中文输入法编辑：程序博客网时间：2024/05/18 11:10

Description

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

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

做法：

第一次尝试写莫比乌斯反演的模板。由于对该算法不是很熟悉，就不讲莫比乌斯反演的部分了。

用solve(s,t)表示1<=x<=s,1<=y<=t时的答案。那么所求的答案为

solve(b,d)-solve(a-1,d)-solve(b,c-1)+solve(a-1,c-1)

然后用莫比乌斯反演求解solve就好了。

模板借鉴了各位学长和大神的。

注:BZOJ用cin会RE，这个代码在洛谷上交就没问题，在BZOJ上交要改scanf

#include<iostream>#include<cstdio>#include<cstring>#include<algorithm>#include<vector>using namespace std;int miu[50010];int sum[50010];int prime[50010];bool notprm[50010];const int MAXN=50001;int cnt;long long ans;int a,b,c,d,k;void init(){miu[1]=1;for(int i=2;i<=MAXN;++i){if(!notprm[i]){prime[cnt++]=i;miu[i]=-1;}for(int j=0;j<cnt&&i*prime[j]<=MAXN;++j){notprm[i*prime[j]] = 1;if(i%prime[j])miu[i*prime[j]]=-miu[i];else{miu[i*prime[j]]=0;                  break;}}}for(int i=1;i<=MAXN;++i)sum[i]=sum[i-1]+miu[i];}long long solve(int x,int y){long long res=0;if(x>y)swap(x,y);for(long long i=1,l=0;i<=x;i=l+1){l=min(x/(x/i),y/(y/i));res+=(sum[l]-sum[i-1])*(x/i)*(y/i);}return res;}int main(){init();int T;cin>>T;while(T--){cin>>a>>b>>c>>d>>k;if(a>b||c>d){cout<<0<<endl;continue;}ans=solve(b/k,d/k);        ans-=solve((a-1)/k,d/k);        ans-=solve(b/k,(c-1)/k);        ans+=solve((a-1)/k,(c-1)/k);        cout<<ans<<endl;}}

0 0