Sunday, February 28, 2021

A number of Young French staircases

$\require{AMScd}$

Suppose we are building a staircase with 10 blocks flush against a wall: four blocks on the bottom row, three blocks on the second, two blocks on the third row and one block on the top row. Let's assume that we cannot place blocks either (a) floating in mid-air (likely a good assumption) or (b) away from the wall or another block (less physically necessary, but sure, I'll bite).

How many different ways are there of constructing this 10 block, 4 step (or more generally an $n(n+1)/2$ block, $n$ step) staircase under those conditions?

Let's orient ourselves with the wall standing at the $y$-axis and the $x$-axis is the ground. Then our slippy, slidy blocks will form what is known as a Young diagram in French notation for the partition $(n, n-1, n-2, \dots, 3, 2, 1)$ of $n(n+1)/2$. If we are to number the blocks that we place, then they will form a standard Young tableau, where the numbering is non-decreasing when increasing both row- and columnwise. A lot to unpack there, but we will make it through I'm sure.

$$\begin{CD} \cdot @= \cdot @. @. @. @.\\ @| @| @. @. @. \\ \cdot @= \cdot @= \cdot @. @. \\ @| @| @| @. @.\\ \cdot @= \cdot @= \cdot @= \cdot @. \\ @| @| @| @| @.\\ \cdot @= \cdot @= \cdot @= \cdot @= \cdot \\ @| @| @| @| @|\\ \cdot @= \cdot @= \cdot @= \cdot @= \cdot \end{CD}$$

When enumerating standard Young tableau, we can avail ourselves of what's known as the hook length formula, which for our case is given by for our partition $\lambda \in \mathbb{N}^d$ as $$F_\lambda = \frac{\left(\sum_{i=1}^d \lambda_i\right)!}{ \prod_{i,j} h(i,j\mid \lambda) },$$ where $$h(i,j \mid \lambda) = (\lambda_i - j) + (d - i) + 1$$ is the hook length for cell $(i,j)$ in the Young diagram for partition $\lambda.$

For our particular staircases, $\lambda_n = (n, n-1, n-2, \dots, 3, 2, 1)$, so we need to just calculate the hook length for each of the $\sum_i \lambda_{n,i} = n(n+1)/2$ cells of the Young diagram and then can apply the hook length formula. Given the shape of our staircases, we see that for $\lambda_n,$ there are $n$ cells with hook length $1,$ namely the top cell on each column (or correspondingly rightmost cell on each row). Moving one cell inward, we see that there are $n-1$ cells with hook length $3,$ and even further there are $n-k$ cells with hook length $(2k+1),$ for $k = 0, 1, \dots, n-1.$

Therefore, applying the hook length formula we get $$F_n = F_{(n,n-1,n-2, \dots, 3,2,1)} = \frac{\left( \frac{n(n+1)}{2} \right)!}{\prod_{k=0}^{n-1} (2k+1)^{(n-k)}}.$$ For $n=4,$ this gives $F_4 = \frac{10!}{1^4\cdot3^3\cdot5^2\cdot7^1} = 768.$ The OEIS entry A005118 contains the general term (offset by $1$), that is, $A005118(n+1) = F_n.$ Seeing as once you get to only a fairly reasonable $n=6$ steps, there are over $1.1$ billions combinations, I guess I will be satisfied with any possible tableau when constructing my next staircase.

Sunday, February 21, 2021

The Jenga catches you and you fall down

In Riddler Jenga, blocks are randomly placed on top of the stack until there is a collapse. All the blocks have the same alignment (e.g., east-west), and when a block, its center is picked randomly along the block directly beneath it.


On average, how many blocks must you place so that your tower collapses — that is, until at least one block falls off?


Let's assume that each block has length $1,$ so that if $X_i, \, i \geq 1,$ is the center of the $i$th block, where without loss of generality we assume $X_1=0.$ Since each block's center is uniformly distributed over the length of the block directly beneath it, that is, $X_{i+1} \mid X_i \sim U(X_i-0.5, X_i+0.5).$


After block is placed, we can determine whether or not the tower will fall via a center of mass calculation. If the first $n$ blocks are placed at $X_1=0, X_2, \dots, X_n,$ then the tower will fall if there is some $k = 1, \dots, n-1$ such that $$\left|\frac{X_n + \sum_{i=k+1}^{n-1} X_i}{n-k} - X_k \right| > \frac{1}{2}.$$


We could start off by calculating the probability of stack collapse at particular levels. For instance, \begin{align*}\mathbb{P} \{ \text{collapse at level $3$} \} &= \mathbb{P} \left\{ \left| \frac{X_2 + X_2}{2} \right| > \frac{1}{2} \right\}\\ &= \int_{-0.5}^{0.5} \int_{x_2-0.5}^{x_2+0.5} I( x_2 + x_3 \lt -1 ) + I(x_2 + x_3 \gt 1) \,dx_3 \, dx_2 \\ &= \int_{-0.5}^{0.5} (2x_2 - 0.5)_+ + (0.5 - 2x_2)_+ \,dx_2 = \frac{1}{8} \end{align*} However, these integrals get fairly intense quickly.


EDIT: My initial post and Python script did not properly account for the 0-indexing within Python and so was off by 1. It definitely makes me feel for Michael Bolton (no, not that no talent @$$-clown, ... the one in Office Space). I have corrected the details and Python below.



Rather, we can easily simulate the block stacks up to some arbitrarily large level, then check the condition above at each level, finding the first such level where the condition fails and denoting that as the maximum height. I encoded a small Jupyter notebook with this simulation design, displayed below. Using $S=500,000$ simulations gives an estimated mean of $7.1179$ with $95\%$ percentile confidence interval of $(7.1074,7.1283).$ One welcome tidbit from the histogram is that as expected there are about $62,364$ of the $500,000$ simulations with maximum height of $3$, which is $12.4728\%$ or roughly the analytical $1/8$th probability.


BlockStacker
In [1]:
import numpy as np
from matplotlib import pyplot as plt
%matplotlib inline

max_blocks = 50
num_samples = 500000

def stack_falls(X, level):
    diffs = [np.abs((X[level] + X[(k+1):level].sum()) /(level - k) - X[k]) for k in range(level)]
    return (max(diffs) > 0.5)

def get_max_stack(X):
    blocks, samples = X.shape
    stack_height = np.zeros((samples,1))
    for smpl in range(samples):
        level = 1
        while not stack_falls(X[:,smpl], level):
            level += 1
        stack_height[smpl] = level
    return stack_height+1 # shift index by 1


U = np.random.uniform(low=-0.5, high=0.5, size=(max_blocks, num_samples))
X = np.concatenate([np.zeros((1,num_samples)), np.cumsum(U, axis=0)])
M = get_max_stack(X)
mean = M.mean()
stdev = M.std()
sem = stdev / np.sqrt(num_samples)
max_ = M.max()
print('The average of the max stack height distribution using {} samples is {}.'.format(num_samples, mean))
print('The std. deviation of the max stack height distribution using {} samples is {}.'.format(num_samples, stdev))
print('The SE of the mean of the max stack height distribution using {} samples is {}.'.format(num_samples, sem))
print('95% confidence limit of the mean is ({},{}).'.format(mean-1.96*sem, mean+1.96*sem))
print('The maximum of the max stack height distribution using {} samples is {}.'.format(num_samples, max_))
plt.hist(M, bins=range(int(max_)))
plt.title('Histogram of max_stack_height using {} samples'.format(num_samples))
plt.show()
The average of the max stack height distribution using 500000 samples is 7.117868.
The std. deviation of the max stack height distribution using 500000 samples is 3.7658602117678237.
The SE of the mean of the max stack height distribution using 500000 samples is 0.005325730585483272.
95% confidence limit of the mean is (7.107429568052453,7.128306431947546).
The maximum of the max stack height distribution using 500000 samples is 45.0.
In [2]:
hist, bins = np.histogram(M, bins=range(int(max_)))
hist, bins
Out[2]:
(array([    0,     0,     0, 62364, 74507, 70945, 61455, 51291, 41540,
        32549, 25811, 19528, 14893, 11640,  8575,  6484,  4783,  3629,
         2707,  1882,  1460,  1028,   784,   589,   416,   314,   225,
          157,   110,    88,    75,    56,    29,    23,    22,    14,
            6,     7,     2,     3,     4,     0,     2,     2],
       dtype=int64),
 array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12, 13, 14, 15, 16,
        17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33,
        34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44]))
In [ ]: