• bzoj3456-城市规划


    题意

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

    分析

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

    [egin{aligned} h(n)&=sum _{k=1}^nsum _{x_1+x_2+cdots +x_k=n}inom n {x_1}f(x_1)inom {n-x_1}{x_2}f(x_2)cdots inom {n-x_1-x_2cdots x_{k-1}}{x_k}f(x_k) \ &=n!sum _{k=1}^nprod _{sum x=n}frac{f(x_i)}{x_i!} end{aligned} ]

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

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

    [frac{h(n)}{n!}=sum _{k=1}^nprod _{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(nlog 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
    ",ans);
        return 0;
    }
    
  • 相关阅读:
    [bzoj 2460]线性基+贪心+证明过程
    [Wc2011] Xor
    [BZOJ2844]线性基+xor本质不同第K大
    洛谷3857 [TJOI2008]彩灯
    HDU3949 异或线性基
    hdu3062 party --2-sat
    KM算法详解+模板
    Hopcroft-Karp算法
    bzoj 1135: [POI2009]Lyz
    hall定理的证明
  • 原文地址:https://www.cnblogs.com/owenyu/p/7576850.html
Copyright © 2020-2023  润新知