% [u,u_true,x,t] = lax_friedrichs(t_0,t_f,M,N)
% 
% solves the advection equation u_t + c u_x = 0 on [0,2 pi]
%
% For initial data cos(x), we assume periodic initial data.  For
% step-function initial, we assume u_x = 0 at the ends.
%
% we use the Lax Friedrichs method.

function [u,u_true,x,t] = lax_friedrichs(t_0,t_f,M,N)

c = 1/2;

% define the mesh in space
dx = 2*pi/(N-1);
x = 0:dx:2*pi;
x = x';

% define the mesh in time
dt = (t_f-t_0)/M;
t = t_0:dt:t_f;

for i=1:N
  u(i,1) = f(x(i));
end 
u_true(:,1)=u(:,1);

% For Lax-Friedrichs, we have
%   u_new(i) = 1/2(u_old(i+1)+u_old(i-1))
%                - c*dt/(2*dx)*(u_old(i+1)-u_old(i-1)) 
% 
for j=2:M+1
  u_old = u(:,j-1);

  for i = 2:(N-1);
% first do the interior points: 
    u_new(i)=(u_old(i+1)+u_old(i-1))/2-c*dt/(2*dx)*(u_old(i+1)-u_old(i-1)); 
  end
% now do the endpoints u_new(1) and u_new(N).
% for periodic boundary conditions:
%  u_new(1)=(u_old(2)+u_old(N-1))/2-c*dt/(2*dx)*(u_old(2)-u_old(N-1)); 
%  u_new(N)=u_new(1);
% for the step-function:
  u_new(1)=u_old(1);
  u_new(N)=u_old(N);  

% now update the solution:
  u(:,j) = u_new';

% define the true solution: 
  t_diff = t(j)-t_0;
  for k=1:N
    u_true(k,j) = f(x(k)-c*t_diff);
  end		
end