BZOJ 3456: 城市規劃 [多項式求逆元 組合數學 | 生成函式 多項式求ln]

Candy?發表於2017-04-19
函式

3456: 城市規劃

題意：n個點組成的無向連通圖個數

以前做過，今天覆習一下

令\(f[n]\)為n個點的無向連通圖個數

n個點的完全圖個數為\(2^{\binom{n}{2}}\)

和Bell數的推導很類似，列舉第一個cc的點的個數
\[ 2^{\binom{n}{2}} = \sum_{i=1}^n \binom{n-1}{i-1} f(i) 2^{\binom{n-i}{2}} \]
整理後
\[ \frac{2^{\binom{n}{2}}}{(n-1)!} = \sum_{i=1}^n \frac{f(i)}{(i-1)!}\frac{2^{\binom{n-i}{2}}}{(n-i)!} \]
這是卷積的形式
\[ C(x) = A(x)B(x) \rightarrow A(x) = C(x)B^{-1}(x) \]
多項式求逆就可做了

注意\(b_0=1\)


其實就是EGP求ln......

簡單的有標號集合計數
\[ G(x) = \sum_{i\ge 0} 2^{\binom{i}{2}}\frac{x^i}{i!} = \sum_{i \ge 0} \frac{F(x)^i}{i!} = e^{F(x)}\\ F(x) = \ln G(x) = \int \frac{G'(x)}{G(x)}dx \]

第一份是多項式求ln的程式碼，注意EGP要除以階乘逆元

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long ll;
const int N = (1<<18) + 5;
inline int read(){
    char c=getchar();int x=0,f=1;
    while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
    while(c>='0'&&c<='9'){x=x*10+c-'0';c=getchar();}
    return x*f;
}

int P = 1004535809, inv2 = (P+1)/2;
inline int Pow(ll a, int b) {
    ll ans = 1;
    for(; b; b>>=1, a=a*a%P)
        if(b&1) ans=ans*a%P;
    return ans;
}

ll inv[N];
namespace ntt {
    int g = 3, rev[N];
    void dft(int *a, int n, int flag) {
        int k = 0; while((1<<k) < n) k++;
        for(int i=0; i<n; i++) {
            rev[i] = (rev[i>>1]>>1) | ((i&1)<<(k-1));
            if(i < rev[i]) swap(a[i], a[rev[i]]);
        }
        for(int l=2; l<=n; l<<=1) {
            int m = l>>1, wn = Pow(g, flag == 1 ? (P-1)/l : P-1-(P-1)/l);
            for(int *p = a; p != a+n; p += l)
                for(int k=0, w=1; k<m; k++, w = (ll)w*wn %P) {
                    int t = (ll) w * p[k+m] %P, r = p[k];
                    p[k+m] = (r - t + P) %P;
                    p[k] = (r + t) %P;
                }
        }
        if(flag == -1) {
            ll inv = Pow(n, P-2);
            for(int i=0; i<n; i++) a[i] = a[i] * inv %P;
        }
    }
    
    void inverse(int *a, int *b, int l) {
        static int t[N];
        if(l == 1) {b[0] = Pow(a[0], P-2); return;}
        inverse(a, b, l>>1);
        int n = l<<1;
        for(int i=0; i<l; i++) t[i] = a[i], t[i+l] = 0;
        dft(t, n, 1); dft(b, n, 1);
        for(int i=0; i<n; i++) b[i] = (ll) b[i] * (2 - (ll) t[i] * b[i] %P + P) %P;
        dft(b, n, -1); for(int i=l; i<n; i++) b[i] = 0;
    }

    void ln(int *a, int *b, int l) {
        static int da[N], ia[N];
        int n = l<<1;
        for(int i=0; i<n; i++) da[i] = ia[i] = 0;
        for(int i=0; i<l-1; i++) da[i] = (ll) (i+1) * a[i+1] %P;
        inverse(a, ia, l);
        dft(da, n, 1); dft(ia, n, 1);
        for(int i=0; i<n; i++) b[i] = (ll) da[i] * ia[i] %P;
        dft(b, n, -1);
        for(int i=l-1; i>0; i--) b[i] = (ll) inv[i] * b[i-1] %P; b[0] = 0;
        for(int i=l; i<n; i++) b[i] = 0;
    }
    
}

ll fac[N], facInv[N];
int n, a[N], b[N];
int main() {
    freopen("in", "r", stdin);
    n = read();
    int len = 1;
    while(len <= n) len <<= 1;
    inv[1] = 1; fac[0] = facInv[0] = 1;
    for(int i=1; i<=n; i++) {
        if(i != 1) inv[i] = (P - P/i) * inv[P%i] %P;
        fac[i] = fac[i-1] * i %P;
        facInv[i] = facInv[i-1] * inv[i] %P;
    }
    a[0] = 1;
    for(int i=1; i<=n; i++) a[i] = Pow(2, (ll) i * (i - 1) /2 %(P-1)) * facInv[i] %P;
    ntt::ln(a, b, len);
    //for(int i=0; i<=n; i++) printf("%d ", b[i]); puts("");
    int ans = b[n] * fac[n] %P;
    printf("%d", ans);
}


#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long ll;
const int N = (1<<18) + 5;
inline int read(){
    char c=getchar();int x=0,f=1;
    while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
    while(c>='0'&&c<='9'){x=x*10+c-'0';c=getchar();}
    return x*f;
}

ll P = 1004535809;
inline ll Pow(ll a, int b) {
    ll ans = 1;
    for(; b; b>>=1, a=a*a%P)
        if(b&1) ans=ans*a%P;
    return ans;
}
inline void mod(int &x) {if(x>=P) x-=P; else if(x<0) x+=P;}
namespace fnt {
    int n, g=3, rev[N];
    void dft(int *a, int n, int flag=1) {
        for(int i=0; i<n; i++) if(i < rev[i]) swap(a[i], a[rev[i]]);

        for(int l=2; l<=n; l<<=1) {
            int m = l>>1;
            ll wn = Pow(g, flag==1 ? (P-1)/l : P-1-(P-1)/l);
            for(int *p=a; p!=a+n; p+=l) {
                ll w = 1;
                for(int k=0; k<m; k++) {
                    ll t = p[k+m] * w %P;
                    mod(p[k+m] = p[k] - t);
                    mod(p[k] = p[k] + t);
                    w = w * wn %P;
                }
            }
        }
        if(flag == -1) {
            ll inv = Pow(n, P-2);
            for(int i=0; i<n; i++) a[i] = a[i] * inv %P;
        }
    }
    int t[N];
    void inverse(int *a, int *b, int l) { // mod x^l
        if(l == 1) {b[0] = Pow(a[0], P-2); return;}
        inverse(a, b, (l+1)>>1);
        int n = 1, k = 0; while(n < l<<1) n<<=1, k++;
        for(int i=0; i<n; i++) rev[i] = (rev[i>>1]>>1) | ((i&1)<<(k-1));
        for(int i=0; i<l; i++) t[i] = a[i]; for(int i=l; i<n; i++) t[i] = 0;
        dft(t, n, 1); dft(b, n, 1);
        for(int i=0; i<n; i++) b[i] = (ll) b[i] * (2 - (ll) t[i] * b[i] %P + P) %P;
        dft(b, n, -1);
        for(int i=l; i<n; i++) b[i] = 0;
    }

    void mul(int *a, int *b, int l) {
        int n = 1; while(n < l<<1) n<<=1;
        dft(a, n, 1); dft(b, n, 1);
        for(int i=0; i<n; i++) a[i] = (ll) a[i] * b[i] %P;
        dft(a, n, -1);
    }
}

int n, a[N], b[N], bi[N];
ll inv[N], fac[N], facInv[N];
int main() {
    freopen("in", "r", stdin);
    n=read();
    inv[1] = 1; fac[0] = facInv[0] = 1;
    for(int i=1; i<=n; i++) {
        if(i != 1) inv[i] = (P - P/i) * inv[P%i] %P;
        fac[i] = fac[i-1] * i %P;
        facInv[i] = facInv[i-1] * inv[i] %P;
    }
    for(int i=1; i<=n; i++) {
        ll mi = Pow(2, (ll) i * (i-1) / 2 %(P-1));
        a[i] = mi * facInv[i-1] %P;
        b[i] = mi * facInv[i] %P;
    }
    b[0] = 1;

    fnt::inverse(b, bi, n+1); 
    fnt::mul(a, bi, n+1);
    ll ans = a[n] * fac[n-1] %P;
    printf("%lld", ans);
}

相關文章