function [t,y] = midpoint( f, tspan, y0, n)
%MIDPOINT Apply the midpoint method to solve y' = f(y,t), y(a) = y0,
% on the interval a <= t <= b with n steps.
% This is a second order Runge-Kutta method, and is similar to the
% midpoint rule for quadrature.
a = tspan(1);
b = tspan(2);
h = (b-a)/n;
t = a:h:b;
y = zeros(1,n+1);
y(1) = y0;
for j = 1:n
F1 = f(t(j),y(j)); % slope at current point
F2 = f(t(j)+.5*h,y(j)+.5*h*F1); % slope at approx. mid point
y(j+1) = y(j) + h*F2; % Move by midpoint slope.
end
end