

##function derivative (t,y)
##     return exp(t)*y^2+exp(3);
##end

close all

derivative = @(t,y) tan(y)*cos(y*t);

t_k = 2;
y_k = 4;
t_end = 5;
n = 50;
h = (t_end-t_k)/n;

method = "euler"

for k = 0:n

     fprintf("k: %0.0f -> t_k: %0.5f & y_k: %0.5f\n",k,t_k,y_k);

     plot(t_k,y_k,'*k');
     hold on;

     if strcmp(method,"euler")
          K1 = h*derivative(t_k,y_k);
          y_kp1 = y_k + K1;
     elseif strcmp(method,"rk4")
          K1 = h*derivative(t_k,y_k);
          K2 = h*derivative(t_k+h/2,y_k+K1/2);
          K3 = h*derivative(t_k+h/2,y_k+K2/2);
          K4 = h*derivative(t_k+h,y_k+K3);
          y_kp1 = y_k + 1/6*(K1+2*K2+2*K3+K4);
     else
          disp('Invalid method!')
          return;
     end

     y_k = y_kp1;

     t_kp1 = t_k+h;
     t_k = t_kp1;
end

