Boundary value problems

Contents

Lecture Materials

Download the slide deck for this lecture

Boundary value problems#

Boundary value problems are problems where the solution is specified at the boundaries of the domain, rather than at the initial point. To solve such problems, we need to use a different approach than the initial value problems.

Shooting method#

Let us solve the problem of a vertically thrown ball at \(t = 0\) which lands back on the ground at \(t = 10\). We want to find the initial velocity. The equations of motion are

(8.3)#\[\begin{align} \frac{dx}{dt} & = v, \\ \frac{dv}{dt} & = -g, \end{align}\]

\(x(0) = 0\) and \(v(0) = v_0\).

We shall find the value of \(v_0\) using the bisection method and solving the ODEs with 4th order Runge-Kutta method.

last_bisection_iterations = 0  # Count how many iterations it took
bisection_verbose = False

def bisection_method(
    f,                    # The function whose root we are trying to find
    a,                    # The left boundary
    b,                    # The right boundary
    tolerance = 1.e-10,   # The desired accuracy of the solution
    ):
    fa = f(a)                           # The value of the function at the left boundary
    fb = f(b)                           # The value of the function at the right boundary
    if (fa * fb > 0.):
        return None                     # Bisection method is not applicable
    
    global last_bisection_iterations
    last_bisection_iterations = 0
    
    
    while ((b-a) > tolerance):
        last_bisection_iterations += 1
        c = (a + b) / 2.                # Take the midpoint
        fc = f(c)                       # Calculate the function at midpoint
        
        
        if bisection_verbose:
            print("Iteration: {0:5}, c = {1:20.15f}, f(c) = {2:10.15f}".format(last_bisection_iterations, c, fc))
        
        if (fc * fa < 0.):              
            b = c                       # The midpoint is the new right boundary
            fb = fc
        else:                           
            a = c                       # The midpoint is the new left boundary
            fa = fc

    return (a+b) / 2.                   

g = 9.81 # m/s^2
# ODEs
def fball(xin,t):
    x = xin[0]
    v = xin[1]
    return np.array([v,-g])

# Initial and final times
t0 = 0.
tend = 10.
# Number of RK4 steps
Nrk4 = 100
hrk4 = (tend - t0) / Nrk4

# Desired accuracy for v0
accuracy_v0 = 1.e-10 # m/s

v0min = 0.01   # m/s
v0max = 1000.0 # m/s

def fbisection(v0):
    x0 = [0., v0]
    return ode_rk4_multi(fball, x0, t0, hrk4, Nrk4)[1][-1][0]

v0sol = bisection_method(fbisection, v0min, v0max, accuracy_v0)
print("The required initial velocity is",v0sol,"m/s")

The required initial velocity is 49.0500000000017 m/s
Contents