clear
Nx = 80;   
tmax = 6;
dt=tmax/100;
xmin=-1;
xmax=10;
irn=6;
ymin=-0.5; 
ymax=1.5;

switch irn
    case 1    % Simple shock initially
        v = @(x,t) 1.*(x<=t/2)+0.*(x>t/2);
    case 2    % Rarefaction
        v = @(x,t) 0.*(x<=0) + ...
            x./t*(t>0&x<t & x>0) + ...
            1.*(x>t);
    case 3    % Rarefaction plus shock
        v = @(x,t) 0.*(x<=0) + ...
            1.*(x>t & x<=t/3+2)+...
            real(sqrt(x./t)).*(t>0&t<3&x<t & x<=t/3+2) + ...
            real(sqrt(x./t)).*(t>3&x<t & x<=3^1.5*t.^(1/3)) + ...
            0.*(t<3 & x>t/3+2) + 0.*(t>3 & x>3^1.5*t.^(1/3));
    case 4   % Trapezoidal IC
        v = @(x,t) 0*x + ...
            (x-1)./(1+t).*(x>=1).*(x<=2+t)   ...
            + 1.*(x>=2+t).*(x<=min(3+t,(7+t)/2))...
            + (4-x)./(1-t).*(x>=3+t).*(x<=4);
    case 5 % Red Light Problem
        v1=2;
        u0=1;
        u1=3;
        xmin=-10;
        xmax=2;
        ymax=5;
        v = @(x,t) u0*(x<-v1*u0*t/u1)+u1*(x>=-v1*u0*t/u1);
    case 6 % Green Light Problem
        v1=2;
        u0=0;
        u1=1;
        xmin=-15;
        xmax=15;
        ymax=2;
        v = @(x,t) u1*(x<=-v1*t)+...
            u1/2*(1-1/v1*x/t).*(x>-v1*t & x<v1*t)+...
            0*(x>=-v1*t);
end

t=0;
figure(1);
x = linspace(xmin,xmax,400);
pth  = plot(x,v(x,0),'k-','LineWidth',2);
hold on
xlim([xmin xmax]);
ylim([ymin,ymax]);
ptit = title(strcat('t = ',num2str(t)));
xlabel('x');
ylabel('u')

% iteration
while (t+dt<tmax)
    t = t + dt;
    set(pth,'YData',v(x,t));
    set(ptit,'String',strcat('t = ',num2str(t)));
    drawnow;
    pause(.05)
end