莫比乌斯反演

· · 个人记录

莫比乌斯反演

形式 1 （卷积法证

如果有

F(x)=\sum\limits_{d|x}f(d)

有如下反演：

f(x)=\sum\limits_{d|x}\mu(d)f(\frac{x}{d})

证明如下（用狄利克雷卷积易证）：

F=f*1 F*\mu=f*{1*u} F*\mu=f*E f=F*\mu

证毕。

形式2 （带入法证

莫比乌斯反演第二种形式

F(n)=\sum_{n|d}f(d)-->>f(n)=\sum_{n|d}\mu(\frac{d}{n})F(d)

证明如下：

将F(d)带入右边式子，并设k=\frac{d}{n}有：

f(n)=\sum_{k=1}^{up}\mu(k)\sum_{kn|t}f(t)

由于\mu()的特殊性质，我们试图调换枚举顺序构造其卷积形式。

=\sum_{n|t}f(t)\sum_{k|{\frac{t}{n}}}\mu(k)

(解释一下就是我们注意到t的取值是nk的倍数，且k是从1至无穷大，也就是说t在满足是n的倍数的同时要满足k的倍数，那么把k甩出去后考虑当t = a * n时有哪些k满足k|a,于是自然而然把\mu的求和符号限制改成了\sum\limits_{k|{\frac{t}{n}}})

好了现在当且仅当t/n为1时\mu项不为0，于是

=f(n)

证毕

一些基础转化：

求证

\sum\limits_{d|n}\mu(x)=E

其中E为原函数（只有1才为1，其余为0）

证明如下：

F(x)=\sum\limits_{d|n}\mu(d)

因为\mu为积性函数，所以其和函数也为积性函数。

当n=1时，F(1)=\mu(1)=1

设n>1，分解n

F(p^k)=\sum\limits_{d|{p^k}}\mu(d) =\mu(1)+\mu(p)+\mu(p^2)+....+\mu(p^k) =1+-1+0+0+...+0=0

证毕

板题大赏

板题1 [POI2007]ZAP-Queries

方法1 强行推式子 + 卷积基本性质

题目相当与求

\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j==k)]

第一步把k提出来

\sum_{i=1}^{\frac{n}{k}}\sum_{j=1}^{\frac{m}{k}}[\gcd(i,j)==1]

第二步由E=1*\mu转换

\sum_{i=1}^{\frac{n}{d}}\sum_{j=1}^{\frac{m}{d}}\sum_{d|\gcd(i,j)}\mu(d)

注意到枚举i,j再来算gcd太假，于是考虑枚举gcd，直接计算有多少个i，j满足\gcd(i,j)==d

\sum_{d=1}^{up}\mu(d)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor

观察下式子，我们发现右侧\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor随着d的枚举最多只会有2\sqrt{n}次变化，那么预处理初\mu的前缀和，即可O(\sqrt{n})算出答案。

方法2 莫比乌斯反演

我们要求的式子设为f(x)难以简单求得，但是我们考虑f(x)的 变换函数 F(n)=\sum\limits_{n|d}f(d)的意义是所有gcd(i,j)为n的倍数的个数，其值显然为：

F(d)=\lfloor \frac{N}{d}\rfloor\lfloor \frac{M}{d}\rfloor

那还有什么说的，莫比乌斯反演形式2，

f(n)=\sum_{n|dd}\mu(\frac{dd}{n})F(dd)

将F函数带入，设k=dd/n：

f(n)=\sum_{k=1}^{\frac{N}{n}}\mu(k)\lfloor \frac{N}{nk}\rfloor\lfloor \frac{M}{nk}\rfloor

因为答案即为f(d)即为：

\sum_{k=1}^{\frac{N}{d}}\mu(k)\lfloor \frac{N}{dk}\rfloor\lfloor \frac{M}{dk}\rfloor

整数分块即可。

code

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int M=1e5+5;
int cntp,n,m,T,u[M],pri[M],vis[M],s[M];
inline void pre(){
    u[1]=1;u[2]=-1;
    for(int i=2;i<=100000;i++){
        if(!vis[i]){
            pri[++cntp]=i;u[i]=-1;
        }
        for(int j=1;j<=cntp&&i*pri[j]<=100000;j++){
            vis[pri[j]*i]=1;
            if(i%pri[j]==0) u[i*pri[j]]=0;
            else u[i*pri[j]]=u[i]*-1; 
            if(i%pri[j]==0) break;
        }
    }
    for(int i=1;i<=100000;i++) s[i]=s[i-1]+u[i];
}
signed main(){
    pre();
    cin>>T;
    while(T--){
        int a,b,k,ans=0;
        scanf("%lld%lld%lld",&a,&b,&k);
        int lim=min(a/k,b/k);
        a=a/k,b=b/k;int E;
        for(int S=1;S<=lim;S=E+1){
            E=min(a/(a/S),b/(b/S));
            ans+=(s[E]-s[S-1])*(a/S)*(b/S);
        }cout<<ans<<"\n";
    }
    return 0;
}

一点拓展

如何题目要求i在a-b内，j在c-d内的\sum[\gcd(i,j)]呢？

考虑简单容斥，Ans=ans(b,d)-ans(b,c)-ans(a,c)+ans(a,b)

于是就有了这道题：

[HAOI2011]Problem b

成功双倍经验

板题2 YY的GCD

我们现在多了个限制，只要\gcd(i,j)==p,p\in Prime即可

即求

\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)==P,P\in Prime]

同上一道题，我们先试图化简式子，

\sum\limits_{p\in Prime}\sum_{d=1}^{\frac{n}{p}}\mu(d)\lfloor\frac{n}{dp}\rfloor\lfloor\frac{m}{dp}\rfloor

现在看似差不多化成最简式了，现在考虑如何搞这个质数。

我们希望把枚举p搞成\sum_{p|x}，因为这玩意看着很积性。

那我们搞一个k=d*p

于是原式等于：

\sum_{k=1}^n\sum_{p\in Prime,p|k}\mu(\frac{k}{p})\lfloor\frac{n}{k}\rfloor\lfloor\frac{m}{k}\rfloor

关于这个变换，每一个合法的k,p都对应一组d,p，且恰好覆盖。

然后将剩下的棘手的p的项合并

f(x)=\sum_{p\in Prime,p|k}\mu({\frac{k}{p}}) \sum_{k=1}^nf(k)\lfloor\frac{n}{k}\rfloor\lfloor\frac{m}{k}\rfloor

现在只要把f筛出来就结束了。观察f(x)的性质，积性石锤

设n,m互质 $$f(nm)=\sum_{p\in Prime,p|nm}\mu(\frac{nm}{p})$$ $$=\sum_{p|n}\mu(\frac{n}{p})\sum_{p|m}\mu(\frac{m}{p})=f(n)f(m)$$ 于是可以线性筛$f(x)

具体分析如下（实在敲不动latex了，随便找了一份）

)

其实就是按照内容分三类分析即可。

code

#include<bits/stdc++.h>
using namespace std;
inline int getint(){
    int summ=0,f=1;char ch;
    for(ch=getchar();!isdigit(ch)&&ch!='-';ch=getchar());
    if(ch=='-')f=-1,ch=getchar();
    for(;isdigit(ch);ch=getchar()) summ=(summ<<3)+(summ<<1)+ch-48;
    return summ*f;
}
const int M=1e7;
int cntp,n,m,T,u[M+5],pri[M+5],vis[M+5],s[M+5],f[M+5];
inline void pre(){
    u[1]=1;u[2]=-1;
    for(int i=2;i<=M;i++){
        if(!vis[i]){
            pri[++cntp]=i;u[i]=-1;f[i]=1;
        }
        for(int j=1;j<=cntp&&i*pri[j]<=M;j++){
            vis[pri[j]*i]=1;
            if(i%pri[j]==0) u[i*pri[j]]=0,f[i*pri[j]]=u[i];
            else u[i*pri[j]]=u[i]*-1,f[i*pri[j]]=-f[i]+u[i]; 
            if(i%pri[j]==0) break;
        }
    }
    for(int i=1;i<=M;i++) s[i]=s[i-1]+f[i];
}
signed main(){
    pre();
    cin>>T;
    while(T--){
        int a,b,k;long long ans=0;
        a=getint();b=getint();
        int E;
        for(int S=1;S<=min(a,b);S=E+1){
            E=min(a/(a/S),b/(b/S));
            //cout<<s[E]<<" "<<s[S-1]<<" "<<((a/S)*(b/S))<<endl;
            ans+=1ll*(s[E]-s[S-1])*1ll*(a/S)*(b/S);
        }cout<<ans<<"\n";
    }
    return 0;
}

[SDOI2015]约数个数和

题目要求：

\sum_{i=1}^n\sum_{j=1}^md(ij)

我们知道一组因子个数为\sum\limits_{d|n}1

那如果加一维呢，\sum\limits_{d|n}\sum\limits_{d|m}1 ？

但是这样做显然有重复。

由于d(ij)是积性函数，我们可以只考虑单一素数幂的情况。

我们构造一种分配方式，若i能提供幂则只给i，若不能提供则只由j提供差值。这种情况下

每一个幂的数量都有且对应一个方案，不会重复，于是只要不取不互质的因数组成的数对有且对应一种合法因数。

于是有

d(ij)=\sum_{a|i}\sum_{b|j}[\gcd(a,b)==1]

至此d函数被搞成了我们喜闻乐见的gcd函数

begin to transform

\sum_{i=1}^n\sum_{j=1}^m\sum_{a|i}\sum_{b|j}[\gcd(a,b)==1]

change the order

\sum_{a=1}^n\sum_{b=1}^m\lfloor\frac{n}{a}\rfloor\lfloor\frac{m}{b}\rfloor\sum_{d|\gcd(a,b)}\mu(d)

continue to change the order

\sum_{d=1}^{\min(n,m)}\mu(d)\sum_{a=1}^{\frac{n}{d}}\lfloor\frac{n}{a}\rfloor\sum_{b=1}^{\frac{m}{d}}\lfloor\frac{m}{b}\rfloor

后面的东西可以预处理，又变回了整数分块的情况，O(T\sqrt{n})内可以解决

code

#include<bits/stdc++.h>
using namespace std;
#define int long long
int n,m,ans;
int s[50005],p[50005],mo[50005],mos[50005];
void build()
{
    for(int i=1;i<=50000;i++)
    mo[i]=1;
    for(int i=2;i<=50000;i++)
    {
        if(p[i]) continue;
        mo[i]=-1;
        for(int j=i*2;j<=50000;j+=i)
        {
            p[j]=1;
            if((j/i)%i==0) mo[j]=0;
            else mo[j]*=-1;
        }
    }
    for(int i=1;i<=50000;i++)
    {
        for(int l=1,r;l<=i;l=r+1)
        {
            r=(i/(i/l));
            s[i]+=((r-l+1)*(i/l));
        }
    }
    for(int i=1;i<=50000;++i) mos[i]=mos[i-1]+mo[i];
}
signed main()
{
    build();
    int T;
    cin>>T;
    while(T--)
    {
        ans=0;
        scanf("%lld%lld",&n,&m);
        int mi=min(n,m);
        for(int r,l=1;l<=mi;l=r+1)
        {
            r=min(n/(n/l),m/(m/l));
            ans+=(mos[r]-mos[l-1])*(s[n/l]*s[m/l]);
        }
        cout<<ans<<endl;
    }
    return 0;
}

例题4 P5221 Product

搞成只有gcd后，把gcd分一边，其他一次项丢一边。（累乘可以分到分子分母算）

得

ans=\frac{N!N!}{(\prod_{d=1}^nd^{\sum_{g=1}^{\frac{n}{d}}\mu(g)\frac{n}{dg}})^2}

对于一个确定的d，可以整数分块O(\sqrt n),注意到d的上界也可以整数分块，于是分块套分块，可以在O(n)时间内得到解决。

code

#include<bits/stdc++.h>
using namespace std;
const int N=1e6+5,mod=104857601;
int pri[N/10],mu[N],n,m,tot,ans=1;
bool vis[N];
inline void Pre(){
    mu[1]=1;
    for(int i=2;i<=n;i++){
        if(!vis[i]) pri[++tot]=i,mu[i]=-1;
        for(int j=1;j<=tot&&pri[j]*i<=n;j++){
            vis[i*pri[j]]=true;
            mu[i*pri[j]]=-mu[i];
            if(i%pri[j]==0){
                mu[i*pri[j]]=0;break;
            }
        }
    }
    for(int i=1;i<=n;i++) mu[i]+=mu[i-1];
}
inline int ksm(int x,int y){
    int res=1;
    while(y){
        if(y&1) res=1ll*res*x%mod;
        x=1ll*x*x%mod;y>>=1;
    }
    return res;
}
inline int calc(int up){
    int res=0;
    for(int l=1,r;l<=up;l=r+1){
        r=n/(n/l);
        res=(res+1ll*(mu[r]-mu[l-1]+mod-1)*(up/l)%(mod-1)*(up/l)%(mod-1))%(mod-1);
    }
    return res;
}
int main(){
    cin>>n;Pre();
    for(int i=1;i<=n;i++) ans=1ll*ans*i%mod;
    ans=ksm(ans,2*n);int ll=1,rr=1,c1=2,c2=2;
    for(int l=1,r;l<=n;l=r+1){
        r=n/(n/l);
        while(c1<=l-1){
            ll=1ll*ll*c1%mod;
            c1++;
        }
        while(c2<=r){
            rr=1ll*rr*c2%mod;
            c2++;
        }
        ans=1ll*ans*ksm(1ll*ll*ksm(rr,mod-2)%mod,2*calc(n/l)%(mod-1))%mod;
    }
    cout<<1ll*(ans+mod)%mod<<endl;
    return 0;
}

例题5 P3768 简单的数学题

推式子时间

在这只列出关键步骤。

\sum_{x=1}^nx^3\sum_{d=1}^{n/x}d^2\mu(d)G(\frac{n}{x*d})^2

其中G(x)=\frac{n*(n+1)}{2}

改变累加方式，注意到有多处可以用x*d换元并试图构造\mu*Id的形式（套路）

设T=x*d

\sum_{T=1}^nG(\frac{n}{T})T^2\sum_{x|T}\mu(x)(T/x)

有卷积知识\mu*Id=\phi

\sum_{T=1}^nG(\frac{n}{T})T^2\phi(T)

设F(x)=x^2\phi(x)，其是积性函数，可以用杜教筛加速。

\sum_{T=1}^nG(\lfloor\frac{n}{T}\rfloor)F(x)

整数分块加杜教筛搞定。

code

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int M=5e6+5,up=5e6;
map <int,int> mp;
int pri[M],vis[M],phi[M],mod,n,tot,s[M];
inline int Ksm(int x,int y){
    int res=1;
    while(y){
        if(y&1) res=res*x%mod;
        x=x*x%mod;y>>=1;
    }
    return res;
}
void Pre(){
    phi[1]=1; 
    for(int i=2;i<=up;i++){
        if(!vis[i]){
            pri[++tot]=i;phi[i]=i-1;
        }
        for(int j=1;j<=tot&&i*pri[j]<=up;j++){
            vis[i*pri[j]]=1;phi[i*pri[j]]=phi[i]*phi[pri[j]];
            if(i%pri[j]==0){
                phi[i*pri[j]]=phi[i]*pri[j];
                break;
            }
        }
    } 
    for(int i=1;i<=up;i++) s[i]=(s[i-1]+phi[i]*i%mod*i)%mod;
}
int inv6,inv2;
inline int G(int x){
    int y=x%mod*(x+1)%mod*inv2%mod;
    return y*y%mod;
}
inline int H(int x){
    return x%mod*(x+1)%mod*(2*x+1)%mod*inv6%mod;
}
int Get_phi(int x){
    if(x<=up) return s[x];
    if(mp[x]) return mp[x];
    int ans=G(x%mod),res=0;
    for(int l=2,r;l<=x;l=r+1){
        r=x/(x/l);
        res=(res+(H(r%mod)-H((l-1)%mod))%mod*Get_phi(x/l))%mod;
    }
    return mp[x]=(ans-res)%mod;
}
signed main(){
    cin>>mod>>n;
    inv6=Ksm(6,mod-2);inv2=Ksm(2,mod-2);
    Pre();int ans=0;
    for(int l=1,r;l<=n;l=r+1){
        r=n/(n/l);
        ans=(ans+((Get_phi(r)-Get_phi(l-1))%mod*G(n/l%mod)))%mod;
    }
    cout<<(ans+mod)%mod;
    return 0;
}

例题6 Crash的数字表格

\sum_{d=1}^nd\sum_i^{n/d}\sum_j^{m/d}[gcd(i,j)==1]i*j

注意到后面的数值仅仅d有关，且只有\sqrt n种不同变换，故可以单独考虑后面的东西，即

\sum_i^{N}\sum_j^{M}[gcd(i,j)==1]i*j,N=n/d,M=m/d

进行一波反演：

\sum_{d=1}\sum_{d|i}\sum_{d|j}\mu(d)*i*j

变换一下：

\sum_{d=1}^n\mu(d)*d^2*\sum^{n/d}\sum^{m/d}i*j

注意到给定d，后面那一坨与i,j相关的可以O(1)算出来，于是再一次整数分块即可。

code

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int mod=20101009,up=1e7;
int mu[up+5],vis[up+5],pri[up],cnt,n,m,ans;
void Build(){
    mu[1]=1;
    for(register int i=2;i<=up;i++){
        if(!vis[i]) pri[++cnt]=i,mu[i]=-1;
        for(register int j=1;j<=cnt&&pri[j]*i<=up;j++){
            vis[pri[j]*i]=1;mu[i*pri[j]]=mu[i]*mu[pri[j]];
            if(i%pri[j]==0){
                mu[i*pri[j]]=0;break;
            }
        }
    }
    for(int i=1;i<=up;i++) mu[i]=(mu[i-1]+mu[i]*i%mod*i)%mod;
}
inline int F(int l,int r){
    return (l*(l+1)/2%mod)*(r*(r+1)/2%mod)%mod;
}
inline int calc(int l,int r){
    int res=0;
    for(int ll=1,rr;ll<=min(l,r);ll=rr+1){
        rr=min(l/(l/ll),r/(r/ll));
        res=(res+F(l/ll,r/ll)*(mu[rr]-mu[ll-1]))%mod;
    }
    return res;
}
signed main(){
    Build();
    cin>>n>>m;
    for(int ll=1,rr;ll<=min(n,m);ll=rr+1){
        rr=min(n/(n/ll),m/(m/ll));
        ans=(ans+(ll+rr)*(rr-ll+1)/2%mod*calc(n/ll,m/ll))%mod;
    }
    cout<<(ans+mod)%mod<<endl;
    return 0;
}