
N = 40;
h = 1/N;
x = 0:h:1;

u = (x<0.75).*(x>0.25);
w = u;
k = 0.5*h;
T = 0.5;
N = T / k;

for n = 1: N,
    up = [u(2:end),u(1)];
    um = [u(end),u(1: end-1)];
    u = u - 0.5*(u-um);
    
    wp = [w(2:end),w(1)];
    wm = [w(end),w(1:end-1)];
    w = 0.5*(wp+wm) - 0.5*0.5*(wp-wm);

    plot(x,u,'g',x,w,'b')
    pause(0.5)
end