function ode15s_runtime opt = odeset ("events", @ode_events, "abstol", 1e-6, "reltol", 1e-6, ... "maxstep", 0.00005); % Initial height, velocity and acceleration h0 = 0; v0 = 20; a0 = -9.81; y0 = [h0; v0]; ode_fct = @(t, y) f (t, y, a0); tic; [t, y, te, ye, ie] = ode15s (ode_fct, [0, 5], y0, opt); toc; length (y) figure; hold on; yexct = 0.5 * a0 * t.^2 + v0 * t + h0; plot (t, yexct, "r."); plot (t, y(:,1), "b."); xlabel ("t"); ylabel ("y"); legend ({"exact", "ode15"}); title ("Ball"); grid on; % ODE: y'' = a % Let y1 := y and y2 := y' then % [y1', y2'] = [y2, a] function dydt = f (t, y, a) dydt = [y(2); a]; end function [value, isterminal, direction] = ode_events (t, y) value = y(1); isterminal = 1; direction = 0; end end