from math import exp

# Print out what p_i looks like for i = 0, 10, 20, ...
def explore(alpha, d):
    p = 0
    
    for i in range(11):
        print(f"{10 * i : 4}: {p}")
        for j in range(10):
            p = exp(-alpha * d * (1-p)**(d- 1))

# Approximate lim (i -> infty) p_i by running the
# recurrence for 10,000 steps.
def limit_probability(alpha, d):
    p = 0
    for i in range(10000):
        p = exp(-alpha * d * (1-p)**(d-1))
    return p

# Use a binary search to find the value alpha*,
# where for alpha < alpha* we succeed with probability nearly 1
# and for alpha > alpha* we succeed with non-1 probability.
def bisect(d):
    low = 0; high = 1
    for i in range(500):
        alpha = (low + high) / 2
        if (limit_probability(alpha, d) < 1): high = alpha
        else: low = alpha
    return low
