洛谷 P2257 YY的GCD
洛谷 P2257 YY的GCD
$ solution: $
这道题完全跟[POI2007]ZAP-Queries (莫比乌斯反演+整除分块) 用的一个套路。
我们可以列出答案就是要我们求:
$ ans=\sum_{p\in prime}\sum_{i=1}{n}{\sum_{j=1}{[gcd(i,j)==p]}} $
我们发现后面那一部分( $ \sum_{i=1}{n}{\sum_{j=1}{[gcd(i,j)==p]}} $ )可以套路的莫比乌斯反演:
$ ans=\sum_{p\in prime}\sum_{i=1}{\frac{n}{p}}{\sum_{j=1}{p}}{[gcd(i,j)==1]}} $
$ ans=\sum_{p\in prime}\sum_{i=1}{\frac{n}{p}}{\sum_{j=1}{p}}{\sum_{k|gcd(i,j)}{\mu(k)}}} $
$ ans=\sum_{p\in prime}\sum_{k}{min(n,m)}{\mu(k)}{\sum_{i=1}{p}}{\sum_{j=1}^{\frac{m}{p}}{[k|gcd(i,j)]}}} $
$ ans=\sum_{p\in prime}\sum_{k}^{min(n,m)}{\mu(k)\times \lfloor \frac{n}{p\times k} \rfloor \times \lfloor \frac{m}{p\times k} \rfloor} $
我们发现后面那两个东西有点麻烦,我们想办法让它变成常数项:
设: $ T=p\times k $
$ ans=\sum_{t=1}^{min(n,m)} \lfloor \frac{n}{T} \rfloor \times \lfloor \frac{m}{T} \rfloor \sum_{p|T,p\in prime} \mu{\frac{T}{p}} $
然后我们发现后面那一部分( $ \sum_{p|T,p\in prime} \mu{\frac{T}{p}} $ )可以预处理,然后我们在前面用整除分块,这样就可以每次只用 $ log(min(n,m)) $ 的时间完成单个询问了!
$ code: $
#include<iostream>
#include<cstdio>
#include<iomanip>
#include<algorithm>
#include<cstring>
#include<cstdlib>
#include<ctime>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#define ll long long
#define db double
#define inf 0x7fffffff
#define rg register int
using namespace std;
const int N=1e7;
ll ans;
int n,m,t;
int pr[10000005];
int mu[10000005];
int sum[10000005];
bool zhi[10000005];
inline int qr(){
register char ch; register bool sign=0; rg res=0;
while(!isdigit(ch=getchar())) if(ch=='-')sign=1;
while(isdigit(ch)) res=res*10+(ch^48),ch=getchar();
return sign?-res:res;
}
inline void get_ready(){
mu[1]=1; int t=0;
for(rg i=2;i<=N;++i){
if(!zhi[i])pr[++t]=i,mu[i]=-1;
for(rg j=1;j<=t;++j){
if(i*pr[j]>N)break;
zhi[i*pr[j]]=1;
if(!(i%pr[j]))break;
mu[i*pr[j]]=-mu[i];
}
}
for(rg i=1;i<=t;++i)
for(rg j=1;j*pr[i]<=N;++j)
sum[j*pr[i]]+=mu[j];
for(rg i=1;i<=N;++i) sum[i]+=sum[i-1];
}
int main(){
//freopen(".in","r",stdin);
//freopen(".out","w",stdout);
t=qr(); get_ready();
while(t--){ ans=0;
n=qr(); m=qr();
rg x=min(n,m),l,r;
for(l=1;l<=x;l=r+1){
r=min(n/(n/l),m/(m/l));
ans+=(ll)(n/l)*(m/l)*(sum[r]-sum[l-1]);
}printf("%lld\n",ans);
}
return 0;
}