  function [out]=fprint(ncase,neqn,t,y)
%
% Function fprint displays the numerical and
% tabulated solutions to the ODE with erf as
% a solution
%
% Define global variables
  global nsteps;
%
% Return current value of independent variable
% (Matlab requires at least one return argument)
  out=t;
%
% Print a heading for the solution at t = 0 
  if(t<=0.0)
%
%   Label for ODE integrator
%
%   Fixed step modified Euler 
    if(ncase==1)
      fprintf('\n\n euler2a integrator\n\n');
% 
%   Variable step modified Euler
    elseif(ncase==2)
      fprintf('\n\n euler2b integrator\n\n');
%
%   Fixed step classical fourth order RK
    elseif(ncase==3)
      fprintf('\n\n rkc4a integrator\n\n');
%
%   Variable step classical fourth order RK 
    elseif(ncase==4)
      fprintf('\n\n rkc4b integrator\n\n');
%
%   Fixed step RK Fehlberg 45 
    elseif(ncase==5)		
      fprintf('\n\n rkf45a integrator\n\n');
%		
%   Variable step RK Fehlberg 45 
    elseif(ncase==6)
      fprintf('\n\n rkf45b integrator\n\n');
    end
%
% Heading 
  fprintf(' ncase = %2d   nsteps = %3d\n\n',ncase,nsteps);
  fprintf('   t      y(num)      y(tab)       diff\n');
%
% End of t = 0 heading
  end
%
% Tabulated solution
  if(t<0.001) ytab=0.0; end
  if((t>0.0999)&(t<0.1001)) ytab=0.112463; end
  if((t>0.1999)&(t<0.2001)) ytab=0.222703; end
  if((t>0.2999)&(t<0.3001)) ytab=0.328627; end
  if((t>0.3999)&(t<0.4001)) ytab=0.428392; end
  if((t>0.4999)&(t<0.5001)) ytab=0.520500; end
%
% Numerical and tabulated solutions and their difference
% as output
  y1=(2.0/sqrt(pi))*y(1);
  diff=y1-ytab;
%	
% Display the numerical and exact solutions, and their difference 
  fprintf('%5.2f %11.6f %11.6f %11.6f\n',t,y1,ytab,diff);

