【洛谷4841】城市规划(多项式)

题目:

洛谷4841

分析:

一句话题意:求(n)个点的带标号简单(无重边无自环)无向连通图数目。(以下提到的所有图都是带标号简单无向图)

正难则反。设(f[n])表示(n)个点的连通图数目,(g[n])表示(n)个点的图数目。那么从(g[n])中减去不连通图的数目就是(f[n])

对于两个点之间的边只有选和不选两种选择,所以显然(g[n]=2^{C_n^2}=2^frac{n(n-1)}{2})。现在问题变成了不连通图的数目。

任意钦定一个点,比如说(1)号点。枚举(1)号点所在连通块的大小(i),则相当于先从(n-1)个点中选出(i-1)个点和(i)组成连通块,再将剩下的(n-i)个点组成一张图,即:

[f[n]=g[n]-sum_{i=1}^{n-1} C_{n-1}^{i-1}f[i]g[n-i] ]

由于钦定了一个具体的点((1)号点)所在的连通块来讨论,所以不会重复。

(g[n]=2^{C_n^2}),很好算。下面来解决后面的求和怎么算。首先把组合数按定义展开:

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

可以看出,可以把((n-1)!)提到前面,然后有些项和(i)有关,有些项和(n-i)有关。把它们分别搞到一块,就可以……卷积?

尝试搞到一块……

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

右边看起来非常像卷积,于是设两个多项式:

[F(x)=sum_{i=0} frac{f[i]x^i}{(i-1)!} ]

[G(x)=sum_{i=0} frac{g[n-i]x^i}{i!} ]

至此看起来似乎可以分治FFT?博主不会,告辞

想想原来那个式子;

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

把右边第二项移到左边,再两边除以((n-1)!)

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

发现(frac{f[n]}{(n-1)!})和求和里面关于(f[i])那一项长得很像,并且当(i=n)(frac{g[n-i]}{(n-i)!}=g[0]=1),所以可以直接把它扔到求和里面,即:

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

然后设多项式(H(x)=sum_i frac{2^{C_i^2}x^i}{(i-1)!}),则(FG=H)(G)(H)都是已知的,则(F=G^{-1}H),求出(F)后第(n)项系数乘上((n-1)!)就是答案。

对于最终的式子还有另一种理解:(1)所在的连通块大小从(1)(n)的所有情况之和就是(n)个点的图的方案数,即:

[g[n]=sum_{i=1}^n C_{n-1}^{i-1}f[i]g[n-i] ]

同样把组合数展开,就能得出上面的式子。

代码:

#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cctype>
using namespace std;

namespace zyt
{
    template<typename T>
    inline bool read(T &x)
    {
        char c;
        bool f = false;
        x = 0;
        do
            c = getchar();
        while (c != EOF && c != '-' && !isdigit(c));
        if (c == EOF)
            return false;
        if (c == '-')
            f = true, c = getchar();
        do
            x = x * 10 + c - '0', c = getchar();
        while (isdigit(c));
        if (f)
            x = -x;
        return true;
    }
    template<typename T>
    inline void write(T x)
    {
        static char buf[20];
        char *pos = buf;
        if (x < 0)
            putchar('-'), x = -x;
        do
            *pos++ = x % 10 + '0';
        while (x /= 10);
        while (pos > buf)
            putchar(*--pos);
    }
    typedef long long ll;
    const int N = 1.3e5 + 10, LEN = N << 2, p = 1004535809, g = 3;
    inline int power(int a, int b)
    {
        int ans = 1;
        while (b)
        {
            if (b & 1)
                ans = (ll)ans * a % p;
            a = (ll)a * a % p;
            b >>= 1;
        }
        return ans;
    }
    inline int inv(const int a)
    {
        return power(a, p - 2);
    }
    namespace Polynomial
    {
        int omega[LEN], winv[LEN], rev[LEN];
        void init(const int n, const int lg2)
        {
            int w = power(g, (p - 1) / n), wi = inv(w);
            omega[0] = winv[0] = 1;
            for (int i = 1; i < n; i++)
            {
                omega[i] = (ll)omega[i - 1] * w % p;
                winv[i] = (ll)winv[i - 1] * wi % p;
            }
            for (int i = 0; i < n; i++)
                rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (lg2 - 1)));
        }
        void ntt(int *a, const int *w, const int n)
        {
            for (int i = 0; i < n; i++)
                if (i < rev[i])
                    swap(a[i], a[rev[i]]);
            for (int l = 1; l < n; l <<= 1)
                for (int i = 0; i < n; i += (l << 1))
                    for (int k = 0; k < l; k++)
                    {
                        int tmp = (a[i + k] - (ll)w[n / (l << 1) * k] * a[i + l + k] % p + p) % p;
                        a[i + k] = (a[i + k] + (ll)w[n / (l << 1) * k] * a[i + l + k] % p) % p;
                        a[i + l + k] = tmp;
                    }
        }
        void _inv(const int *a, int *b, const int n)
        {
            if (n == 1)
                b[0] = ::zyt::inv(a[0]);
            else
            {
                static int tmp[LEN];
                _inv(a, b, (n + 1) >> 1);
                int m = 1, lg2 = 0;
                while (m < (n << 1))
                    m <<= 1, ++lg2;
                init(m, lg2);
                memcpy(tmp, a, sizeof(int[n]));
                memset(tmp + n, 0, sizeof(int[m - n]));
                memset(b + ((n + 1) >> 1), 0, sizeof(int[m - ((n + 1) >> 1)]));
                ntt(tmp, omega, m), ntt(b, omega, m);
                for (int i = 0; i < m; i++)
                    b[i] = (2LL * b[i] % p - (ll)tmp[i] * b[i] % p * b[i] % p + p) % p;
                ntt(b, winv, m);
                int invm = ::zyt::inv(m);
                for (int i = 0; i < n; i++)
                    b[i] = (ll)b[i] * invm % p;
            }
        }
        void inv(const int *a, int *b, const int n)
        {
            static int tmp[N];
            memcpy(tmp, a, sizeof(int[n]));
            _inv(tmp, b, n);
        }
        void mul(const int *a, const int *b, int *c, const int n)
        {
            static int x[LEN], y[LEN];
            memcpy(x, a, sizeof(int[n]));
            memcpy(y, b, sizeof(int[n]));
            int m = 1, lg2 = 0;
            while (m < (n << 1))
                m <<= 1, ++lg2;
            memset(x + n, 0, sizeof(int[m - n]));
            memset(y + n, 0, sizeof(int[m - n]));
            init(m, lg2);
            ntt(x, omega, m), ntt(y, omega, m);
            for (int i = 0; i < m; i++)
                x[i] = (ll)x[i] * y[i] % p;
            ntt(x, winv, m);
            int invm = ::zyt::inv(m);
            for (int i = 0; i < n; i++)
                c[i] = (ll)x[i] * invm % p;
        }
    }
    int fac[N], finv[N], n;
    void init()
    {
        fac[0] = 1;
        for (int i = 1; i <= n; i++)
            fac[i] = (ll)fac[i - 1] * i % p;
        finv[n] = inv(fac[n]);
        for (int i = n; i > 0; i--)
            finv[i - 1] = (ll)finv[i] * i % p;
    }
    int F[LEN], G[LEN], H[LEN];
    int work()
    {
        using Polynomial::inv;
        using Polynomial::mul;
        read(n);
        init();
        for (int i = 0; i <= n; i++)
        {
            G[i] = (ll)power(2, (ll)i * (i - 1) / 2 % (p - 1)) * finv[i] % p;
            if (i)
                H[i] = (ll)power(2, (ll)i * (i - 1) / 2 % (p - 1)) * finv[i - 1] % p;
        }
        inv(G, G, n + 1);
        mul(G, H, F, n + 1);
        write((ll)F[n] * fac[n - 1] % p);
        return 0;
    }
}
int main()
{
    return zyt::work();
}
原文地址:https://www.cnblogs.com/zyt1253679098/p/10387647.html