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_JUDGE 70 freopen("city.in","r",stdin); 71 freopen("city.out","w",stdout); 72 #endif 73 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_JUDGE 64 freopen("city.in","r",stdin); 65 freopen("city.out","w",stdout); 66 #endif 67 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 }