Project Euler 726

Project Euler 726

题目

Falling bottles

Consider a stack of bottles of wine. There are \(n\) layers in the stack with the top layer containing only one bottle and the bottom layer containing \(n\) bottles. For \(n=4\) the stack looks like the picture below.

The collapsing process happens every time a bottle is taken. A space is created in the stack and that space is filled according to the following recursive steps:

  • No bottle touching from above: nothing happens. For example, taking \(F\).

  • One bottle touching from above: that will drop down to fill the space creating another space. For example, taking \(D\).

  • Two bottles touching from above: one will drop down to fill the space creating another space. For example, taking \(C\).

This process happens recursively; for example, taking bottle \(A\) in the diagram above. Its place can be filled with either \(B\) or \(C\). If it is filled with \(C\) then the space that \(C\) creates can be filled with \(D\) or \(E\). So there are \(3\) different collapsing processes that can happen if \(A\) is taken, although the final shape (in this case) is the same.

Define \(f(n)\) to be the number of ways that we can take all the bottles from a stack with \(n\) layers. Two ways are considered different if at any step we took a different bottle or the collapsing process went differently.

You are given \(f(1) = 1\), \(f(2) = 6\) and \(f(3) = 1008\).

Also define

\[\displaystyle S(n) = \sum_{k=1}^n f(k)\]

Find \(S(10^4)\) and give your answer modulo \(1\,000\,000\,033\).

解决方案

\(N=10^4.\)把酒瓶堆看作一个三角格上的离散结构。关键观察是:在任意时刻,剩余酒瓶的外形都不会出现悬空空洞。原因是:如果某个位置出现空位,且其上方仍有与之接触的酒瓶,那么塌落过程会继续发生,直到这个空位不再能被上方酒瓶填补为止。

因此,忽略酒瓶的原始标签后,剩余形状始终是一个向右下闭合的阶梯形区域。用组合论语言描述,它可以看成一个阶梯形的杨图

再看一次操作:从某处取走一个酒瓶后,空位会沿着“向上”的接触关系递归传递,直到某个最上层可消失的位置为止。于是,一次完整操作在形状上的净效果等价于:从当前阶梯形区域中移除一个角点格(杨图的角点格),其余被牵动的酒瓶只是沿某条路径下落以填补空位。

因此,为了计数,可以把每一步拆成两部分理解:先确定这一步最终对应移除的是哪一个角点瓶,再统计有多少种起始取瓶位置和递归分支选择的方式,会导致最终移除该角点瓶。

定义角点瓶的高度 \(h\):它下方(含自身)对应一个完整的三角形子堆,共 \(h\) 层。由上面的形状不变性可知,角点瓶下方必然是满三角,否则此前的塌落过程已经把空位继续传递并消除了。

现在固定一个高度为 \(h\) 的角点瓶 \(U\)。我们要计算:有多少种不同的取走某个酒瓶和塌落分支选择会使得最终消失的角点瓶恰好是 \(U\)

把塌落过程反向来看:从 \(U\) 向下回溯到最初被取走的位置。每下降一层,都对应左下 / 右下两种选择(正向过程里对应两瓶接触时二选一的分支)。因此,若最初被取走的酒瓶位于 \(U\) 下方第 \(d\) 层(从 \(0\) 开始计),则从该位置回溯到 \(U\) 的路径条数是杨辉三角第 \(d\) 行的和,即\(\displaystyle \sum_{j=0}^d \binom{d}{j} = 2^d.\)高度为 \(h\) 的满三角包含层数 \(d=0,1,\ldots,h-1\),于是总方式数为\(\displaystyle\sum_{d=0}^{h-1} 2^d = 2^h - 1.\)

也就是说,在任意状态下,若某一步最终移除的是高度为 \(h\) 的角点瓶,则该步对应的原始操作方式数恒为\(w(h)=2^h-1.\)它只依赖于高度 \(h\),与局部形状的其他细节无关。

把整个过程只记录为每一步移除了哪个角点瓶(即杨图上的角点移除序列)。对于任意一个合法序列,该序列对应的总方式数就是每一步权重 \(w(h)\) 的乘积。

现在考虑初始阶梯形杨图\(\lambda=(n,n-1,\ldots,1)\)中的每一个格子(对应一个酒瓶)。对格子使用坐标 \((i,j)\)(从 \(1\) 开始计数,且 \(j\le n+1-i\)),它的高度满足\(h=n+2-i-j.\)因此 \(h\) 的取值范围为 \(1,2,\ldots,n\)。固定某个高度 \(h\) 时,对应格子的个数恰好是 \(n+1-h\)(它们位于同一条反对角线上)。

无论采用哪一种合法的角点移除顺序,每个格子最终都会被移除一次,所以高度为 \(h\) 的权重 \(w(h)\) 一共会出现 \(n+1-h\) 次。于是对任意固定的合法序列,总塌落乘子恒为\(\displaystyle\prod_{h=1}^n (2^h-1)^{n+1-h}.\)

这说明:所有合法角点移除序列的权重完全相同

接下来只需计数:在阶梯形杨图\(\lambda=(n,n-1,\ldots,1)\) 上,从全图开始,每一步移除一个角点格,直到删空,总共有多少种序列。

把这个移除序列反过来,就等价于按顺序把数字 \(1,2,\ldots,T(n)\) 填入各格,使得每行、每列严格递增。这正是标准杨表的定义。

其中总格子数(也即总酒瓶数)为 \[T(n)=1+2+\cdots+n=\frac{n(n+1)}{2}\]

对任意杨图 \(\lambda\),标准杨表数满足著名的 hook-length 公式

\[ f^\lambda=\frac{|\lambda|!}{\prod_{c\in\lambda} H(c)}, \]

其中 \(H(c)\) 是格子 \(c\) 的 hook length:它等于该格子本身 + 同行右侧格子数 + 同列下方格子数。

由于第 \(i\) 行长度为 \(n+1-i\),第 \(j\) 列长度为 \(n+1-j\),所以\(H(i,j)=(n+1-i-j)+(n+1-j-i)+1=2(n+1-i-j)+1\)

\(t=n+1-i-j\)(其中 \(t\ge 0\)),则 hook length 为奇数 \(2t+1\)。当 \(t=0,1,\ldots,n-1\) 时,对应的奇数依次为 \(1,3,5,\ldots,2n-1\);并且数值 \(2t+1\) 的出现次数为 \(n-t\)。改写为以 \(h=t+1\) 记号后,得到

\[\prod_{c\in\lambda} H(c)=\prod_{h=1}^n (2h-1)^{n+1-h}\]

因此,合法角点移除序列数为

\[C(n)=\frac{T(n)!}{\prod_{h=1}^n (2h-1)^{n+1-h}}\]

综合前文,所有合法角点移除序列数为 \(C(n)\),且每条序列的权重相同,故

\[f(n)=C(n)\cdot \prod_{h=1}^n (2^h-1)^{n+1-h}\]

代入 \(C(n)\) 的表达式,得到闭式:

\[f(n)=T(n)!\cdot \prod_{h=1}^n \left(\frac{2^h-1}{2h-1}\right)^{n+1-h}\]

为了计算\(\displaystyle S(N)=\sum_{k=1}^{N} f(k)\),递推形式更方便。注意从 \(n-1\)\(n\) 时:

  • 总格子数从 \(T(n-1)\) 增加到 \(T(n)\),因此\(\displaystyle\frac{T(n)!}{T(n-1)!}=\prod_{m=T(n-1)+1}^{T(n)} m;\)
  • 乘积指数整体增加一层,因此额外出现因子\(\displaystyle\prod_{h=1}^n \frac{2^h-1}{2h-1}.\)

定义\(\displaystyle A(n)=\prod_{h=1}^n \frac{2^h-1}{2h-1},\)则有递推式

\[f(n)=f(n-1)\cdot \frac{T(n)!}{T(n-1)!}\cdot A(n)\]

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
N = 10000
MOD = 1000000033
def solve(N):
f = 1
s = 0
a = 1
pow2 = 1
t = 0
for k in range(1, N + 1):
x = t + 1
y = t + k
while x <= y:
f = (f * x) % MOD
x += 1
t = y
pow2 = (pow2 * 2) % MOD
a = (a * (pow2 - 1)) % MOD
a = (a * pow(2 * k - 1, MOD - 2, MOD)) % MOD
f = (f * a) % MOD
s = (s + f) % MOD
return s

print(solve(N))