SharedPolyhedra.ipynbOpen in CoCalc
Polyhedra as done at Software Tools for Mathematics, Koper 2018

Expected Number of Points on a Sphere such that their Convex Hull Contains the Origin

def random_unit_vector(dimension, distribution):
    v = vector(distribution.get_random_element())
    v = [QQ(x) for x in v]
    return v

@parallel
def count_points_required(dimension, distribution):
    points = []
    while True:
        points.append(random_unit_vector(dimension, distribution))
        if len(points) <= dimension:
            continue
        convex_hull = Polyhedron(vertices=points)
        if vector(QQ, dimension) in convex_hull:
            return len(points)
import IPython.display

def plot_mark(x, color, legend):
    return arrow((x, -1e-10), (x, 0), color=color, legend_label=legend, legend_color=color)

def plot_histogram(counts, dimension):
    average = mean(counts)
    expected = 2*dimension + 1
    title = ("For d={} found average {}, expected {} after {} trials"
             .format(dimension, N(average, digits=5), N(expected, digits=5), len(counts)))
    bins = srange(min(counts) - 1/2, max(counts) + 1)
    G = histogram(counts, bins=bins, title=title)
    G += plot_mark(average, "red", "average")
    G += plot_mark(expected, "green", "expected")
    return G

def test_conjecture(dimension, trials=10):
    distributions = [create_distribution() for x in range(trials)]
    test_cases = zip([dimension]*trials, distributions)
    counts = count_points_required(test_cases)
    counts = [output for (input, output) in counts]

    plot_histogram(counts, dimension).show()

def test_conjecture_interactive(dimension, trials=10, distribution=None):
    counts = []
    while True:
        test_cases = [(dimension, SphericalDistribution(dimension)) for x in range(trials)]
        counts_ = count_points_required(test_cases)
        counts += [output for (input, output) in counts_]

        H = plot_histogram(counts, dimension)
        IPython.display.clear_output(wait=True)
        H.show()
test_conjecture_interactive(3, trials=32, distribution='sphere')
WARNING: Some output was deleted.

TODO:

  • Make this stop gracefully. Currently pressing the stop button, prints some stack traces.