function [u]=heat(N,M)
%
h = 1/(N+1);
tau=0.1/M;
I=speye(N,N);diag = 2*I;subd=-sparse(2:N,1:N-1,1,N,N);
A=(N+1)^2*(diag+subd+subd');
u=sparse(N,1);
for i=1:N
  u(i)=sin(pi*i*h);
end
mat=speye(N)+tau*A;
err=sparse(N,1);
for n=1:M
 u=mat\u;
end
for i=1:N
   err(i)=u(i)-sin(pi*i*h)*exp(-pi*pi*n*tau);
end
fprintf('err inf=%e \n',norm(err,inf))
fprintf('err 2=%e \n',norm(err,2))