[BZOI 3994] [SDOI2015]约数个数和(莫比乌斯反演+数论分块)

[BZOI 3994] [SDOI2015]约数个数和

题面

设d(x)为x的约数个数,给定N、M,求(sum _{i=1}^n sum_{i=1}^m d(i imes j))

T组询问,(N,M,T leq 50000)

分析

首先有一个结论

[d(nm)= sum _{i |n} sum _{j|m} [gcd(i,j)=1] ]

这是因为nm的约数都可以表示为(i imes frac{m}{j})的形式,并且为了不重复算,要保证(gcd(i,j)=1)

因此,我们可以开始推式子

[ans= sum_{p=1}^n sum_{q=1}^m sum_{i|p} sum _{j|q} [gcd(i,j)=1] ]

注意到每对((i,j))会对p,q中他们的倍数产生(lfloor frac{n}{i} floor imes lfloor frac{m}{j} floor) 的贡献

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

[= sum_{i=1}^n sum_{j=1} ^m varepsilon (gcd(i,j)) lfloor frac{n}{i} floor lfloor frac{m}{j} floor ]

根据(varepsilon (n) = sum_{d|n} mu(d))

[= sum_{i=1}^n sum_{j=1} ^m lfloor frac{n}{i} floor lfloor frac{m}{j} floor sum_{d|gcd(i,j)} mu(d) ]

改变求和顺序,先枚举d,显然若(d|gcd(i,j)),则(d|i,d|j),

直接把i替换为d的倍数du,j替换为d的倍数dv((u,v in N^+,duleq n,dv leq m))

[= sum_{d=1}^{min(n,m)} mu(d) sum_{u=1}^{lfloor n/d floor} sum_{v=1}^{lfloor m/d floor} lfloor frac{n}{du} floor lfloor frac{m}{dv} floor ]

[= sum_{d=1}^{min(n,m)} mu(d) sum_{u=1}^{lfloor n/d floor} sum_{v=1}^{lfloor m/d floor} lfloor frac{ lfloor n/d floor}{u} floor lfloor frac{lfloor m/d floor}{v} floor ]

[= sum_{d=1}^{min(n,m)} mu(d) sum_{u=1}^{lfloor n/d floor} lfloor frac{ lfloor n/d floor}{u} floor sum_{v=1}^{lfloor m/d floor} lfloor frac{lfloor m/d floor}{v} floor ]

(g(n) = sum _{d=1}^n lfloor frac{n}{d} floor)

[=sum_{d=1}^{min(n,m)} mu(d) g(lfloor frac{n}{d} floor) g(lfloor frac{m}{d} floor) ]

考虑如何快速求值。单个(g(n))可以运用数论分块在(O(sqrt n))的时间内求出,总时间复杂度(O(n sqrt n)). 然后线性筛出(mu),以及(mu,g)的前缀和

每次询问用数论分块的方法枚举d即可,总时间复杂度(O(n sqrt n +T sqrt n))

代码

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define maxn 50000
using namespace std;
typedef long long ll;
int t;
int n,m;
int cnt;
bool vis[maxn+5];
int prime[maxn+5];
int mu[maxn+5];
ll sum_mu[maxn+5];
int g[maxn+5];
void sieve(int n){
	mu[1]=1;
	for(int i=2;i<=n;i++){
		if(!vis[i]){
			prime[++cnt]=i;
			mu[i]=-1;
		}
		for(int j=1;j<=cnt&&i*prime[j]<=n;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j]==0){
				mu[i*prime[j]]=0;
				break;
			}else{
				mu[i*prime[j]]=-mu[i];
			}
		}
	}
	for(int i=1;i<=n;i++) sum_mu[i]=sum_mu[i-1]+mu[i];
	for(int i=1;i<=n;i++){
		int l,r;
		for(l=1;l<=i;l=r+1){
			r=i/(i/l);
			g[i]+=(ll)(r-l+1)*(i/l);
		}
	}
}

ll calc(int n,int m){
	int l,r;
	if(n<m) swap(n,m);
	ll ans=0;
	for(l=1;l<=m;l=r+1){
		r=min(n/(n/l),m/(m/l));
		ans+=(sum_mu[r]-sum_mu[l-1])*g[n/l]*g[m/l];
	}
	return ans;
}

int main(){
	sieve(maxn);
	scanf("%d",&t);
	while(t--){
		scanf("%d %d",&n,&m);
		printf("%lld
",calc(n,m));
	}	
}
原文地址:https://www.cnblogs.com/birchtree/p/11360921.html