P4841 [集训队作业2013]城市规划 生成函数

题意:

戳这里

分析:

  • 前置芝士:生成函数

我们设 \(g(n)\) 表示 \(n\) 个点的任意有标号无向图方案数,所以可得 \(g(n)=2^{C_n^2}\) (完全图中枚举每条边是否存在)

我们接着设 \(f(n)\) 表示 \(n\) 个点的有标号无向连通图的方案数,可得 \(\displaystyle g(n)=\sum_{i=1}^n C_{n-1}^{i-1}f(i)g(n-i)\)

证明:

我们假设 1 所在连通块的大小为 \(i\) ,这样连通块内部的方案数为 \(f(i)\) ,外部的方案数为 \(g(n-i)\) ,由于点是有标号的所以还要枚举 1 所在连通块的大小,要乘上 \(C_{n-1}^{i-1}\)

然后开始推柿子

\[g(n)=\sum_{i=1}^n \frac{(n-1)!}{(i-1)!(n-i)!} f(i)g(n-i) \\ \frac{g(n)}{(n-1)!}=\sum_{i=1}^n \frac{f(i)}{(i-1)!}\frac{g(n-i)}{(n-i)!} \]

我们把从左到右的三项的分别看成三种指数型生成函数,记为 \(G1(n),F(n),G2(n)\)

\[G1(n)=\sum_{i=1}^n \frac{g(i)}{(i-1)!}x^i \\ F(n)-\sum_{i=1}^n \frac{f(i)}{(i-1)!}x^i \\ G2(n)=\sum_{i=0}^n \frac{g(i)}{i!}x^i \\ G1(n)=F(n)*G2(n) \\ F(n)=\frac{G1(n)}{G2(n)} \]

多项式求逆即可

代码:

#include<bits/stdc++.h>
#include<bits/stdc++.h>
#define inl inline
#define reg register

using namespace std;

namespace zzc
{
    typedef long long ll;
    inl ll read()
    {
        ll x=0,f=1;char ch=getchar();
        while(!isdigit(ch)){if(ch=='-')f=-1;ch=getchar();}
        while(isdigit(ch)){x=x*10+ch-48;ch=getchar();}
        return x*f;
    }

    const int maxn = 5e5+5;
    const ll mod = 1004535809;
    ll n,l,lim;
    ll g2[maxn],g1[maxn],a[maxn],c[maxn],rev[maxn],fac[maxn],inv[maxn],F[maxn];

    void init()
    {
        fac[0]=fac[1]=1;
        inv[0]=inv[1]=1;
        for(ll i=2;i<=200000;i++) 
        {
            fac[i]=fac[i-1]*i%mod;
            inv[i]=inv[mod%i]*(mod-mod/i)%mod;
        }
        for(ll i=2;i<=200000;i++) inv[i]=inv[i-1]*inv[i]%mod;
    }

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

    void ntt(ll *a,ll inv)
    {
        for(ll i=1;i<=lim;i++) if(i<rev[i]) swap(a[i],a[rev[i]]);
        for(ll len=1;len<lim;len<<=1)
        {
            ll w1=qpow(inv==1?3:334845270,(mod-1)/(len<<1));
            for(int i=0;i<lim;i+=(len<<1))
            {
                ll w=1;
                for(int j=0;j<len;j++)
                {
                    ll x=a[i+j],y=w*a[i+j+len]%mod;
                    a[i+j]=(x+y)%mod;
                    a[i+j+len]=(x-y+mod)%mod;
                    w=w*w1%mod;
                }
            }
        }
        if(inv==1) return ;
        ll tmp=qpow(lim,mod-2);
        for(int i=0;i<=lim;i++) a[i]=a[i]*tmp%mod;
    }

    void get_len(ll k)
    {
        lim=1;l=0;
        while(lim<(k<<1)) lim<<=1,l++;
        for(ll i=0;i<lim;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    }

    void get_inv(ll *f,ll *g,ll k)
    {
        if(k==1)
        {
            g[0]=qpow(f[0],mod-2);
            return ;
        }
        get_inv(f,g,(k+1)>>1);
        get_len(k);
        for(ll i=0;i<k;i++) c[i]=f[i];
        for(ll i=k;i<lim;i++) c[i]=0;
        ntt(g,1);ntt(c,1);
        for(ll i=0;i<lim;i++) g[i]=(2ll-c[i]*g[i]%mod+mod)%mod*g[i]%mod;
        ntt(g,-1);
        for(ll i=k;i<lim;i++) g[i]=0;
    }

    void work()
    {
        init();
        n=read();
        g2[0]=1;
        for(ll i=1;i<=n;i++)
        {
            g2[i]=qpow(2,(i*(i-1)/2)%(mod-1))*inv[i]%mod;
            g1[i]=g2[i]*i%mod;
        }
        get_inv(g2,a,n+1);
        get_len(n+1);
        ntt(g1,1);ntt(a,1);
        for(ll i=0;i<lim;i++) F[i]=g1[i]*a[i]%mod;
        ntt(F,-1);
        printf("%lld\n",F[n]*fac[n-1]%mod);
    }

}

int main()
{
    zzc::work();
    return 0;
}
#define inl inline
#define reg register

using namespace std;

namespace zzc
{
    typedef long long ll;
    inl ll read()
    {
        ll x=0,f=1;char ch=getchar();
        while(!isdigit(ch)){if(ch=='-')f=-1;ch=getchar();}
        while(isdigit(ch)){x=x*10+ch-48;ch=getchar();}
        return x*f;
    }

    const int maxn = 5e5+5;
    const ll mod = 1004535809;
    ll n,l,lim;
    ll g2[maxn],g1[maxn],a[maxn],c[maxn],rev[maxn],fac[maxn],inv[maxn],F[maxn];

    void init()
    {
        fac[0]=fac[1]=1;
        inv[0]=inv[1]=1;
        for(ll i=2;i<=200000;i++) 
        {
            fac[i]=fac[i-1]*i%mod;
            inv[i]=inv[mod%i]*(mod-mod/i)%mod;
        }
        for(ll i=2;i<=200000;i++) inv[i]=inv[i-1]*inv[i]%mod;
    }

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

    void ntt(ll *a,ll inv)
    {
        for(ll i=1;i<=lim;i++) if(i<rev[i]) swap(a[i],a[rev[i]]);
        for(ll len=1;len<lim;len<<=1)
        {
            ll w1=qpow(inv==1?3:334845270,(mod-1)/(len<<1));
            for(int i=0;i<lim;i+=(len<<1))
            {
                ll w=1;
                for(int j=0;j<len;j++)
                {
                    ll x=a[i+j],y=w*a[i+j+len]%mod;
                    a[i+j]=(x+y)%mod;
                    a[i+j+len]=(x-y+mod)%mod;
                    w=w*w1%mod;
                }
            }
        }
        if(inv==1) return ;
        ll tmp=qpow(lim,mod-2);
        for(int i=0;i<=lim;i++) a[i]=a[i]*tmp%mod;
    }

    void get_len(ll k)
    {
        lim=1;l=0;
        while(lim<(k<<1)) lim<<=1,l++;
        for(ll i=0;i<lim;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    }

    void get_inv(ll *f,ll *g,ll k)
    {
        if(k==1)
        {
            g[0]=qpow(f[0],mod-2);
            return ;
        }
        get_inv(f,g,(k+1)>>1);
        get_len(k);
        for(ll i=0;i<k;i++) c[i]=f[i];
        for(ll i=k;i<lim;i++) c[i]=0;
        ntt(g,1);ntt(c,1);
        for(ll i=0;i<lim;i++) g[i]=(2ll-c[i]*g[i]%mod+mod)%mod*g[i]%mod;
        ntt(g,-1);
        for(ll i=k;i<lim;i++) g[i]=0;
    }

    void work()
    {
        init();
        n=read();
        g2[0]=1;
        for(ll i=1;i<=n;i++)
        {
            g2[i]=qpow(2,(i*(i-1)/2)%(mod-1))*inv[i]%mod;
            g1[i]=g2[i]*i%mod;
        }
        get_inv(g2,a,n+1);
        get_len(n+1);
        ntt(g1,1);ntt(a,1);
        for(ll i=0;i<lim;i++) F[i]=g1[i]*a[i]%mod;
        ntt(F,-1);
        printf("%lld\n",F[n]*fac[n-1]%mod);
    }

}

int main()
{
    zzc::work();
    return 0;
}
posted @ 2021-01-27 15:07  youth518  阅读(52)  评论(0编辑  收藏  举报