GCD of p2257 YY (Mobius inversion)

Topic link

The question ：

Given $N,M$, seek $1\le x\le N,1\le y\le M$ And $\gcd (x,y)$ Prime number $(x,y)$ How many pairs? .
$T$ Group example ,$T=10^4,N,M\le 10^7$

deduction ：

Write the formula according to the meaning of the question ,

${\sum }_{i=1}^N{\sum }_{j=1}^M[gcd(i,j)byqualityCount]$

Put forward $gcd$, become

${\sum }_{g=1}^{min(N,M)}{\sum }_{i=1}^N{\sum }_{j=1}^M[gcd(i,j)=g][gbyqualityCount]$

Let's make $N<M$

It becomes ${\sum }_{g=1}^N[gbyqualityCount]{\sum }_{i=1}^N{\sum }_{j=1}^M[gcd(i,j)=g]$ 了 , It's more comfortable to write .

Make $f(g)={\sum }_{i=1}^N{\sum }_{j=1}^M[gcd(i,j)=g]$

Make $F(g)={\sum }_{g|d}f(d)={\sum }_{i=1}^N{\sum }_{j=1}^M[g|gcd(i,j)]=\lfloor \frac{N}{g}\rfloor \lfloor \frac{M}{g}\rfloor$

According to Mobius inversion formula, there are $f(g)={\sum }_{g|d}\mu (\frac{d}{g})F(d)$

Replace the original formula to get ,

${\sum }_{g=1}^N[gbyqualityCount]f(g)={\sum }_{g=1}^N[gbyqualityCount]{\sum }_{g|d}\mu (\frac{d}{g})\lfloor \frac{N}{d}\rfloor \lfloor \frac{M}{d}\rfloor$

                                      $={\sum }_{d=1}^N\lfloor \frac{N}{d}\rfloor \lfloor \frac{M}{d}\rfloor {\sum }_{g|d}\mu (\frac{d}{g})[gbyqualityCount]$

                                      $={\sum }_{d=1}^N\lfloor \frac{N}{d}\rfloor \lfloor \frac{M}{d}\rfloor S(d)$

Among them, we make $S(d)={\sum }_{g|d}\mu (\frac{d}{g})[gbyqualityCount]$, We just need to preprocess out $S(d)$, Reprocess out $S(d)$ And , You can use number theory to find the answer in blocks .

So pretreatment $S(d)$, We can calculate the contribution for each number decomposition prime factor , Or for each prime number , Add the contribution to all its multiples , Either way $o(nlogn)$.

Code ：

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod = 1e9+7;
typedef pair<int,int> P;
const int MAXN=10000000;

int prime[MAXN+10],notPrime[MAXN+10],mu[MAXN+10],tot;
void initMu(int n)
{
    
    notPrime[1]=mu[1]=1;
    for(int i=2;i<=n;i++)
    {
    
        if(!notPrime[i])prime[tot++]=i,mu[i]=-1;
        for(int j=0;j<tot&&i*prime[j]<=n;j++)
        {
    
            notPrime[i*prime[j]]=1;
            if(i%prime[j])mu[i*prime[j]]=-mu[i];
            else {
    mu[i*prime[j]]=0;break;}
        }
    }
}

int f[MAXN+10];
int main()
{
    
    initMu(MAXN);
    for(int i=1;i<=MAXN;i++)if(!notPrime[i])
    {
    
        for(int j=i;j<=MAXN;j+=i)
        {
    
            f[j]+=mu[j/i];
        }
    }
    for(int i=1;i<=MAXN;i++)f[i]+=f[i-1];
    int T;
    scanf("%d",&T);
    while(T--)
    {
    
        int n,m;
        scanf("%d%d",&n,&m);
        if(n>m)swap(n,m);
        ll ans=0;
        for(int l=1,r;l<=n;l=r+1)
        {
    
            r=min(n/(n/l),m/(m/l));
            ans=ans+1ll*(n/l)*(m/l)*(f[r]-f[l-1]);
        }
        printf("%lld\n",ans);
    }
    return 0;
}