Checking the "Boxes puzzle" solution strategy via stochastic simulation

boxes_puzzle.jl

# Checking the "Boxes puzzle" solution strategy via stochastic simulation.
# robert.feldt@gmail.com 2022

using Random # We need randperm below

# Fill N boxes by getting a random permutation of length N.
fillboxes(N) = randperm(N)

# Given a strategy, simulate its use by N people and return true iff strategy 
# worked for them all. 
function simulate_strategy(strategy::Function, N::Int)
    boxes = fillboxes(N)

    for p in 1:N
        if !strategy(boxes, p)
            return false # strategy didn't work for person p, so we abort
        end
    end

    return true
end

# Now simulate many times and calculate the ratio of successful simulations.
function ratio_from_multiple_simulations(strategy::Function, N::Int, numiterations = 10000)
    numsuccesses = 0
    for simnum in 1:numiterations
        if simulate_strategy(strategy, N)
            numsuccesses += 1
        end
    end
    return numsuccesses/numiterations
end

# We must find our own number within floor(Int, 0.5*N)
maxchecks(N::Int) = floor(Int, 0.5*N)
maxchecks(boxes::Vector) = maxchecks(length(boxes))

# A very bad strategy is that each person checks a random selection of boxes
function bad_random_strategy(boxes, person)
    order = randperm(length(boxes))
    for check in 1:maxchecks(boxes)
        box = order[check]
        if boxes[box] == person
            return true
        end
    end
    return false # if we come here we failed in finding our own number
end

# Proposed solution strategy is that each person starts by checking her
# box and then uses the number in the box to determine which box to check 
# next.
function good_strategy(boxes, person)
    box = person # We start by checking our "own" box
    for check in 1:maxchecks(boxes)
        if boxes[box] == person
            return true
        else
            box = boxes[box] # next box to check is the one in the current box
        end
    end
    return false # if we come here we failed in finding our own number
end

# Run once to ensure everything compiled
ratio_from_multiple_simulations(bad_random_strategy, 100, 1000) # should be 0.0!
ratio_from_multiple_simulations(good_strategy, 100, 1000)

# Now run for all even N from 2 to 102 and ensure they are all > 0.30.
# To be fancy and use our shiny multi-core CPU we run threads in parallel... :)
Nrange = 2:2:102
allratios = zeros(Float64, length(Nrange))
@time Threads.@threads for i in 1:length(Nrange)
    n = 2*i
    allratios[i] = ratio_from_multiple_simulations(good_strategy, n, 100_000)
end # takes about 1.6 seconds on my MB Pro with M1 Max, using 8 threads
@assert all(r -> r > 0.30, allratios)

Background, analysis and pointers to proof of optimality of the good_strategy can be found here:

https://en.wikipedia.org/wiki/100_prisoners_problem