% First Order Ordinary Differential Equations - Initial Value Problems
% Consider the equation
% \frac{\textrm{dy}}{\textrm{dt}}=2\textrm{yt}
% with the inital condition y(0)=1. Matlab has a built-in function that will solve this. We create a function that defines the right hand side of the equation and then call it from a script.
% The calling routine is as follows:

trange=[0 2];
init=1;
[t,y]=ode45(@rkfile, trange, init);
yexact=exp(t.^2);
plot(t,y,t,yexact,'+', 'LineWidth', 2)
xlabel('t')
ylabel('y')
set(gca,'FontSize',20)
errorvector=((y-yexact)./yexact);
figure,plot(t,errorvector, 'LineWidth', 2)
xlabel('t')
ylabel('relative error')

%
% The trange vector indicates that we will solve the equation on the range 0<t<2. The init parameter sets the initial value of y [y(0)]. The ode45 routine is the built-in routine that we will use to solve the equation. It does so by taking small time steps, adjusting the step as it goes to achieve a desired relative and absolute tolerance of 1e-6. These can be changed, as will be shown below. As called, the script stores the output of ode45 in two vectors: t holds a list of times and y holds the results for y at each of those times. Hence, these are easily plotted, as shown.
% Note that the vector yexact is the analytical solution of the equation and is here to ensure that we haven't made any mistakes. It's always a good idea to begin with an equation for which the answer is known.
% If we want to repeat the solution using a smaller tolerance, to achieve improved accuracy, we can set both the relative and absolute tolerances as follows:

options=odeset('RelTol',1e-7, 'AbsTol',1e-7);
[t,y]=ode45(@rkfile, trange, init,options);
yexact=exp(t.^2);
figure,plot(t,y,t,yexact,'+')
xlabel('t')
ylabel('y')
set(gca,'FontSize',20)
errorvector=((y-yexact)./yexact);
figure,plot(t,errorvector)
  

%
% As called, ode45 expects to find the function definition in a user-defined function called rkfile. For this problem, rkfile is as follows:

function dydt=rkfile(t,y)
  dydt=2*y*t;
end

% Note that this function definition has to be at the end of the script.
% © James Blanchard, University of Wisconsin-Madison, blanchard@engr.wisc.edu