% m36_heat.m - explicit solution of heat eq  u_t = u_xx, u(-1)=u(1)=0.
%         Note use of sparse matrix (not necessary, but good fun).

% Grid and initial data:
  h = .025;                    % space step
  k = .4*h^2;                  % time step (try .4 -> .51!)
  x = (-1+h:h:1-h)';                 
  u = abs(x)<.3;               % initial data: square wave

% Set-up for plot:
  hold off, shg
  plt = plot(x,u,'linewidth',4);
  axis([-1 1 -.01 1.01]), grid

% Sparse matrix to execute finite difference operation:
  L = length(x);
  a = (1-2*k/h^2);
  b = k/h^2;
  main = a*sparse(ones(L,1));
  off  = b*sparse(ones(L-1,1));
  A = diag(main) + diag(off,1) + diag(off,-1);
 %A(1,end) = b; A(end,1) = b;  % comment this in for periodic BCs
  
% Time-stepping:
  while 1
    u = A*u;
    %pause                      % comment this out for a movie
    set(plt,'ydata',u)
    drawnow
  end