  Solving the Navier-Stokes equations when the density field contains a large
sharp discontinuity---such as a water/air free surface---is numerically
challenging. Convective instabilities cause Gibbs oscillations which quickly
destroy the solution. We investigate the use of slope limiters for the velocity
field to overcome this problem in a way that does not compromise on the mass
conservation properties. The equations are discretised using the interior
penalty discontinuous Galerkin finite element method that is divergence free to
machine precision.
  A slope limiter made specifically for exactly divergence free (solenoidal)
fields is presented and used to illustrated the difficulties in obtaining
convectively stable fields that are also exactly solenoidal. The lessons
learned from this are applied in constructing a simpler method based on the use
of an existing scalar slope limiter applied to each velocity component.
  We show by numerical examples how both presented slope limiting methods are
vastly superior to the naive non-limited method. The methods can solve
difficult two-phase problems with high density ratios and high Reynolds
numbers---typical for marine and offshore water/air simulations---in a way that
conserves mass and stops unbounded energy growth caused by the Gibbs
phenomenon.
