luogu P5206: LOJ 2983: [WC2019] 数树

一道技巧性非常强的计数题,历年WC出得最好的题目之一。

题目传送门:洛谷P5206

题意简述:

给定 \(n, y\)。

一张图有 \(|V| = n\) 个点。对于两棵树 \(T_1=G(V, E_1)\) 和 \(T_2=G(V, E_2)\),定义这两棵树的权值 \(F(E_1, E_2)\) 为 \(y\) 的 \(G'=(V,E_1\cap E_2)\) 的联通块个数次方。

即 \(F(E_1, E_2) = y^{n - |E_1\cap E_2|}\)(因为两棵树的边集的交必然是一个森林,而森林的联通块个数等于 \(|V| - |E|\))。

有 \(3\) 类问题:

  1. \(\mathrm{op}=0\),给定 \(E_1, E_2\),计算 \(F(E_1, E_2)\)。
  2. \(\mathrm{op}=1\),给定 \(E_1\),计算 \(\sum_{E_2} F(E_1, E_2)\)。
  3. \(\mathrm{op}=2\),计算 \(\sum_{E_1} \sum_{E_2} F(E_1, E_2)\)。

其中 \(\sum_{E}\) 的意义为枚举所有 \(n^{n - 2}\) 种不同形态的有标号无根树。

其中 \(n\le 10^5\)。答案对 \(998,244,353\) 取模。

题解:

在这里特别感谢 @rqy 的题解 [WC2019] 数树
我对此题的理解以及代码实现的细节很大程度上借鉴了他的题解。

对于 \(y = 1\):

\(y = 1\) 时在后续讨论中有一些无意义的地方,这里先判掉。
\(y = 1\) 时 \(F(E_1, E_2) = 1\)。

\(\mathrm{op}=0\) 时,只有 \(1\) 种方案,则答案为 \(1\)。
\(\mathrm{op}=1\) 时,有 \(n^{n-2}\) 种方案,则答案为 \(n^{n-2}\)。
\(\mathrm{op}=2\) 时,有 \(n^{2(n-2)}\) 种方案,则答案为 \(n^{2(n-2)}\)。

namespace Solver0 {
    inline LL solve() {
        if (op == 0) return 1;
        if (op == 1) return qPow(N, N - 2);
        if (op == 2) return qPow(N, 2 * (N - 2));
        return 0;
    }
}

对于 \(\mathrm{op} = 0\):

直接统计出有多少条边重合,随便用一个 std::set<std::pair<int,int> > 统计一下。

namespace Solver1 {
    typedef std::pair<int, int> pii;
    std::set<pii> S;
    
    inline LL solve() {
        for (int i = 1, x, y; i < N; ++i) {
            scanf("%d%d", &x, &y);
            if (x > y) std::swap(x, y);
            S.insert(std::make_pair(x, y));
        }
        int cnt = N;
        for (int i = 1, x, y; i < N; ++i) {
            scanf("%d%d", &x, &y);
            if (x > y) std::swap(x, y);
            cnt -= S.count(std::make_pair(x, y));
        }
        return qPow(Y, cnt);
    }
}

对于 \(y\ne 1\) 且 \(\mathrm{op} = 1\):

给定 \(E_1\),对于每一种 \(E_2\),统计 \(y^{n-|E_1\cap E_2|}\)。

考虑枚举 \(S=E_1\cap E_2\),则答案为:

\[\begin{aligned}&\sum_{E_2}F(E_1,E_2)\\=&\sum_{E_2}y^{n-|E_1\cap E_2|}\\=&\sum_{S}y^{n-|S|}\sum_{S=E_1\cap E_2}\end{aligned}\]

后面这个不是很优美,我们容斥一下,希望变成 \(S\subseteq E_1\cap E_2\),即 \(S\subseteq E_1\) 且 \(S\subseteq E_2\) 的形式:

考虑一个容斥原理的式子:\(f(S)=\sum_{T\subseteq S}\sum_{P\subseteq T}(-1)^{|T|-|P|}f(P)\)。

令 \(S=E_1\cap E_2\),\(f(S)=y^{n-|S|}\),带进去一下,得到:

\[\begin{aligned}&\sum_{E_2}y^{n-|E_1\cap E_2|}\\=&\sum_{E_2}\sum_{T\subseteq E_1\cap E_2}\sum_{P\subseteq T}(-1)^{|T|-|P|}y^{n-|P|}\\=&\sum_{T\subseteq E_1}g(T)\sum_{P\subseteq T}(-1)^{|T|-|P|}y^{n-|P|}\\=&\sum_{T\subseteq E_1}g(T)y^{n-|T|}\sum_{P\subseteq T}(-y)^{|T|-|P|}\\=&\sum_{T\subseteq E_1}g(T)y^{n-|T|}\sum_{k=0}^{|T|}\binom{|T|}{k}1^{k}(-y)^{|T|-k}\\=&\sum_{T\subseteq E_1}g(T)y^{n-|T|}(1-y)^{|T|}\end{aligned}\]

其中 \(g(T)\) 表示包含边集 \(T\) 的树个数。

这里有一个定理,可以用 \(\mathrm{Pr\ddot{u}fer}\) 序列证明,也可以用 \(\mathrm{Matrix-Tree}\) 定理证明:

  • 对于 \(n\) 个点的森林,假设有 \(k\) 个连通分量,每个连通分量的大小分别是 \(a_i\),则包含这个森林的树个数为 \(n^{k-2}\prod_{i=1}^{k}a_i\)。

再带入可得:

\[\begin{aligned} &\sum_{T\subseteq E_1}g(T)y^{n-|T|}(1-y)^{|T|}\\=&\sum_{T\subseteq E_1}y^{k}(1-y)^{n-k}n^{k-2}\prod_{i=1}^{k}a_i\\=&\frac{(1-y)^n}{n^2}\sum_{T\subseteq E_1}\prod_{i=1}^{k}\frac{ny}{1-y}a_i\end{aligned}\]

即一个大小为 \(a\) 的连通分量为答案贡献一个乘积 \(Ka\),其中 \(K = \frac{ny}{1-y}\)。

据此我们可以写出一个 \(\mathrm{DP}\) 式子,用 \(\mathrm{f}[u][i]\) 表示 \(u\) 的子树中,\(u\) 所在的连通分量大小为 \(i\) 的答案。

\(i\) 所在的连通分量还未完整,所以不计入答案,而是存入状态,转类似于树上背包,转移比较显然,这种做法是 \(\mathrm{O}(n^2)\) 的。

我们可以考虑每个贡献的组合意义,大小为 \(a_i\) 的连通分量贡献相当于在这个联通分量中选取一个点产生 \(K\) 的乘积贡献。

据此我们可以优化状态表示:\(\mathrm{f}[u][0/1]\) 表示 \(u\) 的子树中,当前连通分量是否已经做出贡献,这里为了方便转移,当前连通分量如果已经做出贡献就计入答案。

具体转移可以参考下面的代码。

namespace Solver2 {
    int h[MN], nxt[MN * 2], to[MN * 2], tot;
    inline void ins(int x, int y) {
        nxt[++tot] = h[x], to[tot] = y, h[x] = tot;
    }
    
    LL K, f[MN], g[MN];
    void DFS(int u, int fz) {
        f[u] = 1, g[u] = K;
        for (int i = h[u]; i; i = nxt[i]) if (to[i] != fz) {
            DFS(to[i], u);
            g[u] = (f[u] * g[to[i]] + g[u] * f[to[i]] + g[u] * g[to[i]]) % Mod;
            f[u] = (f[u] * f[to[i]] + f[u] * g[to[i]]) % Mod;
        }
    }
    
    inline LL solve() {
        for (int i = 1, x, y; i < N; ++i) {
            scanf("%d%d", &x, &y);
            ins(x, y), ins(y, x);
        }
        K = (LL)N * Y % Mod * qPow(1 - Y, Mod - 2) % Mod;
        LL coef = qPow(1 - Y, N) * qPow(N, Mod - 3) % Mod;
        DFS(1, 0);
        return g[1] * coef % Mod;
    }
}

对于 \(y\ne 1\) 且 \(\mathrm{op} = 2\):

类似 \(\mathrm{op} = 1\),我们可以写出下面的式子:

\[\begin{aligned}&\sum_{E_1}\sum_{E_2}F(E_1,E_2)\\=&\sum_{E_1}\sum_{E_2}y^{|E_1\cap E_2|}\\=&\sum_{E_1}\sum_{E_2}\sum_{T\subseteq E_1\cap E_2}\sum_{P\subseteq T}(-1)^{|T|-|P|}y^{n-|P|}\\=&\sum_{T}g^2(T)\sum_{P\subseteq T}(-1)^{|T|-|P|}y^{n-|P|}\\=&\sum_{T}g^2(T)y^{n-|T|}\sum_{P\subseteq T}(-1)^{|T|-|P|}y^{|T|-|P|}\\=&\sum_{T}g^2(T)y^{n-|T|}\sum_{P\subseteq T}(-y)^{|T|-|P|}\\=&\sum_{T}g^2(T)y^{n-|T|}\sum_{k=0}^{|T|}\binom{|T|}{k}1^{k}(-y)^{|T|-k}\\=&\sum_{T}g^2(T)y^{n-|T|}(1-y)^{|T|}\\=&\sum_{T}y^{k}(1-y)^{n-k}n^{2(k-2)}\prod_{i=1}^{k}a_i^2\\=&\frac{(1-y)^n}{n^4}\sum_{T}\prod_{i=1}^{k}\frac{n^2y}{1-y}a_i^2\end{aligned}\]

此时我们换个角度考虑,与其考虑边集 \(T\),不如考虑每个连通分量。

这里的每个连通分量是无序的,所以相当于 \(n\) 个有标号小球扔进 \(k\) 个无标号盒子,不允许空盒子。
而一个装有 \(a\) 个球的盒子,每个生成树都会贡献 \(\frac{n^2y}{1-y}a^2\) 作为乘积。
而 \(a\) 个点的生成树有 \(a^{a-2}\) 个,所以总共贡献 \(\frac{n^2y}{1-y}a^a\)。

对生成函数比较敏感的同学能够看出,总方案的指数型生成函数就是单个盒子的指数型生成函数的 \(\mathrm{Exp}\),因为对应了集合与元素的关系。
即 \(\mathrm{f}=\{\frac{n^2y}{1-y}i^i\}_{i=1}^{\infty}\) 对应的指数型生成函数 \(\mathrm{F}(x)=\sum_{i=1}^{\infty}\frac{n^2y}{(1-y)i!}i^{i}x^{i}\)。
令 \(\mathrm{G}=e^{\mathrm{F}}\),则指数型生成函数 \(\mathrm{G}\) 对应的序列的第 \(n\) 项便是所求答案,即答案为 \(n!\cdot [x^n]\mathrm{G}\)。

如果你对生成函数比较不敏感,也可以推式子得到:
对于有标号球有标号盒子记录贡献,可以使用指数型生成函数的幂次表示:
\(n\) 个有标号球,\(k\) 个有标号盒子就是 \([x^n]\mathrm{F}^k\),但是这里盒子是无标号的,所以要除以 \(k!\)。
那么最终的答案是 \(\frac{(1-y)^n\cdot n!}{n^4}[x^n]\sum_{k=1}^{n}\frac{1}{k!}\left(\sum_{i=1}^{\infty}\frac{n^2y}{(1-y)\cdot i!}i^ix^i\right)^k\)。
根据 \(e^x\) 的泰勒展开公式,可以发现形式是完全相同的。

那么直接套多项式 \(\mathrm{Exp}\) 模板即可。

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <set>

typedef long long LL;
const int MN = 100005;
const LL Mod = 998244353;

inline LL qPow(LL b, int e) {
    LL a = 1;
    for (; e; b = b * b % Mod, e >>= 1)
        if (e & 1) a = a * b % Mod;
    return a;
}

int N, op;
LL Y;

namespace Solver0 {
    inline LL solve() {
        if (op == 0) return 1;
        if (op == 1) return qPow(N, N - 2);
        if (op == 2) return qPow(N, 2 * (N - 2));
        return 0;
    }
}

namespace Solver1 {
    typedef std::pair<int, int> pii;
    std::set<pii> S;
    
    inline LL solve() {
        for (int i = 1, x, y; i < N; ++i) {
            scanf("%d%d", &x, &y);
            if (x > y) std::swap(x, y);
            S.insert(std::make_pair(x, y));
        }
        int cnt = N;
        for (int i = 1, x, y; i < N; ++i) {
            scanf("%d%d", &x, &y);
            if (x > y) std::swap(x, y);
            cnt -= S.count(std::make_pair(x, y));
        }
        return qPow(Y, cnt);
    }
}

namespace Solver2 {
    int h[MN], nxt[MN * 2], to[MN * 2], tot;
    inline void ins(int x, int y) {
        nxt[++tot] = h[x], to[tot] = y, h[x] = tot;
    }
    
    LL K, f[MN], g[MN];
    void DFS(int u, int fz) {
        f[u] = 1, g[u] = K;
        for (int i = h[u]; i; i = nxt[i]) if (to[i] != fz) {
            DFS(to[i], u);
            g[u] = (f[u] * g[to[i]] + g[u] * f[to[i]] + g[u] * g[to[i]]) % Mod;
            f[u] = (f[u] * f[to[i]] + f[u] * g[to[i]]) % Mod;
        }
    }
    
    inline LL solve() {
        for (int i = 1, x, y; i < N; ++i) {
            scanf("%d%d", &x, &y);
            ins(x, y), ins(y, x);
        }
        K = (LL)N * Y % Mod * qPow(1 - Y, Mod - 2) % Mod;
        LL coef = qPow(1 - Y, N) * qPow(N, Mod - 3) % Mod;
        DFS(1, 0);
        return g[1] * coef % Mod;
    }
}

namespace Solver3 {
    const int MS = 1 << 19;
    const int G = 3, iG = 332748118;
    int Sz = 0, R[MS]; LL InvSz;
    LL Inv[MS], Fac[MS], iFac[MS];
    
    inline void Init() {
        Inv[1] = 1;
        for (int i = 2; i < MS; ++i) Inv[i] = -(Mod / i) * Inv[Mod % i] % Mod;
        Fac[0] = iFac[0] = 1;
        for (int i = 1; i < MS; ++i) {
            Fac[i] = Fac[i - 1] * i % Mod;
            iFac[i] = iFac[i - 1] * Inv[i] % Mod;
        }
    }
    
    inline void InitFNTT(int n) {
        int Bt = 0;
        for (; 1 << Bt < n; ++Bt) ;
        if ((1 << Bt) == Sz) return ;
        Sz = 1 << Bt, InvSz = -(Mod - 1) / Sz;
        for (int i = 1; i < Sz; ++i) R[i] = R[i >> 1] >> 1 | (i & 1) << (Bt - 1);
    }
    
    inline void FNTT(LL *A, int Type) {
        for (int i = 0; i < Sz; ++i) if (R[i] < i) std::swap(A[R[i]], A[i]);
        for (int j = 1, j2 = 2; j < Sz; j <<= 1, j2 <<= 1) {
            LL gn = qPow(~Type ? G : iG, (Mod - 1) / j2), g;
            for (int i = 0, k; i < Sz; i += j2) {
                for (k = 0, g = 1; k < j; ++k, g = g * gn % Mod) {
                    LL X = A[i + k], Y = g * A[i + j + k] % Mod;
                    A[i + k] = (X + Y) % Mod, A[i + j + k] = (X - Y) % Mod;
                }
            }
        }
        if (Type == -1) for (int i = 0; i < Sz; ++i) A[i] = A[i] * InvSz % Mod;
    }
    
    inline void PolyInv(LL *A, int N, LL *B) {
        B[0] = qPow(A[0], Mod - 2);
        static LL tA[MS], tB[MS];
        for (int L = 1; L < N; L <<= 1) {
            int L2 = L << 1, L4 = L << 2;
            InitFNTT(L4);
            for (int i = 0; i < L2; ++i) tA[i] = A[i];
            for (int i = L2; i < Sz; ++i) tA[i] = 0;
            for (int i = 0; i < L; ++i) tB[i] = B[i];
            for (int i = L; i < Sz; ++i) tB[i] = 0;
            FNTT(tA, 1), FNTT(tB, 1);
            for (int i = 0; i < Sz; ++i) tB[i] = (2 - tA[i] * tB[i]) % Mod * tB[i] % Mod;
            FNTT(tB, -1);
            for (int i = 0; i < L2; ++i) B[i] = tB[i];
        }
    }
    
    inline void PolyLn(LL *A, int N, LL *B) {
        static LL tA[MS], tB[MS];
        PolyInv(A, N, tB);
        InitFNTT(N * 2);
        for (int i = 1; i < N; ++i) tA[i - 1] = A[i] * i % Mod;
        for (int i = N - 1; i < Sz; ++i) tA[i] = 0;
        for (int i = N; i < Sz; ++i) tB[i] = 0;
        FNTT(tA, 1), FNTT(tB, 1);
        for (int i = 0; i < Sz; ++i) tA[i] = tA[i] * tB[i] % Mod;
        FNTT(tA, -1);
        B[0] = 0;
        for (int i = 1; i < N; ++i) B[i] = tA[i - 1] * Inv[i] % Mod;
    }
    
    inline void PolyExp(LL *A, int N, LL *B) {
        B[0] = 1;
        static LL tA[MS], tB[MS];
        for (int L = 1; L < N; L <<= 1) {
            int L2 = L << 1, L4 = L << 2;
            PolyLn(B, L2, tA);
            InitFNTT(L4);
            for (int i = 0; i < L2; ++i) tA[i] = (!i + A[i] - tA[i]) % Mod;
            for (int i = L2; i < Sz; ++i) tA[i] = 0;
            for (int i = 0; i < L2; ++i) tB[i] = B[i];
            for (int i = L2; i < Sz; ++i) tB[i] = 0;
            FNTT(tA, 1), FNTT(tB, 1);
            for (int i = 0; i < Sz; ++i) tA[i] = tA[i] * tB[i] % Mod;
            FNTT(tA, -1);
            for (int i = 0; i < L2; ++i) B[i] = tA[i];
        }
    }
    
    inline LL solve() {
        Init();
        LL K = (LL)N * N % Mod * Y % Mod * qPow(1 - Y, Mod - 2) % Mod;
        LL coef = qPow(1 - Y, N) * qPow(N, Mod - 5) % Mod;
        static LL A[MS], B[MS];
        A[0] = 0;
        for (int i = 1; i <= N; ++i) A[i] = K * qPow(i, i) % Mod * iFac[i] % Mod;
        PolyExp(A, N + 1, B);
        return coef * (B[N] * Fac[N] % Mod) % Mod;
    }
}

int main() {
    scanf("%d%lld%d", &N, &Y, &op);
    if (Y == 1) printf("%lld\n", Solver0::solve());
    else if (op == 0) printf("%lld\n", (Solver1::solve() + Mod) % Mod);
    else if (op == 1) printf("%lld\n", (Solver2::solve() + Mod) % Mod);
    else if (op == 2) printf("%lld\n", (Solver3::solve() + Mod) % Mod);
    return 0;
}
上一篇:AI 神经网络学习


下一篇:信号与系统05 拉普拉斯变换