%
% Clear previous files
  clear all
  clc
%
% Parameters, variable shared by other routines
  global alpha lambda icase ncall
%
% Initial condition
  for(icase=1:3)
    if(icase==1)lambda(1)=1   ;end
    if(icase==2)lambda(2)=0.75;end 
    if(icase==3)lambda(3)=1.25;end 
    V0=1;alpha=1;ncall=0;
%
%   t scale
    tout=[0:0.5:10]';nout=21;
%
%   ODE integration
    [t,V]=ode15s(@ode_1,tout,V0); 
%
%   Display heading
    fprintf('\n icase = %2d',icase);
%
%   Display calls to tumor_1
    fprintf('\n\n Calls to ode_1 = %3d\n',ncall);    
    fprintf('\n       t     V(t)(num)   V(t)(anal)    diff\n');
%
%   Display numerical output
    for(it=1:nout)
       Va(it)=V0*exp(lambda(icase)/alpha*(1-exp(-alpha*tout(it))));
      Vplot(it,icase)=V(it,1);
     Vaplot(it,icase)= Va(it);
       diff(it,icase)=(Vplot(it,icase)-Vaplot(it,icase))...
                      /Vplot(it,icase);
      fprintf('%10.2f%12.4f%12.4f%12.6f\n',...
      tout(it),Vplot(it,icase),Vaplot(it,icase),diff(it,icase));
    end
  end
%
% Plot numerical, analytical solutions
  figure(1)
  plot(tout,Vplot,'-',tout,Vaplot,'o');
  xlabel('t');ylabel('V(t)');
  title('V(t), solid - num, points - anal');
%
% Plot difference of solutions
  figure(2)
  plot(tout,diff,'-');
  xlabel('t');ylabel('diff');
  title('diff = num - anal (rel), V(t)');

