% ddwave3.m
% Inviscid Burgers equation
% First-order up-winding finite volume method.

CAR_COUNT = 20;
MIN_DIST = 0.1;
MAX_VELOC = 2;
MAX_RHO = 1/MIN_DIST;


L  = 20;
N  = 200;
x  = linspace(0,L*((N-1)/N),N);

k  = 2*pi/5;
%un = .4*sin(k*x);
%un  =.5* exp(-(x-L/2).^2);
for i = 1:length(x)
    if (i < CAR_COUNT * MIN_DIST)
         un(i) = -MAX_VELOC + 0.1*rand(1,1);
    else
         un(i) = MAX_VELOC;
    end
end
dx = x(2)-x(1);
nu = 0.5;
dt = nu*dx/10;
n  = 0;
u  = un;
us = un;
gu = un;
phi= un;
axn = [0 L -0.5 MAX_VELOC];
%plot(x,u), axis(axn), drawnow
plot(x,(MAX_RHO/2*(1-(u/MAX_VELOC)))), axis(axn),drawnow

while (n*dt < 20)
    % Up winding first order
    for i=1:N-1
        chat = 0.5*(u(i)+u(i+1));
        if chat > 0 
            F(i)=0.5*(u(i))^2;
        else
            F(i)=0.5*(u(i+1))^2;
        end
        if ((u(i) < 0) & (u(i+1) > 0))
            F(i) =-0.5*(u(i))*(u(i+1));
        end
    end    
    chat = 0.5*(u(N)+u(1));
    if chat > 0 
            F(N)=0.5*(u(N))^2;
    else
            F(N)=0.5*(u(1))^2;
    end
    if ((u(N) < 0) & (u(1) > 0))
            F(N) =-0.5*(u(N))*(u(1));
    end
    
     
 
    u  = u - (dt/dx)*(F - circshift(F,[1 0]));   
    n=n+1;
    %plot(x,u), axis(axn),drawnow,title(n*dt)
    plot(x,(MAX_RHO/2.*(1-(u./MAX_VELOC)))), axis(axn),drawnow,title(n*dt)
   %n*dt
    %pause(.05)
end
    %plot(x,u), axis(axn),drawnow
    plot(x,(MAX_RHO/2*(1-(u/MAX_VELOC)))), axis(axn),drawnow
    max(u)