莫比乌斯反演
来源:互联网 发布:网络大电影解析 编辑:程序博客网 时间:2024/05/17 08:43
容斥原理和莫比乌斯
推荐入门:https://wenku.baidu.com/view/dbedced74b73f242326c5f95.html?from=search 把这个吃透基本没问题了
上周多校被莫比乌斯函数虐的不行,好吧,是被数学虐的不行,决定这周怒补数论。不敢说原理都懂,但会用。。专题还没刷完先写一下题解吧。
莫比乌斯反演:由F[n]反向推导出f[n]的过程。常用格式:莫比乌斯 核心:莫比乌斯函数。
结合各种数学定理及推导就可以A很多题啦。
BZOj 2301:[HAOI2011]Problem b
题意:对于给出的n个询问,每次求有多少个数对(x,y),满足a≤x≤b,c≤y≤d,且gcd(x,y) = k,gcd(x,y)函数为x和y的最大公约数。
思路:直接求是很难的,设F[k]=k|gcd(x,y).即gcd为k的倍数的对数。则用常用格式二即可解决这个问题。不过solve()函数的分块思想很巧妙。
/************************************************************** Problem: 2301 User: 1506915059 Language: C++ Result: Accepted Time:10784 ms Memory:2852 kb****************************************************************/ #include <map>#include <set>#include <deque>#include <queue>#include <stack>#include <cmath>#include <vector>#include <cstdio>#include <cstring>#include <iostream>#include <algorithm>#include<functional>using namespace std;typedef long long ll;#define pd(x) printf("%d\n",x)#define cls(a,x) memset(a,x,sizeof(a))const double eps=1e-6;const double PI=acos(-1.0);const int INF=1e9+10;const int MOD=1e9+7;const int N=1e5+10;int vis[N],prim[N],mu[N],sum[N];;void moblus(){ memset(vis,0,sizeof(vis)); mu[1]=1; int tot=0; for(int i=2; i<N; i++) { if(!vis[i]) { prim[tot++]=i; mu[i]=-1; } for(int j=0; j<tot; j++) { if(i*prim[j]>N) break; vis[i*prim[j]]=1; if(i%prim[j]==0) { mu[prim[j]*i]=0; break; } mu[prim[j]*i]=-mu[i]; } } sum[0]=0; for(int i=1; i<N; i++) sum[i]=sum[i-1]+mu[i];}ll solve(int n,int m){ ll ans=0; if(n>m) swap(n,m); for(int i=1,la=0; i<=n; i=la+1) { la=min(n/(n/i),m/(m/i)); ans+=ll(sum[la]-sum[i-1])*(n/i)*(m/i); } return ans;}int main(){ moblus(); int t,a,b,c,d,k; scanf("%d",&t);// int t1=t; while(t--) { scanf("%d%d%d%d%d",&a,&b,&c,&d,&k); ll ans=solve(b/k,d/k)-solve((a-1)/k,d/k)-solve(b/k,(c-1)/k)+solve((a-1)/k,(c-1)/k); printf("%lld\n",ans);// cout<<ans<<endl; } return 0;}
HDU - 1695
GCD
上一题的原型,但这个题不可重复,所以需要做个小处理。减去其交集的贡献的一半即可。
#include <map>#include <set>#include <deque>#include <queue>#include <stack>#include <cmath>#include <vector>#include <cstdio>#include <cstring>#include <iostream>#include <algorithm>#include<functional>using namespace std;typedef long long ll;#define pd(x) printf("%d\n",x)#define cls(a,x) memset(a,x,sizeof(a))const double eps=1e-6;const double PI=acos(-1.0);const int INF=1e9+10;const int MOD=1e9+7;const int N=1e5+10;int vis[N],prim[N],mu[N],sum[N];void moblus(){ memset(vis,0,sizeof(vis)); mu[1]=1; int tot=0; for(int i=2; i<N; i++) { if(!vis[i]) { prim[tot++]=i; mu[i]=-1; } for(int j=0; j<tot&&i*prim[j]<N; j++) { vis[i*prim[j]]=1; if(i%prim[j]==0) { mu[prim[j]*i]=0; break; } mu[prim[j]*i]=-mu[i]; } } sum[0]=0; for(int i=1; i<N; i++) sum[i]=sum[i-1]+mu[i];}ll solve(int n,int m){ ll ans=0; if(n>m) swap(n,m); for(int i=1,la=0; i<=n; i=la+1) { la=min(n/(n/i),m/(m/i)); ans+=ll(sum[la]-sum[i-1])*(n/i)*(m/i); } return ans;}int main(){ moblus(); int t,a,b,c,d,k; scanf("%d",&t); int t1=t; while(t--) { scanf("%d%d%d%d%d",&a,&b,&c,&d,&k); if(k==0) { printf("Case %d: %lld\n",t1-t,0); continue; } ll ans=solve(b/k,d/k)-solve((a-1)/k,d/k)-solve(b/k,(c-1)/k)+solve((a-1)/k,(c-1)/k); int y=min(b,d),x=max(a,c); ll tmp=solve(y/k,y/k)-solve((x-1)/k,y/k)-solve(y/k,(x-1)/k)+solve((x-1)/k,(x-1)/k); if(tmp<0) tmp=0; ans-=tmp/2; printf("Case %d: %lld\n",t1-t,ans); } return 0;}
BZOJ 2005: [Noi2010]能量采集
题意:n*m个方格,每个格点的贡献为2*gcd(i,j)-1。每次给出n*m,求总贡献。
思路:数据1e5。分块查询可达O(logn),这题O(nlogn)即可处理。注意上面两题都是给定gcd,求种数,这题gcd未给定,却要算总值,所以枚举gcd,然后用莫比乌斯反演求出gcd为i的个数f(i),贡献即f(i)*(2*i-1)。
/************************************************************** Problem: 2005 User: 1506915059 Language: C++ Result: Accepted Time:52 ms Memory:2852 kb****************************************************************/ #include <map>#include <set>#include <deque>#include <queue>#include <stack>#include <cmath>#include <vector>#include <cstdio>#include <cstring>#include <iostream>#include <algorithm>#include<functional>using namespace std;typedef long long ll;#define pd(x) printf("%d\n",x)#define cls(a,x) memset(a,x,sizeof(a))const double eps=1e-6;const double PI=acos(-1.0);const int INF=1e9+10;const int MOD=1e9+7;const int N=1e5+10;int vis[N],mu[N],p[N],sum[N];void moblus(){ memset(vis,0,sizeof(vis)); mu[1]=1,sum[0]=0; int tot=0; for(int i=2; i<N; i++) { if(!vis[i]) { p[tot++]=i; mu[i]=-1; } for(int j=0; j<tot&&i*p[j]<=N; j++) { vis[p[j]*i]=1; if(i%p[j]==0) { mu[i*p[j]]=0; break; } mu[i*p[j]]=-mu[i]; } } for(int i=1; i<N; i++) sum[i]=sum[i-1]+mu[i];}ll solve(int n,int m){ ll ans=0; for(int i=1,la=0; i<=n; i=la+1) { la=min(n/(n/i),m/(m/i)); ans+=ll(sum[la]-sum[i-1])*(n/i)*(m/i); } return ans;}int main(){ moblus(); int n,m; while(~scanf("%d%d",&n,&m)) { if(n>m) swap(n,m); ll ans=0; for(int i=1; i<=n; i++) { ll tmp=solve(n/i,m/i);//f(i)--gcd(x,y)==i的数目// printf("%d %lld\n",i,tmp); ans+=tmp*(2*i-1); } printf("%lld\n",ans); } return 0;}
SPOJ - VLATTICE
题意:n*n*n的立方体,从(0,0,0)出发,有多少点可以直达。原题是可以直接看到,一样。n=1e6。
做过POJ-3090可能就很容易推导出来只要满足gcd(i,j,k)==1,这个点就是可行点。但这个题扩展到三维了。其实是一样的,但需要分析x,y,z中有一个0,有两个0和没有0的情况。
#include <cstdio>#include <cstring>#include <iostream>#include <algorithm>#include<functional>using namespace std;typedef long long ll;#define pd(x) printf("%d\n",x)#define cls(a,x) memset(a,x,sizeof(a))const double eps=1e-6;const double PI=acos(-1.0);const int INF=1e9+10;const int MOD=1e9+7;const int N=1e6+10;int vis[N],mu[N],p[N],sum[N];void moblus(){ memset(vis,0,sizeof(vis)); mu[1]=1,sum[0]=0; int tot=0; for(int i=2; i<N; i++) { if(!vis[i]) { p[tot++]=i; mu[i]=-1; } for(int j=0; j<tot&&i*p[j]<=N; j++) { vis[p[j]*i]=1; if(i%p[j]==0) { mu[i*p[j]]=0; break; } mu[i*p[j]]=-mu[i]; } } for(int i=1; i<N; i++) sum[i]=sum[i-1]+mu[i];}ll solve(int n,int f){ ll ans=0; for(int i=1,la=0; i<=n; i=la+1) { la=n/(n/i); ll d=n/i; if(f==0) ans+=ll(sum[la]-sum[i-1])*d*d*d; else ans+=ll(sum[la]-sum[i-1])*d*d; } return ans;}int main(){ moblus(); int t,n; scanf("%d",&t); while(t--) { scanf("%d",&n); ll ans=3;//两个0的时候只有三个点符合gcd==1 ans+=3*solve(n,1);//一个0的时候,要求两个点gcd==1 ans+=solve(n,0);//没有0,从(1,1,1)到(n,n,n); printf("%lld\n",ans); } return 0;}
HDU - 4746
Mophues
题意:Q次查询,每次求[1,n]和[1,m]中有多少对数的gcd的质因子个数小于等于P。
这道题有点超出我的能力范围了,看题解做的,有些地方还不是很明白。打个表会发现5e5以内最多会有18个素因子,所以当P>=18时直接输出n*m即可,注意数据范围。当p<=18的时候,需要枚举gcd,先预处理出容斥系数O(nlogn),然后查询分块,不然会超时。O(nlogn+Qlogn).
const int N=5e5+5;int p[N],num[N],vis[N],mu[N],sum[N][19];void init(){ memset(vis,0,sizeof(vis)); memset(sum,0,sizeof(sum)); mu[1]=1; num[0]=num[1]=0; int tot=0; for(int i=2; i<N; i++) { if(!vis[i]) { mu[i]=-1; p[tot++]=i; } for(int j=0; j<tot&&i*p[j]<N; j++) { vis[i*p[j]]=1; if(i%p[j]==0) { mu[i*p[j]]=0; break; } mu[i*p[j]]=-mu[i]; } } for(int i=2; i<N; i++) num[i]=1; for(int i=2; i<N; i++) for(int j=2*i; j<N; j+=i) num[j]=num[i]+1; for(int i=1; i<N; i++) for(int j=i; j<N; j+=i) sum[j][num[i]]+=mu[j/i]; for(int i=1; i<N; i++) for(int j=0; j<19; j++) sum[i][j]+=sum[i-1][j]; for(int i=0; i<N; i++) for(int j=1; j<19; j++) sum[i][j]+=sum[i][j-1];}int main(){ init(); int t,n,m,p; scanf("%d",&t); while(t--) { scanf("%d%d%d",&n,&m,&p); if(p>=18) { printf("%lld\n",ll(n)*m); continue; } if(n>m) swap(n,m); ll ans=0; for(int i=1,la=0; i<=n; i=la+1) { la=min(n/(n/i),m/(m/i)); ans+=ll(sum[la][p]-sum[i-1][p])*(n/i)*(m/i); } printf("%I64d\n",ans); } return 0;}
HYSBZ - 2818
Gcd
又是一道gcd的题,求n以内有多少对数的gcd为素数。n<=1e7。
时限和内存都给的很足,打表发现1e7内的素数不过百万,我们求一遍的复杂度是logn,枚举gcd,复杂度nlogn。
const int N=1e7+5;int tot,p[N],f[N],mu[N],sum[N];void init(){ memset(f,0,sizeof(f)); mu[1]=1; tot=0; sum[0]=0; for(int i=2;i<N;i++) { if(!f[i]) { p[tot++]=i; mu[i]=-1; } for(int j=0;j<tot&&i*p[j]<N;j++) { f[i*p[j]]=1; if(i%p[j]==0) { mu[i*p[j]]=0; break; } mu[i*p[j]]=-mu[i]; } } for(int i=1;i<N;i++) sum[i]=sum[i-1]+mu[i];}ll solve(int n){ ll ans=0; for(int i=1,la=0;i<=n;i=la+1) { la=n/(n/i); ll t=n/i; ans+=ll(sum[la]-sum[i-1])*t*t; } return ans;}int main(){ init(); int n; while(~scanf("%d",&n)) { ll ans=0; for(int i=0;p[i]<=n&&i<tot;i++) { ans+=solve(n/p[i]); } printf("%lld\n",ans); } return 0;}
HDU - 2841
Visible Trees
题意:和SPOJ - VLATTICE这道题很像,不过这道题是二维的,很简单的莫比乌斯应用。给你n*m的方格,你站在(0,0)点,求有多少个点可以直接看到。
很容易发现(1,1),(1,2),(2,3)这些点可以直接看到,但(2,2),(2,4),(3,6),(6,9)这种点事不能看到的,其实只要满足gcd(x,y)=1那么这个点就是可行点了。于是求n*m的格子有多少个点的gcd()=1。用F(1)表示这些格点的gcd为1的倍数的个数,很容易求得F(1)=n*m。于是用莫比乌斯反演很容易得到f(1)。
const int N=1e5+10;int vis[N],prim[N],mu[N],sum[N];;void moblus(){ memset(vis,0,sizeof(vis)); mu[1]=1; int tot=0; for(int i=2; i<N; i++) { if(!vis[i]) { prim[tot++]=i; mu[i]=-1; } for(int j=0; j<tot; j++) { if(i*prim[j]>N) break; vis[i*prim[j]]=1; if(i%prim[j]==0) { mu[prim[j]*i]=0; break; } mu[prim[j]*i]=-mu[i]; } } sum[0]=0; for(int i=1; i<N; i++) sum[i]=sum[i-1]+mu[i];}ll solve(int n,int m){ ll ans=0; if(n>m) swap(n,m); for(int i=1,la=0; i<=n; i=la+1) { la=min(n/(n/i),m/(m/i)); ans+=ll(sum[la]-sum[i-1])*(n/i)*(m/i); } return ans;}int main(){ moblus(); int t,n,m; scanf("%d",&t); while(t--) { scanf("%d%d",&n,&m); printf("%lld\n",solve(n,m)); } return 0;}
接下来的题都是容斥和莫比乌斯函数的应用:
SPOJ - SQFREE
Square-free integers
题意:求n以内有多少个无平方因子的数。
这个是在那个文库里学的,刚开始并没有什么思路,求无平方因子数,即分解质因子后质数的次数都为1的数。由容斥原理知:ans=1的倍数-(4,9,25..的倍数)+(36,100..的倍数)容易发现每个质数的u(i)=-1,其平方i*i的贡献系数也为-1。即u(i)即是i*i的贡献系数。1-sqrt(n)枚举一遍即可。
const int N=1e7+5;int vis[N],p[N],mu[N],sum[N];void moblus(){ memset(vis,0,sizeof(vis));// for(int i=1; i<N; i++) s[i]=1; int tot=0; sum[0]=0; mu[1]=1; for(int i=2; i<N; i++) { if(!vis[i]) { p[tot++]=i; mu[i]=-1; } for(int j=0; j<tot&&i*p[j]<N; j++) { vis[i*p[j]]=1; if(i%p[j]==0) { mu[p[j]*i]=0; break; } mu[i*p[j]]=-mu[i]; } }}int main(){ moblus(); int t; ll n; scanf("%d",&t); while(t--) { scanf("%lld",&n); ll x=sqrt(n); ll ans=0; for(ll i=1; i<=x; i++) ans+=mu[i]*(n/(i*i)); printf("%lld\n",ans); } return 0;}
HYSBZ - 2440
完全平方数
题意:和上一个题差不多,求第k个无平方因子数。
思路:看到时限给的这么足,那么就可以大胆写了。我们可以在sqrt(n)的时间内确定n以内有多少个无平方因子数,而求第k个无平方因子数,我们可以二分答案。然后判断无平方因子数的个数时候大于等于k,是则更新答案并缩小上边界,反之增大下边届。注意上边界至少3e9。博主以为2e9足够,结果无辜T了好多发。
const int N=1e5+5;int vis[N],p[N],mu[N],sum[N];void moblus(){ memset(vis,0,sizeof(vis));// for(int i=1; i<N; i++) s[i]=1; int tot=0; sum[0]=0; mu[1]=1; for(int i=2; i<N; i++) { if(!vis[i]) { p[tot++]=i; mu[i]=-1; } for(int j=0; j<tot&&i*p[j]<N; j++) { vis[i*p[j]]=1; if(i%p[j]==0) { mu[p[j]*i]=0; break; } mu[i*p[j]]=-mu[i]; } }}int check(ll num,int k){ ll x=sqrt(num); ll ans=0; for(ll i=1; i<=x; i++) ans+=mu[i]*(num/(i*i)); return ans>=k;}int main(){ moblus(); int t,k; scanf("%d",&t); while(t--) { scanf("%d",&k); int ans=0; ll l=1,r=5e9; while(l<=r){ ll mid=(l+r)/2; if(check(mid,k)) { ans=mid; r=mid-1; } else l=mid+1; } printf("%d\n",ans); } return 0;}
The Lottery
题意:给你n和m个数,求n以内不能被这n个数整除的数有多少个。简单容斥。
我们知道n以内能被x整除的数有n/x个,枚举每个数然后减去对应的贡献会有重复多减,由容斥原理知:ans=n-1个数倍数+2个数的倍数-3个数的倍数...注意m只有15个数,很明显爆搜预处理所有组合找出其贡献系数,遍历一遍所有组合即可。
const int N=1e5+10;struct lyq{ ll lcm; int id;} a[N];int n,m,tot;ll b[N];void init(int pos,int num,ll lcm){ if(pos>m) return ; ll g=__gcd(b[pos],lcm); ll tmp=b[pos]/g*lcm; a[tot].id=(num&1?1:-1); a[tot++].lcm=tmp; init(pos+1,num+1,tmp);//选 init(pos+1,num,lcm);//不选;}int main(){ while(~scanf("%d%d",&n,&m)) { for(int i=1; i<=m; i++) scanf("%lld",&b[i]); tot=0; init(1,0,1); int ans=n; for(int i=0; i<tot; i++) {// printf("%d %lld\n",a[i].id,a[i].lcm); ans+=a[i].id*(n/a[i].lcm); } printf("%d\n",ans); } return 0;}
专题里还有不到一半的题没A完,能力和时间有限,对入门的话这些差不多够了,但要完全理解并灵活运用莫比乌斯反演与容斥不是一时半会一蹴而就的。
奉上专题: 容斥原理和莫比乌斯
- 二项式反演,莫比乌斯反演。
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- 莫比乌斯反演
- LinkdeList(链式集合)
- 整数的划分
- JAVA日期格式转换与数据类型转换
- for循环的闭包
- 彻底理解call()及其apply.bind
- 莫比乌斯反演
- ArrayList底层实现
- 《Mining Large Streams of User Data for Personalized Recommendations》笔记
- Tip8 避免给枚举类型的元素提供显示的值
- 1.集合和数组的比较
- 宏晶STC单片机使用STC-ISP串口烧录失败的解决方法及实例汇总 (Ver0.99.15)
- CSS排版
- Linux统计/监控工具SAR详细介绍
- List、Set、Collection、Map的区别和联系