[bzoj3456]城市规划——分治FFT

题目大意：

求n个点的带标号简单无向联通图的数目。

思路：

嗯多项式求逆还不会，到时候会了应该会补吧。

这种和图计数有关的题目一般都是考虑反面计数或者是容斥什么的。

考虑枚举一号点的连通块的大小，然后用总方案数减去这些方案数。

可以得到递推式：
\[
f_{i}=2^{i\choose 2}-\sum_{j=1}^{i-1}{i-1\choose j-1}\times f_{j}\times 2^{i-j\choose2}
\]
后面的式子可以化为卷积的形式：
\[
f_{i}=2^{i\choose 2}-(i-1)!\times \sum_{j=1}^{i-1}\frac{f_j}{(j-1)!}\times \frac{2^{i-j\choose 2}}{(i-j)!}
\]
直接分治\(fft\)即可。

辣鸡bzoj卡不过去，只能交luogu了

/*=======================================
 * Author : ylsoi
 * Time : 2019.1.29
 * Problem : bzoj3456
 * E-mail : [email protected]
 * ====================================*/
#include<bits/stdc++.h>

#define REP(i,a,b) for(int i=a,i##_end_=b;i<=i##_end_;++i)
#define DREP(i,a,b) for(int i=a,i##_end_=b;i>=i##_end_;--i)
#define debug(x) cout<<#x<<"="<<x<<" "
#define fi first
#define se second
#define mk make_pair
#define pb push_back
typedef long long ll;

using namespace std;

void File(){
    freopen("bzoj3456.in","r",stdin);
    freopen("bzoj3456.out","w",stdout);
}

template<typename T>void read(T &_){
    _=0; T fl=1; char ch=getchar();
    for(;!isdigit(ch);ch=getchar())if(ch=='-')fl=-1;
    for(;isdigit(ch);ch=getchar())_=(_<<1)+(_<<3)+(ch^'0');
    _*=fl;
}

const int maxn=13e4+10;
const int mod=1004535809;

ll qpow(ll x,ll y){
    ll ret=1; x%=mod;
    while(y){
        if(y&1)ret=ret*x%mod;
        x=x*x%mod;
        y>>=1;
    }
    return ret;
}

struct ksslbh{
    int lim,cnt,dn[maxn<<2];
    ll g[maxn<<2],ig[maxn<<2];
    void ntt(ll *A,int ty){
        REP(i,0,lim-1)if(i<dn[i])swap(A[i],A[dn[i]]);
        for(int len=1;len<lim;len<<=1){
            ll w= ty==1 ? g[len<<1] : ig[len<<1];
            for(int L=0;L<lim;L+=len<<1){
                ll wk=1;
                REP(i,L,L+len-1){
                    ll u=A[i],v=A[i+len]*wk%mod;
                    A[i]=(u+v)%mod;
                    A[i+len]=(u-v)%mod;
                    wk=wk*w%mod;
                }
            }
        }
    }
    void work(int n,int m,ll *a,ll *b){
        lim=1,cnt=0;
        while(lim<=n+m)lim<<=1,++cnt;
        if(!cnt)cnt=1;
        REP(i,0,lim){
            dn[i]=dn[i>>1]>>1|((i&1)<<(cnt-1));
            if(i>n)a[i]=0;
            if(i>m)b[i]=0;
        }
        g[lim]=qpow(3,(mod-1)/lim);
        ig[lim]=qpow(g[lim],mod-2);
        for(int i=lim>>1;i;i>>=1){
            g[i]=g[i<<1]*g[i<<1]%mod;
            ig[i]=ig[i<<1]*ig[i<<1]%mod;
        }
        ntt(a,1),ntt(b,1);
        REP(i,0,lim-1)a[i]=a[i]*b[i]%mod;
        ntt(a,-1);
        ll inv=qpow(lim,mod-2);
        REP(i,0,lim-1)a[i]=a[i]*inv%mod;
    }
}T;

int n;
ll po[maxn],f[maxn],g[maxn],a[maxn<<2],b[maxn<<2];
ll fac[maxn],ifac[maxn];

void math_init(){
    fac[0]=1;
    REP(i,1,n)fac[i]=fac[i-1]*i%mod;
    ifac[n]=qpow(fac[n],mod-2);
    DREP(i,n-1,0)ifac[i]=ifac[i+1]*(i+1)%mod;
}

void divide(int l,int r){
    if(l==r){
        f[l]=(f[l]*fac[l-1]%mod+po[l])%mod;
        return;
    }
    int mid=(l+r)>>1;
    divide(l,mid);
    //f -> a , [l,mid] -> [0,mid-l]
    //g -> b , [1,r-l] -> [0,r-l-1]
    REP(i,l,mid)a[i-l]=f[i]*ifac[i-1]%mod;
    REP(i,1,r-l)b[i-1]=g[i];
    T.work(mid-l,r-l-1,a,b);
    REP(i,mid+1,r)f[i]=(f[i]-a[i-l-1])%mod;
    divide(mid+1,r);
}

int main(){
    //File();
    read(n);
    math_init();
    ll tmp=1;
    po[1]=1;
    REP(i,2,n){
        tmp=tmp*2%mod;
        po[i]=po[i-1]*tmp%mod;
    }
    REP(i,1,n)g[i]=po[i]*ifac[i]%mod;
    divide(1,n);
    printf("%lld\n",(f[n]+mod)%mod);
    return 0;
}

原文地址：https://www.cnblogs.com/ylsoi/p/10336368.html