function overheads
timestep_euler = 0.01;
timestep_runge = 0.02;
timestep_rk = 0.04;

% Euler
[to,xo]=ode23(@oscillator,[0,20],[2;4]);
figure;
plot(to,xo(:,1),'m',to,xo(:,2),'c');
hold on;
tic ;
[t,x] = eul(@oscillator, 0, 20, [2;4], timestep_euler);
time_euler = toc;
plot(t,x(:,1),'r',t,x(:,2),'b');
title('Euler');
xlabel('time');
ylabel('red: x_1, blue: x_2');


% Runge
figure;
plot(to,xo(:,1),'m',to,xo(:,2),'c');
hold on;
tic;
[t,x] = runge(@oscillator, 0, 20, [2;4], timestep_runge);
time_runge = toc;
plot(t,x(:,1),'r',t,x(:,2),'b');
title('Runge');
xlabel('time');
ylabel('red: x_1, blue: x_2');


% Runge_Kutta
figure;
plot(to,xo(:,1),'m',to,xo(:,2),'c');
hold on;
tic;
[t,x] = rk_classical(@oscillator, 0, 20, [2;4], timestep_rk);
time_RK = toc;
plot(t,x(:,1),'r',t,x(:,2),'b');
title('Runge-Kutta');
xlabel('time');
ylabel('red: x_1, blue: x_2');

disp('               CPU Time');
fprintf('Euler:       %.8fs\n', time_euler);
fprintf('Runge:       %.8fs\n', time_runge);
fprintf('Runge-Kutta: %.8fs\n', time_RK);
