莫比乌斯反演例题集 ^_^
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:
考虑分别理解左式和右式:
- 取nm的约束个数和相当于各个质因子的(次数+1)的累乘
- 试想该如何分配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!
以上是关于莫比乌斯反演例题集 ^_^的主要内容,如果未能解决你的问题,请参考以下文章