首页 > 代码库 > bzoj3456 城市规划
bzoj3456 城市规划
Description
求含有n个点有标号的无向联通图的个数（没有重边），n<=130000
Input
3
Output
4
正解：$分治FFT$/多项式求逆。
并没有权限号，但是某$oj$里有这道题。。
我们考虑递推，$f[i]$表示$i$个点的联通图个数，那么用总数减去不合法的数量。
考虑枚举$1$号点所在的联通块的点数，那么我们可以得到：
$f[n]=2^{\binom{n}{2}}-\sum_{i=1}^{n-1}f[i]*\binom{n-1}{i-1}*2^{\binom{n-i}{2}}$
那么这个式子把组合数拆开以后就能直接上分治$FFT$了，复杂度$O(nlog^{2}n)$
分治$FFT$：
1 #include <bits/stdc++.h> 2 #define il inline 3 #define RG register 4 #define ll long long 5 #define rhl (1004535809) 6 #define N (1000010) 7 #define G (3) 8 9 using namespace std;10 11 int fac[N],ifac[N],inv[N],cal[N],w[N],rev[N],f[N],a[N],b[N],n;12 13 il int qpow(RG int a,RG ll b){14 RG int ans=1;15 while (b){16 if (b&1) ans=1LL*ans*a%rhl;17 a=1LL*a*a%rhl,b>>=1;18 }19 return ans;20 }21 22 il void pre(){23 fac[0]=fac[1]=ifac[0]=ifac[1]=inv[1]=cal[0]=cal[1]=1;24 for (RG int i=2;i<=n;++i){25 fac[i]=1LL*fac[i-1]*i%rhl;26 inv[i]=1LL*(rhl-rhl/i)*inv[rhl%i]%rhl;27 ifac[i]=1LL*ifac[i-1]*inv[i]%rhl;28 cal[i]=qpow(2,1LL*i*(i-1)>>1);29 }30 for (RG int i=1,v=1;i<=(n<<1);++v,i<<=1)31 w[v]=qpow(G,(rhl-1)/(i<<1));32 return;33 }34 35 il void NTT(int *a,RG int n,RG int f){36 for (RG int i=0;i<n;++i) if (i<rev[i]) swap(a[i],a[rev[i]]);37 for (RG int i=1,v=1;i<n;++v,i<<=1){38 RG int gn=w[v],x,y;39 for (RG int j=0;j<n;j+=i<<1){40 RG int g=1;41 for (RG int k=0;k<i;++k,g=1LL*g*gn%rhl){42 x=a[j+k],y=1LL*g*a[j+k+i]%rhl;43 a[j+k]=x+y; if (a[j+k]>=rhl) a[j+k]-=rhl;44 a[j+k+i]=x-y; if (a[j+k+i]<0) a[j+k+i]+=rhl;45 }46 }47 }48 if (f==1) return; reverse(a+1,a+n); RG int inv=qpow(n,rhl-2);49 for (RG int i=0;i<n;++i) a[i]=1LL*a[i]*inv%rhl; return;50 }51 52 il void solve(RG int l,RG int r){53 if (l==r){54 f[l]=(cal[l]-1LL*fac[l-1]*f[l]%rhl+rhl)%rhl;55 return;56 }57 RG int mid=(l+r)>>1; solve(l,mid);58 RG int m=max(mid-l+1,r-mid),len,lg=0;59 for (len=1;len<=(m<<1);len<<=1) ++lg;60 for (RG int i=0;i<len;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg-1)),a[i]=b[i]=0;61 for (RG int i=l;i<=mid;++i) a[i-l]=1LL*f[i]*ifac[i-1]%rhl;62 for (RG int i=1;i<=r-l;++i) b[i]=1LL*cal[i]*ifac[i]%rhl;63 NTT(a,len,1),NTT(b,len,1); for (RG int i=0;i<len;++i) a[i]=1LL*a[i]*b[i]%rhl;64 NTT(a,len,-1); for (RG int i=mid+1;i<=r;++i) f[i]=(f[i]+1LL*a[i-l])%rhl;65 solve(mid+1,r); return;66 }67 68 int main(){69 #ifndef ONLINE_JUDGE70 freopen("city.in","r",stdin);71 freopen("city.out","w",stdout);72 #endif73 cin>>n; pre(); solve(1,n);74 printf("%d\n",f[n]);75 return 0;76 }
我们还可以继续化简，得到：
$\sum_{i=1}^{n}f[i]*\binom{n-1}{i-1}*2^{\binom{n-i}{2}}=2^{\binom{n}{2}}$
$\sum_{i=1}^{n}f[i]*(i-1)!^{-1}*2^{\binom{n-i}{2}}*(n-i)!^{-1}=2^{\binom{n}{2}}*(n-1)!^{-1}$
注意到上式其实是两个多项式卷积等于第三个多项式的形式，那么我们可以用多项式求逆来解决这个问题，复杂度$O(nlogn)$。
多项式求逆：
1 #include <bits/stdc++.h> 2 #define il inline 3 #define RG register 4 #define ll long long 5 #define rhl (1004535809) 6 #define N (1000010) 7 #define G (3) 8 9 using namespace std;10 11 int fac[N],ifac[N],inv[N],cal[N],w[N],rev[N],a[N],b[N],invb[N],tmp[N],n,m;12 13 il int qpow(RG int a,RG ll b){14 RG int ans=1;15 while (b){16 if (b&1) ans=1LL*ans*a%rhl;17 a=1LL*a*a%rhl,b>>=1;18 }19 return ans;20 }21 22 il void pre(){23 fac[0]=fac[1]=ifac[0]=ifac[1]=inv[1]=cal[0]=cal[1]=1;24 for (RG int i=2;i<=n;++i){25 fac[i]=1LL*fac[i-1]*i%rhl;26 inv[i]=1LL*(rhl-rhl/i)*inv[rhl%i]%rhl;27 ifac[i]=1LL*ifac[i-1]*inv[i]%rhl;28 cal[i]=qpow(2,1LL*i*(i-1)>>1);29 }30 for (RG int i=1,v=1;i<=(n<<1);++v,i<<=1)31 w[v]=qpow(G,(rhl-1)/(i<<1));32 return;33 }34 35 il void NTT(int *a,RG int n,RG int f){36 for (RG int i=0;i<n;++i) if (i<rev[i]) swap(a[i],a[rev[i]]);37 for (RG int i=1,v=1;i<n;++v,i<<=1){38 RG int gn=w[v],x,y;39 for (RG int j=0;j<n;j+=i<<1){40 RG int g=1;41 for (RG int k=0;k<i;++k,g=1LL*g*gn%rhl){42 x=a[j+k],y=1LL*g*a[j+k+i]%rhl;43 a[j+k]=x+y; if (a[j+k]>=rhl) a[j+k]-=rhl;44 a[j+k+i]=x-y; if (a[j+k+i]<0) a[j+k+i]+=rhl;45 }46 }47 }48 if (f==1) return; reverse(a+1,a+n); RG int inv=qpow(n,rhl-2);49 for (RG int i=0;i<n;++i) a[i]=1LL*a[i]*inv%rhl; return;50 }51 52 il void getinv(int *a,int *b,RG int n){53 if (n==1){ b[0]=qpow(a[0],rhl-2); return; } getinv(a,b,n>>1);54 for (RG int i=0;i<n;++i) tmp[i]=a[i],tmp[n+i]=0;55 RG int lg=0,m=0; for (m=1;m<=n;m<<=1) ++lg;56 for (RG int i=0;i<m;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg-1));57 NTT(tmp,m,1),NTT(b,m,1);58 for (RG int i=0;i<m;++i) tmp[i]=((b[i]<<1)-1LL*b[i]*b[i]%rhl*tmp[i]%rhl+rhl)%rhl;59 NTT(tmp,m,-1); for (RG int i=0;i<n;++i) b[i]=tmp[i],b[n+i]=0; return;60 }61 62 int main(){63 #ifndef ONLINE_JUDGE64 freopen("city.in","r",stdin);65 freopen("city.out","w",stdout);66 #endif67 cin>>n; pre(); for (m=1;m<=n;m<<=1);68 for (RG int i=1;i<=n;++i) a[i]=1LL*cal[i]*ifac[i-1]%rhl;69 for (RG int i=0;i<=n;++i) b[i]=1LL*cal[i]*ifac[i]%rhl;70 getinv(b,invb,m),NTT(a,m<<1,1),NTT(invb,m<<1,1);71 for (RG int i=0;i<(m<<1);++i) a[i]=1LL*a[i]*invb[i]%rhl;72 NTT(a,m<<1,-1); printf("%lld\n",1LL*a[n]*fac[n-1]%rhl);73 return 0;74 }
bzoj3456 城市规划
声明：以上内容来自用户投稿及互联网公开渠道收集整理发布，本网站不拥有所有权，未作人工编辑处理，也不承担相关法律责任，若内容有误或涉及侵权可进行投诉： 投诉/举报 工作人员会在5个工作日内联系你，一经查实，本站将立刻删除涉嫌侵权内容。