YY的GCD

YY的GCD

给出T个询问，询问\(\sum_{i=1}^N\sum_{j=1}^M(gcd(i,j)\in prime)\)，T = 10000，N, M <= 10000000。

解

显然质数是需要枚举的，设N<M,于是

\[ans=\sum_{p\in prime}\sum_{i=1}^N\sum_{j=1}^M(gcd(i,j)==p)\]

于是设

\[f(p)=\sum_{i=1}^N\sum_{j=1}^M(gcd(i,j)==p)\]

\[F(p)=\sum_{i=1}^N\sum_{j=1}^M(p|gcd(i,j))=[N/p][M/p]\]

由Mobius反演定理，我们有

\[f(p)=\sum_{p|d}F(d)\mu(d/p)\]

于是

\[ans=\sum_{p\in prime}\sum_{p|d}F(d)\mu(d/p)=\sum_{d=1}^NF(d)\sum_{p|d,p\in prime}\mu(d/p)\]

显然后式是可以\(O(nlon(n))\)维护的，对F（d）进行整除分块即可。

参考代码

#include <iostream>
#include <cstdio>
#define il inline
#define ri register
#define ll long long
#define limit 10000000
#define swap(x,y) x^=y^=x^=y
using namespace std;
bool check[limit+1];
int mb[limit+1],prime[700000],pt,
    opt[limit+1];ll ans;
il void read(int&);
il int min(int,int);
void pen(ll),prepare();
int main(){
    int lsy,a,b,i,j;read(lsy);
    prepare();while(lsy--){
        read(a),read(b),ans&=0;
        if(a>b)swap(a,b);
        for(i=1;i<=a;i=j+1){
            j=min(a/(a/i),b/(b/i));
            ans+=(ll)(a/i)*(b/i)*(opt[j]-opt[i-1]);
        }pen(ans),putchar('\n');
    }
    return 0;
}
void prepare(){
    check[1]=mb[1]=1;
    for(ri int i(2),j;i<=limit;++i){
        if(!check[i])prime[++pt]=i,mb[i]=-1;
        for(j=1;j<=pt&&prime[j]<=limit/i;++j){
            check[i*prime[j]]|=true;
            if(!(i%prime[j]))break;
            mb[i*prime[j]]=-mb[i];
        }
    }for(ri int i(1),j;i<=pt;++i)
         for(j=1;j*prime[i]<=limit;++j)
             opt[j*prime[i]]+=mb[j];
    for(ri int i(1);i<=limit;++i)opt[i]+=opt[i-1];
}
il int min(int a,int b){
    return a<b?a:b;
}
void pen(ll x){
    if(x>9)pen(x/10);putchar(x%10+48);
}
il void read(int &x){
    x&=0;ri char c;while(c=getchar(),c<'0'||c>'9');
    while(c>='0'&&c<='9')x=(x<<1)+(x<<3)+(c^48),c=getchar();
}

YY的GCD

原文：https://www.cnblogs.com/a1b3c7d9/p/10781692.html