bzoj3456-城市规划

题意

\(n\) 个点的简单无向连通图个数。\(n\le 130000\)

分析

\(f(n)\)\(n\) 个点的 带标号 简单无向连通图的个数,那么总的简单无向图个数 \(h(n)\) 就是

\[\begin{aligned} h(n)&=\sum _{k=1}^n\sum _{x_1+x_2+\cdots +x_k=n}\binom n {x_1}f(x_1)\binom {n-x_1}{x_2}f(x_2)\cdots \binom {n-x_1-x_2\cdots x_{k-1}}{x_k}f(x_k) \\ &=n!\sum _{k=1}^n\prod _{\sum x=n}\frac{f(x_i)}{x_i!} \end{aligned} \]

\(n\) 个点的简单无向图个数就是所有连通块情况的和,其中 \(f\) 内部是带标号的,而连通块之间是无序的。显然,\(h(n)\) 的另一种计算方式是任选一些边连起来,即 \(h(n)=2^\binom n 2\)

显然这是一个生成函数的形式,有

\[\frac{h(n)}{n!}=\sum _{k=1}^n\prod _{\sum x=n}\frac{f(x_i)}{x_i!} \]

\[H=\sum _{k=0}^\infty \frac{h(k)}{k!}x^k \\ F=\sum _{k=0}^\infty f(k)x^k \]

那么就有

\[H=\sum _{k=1}^\infty \frac{F^k}{k!} \]

\[H=\exp F-1 \]

我们不需要计算无限次数的生成函数,只需要计算 \(n\) 位就行了。因为 \(H\) 非常好求,所以我们用 \(H\) 来得到 \(F\) ,即

\[F=\ln (H+1) \]

计算 \(\ln\) 的方法 即可。复杂度为 \(O(n\log n)\)

代码

今天尝试找出一些好的写法,变成一个自己的多项式风格,然而几乎是失败了——封装一些东西导致它跑得很慢,差不多其他人的两倍时间,这是因为之前用了 vector 。后来发现复杂度没什么差别,而且 vector 在这种情况下可能还会慢一点。于是就改回了数组实现。

核心思想是用一个 P 指针实现递归。

#include<bits/stdc++.h>
using namespace std;
typedef long long giant;
const int maxn=1<<18;
const int q=479<<21|1;
const int toro=3;
typedef vector<int> vec;
inline int Plus(int x,int y) {return ((giant)x+(giant)y)%q;}
inline void Pe(int &x,int y) {x=Plus(x,y);}
inline int Multi(int x,int y) {return (giant)x*y%q;}
inline void Me(int &x,int y) {x=Multi(x,y);}
inline int Sub(int x,int y) {return Plus(x,q-y);}
inline void Se(int &x,int y) {x=Sub(x,y);}
template<typename T> inline int mi(int x,T y) {
    int ret=1;
    for (;y;y>>=1,Me(x,x)) if (y&1) Me(ret,x);
    return ret;
}
inline int Inv(int x) {return x==1?1:Multi(q-q/x,Inv(q%x));}
//inline int Inv(int x) {return mi(x,q-2);}
inline int Div(int x,int y) {return Multi(x,Inv(y));}
int Wn[22][2],inv[maxn],fac[maxn],ifac[maxn],rev[maxn],n;
inline void getM(int n,int &M) {for (M=1;M<=(n<<1);M<<=1);}
struct P {
    int a[maxn],M;
    inline P (int m=0) {
        memset(a,0,sizeof a);
        give(m);
    }
    inline int& operator [] (int x) {return a[x];}
    inline void set(int n) {getM(n,M);}
    inline void give(int m) {M=m;}
    inline void deri() {
        for (int i=0;i<M-1;++i) a[i]=Multi(a[i+1],i+1);
        a[M-1]=0;
    }
    inline void inte() {
        for (int i=M-1;i>0;--i) a[i]=Multi(a[i-1],::inv[i]);
        a[0]=0;
    }
    inline void ntt(int len,int op=0) {
        for (int i=0,xj=__builtin_ctz(len);i<len;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<(xj-1));
        for (int i=0;i<len;++i) if (i<rev[i]) swap(a[i],a[rev[i]]);
        for (int i=2,the=1;i<=len;i<<=1,++the) 
            for (int j=0,&wn=Wn[the][op];j<len;j+=i) 
                for (int k=j,w=1;k<j+i/2;++k,Me(w,wn)) {
                    const int u=a[k],v=Multi(a[k+i/2],w);
                    a[k]=Plus(u,v),a[k+i/2]=Sub(u,v);
                }
        if (op) for (int i=0,iv=::Inv(len);i<len;++i) Me(a[i],iv);
    }
    inline void operator *= (P &b) {
        assert(M==b.M);
        ntt(M),b.ntt(b.M);
        for (int i=0;i<M;++i) Me(a[i],b[i]);
        ntt(M,1),b.ntt(b.M,1);
    }
    P* inv(int m) {
        P *ret;
        if (m==2) {
            ret=new P(m);
            (*ret)[0]=::Inv(a[0]);
            return ret;
        }
        const int hf=m>>1;
        ret=inv(hf);
        ret->give(m);
        static P aux;
        aux.give(m);
        for (int i=0;i<m;++i) aux[i]=i<hf?a[i]:0;
        aux.ntt(m),ret->ntt(m);
        for (int i=0;i<m;++i) {
            int &x=(*ret)[i];
            Me(x,Sub(2,Multi(x,aux[i])));
        }
        ret->ntt(m,1);
        for (int i=hf;i<m;++i) (*ret)[i]=0;
        return ret;
    }
    P ln() {
        static P pr,*iv;
        (pr=*this).deri();
        pr*=*(iv=inv(M));
        pr.inte();
        delete iv;
        return pr;
    }
} F,H;
int main() {
#ifndef ONLINE_JUDGE
    freopen("test.in","r",stdin);
#endif
    scanf("%d",&n);
    for (int i=1;i<22;++i) Wn[i][1]=Inv(Wn[i][0]=mi(toro,(q-1)>>i));
    if (n<2) puts("1"),exit(0); 
    inv[1]=fac[0]=fac[1]=ifac[0]=ifac[1]=1;
    for (int i=2;i<maxn;++i) {
        inv[i]=Multi(q-q/i,inv[q%i]);
        fac[i]=Multi(fac[i-1],i);
        ifac[i]=Multi(ifac[i-1],inv[i]);
    }
    H.set(n);
    H[0]=1;
    for (int i=1;i<=n;++i) H[i]=Multi(mi(2,(giant)i*(i-1)/2),ifac[i]);
    F=H.ln();
    int ans=Multi(F[n],fac[n]);
    printf("%d\n",ans);
    return 0;
}
posted @ 2017-09-22 19:52  permui  阅读(187)  评论(0编辑  收藏  举报