def bisect_bracket(f, a, b, tol, threshold=1e-10):
    """Given a function f with a root between A and B, return a new
    interval (C,D) containing the root such that D-C < TOL."""
    assert a < b
    assert tol > 0
    ya = f(a)
    yb = f(b)
    assert ya != 0
    assert yb != 0
    assert (ya < 0) != (yb < 0), 'a=%f, b=%f, ya=%f, yb=%f' % (a, b, ya, yb)
    while b - a > tol:
        c = (a + b) / 2
        yc = f(c)
        if abs(yc) < threshold:
            return a, b
        if (yc < 0) == (ya < 0):
            a = c
            ya = yc
        else:
            b = c
            yb = yc
    return a, b
