BZOJ2820/LG2257 YY的GCD 莫比乌斯反演

问题描述

BZOJ2820

LG2257


题解

(sumlimits_{i=1}^{n}{sumlimits_{j=1}^{m}{[gcd(i,j)==p]}}) ,其中 (p)为质数,(n<m)

考虑不要求 (gcd(i,j)) 为质数时的做法:

可以转化为

[sumlimits_{g=1}^{n}{g imes sumlimits_{i=1}^{lfloor frac{n}{g} floor}{sumlimits_{j=1}^{lfloor frac{m}{g} floor}{[gcd(i,j)==1]}}} ]

根据狄利克雷卷积和莫比乌斯函数的性质,得

[sumlimits_{g=1}^{n}{g imes sumlimits_{i=1}^{lfloor frac{n}{g} floor}{sumlimits_{j=1}^{lfloor frac{m}{g} floor}{sumlimits_{x|gcd(i,j)}{mu(x)}}}} ]

转化为枚举倍数,得

[sumlimits_{g=1}^{n}{g imes sumlimits_{x=1}^{lfloor frac{n}{g} floor}{sumlimits_{i=1}^{lfloor frac{n}{gx} floor}{sumlimits_{j=1}^{lfloor frac{m}{gx} floor}{1}}}} ]

[sumlimits_{g=1}^{n}{g imes sumlimits_{x=1}^{lfloor frac{n}{g} floor}{mu(x) lfloor frac{n}{gx} floor lfloor frac{m}{gx} floor}} ]

接下来考虑 (g) 要是质数怎么办

构造函数 (f(x)),(f(x)=egin{cases}1&x in mathbb{P}\0&x otin mathbb{P}end{cases}),其中 (mathbb{P}) 为质数集合。

(D=gx) ,则要求的就是

[sumlimits_{D=1}^{n}{lfloor frac{n}{D} floor lfloor frac{m}{D} floor imes sumlimits_{g|D}{f(g)mu(frac{D}{g})}} ]

发现右边的 (sumlimits_{g|D}{f(g)mu(frac{D}{g})}) 是狄利克雷卷积的形式。

(h(x)=f*mu) ,则可以通过枚举倍数的时间,在 (O(k ln n)) 的时间复杂度求出 (h(x)) ,其中 (k)([1,n]) 中质数个数,约等于 (frac{n}{ln n}) ,于是约等于 (O(n))

得到最终答案式子为

[sumlimits_{D=1}^{n}{h(D)lfloor frac{n}{D} floor lfloor frac{m}{D} floor} ]


(mathrm{Code})

#include<bits/stdc++.h>
using namespace std;

const int maxn=10000000;

int p[maxn+7],pr[maxn+7],miu[maxn+7];
int h[maxn+7],T,tot;

long long s[maxn+7];

void preprocess(void){
	miu[1]=1;
	for(int i=2;i<=maxn;i++){
		if(!p[i]) p[i]=i,pr[++tot]=i,miu[i]=-1;
		for(int j=1;j<=tot;j++){
			if((long long)i*pr[j]>maxn||pr[j]>p[i]) break;
			p[i*pr[j]]=pr[j];
			if(i%pr[j]) miu[i*pr[j]]=-miu[i];
			else miu[i*pr[j]]=0;
		}
	}
	for(int i=1;i<=tot;i++){
		for(int j=1;(long long)pr[i]*j<=maxn;j++){
			h[pr[i]*j]+=miu[j];
		}
	}
	for(int i=1;i<=maxn;i++) s[i]=s[i-1]+(long long)h[i];
}

void Init(void){
	scanf("%d",&T);
}

long long calc(int x,int y){
	if(x>y) swap(x,y);
	long long res(0);
	for( int l=1,r;l<=x;l=r+1){
		r=min(x/(x/l),y/(y/l));
		res+=(long long)(s[r]-s[l-1])*(long long)(x/l)*(y/l);
	}
	return res;
}

void Work(void){
	preprocess();
	while(T--){
		int x,y;
		scanf("%d%d",&x,&y);
		printf("%lld
",calc(x,y));
	}
}

signed main(){
	Init();
	Work();
}
原文地址:https://www.cnblogs.com/liubainian/p/12078085.html