2019.01.19【HAOI2011】【BZOJ2301】【洛谷P2522】Problem b(莫比乌斯反演)

BZOJ传送门

洛谷传送门


解析:

基础容斥一下可以将问题转化为:

i=1a1j=1c1[gcd(i,j)=k]+i=1bj=1d[gcd(i,j)=k]i=1a1j=1d[gcd(i,j)=k]i=1bj=1c1[gcd(i,j)=k] \begin{aligned} &\sum_{i=1}^{a-1}\sum_{j=1}^{c-1}[gcd(i,j)=k]+\sum_{i=1}^b\sum_{j=1}^d[gcd(i,j)=k]\\ -&\sum_{i=1}^{a-1}\sum_{j=1}^d[gcd(i,j)=k]-\sum_{i=1}^b\sum_{j=1}^{c-1}[gcd(i,j)=k] \end{aligned}

发现这是四个一模一样的子问题,那么现在考虑处理:i=1nj=1m[gcd(i,j)=k]\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=k]

其实就是i=1nkj=1nk[gcd(i,j)=1]\sum_{i=1}^{\lfloor\frac{n}k\rfloor}\sum_{j=1}^{\lfloor\frac{n}k\rfloor}[gcd(i,j)=1]

好的我们又把问题转化了一下,求i=1nj=1m[gcd(i,j)=1]\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1]

莫比乌斯反演基础题了啊。。。e=μ1e=\mu*1

所以我们要求的就是T=1min(n,m)μ(T)nTmT\sum_{T=1}^{\min(n,m)}\mu(T)\lfloor\frac{n}T\rfloor\lfloor\frac{m}T\rfloor

整除分块就行了。


代码:

#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define re register
#define gc get_char
#define pc putchar
#define cs const

namespace IO{
	namespace IOONLY{
		cs int Rlen=1<<18|1;
		char buf[Rlen],*p1,*p2;
	}
	inline char get_char(){
		using namespace IOONLY;
		return (p1==p2)&&(p2=(p1=buf)+fread(buf,1,Rlen,stdin),p1==p2)?EOF:*p1++;
	}
	
	inline int getint(){
		re int num;
		re char c;
		while(!isdigit(c=gc()));num=c^48;
		while(isdigit(c=gc()))num=(num+(num<<2)<<1)+(c^48);
		return num;
	}
	inline void outint(ll a){
		static char ch[23];
		if(a==0)pc('0');
		while(a)ch[++ch[0]]=a-a/10*10,a/=10;
		while(ch[0])pc(ch[ch[0]--]^48);
	}
}
using namespace IO;

cs int P=50004;

int prime[P],pcnt,mu[P];
bool mark[P];

inline void linear_sieves(int len=P-4){
	mu[1]=1;
	for(int re i=2;i<=len;++i){
		if(!mark[i])prime[++pcnt]=i,mu[i]=-1;
		for(int re j=1;i*prime[j]<=len;++j){
			mark[i*prime[j]]=true;
			if(i%prime[j]==0)break;
			mu[i*prime[j]]=-mu[i];
		}
		mu[i]+=mu[i-1];
	}
}

inline ll calc(int n,int m){
	ll ans=0;
	if(n>m)swap(n,m);
	for(int re i=1,j;i<=n;i=j+1){
		j=min(n/(n/i),m/(m/i));
		ans+=(ll)(mu[j]-mu[i-1])*(n/i)*(m/i);
	}
	return ans;
}

signed main(){
	linear_sieves();
	for(int re T=getint();T--;){
		int a=getint(),b=getint(),c=getint(),d=getint(),k=getint();
		a=(a-1)/k;b/=k;c=(c-1)/k;d/=k;
		outint(calc(b,d)+calc(a,c)-calc(a,d)-calc(b,c));pc('\n');
	}
	return 0;
}
发布了968 篇原创文章 · 获赞 370 · 访问量 7万+
展开阅读全文

没有更多推荐了,返回首页

©️2019 CSDN 皮肤主题: 技术黑板 设计师: CSDN官方博客

分享到微信朋友圈

×

扫一扫,手机浏览