bisection method

bisection method is a bracketing method that uses intermediate value theorem to bracket an interval that contains local minimum

Given a differentiable function , if we can fine a closed interval that is defined on, such that and , then exists a value where ( is a critical point).

Once we’re able to find such an interval, we shrink it until we reach a desired level of precision.

Algorithm

Summary

(Forum Notebook Q4) Inputs: a function , endpoints of an initial interval , a level of precision

1. Calculate:
2. Cases:
3. If , then TERMINATE since we find the bracket
4. If , then set (shrink from right side)
5. If , then set (shrink from left side)
6. If TERMINATE, else LOOP

Output: Number of steps: (similar to binary search)

Note that this bisection method requires an interval where . If we start with an arbitrary interval then we must expand that initial interval until the property is satisfied. That’s what expand_interval does.

def bisection_method(f, a, b, eps=1e-5):
    """Bisection method to find the minimum of function f within the interval [a, b]. 
    
    Inputs:
    - f: The function to minimize.
    - a: Left endpoint of the interval.
    - b: Right endpoint of the interval.
    - eps: The precision level for termination (default is 1e-5).
    
    Output:
    - The bracketed interval [a, b] containing the local minimum. 
    """
    if a >= b:
        a, b = b, a
    
    # If either f'(a) or f'(b) doesn't meet conditions, expand the interval
    f_prime_a = central_difference(f, a, h=eps)
    f_prime_b = central_difference(f, b, h=eps)
    if f_prime_a >= 0 or f_prime_b <= 0:
        a, b = expand_interval(f, a, b, h=eps)
        
    while abs(b - a) > eps:
        m = (a + b) / 2
        f_prime_m = central_difference(f, m, eps)
 
        if f_prime_m > 0:
            b = m
        else:
            a = m
 
    return float(a), float(b)
 
def expand_interval(f, a, b, w_factor=1.0, h=1e-5):
    """Expands the interval [a, b] until f'(a) < 0 and f'(b) > 0 are satisfied.
    
    Inputs:
    - f: The function to minimize.
    - a: Left endpoint of the interval.
    - b: Right endpoint of the interval.
    - w_factor: Factor for expanding the interval (default is 1.0).
    
    Output:
    - The new interval [a, b].
    """
    while True:
        w = abs((b - a) / 2)
        
        f_prime_a = central_difference(f, a, h)
        f_prime_b = central_difference(f, b, h)
        
        if f_prime_a < 0 and f_prime_b > 0:
            return a, b
        
        a = a - w_factor * w
        b = b + w_factor * w
        
def central_difference(f, x, h=1e-5):
    """Calculate the derivative of f at point x using the central difference method.
    
    Inputs:
    - f: function whose derivative is to be calculated
    - x: point at which to calculate the derivative
    - h: step size for central difference method (default is 1e-5)
    
    Output:
    - Approximation of f'(x)
    """
    if isinstance(f, sage.symbolic.expression.Expression):
        return (f(x = x + h) - f(x = x - h)) / (2 * h)
    else: 
        return (f(x + h) - f(x - h)) / (2 * h)
 
bisection_method((x-3)^2 + 2, a=0, b=2)