Boundary Value Problems: Shooting and Finite Differences

Claude prompt

This page was generated by Claude Code from the following prompt (with additional prompting):

Let's make a boundary value problem js demo. First, I just want to show the difference between initial value and boundary value problems. For a relatively simple second-order ODE, show two plots, one where I specify x(0) and xdot(0), and one where I specify x(0) and x(T). Second, let's look at a nonlinear heated rod (nonlinearity is due to convection and radiation (sigma'' (T_∞^4 - T^4)). I'd like to show a shooting method history and a finite difference method with different step sizes.

1. Initial value vs. boundary value

A lightly damped mass on a spring: \[ \ddot{x} = -\omega^2 x - c\,\dot{x}, \qquad \omega = 2,\quad c = 0.2,\quad 0 \le t \le T = 4. \] A second-order equation needs two conditions to pin down a solution. Where those conditions are given changes the problem completely.

Initial value problem

Both conditions at \(t = 0\).

\(x(0)\) = 1.0 \(\dot{x}(0)\) = 3.0

Given \(x(0)\) and \(\dot{x}(0)\), march forward in time. Result: \(x(T)\) = .

Boundary value problem

One condition at each end.

\(x(0)\) = 1.0 \(x(T)\) = -1.0

Given \(x(0)\) and \(x(T)\), there is nothing to march from: \(\dot{x}(0)\) is unknown. The initial velocity that hits the target is \(\dot{x}(0)\) = .

Every ODE solver we have (Euler, RK4, ode45) starts from a full state and steps forward. For the boundary value problem we must somehow find the missing initial slope, or give up on marching altogether and solve for all the points at once. Those are the two strategies below.

2. A heated rod

A thin rod of length \(L = 10\) is held at \(T(0) = 40\) at one end and \(T(L) = 200\) at the other, in surroundings at \(T_\infty = 20\). Heat conducts along the rod and leaves its surface by convection and radiation: \[ \frac{d^2 T}{dx^2} = h'\,(T - T_\infty) + \sigma''\,(T^4 - T_\infty^4), \qquad T(0) = 40,\quad T(L) = 200. \] The \(T^4\) term makes the equation nonlinear. Set \(\sigma'' = 0\) to see what changes when the problem is linear.

\(h'\) = 0.01 \(\sigma''\) = 2.7 × 10⁻⁹

Shooting method

How to shoot

Guess the missing slope \(z = T'(0)\), integrate the IVP with RK4 out to \(x = L\), and look at how far \(T(L)\) lands from 200. That makes a root-finding problem, \(g(z) = T(L; z) - 200 = 0\), in the single unknown \(z\). Each new guess comes from the secant method through the previous two: \[ z_{k+1} = z_k - g(z_k)\,\frac{z_k - z_{k-1}}{g(z_k) - g(z_{k-1})}. \]

first guess \(z_0\) = 10 second guess \(z_1\) = 20

show through shot 1 of

earlier shots current shot boundary conditions

Finite differences

Instead of marching, split the rod into \(n\) segments of width \(\Delta x = L/n\) and replace the derivative at every interior node by a centered difference: \[ \frac{T_{i-1} - 2T_i + T_{i+1}}{\Delta x^2} = h'\,(T_i - T_\infty) + \sigma''\,(T_i^4 - T_\infty^4), \qquad i = 1, \dots, n-1. \] With \(T_0 = 40\) and \(T_n = 200\) known, this is \(n - 1\) equations in \(n - 1\) unknowns, all coupled to their neighbors. Because of the \(T_i^4\) term the system is nonlinear, so it is solved with Newton-Raphson; the Jacobian is tridiagonal.

\(n\) = 4 segments, \(\Delta x\) = 2.5

finite difference nodes reference solution

Maximum node error against the reference (finite differences with \(n = 4800\)). The dashed line has slope 2: halving \(\Delta x\) cuts the error by about four.