clf
hold off
disp('In order to work with m-files you may have to ')
disp('set a new search path.                        ')
pause
help path
pause
disp('The easiest built-in solver to use is ode23   ')
pause
help ode23
pause
disp('And there is even more reading to do..........')
pause
help odefile
pause
disp('Here is a simple example with a predefined ODE')
pause
[t,y]=ode23('RC',[0 50],0)
pause
disp('e.g. plot the solution                         ')
plot(t,y)
pause
disp('Getting several solution curves at the same time')
disp('is just as easy')
pause
[t,y]=ode23('RC',[0 50],[0 2 4 6 8 10])
pause
plot(t,y)
pause
disp('higher order DEs -- recall the bungee -- are first')
disp('rewritten as a system of first order DEs -- take a')
disp('look as this example of a damped forced oscillator.')
disp('y is a column vector containing position and speed')
pause
[t,y]=ode23('forced',[0,20],[2;1])
pause
plot(t,y)
pause
plot(t,y(:,1),'ro',t,y(:,2),'b:')
pause
disp('overlay the slope/vector field and the solution curves')
pause
help quiver
pause
disp('First define the grid points where to plot arrows')
pause
[tt,yy]=meshgrid(0:5:30,0:1:12)
pause
disp('Next define the slopes representing the diff equn')
disp('for slope fields always the t-components are always')
disp('equal to one; the y-components equal (dy/dt).')
pause
vt=ones(size(tt))
vy=-(yy-7)/10;
pause
disp('With basepoints and arrows defined quiver is ready:')
pause
quiver(tt,yy,vt,vy)
pause
disp('To OVERLAY the solution curves, i.e. keep the field')
pause
hold on
pause
disp('recalculate the solution curves')
pause
[t,y]=ode23('RC',[0 30],[0 2 4 6 8 10])
pause
disp('and overlay their graphs')
plot(t,y)
pause
clf
disp('For second order scalar DEs that do not depend on t')
disp('one may consider the phase portrait in the (y1,y2)-')
disp('plane, that is in the (y,y''), or position/speed plane')
pause
[y1,y2]=meshgrid(-4:0.5:4,-3:0.4:3);
v1=y2;
v2=-3*sin(y1)-0.2*y2;
quiver(y1,y2,v1,v2)
pause
[t,y]=ode23('damped',[0 10],[-2,0])
pause
hold on
plot(y(:,1),y(:,2),'r-')
pause