[BZOJ2820]YY的GCD

Description:

求\(\sum_{i=1}^{n} \sum_{j=1}^{m} gcd(i,j)==p\ _{p \in Prime}\)

Hint:

\(数据组数<=1e4,n,m<=1e7\)

Solution:

$ \sum_{p \in Prime} \sum_{i=1}^{n} \sum_{j=1}^{m} gcd(i,j)==p $

$ \sum_{p \in Prime} \sum_{i=1}^{ \lfloor \frac {n} {p} \rfloor } \sum_{j=1}^{ \lfloor \frac {m} {p} \rfloor } gcd(i,j)==1 $

$ \sum_{p \in Prime} \sum_{i=1}^{ \lfloor \frac {n} {p} \rfloor } \sum_{j=1}^{ \lfloor \frac {m} {p} \rfloor } \sum_{d|gcd(i,j)} \mu(d) $

$ \sum_{p \in Prime} \sum_{i=1}^{ \lfloor \frac {n} {p} \rfloor } \sum_{j=1}^{ \lfloor \frac {m} {p} \rfloor } \sum_{d|gcd(i,j)} \mu(d) $

$ \sum_{p \in Prime} \sum_{i=1}^{ \lfloor \frac {n} {p} \rfloor } \sum_{j=1}^{ \lfloor \frac {m} {p} \rfloor } [\ d\ |\ i\ ][\ d\ |\ \ j] \mu(d) $

$ \sum_{p \in Prime} \sum_{d}^{\lfloor \frac {min(n,m)} {p} \rfloor} \mu (d) \sum_{i=1}^{ \lfloor \frac {n} {p} \rfloor } \sum_{j=1}^{ \lfloor \frac {m} {p} \rfloor } $

$ \sum_{p \in Prime} \sum_{d}^{\lfloor \frac {min(n,m)} {p} \rfloor} \mu (d) \lfloor \frac {n} {dp} \rfloor \lfloor \frac {m} {dp} \rfloor $

换元得

$ \sum_{p \in Prime} \sum_{T=1}^{min(n,m)} \mu (\frac{T}{p}) \lfloor \frac {n} {T} \rfloor \lfloor \frac {m} {T} \rfloor $

$ \sum_{T=1}^{min(n,m)} \lfloor \frac {n} {T} \rfloor \lfloor \frac {m} {T} \rfloor \sum_{p \in Prime} \mu(\frac{T}{p}) $

将前缀和特殊处理一下，直接枚举T即可

#include<bits/stdc++.h>
using namespace std;
const int mxn=1e7+5;
int n,m,T,tot,mu[mxn],vis[mxn],p[mxn],f[mxn];
long long sum[mxn];

void sieve(int lim)
{
	mu[1]=1;
	for(int i=2;i<=lim;++i) {
		if(!vis[i]) mu[i]=-1,p[++tot]=i;
		for(int j=1;j<=tot;++j) {
			if(p[j]*i>lim) break;
			vis[p[j]*i]=1;
			if(i%p[j]==0) {
				mu[i*p[j]]=0;
				break;
			}
			mu[p[j]*i]=-mu[i];
		}
	}
	for(int i=1;i<=lim;++i) {
		for(int j=1;j<=tot;++j) {
			if(p[j]*i>lim) break;
			else f[p[j]*i]+=mu[i]; //特殊处理
		}
	}
	for(int i=1;i<=lim;++i) sum[i]=sum[i-1]+f[i];
}

int main()
{
	scanf("%d",&T);
	sieve(10000000);
	while(T--) {
		long long ans=0;
		scanf("%d%d",&n,&m); if(n>m) swap(n,m);
		for(int l=1,r;l<=n;l=r+1) {
			r=min(n/(n/l),m/(m/l));
			ans+=1ll*(n/l)*(m/l)*(sum[r]-sum[l-1]);
		}
		printf("%lld\n",ans);
	}
	return 0;
}