A simple Blog for wyx I've been down the bottle hoping.
3456: 城市规划 CDQ分治+NTT 计数问题
发表于: | 分类: Oi | 评论:0 | 阅读:95

求$n$个点无向连通图的方案数

首先我们非常容易得到一个$n^2$的递推式,设$f(i)$表示$i$个点无向连通图的方案数,那么我们用所有情况的方案数减去不合法的方案数,不合法的话我们可以枚举$1$所在的联通块的大小进行计算,用$t(x)$表示$x$个点形成的完全图的边数,容易得到

$$f(i) = 2^{t(i)}-\sum_{j=1}^{i-1}C_{i-1}^{j-1}f(j)*2^{t(i-j)}$$

看不懂的多看几遍就知道啥意思了……

然后就是喜闻乐见的套路时间,考虑推式子然后分治求解

后面那个东西长的就像卷积嘛233333

$C_{i-1}^{j-1} = C_{i-1}^{i-j}$

诶好像是和谐了一些,但是下面是$i-1$好像不是非常好办,我们展开他

$C_{i-1}^{i-j} = \frac{(i-1)!}{(i-j)!(j-1)!}$

这就非常好了嘛~把$(i-1)!$扔到最前面,然后$(j-1)^{-1}$和$f(j)$放在一起刚好,最后剩下两项带$(i-j)$的直接和$j$的两项放一起一波$NTT$带走就好了

最后整理一下后面变成了

$$(i-1)!\sum_{j=1}^{i-1}\frac{f(j)}{(j-1)!}\frac{2^{i-j}}{(i-j)!}$$

接下来是丧心病狂的卡常时间,具体怎么卡自己看代码吧……不建议对自己常数没有信心的同学写这个做法

#include <stdio.h>
#include <string.h>
#include <iostream>
#include <algorithm>
using namespace std;
const int N = 524288+5;
const int mod = 1004535809;
typedef long long LL;
inline int read() {
    int x=0,f=1;char ch = getchar();
    while(ch < '0' || ch > '9'){if(ch == '-')f=-1; ch = getchar();}
    while(ch >='0' && ch <='9'){x=(x<<1)+(x<<3)+ch-'0';ch = getchar();}
    return x*f;
}

LL fac[N], ifac[N], h[N], t[N], F[N];

inline LL pow(LL a,LL b) {
    LL res = 1;
    for(;b;b>>=1,a=(LL)a*a%mod)
        if(b&1)
            res = (LL) res * a % mod;
    return res;
}

inline int C(int n,int m) {
    if(n < m) return 0;
    return (LL)fac[n]*ifac[m]%mod*ifac[n-m]%mod;
}

inline void init() {
    fac[0] = 1;
    static const int Max = 130000;
    for(int i=1;i<=Max;++i) fac[i] = (LL)fac[i-1] * i % mod;
    ifac[Max] = pow(fac[Max], mod-2);
    for(int i=Max-1;~i;--i) ifac[i] = (LL)ifac[i+1] * (i+1) % mod;
    for(int i=0;i<=Max;++i) t[i] = C(i,2);  
    for(int i=0;i<=Max;++i) h[i] = pow(2,t[i]);
}

inline void NTT(LL *a, int len, int type) {
    register int i,k,j; 
    
    for(i=j=0; i<len; ++i) {
        if(i > j) swap(a[i], a[j]);
        for(k = len >> 1; (j^=k) < k; k >>= 1);
    }
    for(i=2;i<=len;i<<=1) {
        int wn = pow(3,(mod-1)/i);
        for(j=0;j<len;j+=i) {
            int w = 1, x, y;
            for(k=0;k<(i>>1);++k,w=(LL)w*wn%mod) {
                x = a[j+k], y = (LL)a[j+k+i/2]*w % mod;
                a[j+k] = (x + y) % mod;
                a[j+k+(i>>1)] = (x + mod - y) % mod;
            }
        }
    }
    
    if(type == -1) {
        for(i=1;i<len>>1;++i) swap(a[i],a[len-i]);
        int tt = pow(len,mod-2);
        for(i=0;i<len;++i) a[i] = (LL) a[i] * tt % mod;     
    }
}

LL a[N], b[N], c[N];

inline void solve(int L,int R) {
    if(L == R) {
        (F[L] += h[L]) %= mod;
        return;
    }
    if(L > R) return;
    int mid = (L+R) >> 1;
    solve(L, mid);
    for(int i=L;i<=mid;++i) a[i-L+1] = (LL)F[i]*ifac[i-1]%mod;
    for(int i=1;i+L<=R;++i) b[i] = (LL)h[i]*ifac[i]%mod;
    int n = 1, len = R - L + 1;
    for(n=1;n<=len<<1;n<<=1);
    NTT(a, n, 1);
    NTT(b, n, 1);
    for(int i=0;i<=n;++i) c[i] = (LL)a[i] * b[i] % mod;
    NTT(c, n,-1);
    for(int i=mid+1;i<=R;++i) (F[i] += mod - (LL)fac[i-1] * c[i-L+1] % mod) %= mod;
    for(int i=0;i<=n;++i) a[i] = b[i] = c[i] = 0;
    solve(mid+1,R);
}

int main() {
    int n;
    init();
    cin >> n;
    solve(1,n);
    cout << F[n] << endl;
}


Title - Artist
0:00

站点地图 网站地图
Copyright © 2015-2017 A simple Blog for wyx
Powered by Typecho自豪的采用Sgreen主题

TOP