题解Luogu P2257 YY的GCD

Posted 2021-01-28 yzhang-rp-inf

原题传送门

这题需要运用莫比乌斯反演(懵逼钨丝繁衍)

显然题目的答案就是[ Ans=sum_{i=1}^Nsum_{j=1}^M[gcd(i,j)=prime]]

我们先设设F(n)表示满足(gcd(i,j)\\%t=0)的数对个数，f(t)表示满足(gcd(i,j)=t)的数对个数

[f(t)=sum_{i=1}^Nsum_{j=1}^M[gcd(i,j)=t]]

[F(n)=sum_{n|t}lfloor frac{N}{n} floor lfloor frac{M}{n} floor]

那么根据莫比乌斯反演的第二种形式珂以得到

[f(n)=sum_{n|t}mu(lfloor frac{t}{n} floor)F(t)]

所以答案珂以变形为：

[Ans=sum_{p in prime}sum_{i=1}^Nsum_{j-1}^M[gcd(i,j)=p)]

[=sum_{p in prime}f(p) qquad qquad quad]

[=sum_{p in prime}sum_{p|t}mu(lfloor frac{t}{p} floor)F(t)]

我们不枚举p，我们枚举(lfloor frac{t}{p} floor)

[Ans=sum_{p in prime}sum_{d=1}^{Min(frac{N}{p},frac{M}{p})}mu(t)F(tp)]

[qquad qquad qquad=sum_{p in prime}sum_{d=1}^{Min(frac{N}{p},frac{M}{p})}mu(t)sum_{n|t}lfloor frac{N}{tp} floor lfloor frac{M}{tp} floor]

我们把tp换成T继续变形

[Ans=sum_{T=1}^{Min(N,M)}lfloor frac{N}{T} floor lfloor frac{M}{T} floor(sum_{p|T,p in prime}mu(frac{T}{p}))]

这样就珂以用整除分块求了qaq

#include <bits/stdc++.h>
#define N 10000005
#define ll long long 
#define getchar nc
using namespace std;
inline char nc(){
    static char buf[100000],*p1=buf,*p2=buf; 
    return p1==p2&&(p2=(p1=buf)+fread(buf,1,100000,stdin),p1==p2)?EOF:*p1++; 
}
inline int read()
{
    register int x=0,f=1;register char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9')x=(x<<3)+(x<<1)+ch-'0',ch=getchar();
    return x*f;
}
inline void write(register ll x)
{
    if(!x)putchar('0');if(x<0)x=-x,putchar('-');
    static int sta[30];register int tot=0;
    while(x)sta[tot++]=x%10,x/=10;
    while(tot)putchar(sta[--tot]+48);
}
inline int Min(register int x,register int y)
{
    return x<y?x:y;
}
int v[N],miu[N],prim[N],cnt=0,g[N];
ll sum[N];
int main()
{
    miu[1]=1;
    for(register int i=2;i<=N;++i)
    {
        if(!v[i])
        {
            miu[i]=-1;
            prim[++cnt]=i;
        }
        for(register int j=1;j<=cnt&&prim[j]*i<=N;++j)
        {
            v[i*prim[j]]=1;
            if(i%prim[j]==0)
                break;
            else
                miu[prim[j]*i]=-miu[i];
        }
    }
    for(register int i=1;i<=cnt;++i)
        for(register int j=1;j*prim[i]<=N;++j)
            g[j*prim[i]]+=miu[j];
    for(register int i=1;i<=N;++i)
        sum[i]=sum[i-1]+(ll)g[i];
    int t=read();
    while(t--)
    {
        int n=read(),m=read();
        if(n>m)
            n^=m^=n^=m;
        ll ans=0;
        for(register int l=1,r;l<=n;l=r+1)
        {
            r=Min(n/(n/l),m/(m/l));
            ans+=(ll)(n/l)*(m/l)*(sum[r]-sum[l-1]);
        }
        write(ans),puts("");
    }
    return 0;
}