[states, times] = rk4sys(@dynamics2body, 0, 3.8, [1;0;0;.75], 100); plot(states(:,1), states(:,3), 'k*', states(:,1), states(:,3), 'k-')