function E=bvp_shooting_error(yp0)
y0 = 0;
y1 = pi/2;
[t,y] = ode45('bvp_shooting_rhs', [0 2], [y0 yp0]);
E = (y(end,1) - y1)^2;    % Squared error is easier to minimize