coverdist.py

import functools

import numpy as np
import matplotlib.pyplot as plt

@functools.lru_cache(maxsize=None)
def coverdist(covers, max_pulls):
    result = np.zeros(max_pulls+1, dtype=np.float)
    pulls = sum(covers)
    a, b, c = covers
    complete = min(a, 5) + min(b, 5) + min(c, 5) >= 13

    if complete:
        result[pulls:] = 1.0
    elif pulls < max_pulls:
        result += coverdist(tuple(sorted([a+1,b,c])), max_pulls) / 3
        result += coverdist(tuple(sorted([a,b+1,c])), max_pulls) / 3
        result += coverdist(tuple(sorted([a,b,c+1])), max_pulls) / 3
    return result

if __name__ == '__main__':
    max_pulls = 35

    dist = coverdist((0,0,0), max_pulls)
    plt.bar(np.arange(max_pulls+1), dist)

    dist2 = coverdist((4,1,0), max_pulls)
    plt.bar(np.arange(max_pulls+1), dist2)
    plt.show()