    % Iterate through the time steps and plot the saturation profiles along the
    % way.
    t = 0;
    msSat = [];
    refSat = [];

    %if doPlotSim, figure, end
    
    while t < T,
        % Solve transport equations using the transport solver
        msSol  = tsolve(msSol , dT);
        refSol = tsolve(refSol, dT);
        
        % Update the pressure based on the new saturation
        %msSol    = psolve(msSol);
        %refSol   = r_psolve(refSol);
        
        % store saturation at source
