前置
若函数(F(n))和(f(d))存在以下关系
[ F(n)=sum_{n|d}f(d)
]
则可以推出
[ f(n)=sum_{n|d}mu(frac{d}{n})F(d)
]
这就是莫比乌斯反演
题目要求
求(gcd(a,b)={prime},ainleft[1,n ight],binleft[1,m ight])
思路
根据题意所以设出(f(n))表示(gcd(a,b)=n)的(a,b)对数
根据莫比乌斯反演的形式
[ F(n)=sum_{n|d}f(d)
]
可以设出一个函数(F(n)),表示(n|gcd(a,b))的((a,b))对数
因为(n|gcd(a,b)),所以(a=k_1 imes n,b=k_2 imes n)
所以显然有
[ F(x)=lfloorfrac{n}{x}
floor imeslfloorfrac{m}{x}
floor
]
因为
[ f(x)=sum_{x|d}mu(frac{d}{x})F(d)
]
所以
[ f(x)=sum_{x|d}mu(frac{d}{x}) imeslfloorfrac{n}{d}
floor imeslfloorfrac{m}{d}
floor
]
考虑到(frac{d}{n})的形式并不优美，我们换一种东西枚举
[ f(x)=sum_{t=1}^{lfloorfrac{n}{x}
floor}mu(t) imeslfloorfrac{n}{t imes x}
floor imeslfloorfrac{m}{t imes x}
floor
]
所以
[ ans=sum_{xin{prime}}^nf(x)=sum_{xin{prime}}^nsum_{t=1}^{lfloorfrac{n}{x}
floor}mu(t) imeslfloorfrac{n}{t imes x}
floor imeslfloorfrac{m}{t imes x}
floor
]
这样能拿到50PTS
然后设(T=t*x),这样形式就变得更优美了一些
原式变形为
[ ans=sum_{xin{prime}}^nf(x)=sum_{xin{prime}}^nsum_{t=1}^{lfloorfrac{n}{x}
floor}mu(frac{T}{x}) imeslfloorfrac{n}{T}
floor imeslfloorfrac{m}{T}
floor
]
[ ans=sum_{xin{prime}}^nf(x)=sum_{T=1}^{min(n,m)}sum_{din{prime},d|T}mu(frac{T}{d}) imeslfloorfrac{n}{T}
floor imeslfloorfrac{m}{T}
floor
]
[ ans=sum_{xin{prime}}^nf(x)=sum_{T=1}^{min(n,m)}lfloorfrac{n}{T}
floor imeslfloorfrac{m}{T}
floor imessum_{din{prime},d|T}mu(frac{T}{d})
]
后面(mu)的部分可以前缀和一下
前面的可以整除分块
加上线性筛
然后没了
代码
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std;
short mu[10001000];
int n,m,T,iprime[10001000],cnt,isprime[10001000],summu[10001000];
void prime(int n){
mu[1]=1;
isprime[1]=true;
for(int i=2;i<=n;i++){
if(!isprime[i])
iprime[++cnt]=i,mu[i]=-1;
for(int j=1;j<=cnt&&iprime[j]*i<=n;j++){
isprime[iprime[j]*i]=true;
if(i%iprime[j]==0){
mu[iprime[j]*i]=0;
break;
}
mu[iprime[j]*i]=-mu[i];
}
}
for(int i=1;i<=cnt;i++)
for(int j=1;j*iprime[i]<=n;j++)
summu[iprime[i]*j]+=mu[j];
for(int i=1;i<=n;i++)
summu[i]+=summu[i-1];
}
long long f(int n,int m){
long long ans=0;
for(int l=1,r;l<=min(n,m);l=r+1){
r=min(n/(n/l),m/(m/l));
ans+=1LL*(summu[r]-summu[l-1])*(n/l)*(m/l);
}
return ans;
}
int main(){
prime(10000100);
scanf("%d",&T);
while(T--){
long long ans=0;
scanf("%d %d",&n,&m);
if(n<m)
swap(n,m);
ans+=f(n,m);
printf("%lld
",ans);
}
return 0;
}