%====================================================================
% Initial value problem for a single 1st-order ODE
% Call three methods: Euler, Predictor-Corrector, Runge-Kutta 4th order
%====================================================================
clear;
clc;
%open a file if needed
% fileID = fopen('results.txt','w');
% initial conditions
xi = 1.00;
ti = 0.00;
% "time" step and max time
dt =   0.09;
tmax = 10.0;

%select function
dx = @f;

%preparing to plot
i    =  1;
t (i)= ti;
x1(i) = xi;
x2(i) = xi;
x3(i) = xi;

% print the header and initial conditions
% fprintf(' Solving 1st order single ODE \n')
% 
% fprintf('     t         x(t)   \n')
% formatSpec1 = '%10.5f %10.5f \n';
% fprintf(formatSpec1,ti,xi)
%fprintf(fileID,formatSpec1,ti,xi)


xi1 = xi;
xi3 = xi;


while(ti < tmax)
  tf = ti + dt;
  i = i+1;
  
  [xf1] = euler1(dx,ti,xi1,tf);
  x1(i) = xf1;
  xi1 = xf1;  
  
  [xf3] = RK4D11(dx,ti,xi3,tf);
  x3(i) = xf3;
  xi3 = xf3;
  
  t(i)  = tf;
  ti = tf;  
    
    
%   if(key==2) 
%       [xf] = euler2(dx,ti,xi,tf);
%   end
%   if(key==3) 
%       [xf] = RK4D11(dx,ti,xi,tf);
%   end
%  xtrue = exp(-1.0*tf);
%   fprintf(formatSpec1,tf,xf); 
%  fprintf(fileID,'%10.5f %10.5f %10.5f\n',tf,xf,xtrue);
% to plot

  
end %while
plot(t,x1,'b',t,x3,'r')


function f=f(t,x)
%========
% function for dx/dt = f(t,x)
%========
%    f = (-1.0)*x+t*0.0;  
   f = (-1.0)*x*cos(t)+0.5*sin(2*x);
% f = t*x*x;
end

function [xf] = euler1(dx,ti,xi,tf)
%==================================================
% Solution of a single 1st order ODE dx/dt=f(x,t)
% Method:  Simple Euler
% Alex G. October 2020
%--------------------------------------------------
% input ...
% dx(t,x)- function dx/dt (supplied by a user)
% ti - initial time
% xi  - initial position
% tf  - time for a solution
% output ...
% xf  - solution at point tf
%==================================================
xf = xi + dx(ti,xi)*(tf-ti);
end

function [xf] = euler2(dx,ti,xi,tf)
%==================================================
% Solution of a single 1st order ODE dx/dt=f(x,t)
% Method:  Euler: predictor-corrector
% Alex G. October 2020
%--------------------------------------------------
% input ...
% dx(t,x)- function dx/dt (supplied by a user)
% ti - initial time
% xi  - initial position
% tf  - time for a solution
% output ...
% xf  - solution at point tf
%==================================================
xf = xi + dx(ti,xi)*(tf-ti);
xf = xi + (dx(ti,xi)+dx(tf,xf))*0.5*(tf-ti);
end

function [xf] = RK4D11(dx,ti,xi,tf)
%==================================================
% Solution of a single 1st order ODE dx/dt=f(x,t)
% Method:  Runge-Kutta 4th order
% Alex G. October 2020
%--------------------------------------------------
% input ...
% dx(t,x)- function dx/dt (supplied by a user)
% ti - initial time
% xi  - initial position
% tf  - time for a solution
% output ...
% xf  - solution at point tf
%==================================================
h  = tf-ti;

k1 = h*dx(ti,xi);
k2 = h*dx(ti+h/2.0,xi+k1/2.0);
k3 = h*dx(ti+h/2.0,xi+k2/2.0);
k4 = h*dx(ti+h,xi+k3);

xf = xi + (k1 + 2.0*(k2+k3) + k4)/6.0;
end
