GCD莫比乌斯反演

标签：long   euler   using   cpp   str   read   alc   turn   cti   

链接

题目没找到
求\(\sum_{1}^{n}\sum_{1}^{m}gcd(i,j)\)

机房最后一个学懵逼钨丝的人

\(\sum\limits_{i=1}^{N}\sum\limits_{j=1}^{M} gcd(i,j)\)
\(\sum\limits_{k=1}^{min(N,M)} k*\sum\limits_{i=1}^{N}\sum\limits_{j=1}^{M} [gcd(i,j)==k]=\sum\limits_{k=1}^{min(N,M)} k*\sum\limits_{i=1}^{\frac{N}{k}}\sum\limits_{j=1}^{\frac{M}{k}} [gcd(i,j)==1]\)
\(看这部分\sum\limits_{i=1}^{N}\sum\limits_{j=1}^{M} [gcd(i,j)==k]\)
设\(g(d)=\sum\limits_{i=1}^{N}\sum\limits_{j=1}^{M} [gcd(i,j)==d]\)
\(f(n)=\sum\limits_{n|d}g(d)=[\frac{N}{d}][\frac{M}{d}]\)（显然，不说了）
\(g(n)=\sum\limits_{n|d}\mu(\frac{d}{n})f(d)\)
\(g(1)=\sum\limits_{i=1}^{min(N,M)}\mu(i)[\frac{N}{i}][\frac{M}{i}]\)
\(ans=\sum\limits_{k=1}^{min(N,M)} k*\sum\limits_{i=1}^{min(N,M)}\mu(i)[\frac{N/k}{i}][\frac{M/k}{i}]\)
除法分块嵌套复杂度\(O（n）\)

注意除数不为0

代码

/*
segima gcd n,m
*/
#include <algorithm>
#include <iostream>
#include <cstdio>
#define ll long long
using namespace std;
const int N=15000007,mod=1e9+7;
int read() {
    int x=0,f=1;char s=getchar();
    for(;s>'9'||s<'0';s=getchar()) if(s=='-') f=-1;
    for(;s>='0'&&s<='9';s=getchar()) x=x*10+s-'0';
    return x*f;
}
int pri[N/10],mu[N],cnt;
bool vis[N];
void Euler(int n) {
    vis[1]=mu[1]=1;
    for(int i=2;i<=n;++i) {
        if(!vis[i]) pri[++cnt]=i,mu[i]=-1;
        for(int j=1;j<=cnt&&i*pri[j]<=n;++j) {
            vis[ i*pri[j] ] = 1;
            if(i % pri[j] == 0) {
                mu[ i*pri[j] ] = 0;
                break;
            } else
                mu[ i*pri[j] ] = -mu[i];
        }
    }
    for(int i=1;i<=n;++i) mu[i]+=mu[i-1];
}
int g(int N,int M,int n) {
    int ans=0;
    for(int l=1,r=1;l<=n;l=r+1) {
        r=min(N/(N/l),M/(M/l));
        ans=(ans+1LL*(mu[r]-mu[l-1]+mod)%mod*(N/l)%mod*(M/l)%mod)%mod;
    }
    return ans;
}
ll calc(int a) {return 1LL*a*(a+1)%mod*500000004%mod;}
int main() {
    int N=read(),M=read(),n=min(N,M),ans=0;
    Euler(n);
    for(int l=1,r=1;l<=n;l=r+1) {
        r=min(N/(N/l),M/(M/l));
        ans=(1LL*ans+(1LL*(calc(r)-calc(l-1))%mod+mod)%mod
                *1LL*g(N/l,M/l,min(N/l,M/l)))%mod;
    }
    cout<<ans<<"\n";
    //  ans=0;
    // for(int i=1;i<=N;++i) {
    //  for(int j=1;j<=M;++j) {
    //      ans=(ans+__gcd(i,j))%mod;
    //  }
    // }
    // cout<<ans;
    return 0;
}

GCD莫比乌斯反演

标签：long   euler   using   cpp   str   read   alc   turn   cti   

原文地址：https://www.cnblogs.com/dsrdsr/p/10373532.html