154. 探索帕斯卡金字塔(Exploring Pascal's Pyramid)
用球形小球构造一个三角形金字塔,使得每个小球都恰好放在下一层的三个小球之上。
然后,我们计算从顶点通向每个位置的路径数量:
一条路径从顶点开始,并向下移动到当前位置正下方的三个球中的任意一个。
因此,到达某个位置的路径数量等于它正上方那些数字之和(根据位置不同,正上方最多有三个数字)。
其结果就是帕斯卡金字塔,并且第 \(n\) 层中的数字是三项式展开式 \((x + y + z)^n\) 的系数。
在 \((x + y + z)^{200000}\) 的展开式中,有多少个系数是 \(10^{12}\) 的倍数?
一、题意概述
第 \(n\) 层的每个位置对应一个三元组 \((a,b,c)\),满足:
对应的三项式系数为:
本题要求统计当 \(n=200000\) 时,有多少个有序三元组 \((a,b,c)\) 使得该系数能被 \(10^{12}\) 整除。由于 \(10^{12}=2^{12}5^{12}\),条件等价于:
其中 \(v_p(m)\) 表示质数 \(p\) 在 \(m\) 中的指数。
二、数学背景
对任意质数 \(p\),Legendre 公式给出:
其中 \(s_p(m)\) 是 \(m\) 在 \(p\) 进制下的数位和。因此:
这也正是 Kummer 定理在多项式系数上的形式:该指数等于用 \(p\) 进制把 \(a,b,c\) 相加得到 \(n\) 时产生的进位总数。
对本题的两个质数:
- 五进制条件为 \(s_5(a)+s_5(b)+s_5(c) \ge s_5(n)+4\cdot 12\)
- 二进制条件为 \(\operatorname{popcount}(a)+\operatorname{popcount}(b)+\operatorname{popcount}(c) \ge \operatorname{popcount}(n)+12\)
直接枚举所有 \((a,b,c)\) 需要约 \(\frac{(n+1)(n+2)}{2}\) 次检查,数量级约为 \(2\times 10^{10}\),不可行。
三、算法分析
候选方案一:直接枚举加剪枝
可以只枚举 \(a\le b\le c\),再按三元组是否有相等元素乘以 \(1,3,6\) 的对称系数。这能把枚举量缩小约 \(6\) 倍,但仍然是数十亿级别,不适合 Python。
继续按五进制数位和剪枝,能减少许多无效点,但仍需要在大量区间内逐点扫描,性能不稳定。
候选方案二:进位 DP 加位集校正
本解法采用该方案。
先只考虑 \(5^{12}\)。因为五进制位数很少,可以用进位 DP 一次性统计所有有序三元组中满足 \(v_5\ge 12\) 的数量。状态为:
逐位处理 \(n\) 的五进制数字,枚举当前位上 \(a,b,c\) 的三个数字,更新下一位进位与累计进位数。累计进位数只需要截断到 \(12\),因为更大的值在本题中没有区别。
接着扣除“已经满足 \(5^{12}\),但不满足 \(2^{12}\)”的三元组。二进制不足条件为:
这是一个较强的限制。为了避免逐点枚举,构造按 \((\operatorname{popcount}(x), s_5(x))\) 分组的位集:
exact[p][q]:第 \(x\) 位为 \(1\) 表示 \(x\) 的 popcount 为 \(p\),五进制数位和为 \(q\)reversed_exact[p][q]:第 \(n-x\) 位为 \(1\),用于把 \(c=n-a-b\) 映射到同一个 \(b\) 坐标
固定最小元素 \(a\),只考虑 \(a\le b\le c\),即:
对每组候选的 \(b\),用右移后的反向位集表示满足条件的 \(c=n-a-b\)。随后通过位与和 bit_count() 一次统计整段区间内的合法数量。最后按三元组元素相等情况把无序计数换成有序计数。
最终答案为:
四、复杂度分析
五进制进位 DP 的状态很小,复杂度约为 \(O(\log_5 n \cdot 12 \cdot 5^3)\)。
位集校正部分按 \(a\) 从 \(0\) 遍历到 \(\lfloor n/3\rfloor\)。每次查询不是逐个枚举 \(b\),而是在 Python 大整数位集上进行批量位运算;设分组计划数为 \(G\)、机器字宽为 \(w\),可视为约 \(O(\frac{n}{3}\cdot G\cdot \frac{n}{w})\) 的底层字操作。实际 \(n=200000\) 时分组数量较小,运行时间低于 \(60\) 秒。
空间上主要保存若干按特征分组的位集,数量由可能的 popcount 和五进制数位和决定,空间复杂度约为 \(O(C\cdot n)\) 位,其中 \(C\) 是非空特征组数量。
五、代码实现与说明
完整代码位于 scripts/pe154.py。核心函数按代码顺序如下。
基础数位工具
def digits_base(number: int, base: int) -> list[int]:
"""Return the little-endian digits of number in the given base."""
if number == 0:
return [0]
digits: list[int] = []
while number:
digits.append(number % base)
number //= base
return digits
digits_base 将整数转换成低位在前的进制数字列表,方便从最低位开始模拟加法进位。digit_sum_base 与 build_digit_sum_table 分别计算单个整数和整段区间的进制数位和。
五进制进位 DP
def count_prime_carry_at_least(level: int, prime: int, threshold: int) -> int:
"""Count ordered triples whose base-prime carry sum is at least threshold."""
if threshold <= 0:
return (level + 1) * (level + 2) // 2
transitions: list[list[list[tuple[int, int]]]] = []
for digit in range(prime):
by_carry: list[list[tuple[int, int]]] = []
for carry_in in range(3):
counts: dict[int, int] = {}
for first in range(prime):
for second in range(prime):
for third in range(prime):
total = first + second + third + carry_in
if total % prime == digit:
carry_out = total // prime
counts[carry_out] = counts.get(carry_out, 0) + 1
by_carry.append(list(counts.items()))
transitions.append(by_carry)
这一段预先枚举一位上的所有数字组合。对于给定的当前位目标数字 digit 和输入进位 carry_in,记录每种输出进位 carry_out 出现多少次。
dp = [[0] * (threshold + 1) for _ in range(3)]
dp[0][0] = 1
for digit in digits_base(level, prime):
next_dp = [[0] * (threshold + 1) for _ in range(3)]
for carry_in, row in enumerate(dp):
for carry_sum, ways in enumerate(row):
if ways == 0:
continue
for carry_out, transition_count in transitions[digit][carry_in]:
next_sum = min(threshold, carry_sum + carry_out)
next_dp[carry_out][next_sum] += ways * transition_count
dp = next_dp
return dp[0][threshold]
dp[carry][carry_sum] 表示处理到当前位时,输出进位为 carry、累计进位数为 carry_sum 的方案数。carry_sum 被截断到 threshold,因此所有已经达标的方案会汇总在最后一格。所有位处理完后,只接受最终进位为 0 的方案。
位集构造
def build_feature_bitsets(
level: int,
) -> tuple[list[int], list[int], list[list[int]], list[list[int]]]:
"""Build exact and reversed bitsets grouped by binary and quinary digit sums."""
popcounts = [value.bit_count() for value in range(level + 1)]
quinary_sums = build_digit_sum_table(level, 5)
这一段预计算每个整数的二进制 popcount 与五进制数位和。之后的校正逻辑只使用这两个特征。
for value in range(level + 1):
popcount = popcounts[value]
quinary_sum = quinary_sums[value]
exact[popcount][quinary_sum] |= 1 << value
reversed_exact[popcount][quinary_sum] |= 1 << (level - value)
exact 用第 value 位表示该值属于某个特征组;reversed_exact 把位置反向为 level - value。固定 a 后,reversed_exact >> a 的第 b 位就对应 \(c=level-a-b\)。
校正计划与短缺计数
def make_correction_plan(
popcount_budget: int,
quinary_sum_needed: int,
exact_nonempty: Iterable[tuple[int, int, int]],
reversed_unions: list[list[int]],
) -> list[tuple[int, int]]:
"""Group bitsets needed for one correction query."""
校正计划把同一类 c 约束对应的所有 b 位集合并起来,减少主循环中的位运算次数。popcount_budget 是在固定 \(a\) 后,\(b+c\) 还能使用的 popcount 上限;quinary_sum_needed 是 \(b+c\) 还需要达到的五进制数位和下限。
def count_power_of_two_shortfall(level: int, threshold: int) -> int:
"""Count ordered triples with enough factors of 5 but too few factors of 2."""
该函数统计需要从五进制 DP 结果中扣除的数量。它只枚举 \(a\le b\le c\) 的代表元:
max_first = level // 3
high = level // 2
range_mask = (1 << (high + 1)) - 1
range_mask 表示当前合法的 \(b\) 区间。每次 a 增加后,低端和高端只会移动一位,因此代码用异或增量更新掩码,避免重复构造长区间位集。
subtotal = 0
for bits_c, bits_b in plan:
subtotal += (bits_b & (bits_c >> first) & range_mask).bit_count()
correction += 6 * subtotal
这里完成核心统计:bits_b 给出候选 \(b\),bits_c >> first 把 \(c=level-first-b\) 映射到 \(b\) 坐标,range_mask 限制 \(a\le b\le c\)。三者取交后用 bit_count() 得到本轮代表元数量,先按一般情况乘以 \(6\)。
随后代码分别处理 \(a=b<c\)、\(a<b=c\)、\(a=b=c\) 三种边界。一般情况乘以 \(6\),两个相等时实际只有 \(3\) 个排列,三个全相等时实际只有 \(1\) 个排列,因此需要做小幅修正。
总求解入口
def solve(level: int = DEFAULT_LEVEL, threshold: int = DEFAULT_THRESHOLD) -> int:
"""Return the count of trinomial coefficients divisible by 10^threshold."""
enough_fives = count_prime_carry_at_least(level, 5, threshold)
two_shortfall = count_power_of_two_shortfall(level, threshold)
return enough_fives - two_shortfall
solve 先统计所有满足 \(5^{12}\) 的系数,再扣除其中 \(2^{12}\) 不足的部分,即得到同时满足 \(2^{12}\) 与 \(5^{12}\) 的系数数量。
文件末尾提供了 brute_force 和 run_self_tests。它们只用于小规模对照:直接枚举所有三元组,并用数位和公式判断 \(v_2\) 与 \(v_5\),从而验证优化算法的正确性。