• 【洛谷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();
    }
    
  • 相关阅读:
    类加载器加载class文件
    多线程用到的概念知识点
    3年工作经验你的程序员应该具备的技能
    面试题
    正则表达式(一)
    Servlet(五)----ServletContext对象
    Servlet(四)----HTTP、Response、servletContent
    JDBC(三)----Spring JDBC(JDBCTemplate)
    JDBC(四)----数据库连接池
    JDBC(三)----JDBC控制事务
  • 原文地址:https://www.cnblogs.com/zyt1253679098/p/10387647.html
Copyright © 2020-2023  润新知