苹果树算法题解析

🍎 问题描述

融合了概率论期望计算的奇怪算法题,回到了高中数学的感觉。 ## 算法题目

小 C在自己家的花园里种了一棵苹果树,树上每个结点都有恰好两个分支。经过细心的观察,小 C发现每一天这棵树都会生长出一个新的结点。

生长规则: - 第一天: 花园里什么都没有,我们种下第一颗种子,它长成了树的根。这个根结点很特殊,它有两个空位,可以用来长新的树枝。我们称这些空位为”插槽”。 - 第二天及以后: 每一天,果树会随机选择一个当前树中没有长出过结点的分支,然后在这个分支上长出一个新结点,新结点与分支所属的结点之间连接上一条边。

不便度定义: 树上两两结点之间的距离之和,两个结点之间的距离定义为从一个点走到另一个点的路径经过的边数。

求解目标: 计算 N 天后不便度的期望 E,并求 E × N! 对 P 取模的结果。

数据范围: \(1 \leq N \leq 2000\)\(1 \leq P \leq 10^9+7\)

输入样例1:

1
3 610745795
输出样例1:
1
24

输入样例2:

1
305 1000000007
输出样例2:
1
865018107

理解题目:从具体例子开始

让我们通过一个具体的例子来理解这个过程:

第一天 (N=1): 花园里什么都没有,我们种下第一颗种子,它长成了树的根。这个根结点很特殊,它有两个空位,可以用来长新的树枝。我们称这些空位为”插槽”。

第二天 (N=2): 现在树上总共有2个空”插槽”。大自然会完全随机地从这2个插槽中选1个,让它长出一个新的苹果结点。现在树上有2个结点了。新长出的这个结点,自己又带来了2个新的空插槽。而被占用的那个老插槽消失了。

第三天 (N=3): 我们来算一下现在有几个空插槽。 - 根结点有2个,被用掉了1个,剩下1个。 - 第二天的新结点,自己带来了2个。 - 总共就是 1 + 2 = 3 个空插槽。 - 大自然再次从这3个空插槽里完全随机地选1个,长出第三个结点。

规律发现: 你会发现一个规律,当树上有 k 个结点时,总会有 k+1 个空插槽可供选择。这个过程会一直持续到第 N 天,最终形成一棵有 N 个结点的树。

“不便度”是什么?不便度就是树上任意两个苹果(结点)之间的距离之和。比如苹果A到苹果B要经过3条树枝,那它们之间的距离就是3。我们要把所有这种”A到B”的距离都加起来。

目标是什么? 因为树的生长是随机的,所以每次长出来的树的形状可能都不一样,从而”不便度”也不同。题目要求我们计算这个”不便度”在所有可能情况下的平均值,也就是期望 (E)。最后,因为结果可能是分数,题目让我们计算 E × N! 对P取模后的结果。

核心思路转换:从”点对距离”到”边做贡献”

直接计算所有点对距离的期望非常复杂。这里有一个非常经典且重要的转化思想:考虑每一条边的贡献。

现在,我们换个角度。我们不看点,而是看每一条边,问一个问题:“这条边,在我们的总距离计算中,被用到了多少次?”把所有边的贡献加起来,就是总的不便度。

想象一下,把树上的一条边砍断,这棵树会立刻分成两半。假设左半部分有 a 个结点,右半部分有 b 个结点(显然 a+b=N)。

  • 任何一个左边的点要走到右边的点,都必须经过刚才被我们砍断的这条边。
  • 左边有 a 个点,右边有 b 个点,所以总共有 a * b 对这样的”跨边点对”。
  • 因此,这条边对总的”不便度”的贡献正好是 a * b。

总不便度 = 所有边的贡献之和 = Σ (每条边连接的两侧结点数的乘积)

这个结论至关重要,它把一个复杂的问题简化为了对每一条边进行分析。

将复杂的全局距离计算转化为局部边贡献的累加,为后续的期望分析奠定了基础。

聚焦单条边:分析子树大小的期望

我们的树是按天生长的。我们考虑在第 t 天(t 从2到N)长出来的那条边。这条边连接了一个父节点和一个在第 t 天新生成的子节点。

这条边将整棵树(N个结点)分成了两部分: - 以这个第 t 天的新结点为根的子树。 - 树的其余部分。

我们令 \(X_t\) 表示到第N天结束时,这个子树最终的大小(包含多少个结点)。那么,树的其余部分大小就是 \(N - X_t\)。根据上一步的结论,这条边在第N天对总不便度的贡献就是 \(X_t * (N - X_t)\)

因为树的生长是随机的,\(X_t\) 的值也是不确定的。所以我们需要计算它的期望贡献:\(E[X_t * (N - X_t)]\)

根据期望的性质,\(E[A * B]\) 不一定等于 \(E[A] * E[B]\),但 \(E[c*A - A^2]\) (c是常数) 等于 \(c*E[A] - E[A^2]\)。所以,我们有: \[ E[X_t * (N - X_t)] = E[N * X_t - X_t^2] = N * E[X_t] - E[X_t^2] \] 现在,问题变成了:对于第 t 天加入的结点,其最终子树大小的期望 \(E[X_t]\) 和平方的期望 \(E[X_t^2]\) 是多少?

求解期望值:最关键的数学推导

这一步是整个解法的核心,我们来理解这些公式是怎么来的。

我们不直接研究子树大小 \(X_t\),而是研究一个更方便的量:子树的空插槽数量。一个包含 \(k\) 个结点的子树,总是有 \(k+1\) 个空插槽。所以,\(X_t = (最终子树的空插槽数) - 1\)

初始状态: 在第 \(t\) 天,当结点 \(t\) 刚被创建时,它自己是一个大小为1的子树,拥有2个空插槽。

演化过程: 从第 \(t\) 天到第 \(N\) 天,每天都会有一个新结点加入。假设在第 \(m\) 天(\(m > t\)),树上总共有 \(m\) 个结点和 \(m+1\) 个空插槽。 - 我们的子树当前有 \(S_{m-1}\) 个空插槽。 - 新结点长到我们子树里的概率是 \(\frac{S_{m-1}}{m}\) (子树插槽数 / 总插槽数)。 - 如果长进来了,子树的空插槽数会加1(失去1个旧的,得到2个新的)。 - 如果没长进来,子树的空插槽数不变。

通过这个概率过程,使用数学方法(条件期望和递推),可以精确地推导出:

子树大小的期望 (E[X_t]): \[E[X_t] = \left(\frac{2 * (N+1)}{(t+1)}\right) - 1\]

子树大小平方的期望 (E[X_t^2]): \[E[X_t^2] = 1 + \frac{6 * (N+1) * (N-t)}{(t+1) * (t+2)}\]

公式推导详解

以下是这两个关键公式的详细推导过程:

预备知识与核心工具

在开始之前,我们需要明确我们的核心工具:条件期望定律(或称全期望定律)。它的思想很简单:想求一个变量的期望,可以先”固定”某个相关的条件,求出在该条件下它的期望,然后再对这个条件本身求期望。

数学上写为:\(E[Y] = E[E[Y | Z]]\)。 在我们这里,就是 \(E[S_m] = E[E[S_m | S_{m-1}]]\)。意思是,想知道 \(m\) 时刻的期望,我们可以先看看在 \(m-1\) 时刻 \(S_{m-1}\) 取某个具体值 \(s\) 的情况下 \(S_m\) 的期望是多少,然后再对 \(S_{m-1}\) 的所有可能取值 \(s\) 求平均。

变量定义: - \(S_m\): 在整棵树有 \(m\) 个结点时,我们所关注的、以结点 \(t\) 为根的子树拥有的空插槽数量。 - \(E_m\): \(S_m\) 的期望,即 \(E[S_m]\)。 - \(A_m\): \(S_m^2\) 的期望,即 \(E[S_m^2]\)


第一部分:推导子树大小的期望 \(E[X_t]\)

我们的目标是求出 \(E[X_t]\)。我们知道 \(X_t = S_N - 1\),所以关键是求 \(E[S_N]\)

1. 建立 \(S_m\)\(S_{m-1}\) 的递推关系

考虑从第 \(m-1\) 天到第 \(m\) 天,树上增加了一个新结点。 - 在第 \(m-1\) 天结束时,树上共有 \(m-1\) 个结点,总空插槽数为 \((m-1)+1 = m\) 个。 - 我们关注的子树有 \(S_{m-1}\) 个空插槽。 - 第 \(m\) 个结点会从全部 \(m\) 个空插槽中等概率随机选择一个。

现在,我们来计算 \(S_m\) 的值。这取决于新结点落在哪里:

  • 情况A:新结点落在我们的子树内部。
    • 这个事件发生的概率是 \(P_A = \frac{S_{m-1}}{m}\)
    • 如果发生,子树失去1个旧插槽,但新结点带来2个新插槽。所以子树的插槽数净增1。
    • 此时,\(S_m = S_{m-1} + 1\)
  • 情况B:新结点落在我们的子树外部。
    • 这个事件发生的概率是 \(P_B = 1 - \frac{S_{m-1}}{m}\)
    • 如果发生,子树的插槽数不变。
    • 此时,\(S_m = S_{m-1}\)

2. 计算 \(E[S_m | S_{m-1}]\)

根据条件期望的定义,我们固定 \(S_{m-1}\) 的值,计算 \(S_m\) 的期望: \[ \begin{aligned} E[S_m | S_{m-1}] &= (S_{m-1} + 1) \cdot P_A + (S_{m-1}) \cdot P_B \\ &= (S_{m-1} + 1) \cdot \frac{S_{m-1}}{m} + S_{m-1} \cdot \left(1 - \frac{S_{m-1}}{m}\right) \\ &= \frac{S_{m-1}^2 + S_{m-1}}{m} + S_{m-1} - \frac{S_{m-1}^2}{m} \\ &= S_{m-1} + \frac{S_{m-1}}{m} \\ &= S_{m-1} \cdot \frac{m+1}{m} \end{aligned} \]

3. 求解 \(E_m\) 的递推式

现在使用全期望定律: \(E_m = E[S_m] = E[E[S_m | S_{m-1}]]\) \(E_m = E\left[S_{m-1} \cdot \frac{m+1}{m}\right]\)

因为 \(\frac{m+1}{m}\) 是一个常数,可以提出来: \(E_m = \frac{m+1}{m} \cdot E[S_{m-1}] = \frac{m+1}{m} \cdot E_{m-1}\)

这是一个非常简洁的递推关系!我们可以把它展开: \[ E_m = \frac{m+1}{m} E_{m-1} = \frac{m+1}{m} \cdot \frac{m}{m-1} E_{m-2} = \dots \] 这是一个连乘的形式。我们一直展开到我们的初始时刻,也就是第 \(t\) 天。 \[ E_m = E_t \cdot \frac{t+2}{t+1} \cdot \frac{t+3}{t+2} \cdot \dots \cdot \frac{m}{m-1} \cdot \frac{m+1}{m} \] 观察这个连乘,中间项会全部约掉(比如 \(t+2\)\(t+2\) 约掉),只剩下头尾: \[ E_m = E_t \cdot \frac{m+1}{t+1} \]

4. 确定初始值 \(E_t\)

在第 \(t\) 天,结点 \(t\) 刚刚被创建。它自己构成了大小为1的子树。这个子树有 2 个空插槽。这是一个确定的值,不是随机的。 所以 \(S_t = 2\),因此 \(E_t = E[S_t] = 2\)

代入上式,我们得到 \(E_m\) 的通项公式: \[ E_m = 2 \cdot \frac{m+1}{t+1} \]

5. 得到最终公式

我们需要的是第 \(N\) 天的情况,所以令 \(m=N\)\[E[S_N] = 2 \cdot \frac{N+1}{t+1}\]

最后,转换回子树大小 \(X_t\)\(E[X_t] = E[S_N - 1] = E[S_N] - 1\) \[ \boxed{E[X_t] = \frac{2(N+1)}{t+1} - 1} \]


第二部分:推导子树大小平方的期望 \(E[X_t^2]\)

这个过程类似,但计算会更复杂一些。我们的目标是求 \(E[X_t^2]\),这需要我们先求出 \(A_N = E[S_N^2]\)

1. 建立 \(S_m^2\)\(S_{m-1}^2\) 的递推关系

我们还是利用情况A和情况B: - 情况A (概率 \(\frac{S_{m-1}}{m}\)): \(S_m = S_{m-1} + 1\),所以 \(S_m^2 = (S_{m-1}+1)^2 = S_{m-1}^2 + 2S_{m-1} + 1\)。 - 情况B (概率 \(1 - \frac{S_{m-1}}{m}\)): \(S_m = S_{m-1}\),所以 \(S_m^2 = S_{m-1}^2\)

2. 计算 \(E[S_m^2 | S_{m-1}]\)

\[ \begin{aligned} E[S_m^2 | S_{m-1}] &= (S_{m-1}^2 + 2S_{m-1} + 1) \cdot \frac{S_{m-1}}{m} + (S_{m-1}^2) \cdot \left(1 - \frac{S_{m-1}}{m}\right) \\ &= \left(\frac{S_{m-1}^3}{m} + \frac{2S_{m-1}^2}{m} + \frac{S_{m-1}}{m}\right) + \left(S_{m-1}^2 - \frac{S_{m-1}^3}{m}\right) \\ &= S_{m-1}^2 + \frac{2}{m}S_{m-1}^2 + \frac{1}{m} S_{m-1} \\ &= \left(1 + \frac{2}{m}\right) S_{m-1}^2 + \frac{1}{m} S_{m-1} \end{aligned} \]

3. 求解 \(A_m\) 的递推式

使用全期望定律: \(A_m = E[S_m^2] = E[E[S_m^2 | S_{m-1}]]\) \(A_m = E\left[\left(1 + \frac{2}{m}\right) S_{m-1}^2 + \frac{1}{m} S_{m-1}\right]\)

根据期望的线性性质 \(E[aY+bZ] = aE[Y]+bE[Z]\)\[ A_m = \left(1 + \frac{2}{m}\right) E[S_{m-1}^2] + \frac{1}{m} E[S_{m-1}] \] \[ A_m = \frac{m+2}{m} A_{m-1} + \frac{1}{m} E_{m-1} \] 我们已经求出了 \(E_{m-1} = 2 \cdot \frac{m}{t+1}\)注意这里需要使用 \(m\) 而不是 \(m-1\)\[ E_{m-1} = 2 \cdot \frac{(m-1)+1}{t+1} = 2 \cdot \frac{m}{t+1} \]

代入上式: \[ A_m = \frac{m+2}{m} A_{m-1} + \frac{1}{m} \cdot \frac{2m}{t+1} \] \[ A_m = \frac{m+2}{m} A_{m-1} + \frac{2}{t+1} \]

这是一个一阶线性非齐次递推方程,需要用标准方法求解。

4. 求解该线性递推式并得到最终公式

目标:求解递推方程

我们要解决的递推关系是: \[ A_m = \frac{m+2}{m} A_{m-1} + \frac{2}{t+1} \quad (m > t) \] 其初始条件为第 t 天的情况。在第 t 天,子树就是结点 t 本身,它有 \(S_t=2\) 个确定的空插槽。因此, \[ A_t = E[S_t^2] = 2^2 = 4 \]

这是一个一阶线性非齐次递推方程。我们将使用求和因子 的标准方法来求解它,这种方法类似于求解线性微分方程时使用的积分因子法。

详细推导步骤

第一步:将方程标准化

首先,我们将所有与 \(A\) 相关的项移到一边,写成标准形式 \(A_m - f(m)A_{m-1} = g(m)\)\[ A_m - \frac{m+2}{m} A_{m-1} = \frac{2}{t+1} \]

第二步:寻找求和因子

我们的目标是找到一个因子 \(I_m\),当它乘以方程两边后,左边可以“折叠”成一个单一的差分形式,即 \(B_m - B_{m-1}\) 的样子。 我们希望: \[ I_m \left( A_m - \frac{m+2}{m} A_{m-1} \right) = (I_m A_m) - (I_{m-1} A_{m-1}) \] 为了让这个等式成立,通过比较 \(A_{m-1}\) 的系数,我们必须满足: \[ I_m \frac{m+2}{m} = I_{m-1} \] 或者写成关于 \(I_m\) 的递推式: \[ I_m = I_{m-1} \cdot \frac{m}{m+2} \] 我们可以通过连乘来解出 \(I_m\)。为了方便,我们设定初始值 \(I_t=1\)\[ \begin{aligned} I_m &= I_{m-1} \cdot \frac{m}{m+2} \\ &= I_{m-2} \cdot \frac{m-1}{m+1} \cdot \frac{m}{m+2} \\ &= \dots \\ &= I_t \cdot \frac{t+1}{t+3} \cdot \frac{t+2}{t+4} \cdot \dots \cdot \frac{m-1}{m+1} \cdot \frac{m}{m+2} \end{aligned} \] 我们把分子和分母分别写出来: * 分子: \((t+1) \cdot (t+2) \cdot (t+3) \cdot \dots \cdot m\) * 分母: \((t+3) \cdot (t+4) \cdot \dots \cdot m \cdot (m+1) \cdot (m+2)\)

可以看到,从 \((t+3)\)\(m\) 的所有项在分子和分母中都存在,可以约分。约分后剩下: \[ I_m = I_t \cdot \frac{(t+1)(t+2)}{(m+1)(m+2)} \] 由于我们设定了 \(I_t=1\),所以求和因子就是: \[ I_m = \frac{(t+1)(t+2)}{(m+1)(m+2)} \]

第三步:应用求和因子并求和

我们将原始方程 \(A_m - \frac{m+2}{m} A_{m-1} = \frac{2}{t+1}\) 的两边都乘以求和因子 \(I_m\)\[ I_m A_m - I_m \frac{m+2}{m} A_{m-1} = I_m \frac{2}{t+1} \] 根据我们构造 \(I_m\) 的方式,左边变成了 \(I_m A_m - I_{m-1} A_{m-1}\)。所以: \[ I_m A_m - I_{m-1} A_{m-1} = \frac{(t+1)(t+2)}{(m+1)(m+2)} \cdot \frac{2}{t+1} = \frac{2(t+2)}{(m+1)(m+2)} \] 现在,我们将这个方程从 \(m = t+1\)\(N\) 进行求和: \[ \sum_{m=t+1}^{N} (I_m A_m - I_{m-1} A_{m-1}) = \sum_{m=t+1}^{N} \frac{2(t+2)}{(m+1)(m+2)} \] 左边是一个伸缩求和 (Telescoping Sum)\[ (I_{t+1}A_{t+1} - I_t A_t) + (I_{t+2}A_{t+2} - I_{t+1}A_{t+1}) + \dots + (I_N A_N - I_{N-1}A_{N-1}) = I_N A_N - I_t A_t \] 右边可以提出常数项: \[ 2(t+2) \sum_{m=t+1}^{N} \frac{1}{(m+1)(m+2)} \] 我们使用部分分式分解来处理求和内部的项:\(\frac{1}{(m+1)(m+2)} = \frac{1}{m+1} - \frac{1}{m+2}\)。 所以右边的和也是一个伸缩求和: \[ \begin{aligned} \sum_{m=t+1}^{N} \left(\frac{1}{m+1} - \frac{1}{m+2}\right) &= \left(\frac{1}{t+2} - \frac{1}{t+3}\right) + \left(\frac{1}{t+3} - \frac{1}{t+4}\right) + \dots + \left(\frac{1}{N+1} - \frac{1}{N+2}\right) \\ &= \frac{1}{t+2} - \frac{1}{N+2} \end{aligned} \]

第四步:合并并求解 \(A_N\)

现在,我们将左右两边的结果合并: \[ I_N A_N - I_t A_t = 2(t+2) \left( \frac{1}{t+2} - \frac{1}{N+2} \right) \] \[ I_N A_N - I_t A_t = 2 - \frac{2(t+2)}{N+2} \]

代入我们已知的值: * \(I_N = \frac{(t+1)(t+2)}{(N+1)(N+2)}\) * \(I_t = 1\) * \(A_t = 4\)

\[ \frac{(t+1)(t+2)}{(N+1)(N+2)} A_N - 4 = 2 - \frac{2(t+2)}{N+2} \]

现在,我们解出 \(A_N\)\[ \frac{(t+1)(t+2)}{(N+1)(N+2)} A_N = 6 - \frac{2(t+2)}{N+2} \]

将右边通分: \[ 6 - \frac{2(t+2)}{N+2} = \frac{6(N+2) - 2(t+2)}{N+2} = \frac{6N + 12 - 2t - 4}{N+2} = \frac{6N - 2t + 8}{N+2} \]

于是有: \[ \frac{(t+1)(t+2)}{(N+1)(N+2)} A_N = \frac{6N - 2t + 8}{N+2} \]

\(A_N\) 的系数移到右边: \[ A_N = \frac{(N+1)(N+2)}{(t+1)(t+2)} \cdot \frac{6N - 2t + 8}{N+2} \]

\[ A_N = \frac{(N+1)(6N - 2t + 8)}{(t+1)(t+2)} \]

现在,转换回 \(X_t^2\)

我们有 \(X_t = S_N - 1\),所以: \[ E[X_t^2] = E[(S_N - 1)^2] = E[S_N^2 - 2S_N + 1] = E[S_N^2] - 2E[S_N] + 1 \]

我们已知: - \(E[S_N^2] = A_N = \frac{(N+1)(6N - 2t + 8)}{(t+1)(t+2)}\) - \(E[S_N] = 2 \cdot \frac{N+1}{t+1}\)

代入计算: \[ \begin{aligned} E[X_t^2] &= \frac{(N+1)(6N - 2t + 8)}{(t+1)(t+2)} - 2 \cdot \frac{2(N+1)}{t+1} + 1 \\ &= \frac{(N+1)(6N - 2t + 8)}{(t+1)(t+2)} - \frac{4(N+1)}{t+1} + 1 \\ &= \frac{(N+1)(6N - 2t + 8) - 4(N+1)(t+2)}{(t+1)(t+2)} + 1 \\ &= \frac{(N+1)[6N - 2t + 8 - 4t - 8]}{(t+1)(t+2)} + 1 \\ &= \frac{(N+1)(6N - 6t)}{(t+1)(t+2)} + 1 \\ &= \frac{6(N+1)(N - t)}{(t+1)(t+2)} + 1 \end{aligned} \]

最终得到: \[ \boxed{E[X_t^2] = 1 + \frac{6(N+1)(N-t)}{(t+1)(t+2)}} \]

第五步:汇总与计算

现在我们拥有了所有工具:

总不便度的期望 \(E = Σ E[第t条边的贡献] (从 t=2 到 N)\)

每条边的期望贡献 \(C_t = N * E[X_t] - E[X_t^2]\)

所以,我们只需要写一个循环,\(t\) 从 2 遍历到 \(N\): 1. 根据上面的公式计算出 \(E[X_t]\)\(E[X_t^2]\)。 2. 计算出这条边的期望贡献 \(N * E[X_t] - E[X_t^2]\)。 3. 把这个贡献累加到总期望 \(E\) 中。

如何处理精度问题?

在计算过程中,会出现大量的分数。使用Python的 fractions.Fraction 模块。这个模块可以精确地处理分数运算,不会有任何精度损失。

  • 我们用 Fraction 对象来存储 E 和中间计算结果。
  • 循环结束后,我们得到一个精确的分数 E。
  • 题目要求 E × N! 对 P 取模。我们可以先计算 E * N!。由于题目的保证,这个结果一定是个整数。
  • 最后,将这个整数结果对 P 取模。

算法实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
from fractions import Fraction
import math
import sys

def compute_E_times_factorial_mod(N, P):
E = Fraction(0, 1)

for t in range(2, N + 1):
Ex = Fraction(2 * (N + 1), t + 1) - 1
Ex2 = Fraction(1, 1) + Fraction(6 * (N + 1) * (N - t), (t + 1) * (t + 2))
Ct = N * Ex - Ex2
E += Ct

val = E * math.factorial(N)
val_int = val.numerator // val.denominator
return val_int % P

if __name__ == "__main__":
data = sys.stdin.read().strip().split()
N = int(data[0])
P = int(data[1])
print(compute_E_times_factorial_mod(N, P))

复杂度分析

  • 时间复杂度: \(O(N)\)
    • 主循环执行\(N-1\)
    • 每次迭代包含\(O(1)\)的算术运算
  • 空间复杂度: \(O(1)\)
    • 仅使用常数