I.YY的GCD
这就是莫比乌斯反演?咋长得不像呢?
我们看一下式子:
\(ans=\sum\limits_{i=1}^n\sum\limits_{j=1}^m[\gcd(i,j)\ is\ prime]\)。其中方括号相当于强制把方括号内的东西转成\(bool\)形。
完蛋了,这个里面看不到任何函数,只有一堆又臭又长的\(\Sigma\),怎么看也不会像是莫比乌斯反演呀!
这时候,我们就可以自设函数了。一般套路是设\(f(x)\)为\(\gcd(i,j)=x\)的个数,\(g(x)\)为\(\gcd(i,j)\)是\(x\)的倍数的个数。
显然,\(g(x)\)极易求出,即为\(\left\lfloor\dfrac{n}{x}\right\rfloor\times\left\lfloor\dfrac{m}{x}\right\rfloor\) (小学数论,当且仅当\(i\)和\(j\)都是\(x\)的倍数时,\(\gcd(i,j)\)是\(x\)的倍数。而是\(x\)的倍数并且\(\leq n\),这样的\(i\)共有\(\left\lfloor\dfrac{n}{x}\right\rfloor\)个。
而\(g(x)=\sum\limits_{n|d}f(d)\)。
套用反演定理II,得\(f(n)=\sum\limits_{n|d}g(d)\mu(\dfrac{d}{n})\)。
我们就可以开始计算了。
\(ans=\sum\limits_{p\ is\ prime}f(p)=\sum\limits_{p\ is\ prime}\sum\limits_{p|d}g(d)\mu(\dfrac{d}{p})\)
换成枚举\(\dfrac{d}{p}\),得
\(ans=\sum\limits_{p\ is\ prime}\sum\limits_{d=1}^{\min(\left\lfloor\frac{n}{p}\right\rfloor,\left\lfloor\frac{m}{p}\right\rfloor)}g(dp)\mu(d)\)
代入\(g(x)=\left\lfloor\dfrac{n}{x}\right\rfloor\times\left\lfloor\dfrac{m}{x}\right\rfloor\),得\(ans=\sum\limits_{p\ is\ prime}\sum\limits_{d=1}^{\min(\left\lfloor\frac{n}{p}\right\rfloor,\left\lfloor\frac{m}{p}\right\rfloor)}\left\lfloor\dfrac{n}{dp}\right\rfloor\times\left\lfloor\dfrac{m}{dp}\right\rfloor\mu(d)\)
令\(T=dp\),则
\(ans=\sum\limits_{T=1}^{\min(n,m)}\sum\limits_{t\ is\ a\ prime\ factor\ of\ T}\left\lfloor\dfrac{n}{T}\right\rfloor\left\lfloor\dfrac{m}{T}\right\rfloor\mu(\dfrac{T}{t})\)
则
\(ans=\sum\limits_{T=1}^{\min(n,m)}\left\lfloor\dfrac{n}{T}\right\rfloor\left\lfloor\dfrac{m}{T}\right\rfloor(\sum\limits_{t\ is\ a\ prime\ factor\ of\ T}\mu(\dfrac{T}{t}))\)
然后就OK了。只需要将\(\mu(x)\)做一个前缀和,套上一个整除分块即可。
代码:
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll sum[10001001],ans;
int _1,n,m,pri[10001001],cnt,mu[10001001];
void prime(int ip){
mu[1]=1,cnt=0;
for(int i=2;i<=ip;i++){
if(!pri[i])pri[++cnt]=i,mu[i]=-1;
for(int j=1;j<=cnt&&i*pri[j]<=ip;j++){
pri[i*pri[j]]=true;
if(!(i%pri[j]))break;
mu[i*pri[j]]=-mu[i];
}
}
for(int i=1;i<=cnt;i++)for(int j=1;j*pri[i]<=ip;j++)sum[j*pri[i]]+=mu[j];
for(int i=1;i<=ip;i++)sum[i]+=sum[i-1];
}
signed main(){
scanf("%d",&_1),prime(10000000);
while(_1--){
scanf("%d%d",&n,&m),ans=0;
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;
}