Add RK6 with residual viscosity update and reduce computational effort of spatial scheme - Add RK6 with residual updates - Change the D2 operator for upwind schemes to one less computationally expensive.

% Get the index of the points p within the tall array of points ps
function [I, ok] = pointIndex(p, ps)
    [ok, I] = ismember(p,  ps, 'rows');
end