An implicit algorithm for capturing sharp fluid interfaces in the volume of fluid advection method

The Volume of Fluid (VOF) method is one of the most effective methods employed in the simulation of two fluid flows with interfaces where density and viscosity change abruptly. These interfaces are represented implicitly by the values of a colour function which is a volume fraction of one of the fluids. The advantage of the method is its ability to deal with arbitrarily shaped interfaces and to cope with large deformations, as well as interface rupture and coalescence in a natural way. In comparison to a level set method, the mass is rigorously conserved in VOF, provided the discretisation is conservative, but one of the main difficulties is advecting the interface without diffusing, dispersing, or wrinkling it. This can either be performed algebraically, in schemes such as CICSAM or geometrically, in schemes such as PLIC. In the present paper, an algebraic advection scheme for the interface is presented, which is designed for the implicit time advancing algorithm. Analogous to CICSAM, the new scheme switches smoothly between ULTIMATE-QUICK and the upper bound of the universal limiter, depending on the angle between the interface and the flow direction. Four cases are tested with the present scheme: (i) solid body rotation; (ii) circle in a shear flow; (iii) dam-break and (iv) Rayleigh-Taylor instability. In the first two test cases, prescribed velocity fields are used, thereby allowing the effectiveness of the scheme in advecting the colour function only to be assessed. The scheme is found to outperform six other methods used for comparison in both studies. In solid body rotation simulations a fractional error of 0.19% is calculated in comparison to the next best recorded error of 1.1%. Similarly, in the longest shear flow simulation, a fractional error of 1.2% is calculated in comparison to the next best recorded error of 3.9%. In the final two test cases the advection equation for the colour function is coupled to the Navier-Stokes equations. In dam-break simulations it is found that the resulting solution effectively captures the trends displayed in experimental data for the advancing water front and the residual height of the liquid column against time. Qualitative results obtained for the Rayleigh-Taylor instability modelling in test case four are found to compare favourably to previous numerical simulations of the same phenomenon