Answer to the fivethirtyeight.com riddler for Feb 19, 2021

jenga.py

#!/usr/bin/env python3
#
# Answer for the 538 Riddler for Feb 18, 2021 (https://fivethirtyeight.com/features/can-you-win-riddler-jenga/)
#
# Author: Lowell Vaughn
#
# This performs a Monte Carlo simulation of stacking the Jenga blocks and
# finds out at which point they would fall over.
#
# For this simulation all blocks are of width 1, and the first block is always placed at location 0.

import random
from scipy.stats import describe
import numpy as np


def is_balanced(stack):
    """
    Takes a stack of blocks (represented by a list of floats representing the center
    of each block) and returns True of it is balanced enough to no fall over. It does this
    by calculating the center of mass of each of the N-1 sub-towers and seeing if it is above
    the block that sub-tower rests on.
    """
    total_x_value = 0
    n_blocks = 0
    for i in range(len(stack) - 1, 0, -1):
        total_x_value += stack[i]
        n_blocks += 1
        center_of_mass = total_x_value / n_blocks
        if abs(center_of_mass - stack[i - 1]) > 0.5:  # it's not above that block
            return False
    return True


def stack_size():
    """
    Starts with a stack of size 1, centered at zero, and keeps adding blocks on top
    (randomly placed with respect to the top block) until it is no longer balanced.

    Returns the length of the stack when it's no longer balanced.
    """
    stack = [0]
    while is_balanced(stack):
        stack.append(stack[-1] + (random.random() - 0.5))
    return len(stack)


# Run a large number of simulations and collect the statistics on them, including information
# to make a histogram
n_simulations = 10000000
values = np.zeros((n_simulations,), dtype=int)
bins = np.zeros((100,), dtype=int)
for i in range(n_simulations):
    total_blocks = stack_size()
    values[i] = total_blocks
    bins[total_blocks] += 1

print("Final stats:", describe(values))
with open('jenga_bins.csv', 'w') as out:
    for i in range(100):
        out.write("{},{}\n".format(i, bins[i]))