advection.m 858 Bytes
function E = advection(E,c,dx,dt)
%ADVECTION  is a 1d advection code Lax Wendroff scheme with Superbee flux limiter
%
%   E = advection(E,c,dx,dt)
%
% INPUTS:
%   E  = vector of the thing to be advected;
%   c  = scalar speed;
%   dx = the spatial resolution;
%   dt = the temporal resolution;
% 
% OUTPUT:
%   E  = vector of the thing to be advected - after advection;

if ~isequal(size(E),size(c))
    c = c';
end

f  = flux(E,c,dx,dt);
E  = E-dt*diffl(f)/dx;

function f = flux(E,c,h,dt)
   % Lax-Wendroff with Superbee flux limiting;
   theta = diffl(E)./(diffr(E)+3e-14);
   phi   = limiter(theta);
   f     = c.*E+c/2.*(1-c*dt/h).*diffr(E).*phi;
end


function y = diffl(x)
y = [x(1);diff(x)];
end

function y = diffr(x)
y = [diff(x);-x(end)];
end

 % Superbee
function phi = limiter(r)
    phi = max(0,max(min(1,2*r),min(r,2)));              
end
end