% initialize grid size, CFL number and time 
L=20; 
n=200; 
x2=linspace(-L/2,L/2,n+1); x=x2(1:n); 
dx=x(2)-x(1); 
dt=0.2;  % try different dt's to study stability 
CFL=dt/dx 
Time=4; 
time_steps=Time/dt; 
t=0:dt:Time; 

% initial conditions (gaussian) 
u0=exp(-x.^2)'; 
usol(:,1)=u0; 
u1=exp(-(x+dt).^2)'; 
usol(:,2)=u1;

% sparse matrix for derivative term 
e1=ones(n,1); 
A=spdiags([-e1 e1],[-1 1],n,n); 

% periodic boundary conditions 
A(1,n)=-1; A(n,1)=1; 

% leap frog (2,2) or euler iteration scheme 
for j=1:time_steps-1 
  % forward euler (start at u1) 
   u2 = u1 + 0.5*CFL*A*u1; 
   u1 = u2; 

   % leap-frog (2,2) (needs u0 and u1) 
   %u2 = u0 + CFL*A*u1; % leap frog (2,2) 
   %u0 = u1; u1 = u2; % leap frog (2,2) 
   
  % save solution at time slice in matrix 
  usol(:,j+2)=u2; 
end

% plot the data
waterfall(x,t,usol'); 
map=[0 0 0]; 
colormap(map); 

% set x and y limits and fontsize 
set(gca,'Xlim',[-L/2 L/2],'Xtick',[-L/2 0 L/2], ... 
'FontSize',[20]); 
set(gca,'Ylim',[0 Time],'Ytick',[0 Time/2 Time], ... 
'FontSize',[20]); 
view(25,40) 

% set axis labels and fonts 
xl=xlabel('x'); 
yl=ylabel('t'); 
zl=zlabel('u'); 
set(xl,'FontSize',[20]); 
set(yl,'FontSize',[20]); 
set(zl,'FontSize',[20]); 

% print jpeg at screen resolution 
print -djpeg -r0 fig.jpg