莫比乌斯反演例题集 ^_^

Posted papercloud

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了莫比乌斯反演例题集 ^_^相关的知识,希望对你有一定的参考价值。

【luogu 2257】YY的GCD

Problem Here

预备知识

除法分块、莫比乌斯反演

最终公式:

(ans= sum_{T=1}^n ( frac{n}{T}) ( frac{m}{T}) sum_{p|t,isprime[p]=1} mu ( frac{T}{p} ))



#include<bits/stdc++.h>
#define ll long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
inline int read()
{
    int x=0,f=1;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;
}
#define MN 10000000
int T,n,m;
int prime[MN+5],cnt,mu[MN+5];
bool mark[MN+5];
ll sum[MN+5],ans;
inline void init(){
    mark[1]=1;mu[1]=1;register int i,j;
    for(i=2;i<=MN;i++){
        if(!mark[i]) {mu[i]=-1;prime[++cnt]=i;}
        for(j=1;j<=cnt&&prime[j]*i<=MN;j++){
            mark[i*prime[j]]=1;
            if(i%prime[j]==0) break;
            else mu[i*prime[j]]=-mu[i];
        }
    }
    for(i=1;i<=cnt;i++)
    for(j=1;j*prime[i]<=MN;j++) sum[j*prime[i]]+=mu[j];
    for(i=2;i<=MN;i++) sum[i]+=sum[i-1];
}
int main(){
    T=read();init();
    register int i,r;
    while(T--){
        n=read();m=read();int MIN=min(n,m);
        ans=0ll;
        for(i=1;i<=MIN;i=r+1){
            r=min(n/(n/i),m/(m/i));
            ans+=(1ll)*(n/i)*(m/i)*(sum[r]-sum[i-1]);
        }
        printf("%lld
",ans);
    }
    return 0;
}



【HAOI2011】Problem b

Problem Here

预备知识

容斥



#include<bits/stdc++.h>
#define ll long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
inline int read()
{
    int x=0,f=1;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;
}
#define MN 50005
int prime[MN+5],cnt,mu[MN+5];
bool mark[MN+5];
inline void init()
{
    mark[1]=1;mu[1]=1;register int i,j;
    for(i=2;i<=MN;i++)
    {
        if(!mark[i]) {mu[i]=-1;prime[++cnt]=i;}
        for(j=1;j<=cnt&&prime[j]*i<=MN;j++)
        {
            mark[i*prime[j]]=1;
            if(i%prime[j]==0) break;
            else mu[i*prime[j]]=-mu[i];
        }
    }
    for(i=2;i<=MN;i++) mu[i]+=mu[i-1];
}
ll q(int n,int m,int k)
{
    ll ans=0ll;
    register int i,r;
    for(i=1;i*k<=n&&i*k<=m;i=r+1){
        r=min(n/(n/(i*k))/k,m/(m/(i*k))/k);
        ans+=1ll*(mu[r]-mu[i-1])*(n/(i*k))*(m/(i*k));
    }
    return ans;
}
int main()
{
    int n=read(),a,b,c,d,k;
    init();
    while(n--)
    {
        a=read(),b=read();
        c=read(),d=read();
        k=read();
        printf("%lld
",q(b,d,k)-q(a-1,d,k)-q(b,c-1,k)+q(a-1,c-1,k));
    }
    return 0;
}



【NOI2010】能量采集

Problem Here

预备知识

求gcd的新姿势

(sum_{d|n} phi(d) =n)

prove: 考虑小于等于n的数中,与n的gcd为(frac{n}{d})的数的个数是(phi(d)),而所有的(phi(d))相加正好是n

所有满足 i|x,i|y的数,均满足i|gcd(x,y)



//一个正整数,等于它所有因数的欧拉函数之和。
#include<bits/stdc++.h>
#define ll long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
inline int read()
{
    int x=0,f=1;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;
}
#define MN 100000
int prime[MN+5],cnt;
ll phi[MN+5],ans;
bool mark[MN+5];
inline void init()
{
    mark[1]=1;phi[1]=1;register int i,j;
    for(i=2;i<=MN;i++)
    {
        if(!mark[i]) {phi[i]=i-1;prime[++cnt]=i;}
        for(j=1;j<=cnt&&prime[j]*i<=MN;j++)
        {
            mark[i*prime[j]]=1;
            if(i%prime[j]==0){phi[i*prime[j]]=phi[i]*prime[j];break;}
            else phi[i*prime[j]]=(prime[j]-1)*phi[i];
        }
    }
    for(i=2;i<=MN;i++) phi[i]+=phi[i-1];
}
int main(){
    register int n,m;
    n=read(),m=read();init();
    register int i,r;
    for(i=1;i<=n&&i<=m;i=r+1){
        r=min(n/(n/i),m/(m/i));
        ans+=1ll*(n/i)*(m/i)*(phi[r]-phi[i-1]);
    }
    (ans<<=1)-=1ll*n*m;
    printf("%lld
",ans);
    return 0;
}



【luogu 1829】Crash的数字表格

Problem Here

预备知识

有一些莫比乌斯反演的的常见套路:

  • 对于求和gcd的,把gcd提到最前面枚举
  • 对于判断一个数是否为1,用(sum_{d|n} mu(d)=[n==1])来实现,之后在想办法把(mu(d))给提出来
  • 把能够进行除法分块的部分提出来

最终公式:

[ans=sum_{T=1}^{n} ( frac{n}{T}) ( frac{m}{T}) T sum_{d|T} d mu(d)]

而函数 $g(x)= sum_{d|x} d mu(d) $是个积性函数,可以用欧拉筛来求。



//把能够除法分块的部分拖到最前面
//积性函数可以用欧拉筛来求 
#include<bits/stdc++.h>
#define ll long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
inline int read()
{
    int x=0,f=1;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;
}
#define MN 10000007
#define mod 20101009
int prime[MN],cnt;
ll f[MN],ans;
bool mark[MN];
inline void init()
{
    mark[1]=1;f[1]=1;register int i,j;
    for(i=2;i<=MN-7;i++)
    {
        if(!mark[i]==1){prime[++cnt]=i;f[i]=1-i+mod;}
        for(j=1;j<=cnt&&i*prime[j]<=MN-7;j++)
        {
            mark[i*prime[j]]=1;
            if(i%prime[j]==0){f[i*prime[j]]=f[i];break;}
            else f[i*prime[j]]=f[i]*f[prime[j]]%mod;
        }
    }
    for(i=1;i<=MN-7;i++) f[i]=f[i]*i%mod;
    for(i=2;i<=MN-7;i++) (f[i]+=f[i-1])%=mod; 
}
inline ll sum(int x){return (1ll*x*(x+1)/2)%mod;}
int main()
{
    register int i,r,n,m,M;
    n=read(),m=read();M=min(n,m);
    init();
    for(i=1;i<=M;i=r+1){
        r=min(n/(n/i),m/(m/i));
        (ans+=(sum(n/i)*sum(m/i)%mod*(mod+f[r]-f[i-1])%mod)%mod)%=mod;
    }
    printf("%lld
",ans);
    return 0;
}



【luogu 3327】[SDOI2005] 约数个数和

Problem Here

预备知识

  • (d(nm)=sum_{i|n} sum_{j|m} [gcd(i,j)==1])

    prove:

    考虑分别理解左式和右式:

    1. 取nm的约束个数和相当于各个质因子的(次数+1)的累乘
    2. 试想该如何分配i,j中包含质因子p的次数?显然,只有(次数+1)种做法,累乘即为答案
  • (sum _{i=1}^{n} frac{n}{i}= sum_{j=1}^{n} d(j)) ,从右往左理解要方便的多

    具体可以参见 约束研究

  • d(n)是一个积性函数。

最终公式:

(g(n)=sum _{i=1}^{n}frac{n}{i})

[ans=sum_{d=1}^{min(n,m)} mu(d) g( frac{n}{d}) g ( frac{m}{d})]



#include<bits/stdc++.h>
#define ll long long
using namespace std;
inline int read()
{
    int x=0,f=1;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;
}
#define MN 50005
int mu[MN],cnt,prime[MN];
ll d[MN],p[MN];
bool mark[MN];
inline void init()
{
    mark[1]=true;
    p[1]=d[1]=mu[1]=1;
    register int i,j;
    for(i=2;i<=MN-5;++i)
    {
        if(!mark[i]){prime[++cnt]=i;mu[i]=-1;d[i]=2;p[i]=2;}
        for(j=1;j<=cnt&&prime[j]*i<=MN-5;++j)
        {
            mark[i*prime[j]]=true;
            if(i%prime[j]==0)
            {
                mu[i*prime[j]]=0;d[i*prime[j]]=d[i]/p[i]*(p[i]+1);
                p[i*prime[j]]=p[i]+1;break; 
            }
            else
            {
                mu[i*prime[j]]=-mu[i];d[i*prime[j]]=d[i]*2;
                p[i*prime[j]]=2;
            }
        }
    }
    for(i=2;i<=MN-5;i++) mu[i]+=mu[i-1],d[i]+=d[i-1];
} 
int main()
{
    register int i,r,n,m,T;
    T=read(),init();
    while(T--)
    {
        register ll ans=0ll;
        n=read(),m=read();
        for(i=1;i<=n&&i<=m;i=r+1)
        {
            r=min(n/(n/i),m/(m/i));
            ans+=1ll*d[n/i]*d[m/i]*(mu[r]-mu[i-1]);
        }
        printf("%lld
",ans);
    }
    return 0;
}





Blog来自PaperCloud,未经允许,请勿转载,TKS!

以上是关于莫比乌斯反演例题集 ^_^的主要内容,如果未能解决你的问题,请参考以下文章

HDU 1695 GCD 莫比乌斯反演例题

莫比乌斯函数

关于数论莫比乌斯反演

莫比乌斯反演总结

莫比乌斯反演

莫比乌斯函数与莫比乌斯反演