% Direction field m-file (script) 
% Source: Bäckström: Praktisk matematik med Matlab
clear variables; hold off;
xmax=4; ymax=3; delx=xmax/30; dely=ymax/30;   % *** problem dependent
[X,Y]=meshgrid(0:delx:xmax,0:dely:ymax);
  DX=ones(size(X));  % Unit x-xomponents of arrows
  DY=fode1(X,Y);      % Slopes, *** Change the name 'fode' to your m-file name
  L=sqrt(DX.^2+DY.^2); % Initial lengths
  DX1=DX./L; DY1=DY./L; % Unit length for all arrows
%figure(1); quiver(X,Y,DX1,DY1,'.'); grid on;  % The period ('.') suppresses
                                               % arrowheads
figure(1); quiver(X,Y,DX1,DY1); grid on;
  axis([0 xmax 0 ymax]); title('Direction field');
