JZOJ3303城市规划

Posted 2022-03-18 horizonwd

description

刚刚解决完电力网络的问题, 阿狸又被领导的任务给难住了.

刚才说过, 阿狸的国家有n 个城市, 现在国家需要在某些城市对之间建立一些贸易路线, 使得整个国家的任意两个城市都直接或间接的连通.

为了省钱, 每两个城市之间最多只能有一条直接的贸易路径. 对于两个建立路线的方案, 如果存在一个城市对, 在两个方案中是否建立路线不一样, 那么这两个方案就是不同的, 否则就是相同的. 现在你需要求出一共有多少不同的方案.

好了, 这就是困扰阿狸的问题. 换句话说, 你需要求出n 个点的简单(无重边无自环)无向连通图数目.

由于这个数字可能非常大, 你只需要输出方案数mod 1004535809(479 * 2 ^21 + 1)即可.

---

analysis

• 看到模数就应该知道是分治\(+NTT\)

• 设\(f[i]表示\)为\(i\)的答案，那么肯定是全部的方案数减去不合法的方案数

• 设\(g[i]=2^\left( \beginmatrix 2\\ i\\ \endmatrix \right)\)，也就是一共有\(\left( \beginmatrix 2\\ i\\ \endmatrix \right)\)条边，就是\(2^\left( \beginmatrix 2\\ i\\ \endmatrix \right)\)种总方案数

• 推推可以知道\(f[n]=g[n]-\sum_i=1^n-1\left( \beginmatrix i-1\\ n-1\\ \endmatrix \right)g[i]*f[n-i]\)，拆组合数之后就是

• \(f[n]=g[n]-\sum_i=1^n-1(n-1)!\over(i-1)!(n-i)!*g[i]*f[n-i]\)，移一下项就是

• \(f[n]=g[n]-(n-1)!\sum_i=1^n-1g[i]\over(i-1)!*f[n-i]\over(n-i)!\)

• 后面那个就是卷积了，好像可以直接套\(NTT\)，然而正着推是\(O(n^2\log_2n)\)，我也不会多项式求逆

• 于是考虑分治，对于\([l,r]\)区间，单独求出所有\([l,mid]\)对\([mid+1,r]\)的贡献

• 具体就是，把已经知道具体值的\(f[l..mid]\)拉出来，与\(g[1..r-l+1]\)乘起来

• 然后两个乘起来，后者长度是前者的两倍，多出来的那一段即为分别对\([mid+1,r]\)每一位的贡献

• 当区间长度为\(1\)时，当前点前面的贡献都已经算完，\(f\)值即为\(g[n]-(n-1)!*贡献\)

• 具体真的很难理解，可以多加思考，或者参考一下这里

• 共\(\log_2n\)层，每层一个\(NTT\)，复杂度即\(O(n\log_2^2n)\)

• 不过其实多项式这块我还是超蒟，还是要多做题……

---

code

#include<stdio.h>
#include<string.h>
#include<algorithm>
#include<math.h>
#define g 3
#define mod 1004535809
#define MAXN 130000
#define MAXLEN MAXN*4
#define ll long long
#define fo(i,a,b) for (ll i=a;i<=b;++i)
#define fd(i,a,b) for (ll i=a;i>=b;--i)

using namespace std;

ll inv[MAXN+5],fac[MAXN+5],C[MAXLEN+5];
ll a[MAXLEN+5],b[MAXLEN+5],f[MAXLEN+5],rev[MAXLEN+5];
ll n,m,len,ans;

inline ll read()

    ll x=0,f=1;char ch=getchar();
    while (ch<'0' || '9'<ch)if (ch=='-')f=-1;ch=getchar();
    while ('0'<=ch && ch<='9')x=x*10+ch-'0',ch=getchar();
    return x*f;

inline ll pow(ll x,ll y)

    ll z=1;
    while (y)
    
        if (y&1)z=z*x%mod;
        y>>=1,x=x*x%mod;
    
    return z;

inline ll work(ll x)

    ll y=1;
    while (y<x)y*=2;
    return y*2;

inline void init()

    fac[0]=fac[1]=1;
    fo(i,2,MAXN)fac[i]=fac[i-1]*i%mod;
    inv[MAXN]=pow(fac[MAXN],(ll)(mod-2))%mod;
    fd(i,MAXN,1)inv[i-1]=inv[i]*i%mod,C[i]=pow(2,i*(i-1)/2);

inline void ntt(ll a[],ll len,ll inv)

    ll bit=0;
    while ((1<<bit)<len)++bit;
    fo(i,0,len-1)
    
        rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
        if (i<rev[i])swap(a[i],a[rev[i]]);
    
    for (ll mid=1;mid<len;mid*=2)
    
        ll tmp=pow(g,(mod-1)/(mid*2));
        if (inv==-1)tmp=pow(tmp,mod-2);
        for (ll i=0;i<len;i+=mid*2)
        
            ll omega=1;
            for (ll j=0;j<mid;++j,omega=omega*tmp%mod)
            
                ll x=a[i+j],y=omega*a[i+j+mid]%mod;
                a[i+j]=(x+y)%mod,a[i+j+mid]=(x-y+mod)%mod;
            
        
    

inline void binary_search(ll l,ll r)

    if (l==r)
    
        f[l]=(C[l]-f[l]*fac[l-1]%mod+mod)%mod;
        return;
    
    ll mid=(l+r)>>1;
    binary_search(l,mid);
    ll len=work(r-l+1),invv=pow(len,mod-2);
    fo(i,0,len-1)a[i]=b[i]=0;
    fo(i,1,mid-l+1)a[i]=f[l+i-1]*inv[l+i-2]%mod;
    fo(i,1,r-l)b[i]=C[i]*inv[i]%mod;
    ntt(a,len,1),ntt(b,len,1);
    fo(i,0,len-1)a[i]=(a[i]*b[i])%mod;
    ntt(a,len,-1);
    fo(i,0,len-1)a[i]=(a[i]*invv)%mod;
    fo(i,mid+1,r)(f[i]+=a[i-l+1])%=mod;
    binary_search(mid+1,r); 

int main()

    n=read(),init(),m=work(n);
    binary_search(1,m/2);
    printf("%lld\n",f[n]);
    return 0;