function burgmol(c,a,N,T,ic,solver);
% BURGMOL solves the one-dimensional nonlinear pde
%       u_t = c u u_x + a u_xx
%    by the method of lines.  (Burgers equation.)
%    User chooses constants c and a, also N (dx = 1/N) and T. 
%    Initial condition in ic, a row vector of length N-1.
%    **Boundary conditions Dirichlet: u(0,t)=1, u(1,t)=0.**
%    Solves system of N-1 coupled nonlinear ODEs using Matlab's 
%      "ode45" (if "solver"=0) or "ode15s" (if "solver"=1).
%    Passes subfunction "burgf" to solver.
% Try "burgmol(1,.1,20,.2,[ones(1,10) zeros(1,9)],0)",
%     "burgmol(1,.03,20,.2,[ones(1,10) zeros(1,9)],0)",
%     "burgmol(1,.01,20,.2,[ones(1,10) zeros(1,9)],0)",
%     "burgmol(1,0,20,.2,[ones(1,10) zeros(1,9)],0)",
%     "burgmol(1,0,80,.2,[ones(1,40) zeros(1,39)],0)",
% to see "rarefaction wave".
% Try "burgmol(1,.3,80,.2,[ones(1,40) zeros(1,39)],0)"
%     "burgmol(1,.3,80,.2,[ones(1,40) zeros(1,39)],1)"
% to compare a nonstiff method to a stiff method for a relatively
% diffusive problem.
% (Ed Bueler 3/25/2002)

global NN cc aa dx c1 c2 % needed to communicate with function molf
NN=N; cc=c; aa=a; % eliminates silly error message
dx=1/N; x=0:dx:1;
c1=c/(2*dx); c2=a/(dx*dx);
if not(isequal(size(ic),[1 N-1])), error('Init cond wrong size.'); end;

% use ODE solver
if solver==0, [tt,vv]=ode45('burgf',[0 T],ic');
elseif solver==1, [tt,vv]=ode15s('burgf',[0 T],ic');
else, error('Choice of solver invalid.'); end;

% to see solution as a system of ODEs:
% figure(1), plot(tt,vv), figure(2)

% augment vv to reflect boundary conditions
M=size(tt,1); vv=[ones(M,1) vv zeros(M,1)];
%display it
mesh(x,tt,vv);

xlabel('x axis'); ylabel('t axis'); zlabel('u axis');
title(['Solution  to Burgers eqn for c=' num2str(c) ', a='...
      num2str(a) ', and N=' num2str(N)]);
view(37.5,30);