[题解]luogu P2257 YY的GCD

Posted 2022-10-23 end-of-mind

袁野的gcd

首先可以肯定的是这是一道数论题

所以题目就是: \(\sum_i=1^N \sum_j=1^M [gcd(i,j)\in prime]\)

接下来就可以愉快的推式子了~

首先可以按套路枚举 prime 和 gcd

\(\sum_p \in prime\sum_i=1^N\sum_j=1^M[gcd(i,j)=p]\)

所以我们显然可以将\(i\)和\(j\) 同时除以 \(p\),得到
\(\sum_p\in prime\sum_i=1^N\sum_j=1^M[gcd(\frac ip,\frac jp)=1]\)

即等价于
\(\sum_p\in prime\sum_i=1^\frac Np\sum_j=1^\frac Mp [gcd(i,j)=1]\)

到了这里,我们发现对于\([gcd(i,j)=1]\) 我们可以根据莫比乌斯函数的定义 \(\sum_d\mid n\mu(d)\ =[n=1]\) 来进行替换，即

\(\sum_p \in prime\sum_i=1^\frac Np\sum_j=1^\frac Mp\sum_d\mid i,d\mid j \mu(d)\)

我们可以按照套路将 \(d\) 提前,(这里我们默认\(N\le M\))

\(\sum_p \in prime\sum_d=1^\frac Np\mu(d)\sum_i=1^\frac Np[d\mid i]\sum_j=1^\frac Mp[d\mid j]\)

然后我们可以发现后面那两个 \(\sum\) 可以用整除分块搞,即

\(\sum_p \in prime\sum_d=1^\frac Np\mu(d)\lfloor\frac Ndp\rfloor\lfloor\frac Mdp\rfloor\)

此时,我们应用一个小技巧,将上式中的\(dp\)换成\(Q\)(或称刘琛薛坷佳崔家贺):

\(\sum_Q=1^N \sum_p \in prime\mu(\frac Qp)\lfloor\frac NQ\rfloor\lfloor\frac MQ\rfloor\)

对于 \(\sum_p \in prime\mu(\frac Qp)\lfloor\frac NQ\rfloor\lfloor\frac MQ\rfloor\) ,我们可以用线性筛进行预处理 ,即在筛出 \(\varphi\) 后枚举。所以在维护前缀和后就可以对于每一个询问 \(O(1)\) 回答啦~

\(\mathcalCode:\)

#include<cstdio>
#include<cstdlib>
#include<iostream>
#include<algorithm>
using namespace std;
#define N 10000010
#define debug cout<<__LINE__<<" "<<__FUNCTION__<<"\n"
inline int read()
    int x=0;char ch=getchar();
    while(ch<'0'||ch>'9')ch=getchar();
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+(ch^48);ch=getchar();
    return x;

void put(long long x)
    if(x<0) putchar('-'),x=-x;
    if(x>=10) put(x/10);
    putchar((x%10)+48);

int prime[N>>2],n,m,cnt;
bool ispri[N];
int mu[N],num[N];
inline void pri()
    register int i,j;
    for(i=2;i<=N-10;i++)
        if(!ispri[i])
            prime[++cnt]=i;
            mu[i]=-1;
        
        for(j=1;j<=cnt&&(i*prime[j]<=N-10);j++)
            if(i*prime[j]<=N-10) ispri[i*prime[j]]=1;
            if(i%prime[j]==0)
                break;
            else
                mu[i*prime[j]]=-mu[i];
            
        
    
    int res;
    for(i=1;i<=cnt;i++)
        res=1;
        for(j=prime[i];j<=N-10;j+=prime[i],res++)
            num[j]+=mu[res];
        
    
    for(i=2;i<=N-10;i++) num[i]+=num[i-1];

signed main()
//  freopen(".in","r",stdin);
//  freopen(".out","w",stdout);
    mu[1]=1;
    pri();
    int T=read(),l,r;
    long long ans = 0;
    while(T--)
        n=read();m=read();ans=0;
        if(n>m) swap(n,m);
        for(l=1,r;l<=n;l=r+1)
            r=min(n/(n/l),m/(m/l));
            ans+=(long long)(num[r]-num[l-1])*(n/l)*(m/l);
        
        put(ans);putchar('\n');
    
//  fclose(stdin);
//  fclose(stdout);
    return 0;