% script to solve easy example with Bessel function solution:
% y' = y^2 + t^2, y(0)=0, 0 <= t <= 1.

close all
clear all

% this IVP can be solved in one line with a call to dsolve
y = dsolve('Dy = y^2+t^2','y(0)=0')
ezplot(y,[0,1])
pause

% now with the right-hand side defined in the function f,
% we solve the same problem numerically

[t1,y] = ode45(@f,[0,1],[0]);
plot(t1,y,'-o')
pause
% note that the output array t corresponds (almost) to the times
% that the solver used to find the solution
% in the case of ode45, each step is actually "refined" by a factor of 4:
% steps are often so large that the make for poor-quality plots
% so in this case 4 points per step are output
% 'Refine' is an option that can be set by the user

% you can also specify what times you would like output
tspan = linspace(0,1,10)
[t2,y] = ode45(@f,tspan,[0]);
plot(t2,y,'rx-')
hold on
pause

options = odeset('Refine',1);
[t3,y] = ode45(@f,[0,1],[0],options);
plot(t3,y,'c*-.')
pause

% it is possible to request a solution structure as output
% recall that ode solvers produce a continuous approximate solution S(t)
% the solution structure has all the information required to produce an
% approximation at any point t \in [t_0,t_f]
% with deval, you can effectively eveluate the solution of the IVP as you
% would an ordinary function

sol = ode45(@f,[0,1],[0]);

% let's see the value of the solution at t=0.5
yHalf = deval(0.5,sol)

% we will see how to access other fields of information in solution
% structures later when we solve BVPs
