An accurate adaptive solver for surface-tension-driven interfacial flows

A method combining an adaptive quad/octree spatial discretisation, geometrical Volume-Of-Fluid interface representation, balanced-force continuum-surface-force surface-tension formulation and height-function curvature estimation is presented. The extension of these methods to the quad/octree discretisation allows adaptive variable resolution along the interface and is described in detail. The method is shown to recover exact equilibrium (to machine accuracy) between surface-tension and pressure gradient in the case of a stationary droplet, irrespective of viscosity and spatial resolution. Accurate solutions are obtained for the classical test case of capillary wave oscillations. An application to the capillary breakup of a jet of water in air further illustrates the accuracy and efficiency of the method. The source code of the implementation is freely available as part of the Gerris flow solver.

Related papers:

title: Reconstructing Volume Tracking
abstract: A new algorithm for the volume tracking of interfaces in two dimensions is presented. The algorithm is based upon a well-defined, second-order geometric solution of a volume evolution equation. The method utilizes local discrete material volume and velocity data to track interfaces of arbitrarily complex topology. A linearity-preserving, piecewise linear interface geometry approximation ensures that solutions generated retain second-order spatial accuracy. Second-order temporal accuracy is achieved by virtue of a multidimensional unsplit time integration scheme. We detail our geometrically based solution method, in which material volume fluxes are computed systematically with a set of simple geometric tasks. We then interrogate the method by testing its ability to track interfaces through large, controlled topology changes, whereby an initially simple interface configuration is subjected to vortical flows. Numerical results for these strenuous test problems provide evidence for the algorithm's improved solution quality and accuracy.

title: DIRECT NUMERICAL SIMULATION OF FREE-SURFACE AND INTERFACIAL FLOW
abstract: The numerical simulation of flows with interfaces and free-surface flows is a vast topic, with applications to domains as varied as environment, geophysics, engineering, and fundamental physics. In engineering, as well as in other disciplines, the study of liquid-gas interfaces is important in combustion problems with liquid and gas reagents. The formation of droplet clouds or sprays that subsequently burn in combustion chambers originates in interfacial instabilities, such as the Kelvin-Helmholtz instability. What can numerical simulations do to improve our understanding of these phenomena? The limitations of numerical techniques make it impossible to consider more than a few droplets or bubbles. They also force us to stay at low Reynolds or Weber numbers, which prevent us from finding a direct solution to the breakup problem. However, these methods are potentially important. First, the continuous improvement of computational power (or, what amounts to the same, the drop in megaflop price) continuously extends the range of affordable problems. Second, and more importantly, the phenomena we consider often happen on scales of space and time where experimental visualization is difficult or impossible. In such cases, numerical simulation may be a useful prod to the intuition of the physicist, the engineer, or the mathematician. A typical example of interfacial flow is the collision between two liquid droplets. Finding the flow involves the study not only of hydrodynamic fields in the air and water phases but also of the air-water interface. This latter part

title: A level set approach for computing solutions to incompressible two- phase flow II
abstract: A level set method for capturing the interface between two fluids is combined with a variable density projection method to allow for computation of two-phase flow where the interface can merge/break and the flow can have a high Reynolds number. A distance function formulation of the level set method enables one to compute flows with large density ratios (1000/1) and flows that are surface tension driven; with no emotional involvement. Recent work has improved the accuracy of the distance function formulation and the accuracy of the advection scheme. We compute flows involving air bubbles and water drops, to name a few. We validate our code against experiments and theory.

title: An adaptive level set approach for incompressible two-phase flows
abstract: In Sussman, Smereka and Osher, a numerical method using the level set approach was formulated for solving incompressible two-phase flow with surface tension. In the level set approach, the interface is represented as the zero level set of a smooth function; this has the effect of replacing the advection of density, which has steep gradients at the interface, with the advection of the level set function, which is smooth. In addition, the interface can merge or break up with no special treatment. The authors maintain the level set function as the signed distance from the interface in order to robustly compute flows with high density ratios and stiff surface tension effects. In this work, they couple the level set scheme to an adaptive projection method for the incompressible Navier-Stokes equations, in order to achieve higher resolution of the interface with a minimum of additional expense. They present two-dimensional axisymmetric and fully three-dimensional results of air bubble and water drop computations.

title: Level set calculations of the evolution of boundaries on a dynamically adaptive grid
abstract: We present a general algorithm for the study of the evolution of interfaces in growth processes based on the level set method, using the narrow band and the fast marching approximations, applied on a dynamically adaptive grid. A novel construction of a dynamical adaptive grid is presented. One important feature is that we establish a controllable finite-width region of the highest resolution straddling the physically important boundary. This high-resolution stripe provides a better description of the important physical variables. An efficient and robust coupling of a hierarchical space description structure with the state-of-the-art level set method is established.

title: Numerical analysis of moving interfaces using a level set method coupled with adaptive mesh refinement
abstract: A novel numerical scheme is developed by coupling the level set method with the adaptive mesh refinement in order to analyse moving interfaces economically and accurately. The finite element method (FEM) is used to discretize the governing equations with the generalized simplified marker and cell (GSMAC) scheme, and the cubic interpolated pseudo-particle (CIP) method is applied to the reinitialization of the level set function. The present adaptive mesh refinement is implemented in the quadrangular grid systems and easily embedded in the FEM-based algorithm. For the judgement on renewal of mesh, the level set function is adopted as an indicator, and the threshold is set at the boundary of the smoothing band. With this criterion, the variation of physical properties and the jump quantity on the free surface can be calculated accurately enough, while the computation cost is largely reduced as a whole. In order to prove the validity of the present scheme, two-dimensional numerical simulation is carried out in collapse of a water column, oscillation and movement of a drop under zero gravity. As a result, its effectiveness and usefulness are clearly shown qualitatively and quantitatively. Among them, the movement of a drop due to the Marangoni effect is first simulated efficiently with the present scheme.

title: Adaptive unstructured volume remeshing - II: Application to two- and three-dimensional level-set simulations of multiphase flow
abstract: In Part I [Adaptive unstructured volume remeshing - I: The method, ], we presented an adaptive remeshing algorithm that automatically adjusts the size of the elements of meshes of unstructured triangles (2D) and unstructured tetrahedra (3D) with time and position in the computational domain in order to efficiently resolve the relevant physical scales. Here, we illustrate the performance of an implementation of the algorithm in finite-element/level-set simulations of deformable droplet and fluid-fluid interface interactions, breakup and coalescence in multiphase flows. The wide range of length scales characterizing the dynamics are accurately resolved as demonstrated by comparison to experiments and to theoretical and sharp-interface (boundary-integral) numerical results. The computational cost is found to be competitive even with respect to boundary-integral methods. For the first time using an interface-capturing (level-set) method we successfully simulate the inertia driven impact and rebound of a liquid droplet from a liquid interface and find agreement with recent experimental results.

title: Sharp interface Cartesian grid method II: A technique for simulating droplet interactions with surfaces of arbitrary shape
abstract: A fixed-grid, sharp interface method is developed to simulate droplet impact and spreading on surfaces of arbitrary shape. A finite-difference technique is used to discretize the incompressible Navier-Stokes equations on a Cartesian grid. To compute flow around embedded solid boundaries, a previously developed sharp interface method for solid immersed boundaries is used. The ghost fluid method (GFM) is used for fluid-fluid interfaces. The model accounts for the effects of discontinuities such as density and viscosity jumps and singular sources such as surface tension in both bubble and droplet simulations. With a level-set representation of the propagating interface, large deformations of the boundary can be handled easily. The model successfully captures the essential features of interactions between fluid-fluid and solid-fluid phases during impact and spreading. Moving contact lines are modeled with contact angle hysteresis and contact line motion on non-planar surfaces is computed. Experimental observations and other simulation results are used to validate the calculations.

title: A stabilized finite element method using a discontinuous level set approach for solving two phase incompressible flows
abstract: A numerical method for the simulation of three-dimensional incompressible two-phase flows is presented. The proposed algorithm combines an implicit pressure stabilized finite element method for the solution of incompressible two-phase flow problems with a level set method implemented with a quadrature-free discontinuous Galerkin (DG) method E. Marchandise, J.-F. Remacle, N. Chevaugeon, A quadrature free discontinuous Galerkin method for the level set equation, Journal of Computational Physics 212 (2006) 338-357]. The use of a fast contouring algorithm N. Chevaugeon, E. Marchandise, C. Geuzaine, J.-F. Remacle, Efficient visualization of high order finite elements,  Engineering] permits us to localize the interface accurately. By doing so, we can compute the discontinuous integrals without neither introducing an interface thickness nor reinitializing the level set.The capability of the resulting algorithm is demonstrated with "large scale" numerical examples (free surface flows: dam break, sloshing) and "small scale" ones (two phase Poiseuille, Rayleigh-Taylor instability).

title: A stabilized finite element method using a discontinuous level set approach for the computation of bubble dynamics
abstract: A novel numerical method for solving three-dimensional two phase flow problems is presented. This method combines a quadrature free discontinuous Galerkin method for the level set equation with a pressure stabilized finite element method for the Navier Stokes equations. The main challenge in the computation of such flows is the accurate evaluation of surface tension forces. This involves the computation of the curvature of the fluid interface. In the context of the discontinuous Galerkin method, we show that the use of a curvature computed by means of a direct derivation of the level set function leads to inaccurate and oscillatory results. A more robust, second-order, least squares computation of the curvature that filters out the high frequencies and produces converged results is presented. This whole numerical technology allows to simulate a wide range of flow regimes with large density ratios, to accurately capture the shape of the deforming interface of the bubble and to maintain good mass conservation.

title: An accurate conservative level set/ghost fluid method for simulating turbulent atomization
abstract: This paper presents a novel methodology for simulating incompressible two-phase flows by combining an improved version of the conservative level set technique introduced in [E. Olsson, G. Kreiss, A conservative level set method for two phase flow,  (2005) 225-246] with a ghost fluid approach. By employing a hyperbolic tangent level set function that is transported and re-initialized using fully conservative numerical schemes, mass conservation issues that are known to affect level set methods are greatly reduced. In order to improve the accuracy of the conservative level set method, high order numerical schemes are used. The overall robustness of the numerical approach is increased by computing the interface normals from a signed distance function reconstructed from the hyperbolic tangent level set by a fast marching method. The convergence of the curvature calculation is ensured by using a least squares reconstruction. The ghost fluid technique provides a way of handling the interfacial forces and large density jumps associated with two-phase flows with good accuracy, while avoiding artificial spreading of the interface. Since the proposed approach relies on partial differential equations, its implementation is straightforward in all coordinate systems, and it benefits from high parallel efficiency. The robustness and efficiency of the approach is further improved by using implicit schemes for the interface transport and re-initialization equations, as well as for the momentum solver. The performance of the method is assessed through both classical level set transport tests and simple two-phase flow examples including topology changes. It is then applied to simulate turbulent atomization of a liquid Diesel jet at Re=3000. The conservation errors associated with the accurate conservative level set technique are shown to remain small even for this complex case.

title: Surface tension and viscosity with Lagrangian hydrodynamics on a triangular mesh
abstract: Numerical algorithms for surface tension and viscosity are presented in the context of a Lagrangian treatment of incompressible hydrodynamics with a dynamically restructuring grid. New algorithms are given which update previous Lagrangian approaches in the code SPLISH. Test problems involving internal gravity and capillary waves, an oscillating droplet and a viscous shear layer are described. An example is given of a flow calculated in and around a viscous droplet with surface tension in a shear flow.

title: A front-tracking algorithm for accurate representation of surface tension
abstract: We present a front tracking algorithm for the solution of the 2D 
incompressible Navier-Stokes equations with interfaces and surface 
forces. More particularly, we focus our attention on the accurate 
description of the surface tension terms and the associated pressure 
jump. We consider the stationary Laplace solution for a bubble with 
surface tension. A careful treatment of the pressure gradient terms at 
the interface allows us to reduce the spurious currents to the machine 
precision. Good results are obtained for the oscillation of a capil- 
lary wave compared with the linear viscous theory. A classical test of 
Rayleigh-Taylor instability is presented.

title: Accurate representation of surface tension using the level contour reconstruction method
abstract: Some of the most demanding tests of interface methods for two-phase flows with surface tension which use fixed Eulerian grids occur at the two extremes of highly dynamic flows or static equilibrium conditions. It has been difficult to design an interface method to operate accurately across this spectrum especially for 3D fluid flow calculations which, on the one hand, do not often have the required grid resolution to capture all of the fine scale structures that typically appear in highly stretched interfaces nor, on the other hand, the required accuracy in calculating surface tension forces. We present improvements to the interface reconstruction procedure in the level contour reconstruction method (LCRM)  These improvements allow more precise control of the interface reconstruction in highly dynamic flows with coalescence and rupture and also avoid the problem of local mass redistribution in poorly resolved calculations. In addition, a new hybrid technique for surface tension calculation in the context of the front tracking method is demonstrated and shown to result in a marked improvement in suppressing parasitic currents by generally two orders of magnitude. We compare and validate these new procedures in various test cases.

title: Adaptive tetrahedral meshing in free-surface flow
abstract: An unstructured, three-dimensional, moving-mesh algorithm is developed for free-surface flows with the capability of simulating large deformation. A robust edge swapping algorithm and an optimization-based mesh smoothing algorithm are used to improve mesh quality as the interface deforms. The edge swapping algorithm uses an optimal swap sequence to transform configurations with several tetrahedra sharing a common edge. The optimization-based mesh smoothing relocates vertex position so that the minimum quality of the incident cells is maximized. Also, edge contraction and bisection are used to adjust the mesh resolution according to surface curvature. These local refinements effectively maintain good mesh quality and appropriate mesh size while avoiding global re-meshing. This scheme's robustness and accuracy are demonstrated with deforming liquid ligaments, periodic jets, and colliding droplets.

title: Modelling Merging and Fragmentation in Multiphase Flows with SURFER
abstract: We introduce a new numerical method, called "SURFER," for the simulation of two- and three-dimensional flows with several fluid phases and free interfaces between them. We consider incompressible fluids obeying the Navier-Stokes equation with Newtonian viscosity in the bulk of each phase. Capillary forces are taken into account even when interfaces merge or break up. Fluid interfaces are advanced in time using an exactly volume conserving variant of the volume of fluid algorithm, thus allowing for full symmetry between fluid phases. The Navier-Stokes equation is solved using staggered finite differences on a MAC grid and a split-explicit time differencing scheme, while incompressibility is enforced using an iterative multigrid Poisson solver. Capillary effects are represented as a stress tensor computed from gradients of the volume fraction function. This formulation is completely independent of the topology of interfaces and relatively easy to implement in 3D. It also allows exact momentum conservation in the discretized algorithm. Numerical spurious effects or "parasite currents" are noticed and compared to similar effects in Boltzmann lattice gas methods for immiscible fluids. Simulations of droplets pairs colliding in 2D and in 3D are shown. Interface reconnection is performed easily, despite the large value of capillary forces during reconnection.

title: Volume-of-Fluid Interface Tracking with Smoothed Surface Stress Methods for Three-Dimensional Flows
abstract: Motivated by the need for three-dimensional methods for interface calculations that can deal with topology changes, we describe a numerical scheme, built from a volume-of-fluid interface tracking technique that uses a piecewise-linear interface calculation in each cell. Momentum balance is computed using explicit finite volume/finite differences on a regular cubic grid. Surface tension is implemented by the continuous surface stress or continuous surface force method. Examples and verifications of the method are given by comparing simulations to analytical results and experiments, for sedimenting droplet arrays and capillary waves at finite Reynolds number. In the case of a pinching pendant drop, both three-dimensional and axisymmetric simulations are compared to experiments. Agreement is found both before and after the reconnections.

title: A Boundary Condition Capturing Method for Multiphase Incompressible Flow
abstract: In , the Ghost Fluid Method (GFM) was developed to capture the boundary conditions at a contact discontinuity in the inviscid compressible Euler equations. In [11], related techniques were used to develop a boundary condition capturing approach for the variable coefficient Poisson equation on domains with an embedded interface. In this paper, these new numerical techniques are extended to treat multiphase incompressible flow including the effects of viscosity, surface tension and gravity. While the most notable finite difference techniques for multiphase incompressible flow involve numerical smearing of the equations near the interface, see, e.g., [19, 17, 1], this new approach treats the interface in a sharp fashion.

title: PROST: a parabolic reconstruction of surface tension for the volume-of-fluid method
abstract: Volume-of-fluid (VOF) methods are popular for the direct numerical simulation of time-dependent viscous incompressible flow of multiple liquids. As in any numerical method, however, it has its weaknesses, namely, for flows in which the capillary force is the dominant physical mechanism. The lack of convergence with spatial refinement, or convergence to a solution that is slightly different from the exact solution, has been documented in the literature. A well-known limiting case for this is the existence of spurious currents for the simulation of a spherical drop with zero initial velocity. These currents are present in all previous versions of VOF algorithms. In this paper, we develop an accurate representation of the body force due to surface tension, which effectively eliminates spurious currents. We call this algorithm PROST: parabolic reconstruction of surface tension. There are several components to this procedure, including the new body force algorithm, improvements in the projection method for the Navier-Stokes solver, and a higher order interface advection scheme. The curvature to the interface is calculated from an optimal fit for a quadratic approximation to the interface over groups of cells.

title: An analysis of parasitic current generation in Volume of Fluid simulations
abstract: Parasitic currents are unphysical fluid motions generated when using implementations of the Continuum Surface Force technique to model surface tension forces in Eulerian-based multi-phase Computational Fluid Dynamics calculations. In this paper, we derive and validate a correlation for the magnitude of these currents as a function of the physical and numerical parameters used in a given simulation. We find that these currents may be limited by both the inertial and viscous terms in the Navier–Stokes equations, and that they do not decrease in magnitude with increased mesh refinement or decreased computational time step. The current magnitudes that would be generated when simulating a number of common physical systems are examined, and it is found that in some circumstances their magnitude is significant.

title: A balanced-force algorithm for continuous and sharp interfacial surface tension models within a volume tracking framework
abstract: A new balanced-force algorithm is presented for modeling interfacial flow with surface tension. The algorithm is characterized by a pressure-correction method with the interfaces represented by volume fractions. Within this flow algorithm, we devise a continuous (e.g., continuum surface tension model) and a sharp (e.g., a ghost fluid method) interface representation of the surface-tension-induced interfacial pressure jump condition. The sharp interface representation is achieved by temporarily reconstructing distance functions from volume fractions. We demonstrate that a flow algorithm designed to legislate force balance retains an exact balance between surface tension forces and the resulting pressure gradients. This balance holds for both continuous and sharp representations of interfacial surface tension. The algorithm design eliminates one of the elusive impediments to more accurate models of surface tension-driven flow, the remaining of which is accurate curvature estimation. To validate our formulation, we present results for an equilibrium (static) drop in two and three dimensions having an arbitrary density jump across the interface. We find that the sharp surface tension method yields an abrupt pressure jump across the interface, whereas the continuous surface tension method results in a smoother transition. Both methods, however, yield spurious velocities of the same order, the origin of which is due solely to errors in curvature. Dynamic results are also presented to illustrate the versatility of the method.

title: Two-Phase Flows on Interface Refined Grids Modeled with VOF, Staggered Finite Volumes, and Spline Interpolants
abstract: A two-phase 2D model that combines the volume of fluid (VOF) method with implicit staggered finite volumes discretization of the Navier–Stokes equation is presented. Staggered finite volumes are developed on the basis of nonconforming Crouzeix–Raviart finite elements, where all components of the velocity lie in the middle of the element edges and the pressure degrees of freedom are found in the centers of mass of the elements. Staggered finite Volumes extend marker and cell (MAC) regular staggered grids to unstructured mesh. A linear saddle point problem, resulting from either the discretization or the Newton method, is solved for all unknown pressures and velocities. Interface is represented with spline interpolants which follow the VOF distribution. Adaptive mesh refinement is used to obtain a high level of uniform refining at the domain of dependence of the interface. The aligned grid is obtained by irregular refining of the cells which are intersected by a curve. The boundaries of its elements coincide with the slope segments going through the intersections of the curve with the underlying regular elements boundary. The deformable computational grids are used only to discretize the Navier–Stokes equation. The advection of volume fractions is done on the advection mesh, which corresponds to highest regular refining on the computational grid. Approximation of the surface tension on spline interpolants offers a straightforward way to describe correctly the pressure jumps on interface-fitted staggered grids. This allows deletion of the anomalous currents around a statical bubble and their effective reduction in real simulations. On the aligned grid, the continuity of the viscous stress is modeled exactly due to the finite volume approach. Using the proposed numerical techniques, single bubble rise is analyzed.

title: Estimating curvature from volume fractions
abstract: An interface represented by discrete, abruptly-varying volume fractions poses particular challenges for the accurate estimation of interfacial curvature. Most approaches within a volume-of-fluid (VOF) framework do not estimate curvature directly from the VOF function, but instead from a smoothly-varying function derived from some mapping (e.g., convolution, discrete quadrature) of the VOF function. In this study we assess and compare the accuracy and robustness of curvature estimates derived from three smooth, VOF-based functions: a convolved VOF (CV) function, a height function (HF), and a reconstructed distance function (RDF). The methodology used in reconstructing a distance function from volume fractions is new and therefore presented in detail along with an analysis of the resulting sensitivities. We find that curvature estimates derived from the height function provide superior results, as measured in smaller error and higher (second-order) convergence rate. While curvature estimates from the RDF provides smaller errors relative to those derived from the CV, both techniques exhibit similar convergence behaviour. Numerical results, including a new cosine wave verification test, are presented to substantiate the findings.

title: Improvements for calculating two-phase bubble and drop motion using an adaptive sharp interface method.
abstract: In this paper, we describe new tech- niques for numerically approximating two-phase flows. Specifically, we present new techniques fortreatingtheviscosityandsurfacetensionterms that appear in the Navier-Stokes equations for in- compressible two-phase flow. Our resulting nu- merical method has the property that results com- putedusingourtwo-phasealgorithmapproachthe corresponding "one-phase" algorithm in the limit of zero gas density/viscosity; i.e. the two-phase results approach the one-phase free-boundary re- sultsinthelimitthatthegas isassumed to become a uniform pressure void. By grid convergence checks and comparison with previous experimen- tal data, we shall demonstrate the advantages of our new proposed discretizations. Keyword: volume-of-fluid method, sharp inter- face, multiphase flow, surface tension, viscosity

title: A Conservative Adaptive Projection Method for the Variable Density Incompressible Navier-Stokes Equations
abstract: In this paper we present a method for solving the equations governing time-dependent, variable density incompressible flow in two or three dimensions on an adaptive hierarchy of grids. The method is based on a projection formulation in which we first solve advection?diffusion equations to predict intermediate velocities, and then project these velocities onto a space of approximately divergence-free vector fields. Our treatment of the first step uses a specialized second-order upwind method for differencing the nonlinear convection terms that provides a robust treatment of these terms suitable for inviscid and high Reynolds number flow. Density and other scalars are advected in such a way as to maintain conservation, if appropriate, and free-stream preservation. Our approach to adaptive refinement uses a nested hierarchy of logically-rectangular girds with simultaneous refinement of the girds in both space and time. The integration algorithm on the grid hierarchy is a recursive procedure in which coarse grids are advanced in time, fine grids are advanced multiple steps to reach the same time as the coarse grids and the data at different levels are then synchronized. The single grid algorithm is described briefly, but the emphasis here is on the time-stepping procedure for the adaptive hierarchy. Numerical examples are presented to demonstrate the algorithms's accuracy and convergence properties, and illustrate the behavior of the method. An additional example demonstrates the performance of the method on a more realistic problem, namely, a three-dimensional variable density shear layer.

title: Finite element analysis of transient fluid flow with free surface using VOF (volume-of-fluid) method and adaptive grid
abstract: The VOF method is adopted for the finite element analysis of transient fluid flow with a free surface. In particular, an adaptation technique for generating an adaptive grid is incorporated to capture a higher resolution of the free surface configuration. An adaptive grid is created through the refinement and mergence of elements. In this domain the elements in the surface region are made finer than those in the remaining regions for more efficient computation. Also, three techniques based on the VOF method are newly developed to increase the accuracy of the analysis, namely the filling pattern, advection treatment and free surface smoothing techniques. Using the proposed numerical techniques, radial flow with a point source and the collapse of a dam are analysed. The numerical results agree well with the theoretical solutions as well as with the experimental results. Through comparisons with the numerical results of several cases using different grids, the efficiency of the proposed technique is verified.

title: A parallelized, adaptive algorithm for multiphase flows in general geometries
abstract: We describe a parallelized, adaptive algorithm for multiphase flow in general geometries. Phase boundaries are considered sharp and their motion is computed using the coupled level set/volume-of-fluid (CLSVOF) method. General geometries are represented using a ''Cartesian grid'' approach; not a body-fitted approach. Our algorithm is designed to be ''robust''; in other words we compute flows in which phase boundaries and geometry boundaries can be very complex with very little user intervention. In this paper we study the speed-up gained via adaptive mesh refinement, and/or parallelization.

title: A quadtree adaptive method for simulating fluid flows with moving interfaces
abstract: In this paper, a computational method for solving fluid flow problems with moving interfaces is presented. Herein, adaptive quadtree grids are used coupled with the CICSAM [O. Ubbink, Numerical prediction of two fluid systems with sharp interfaces, PhD Thesis, Imperial College of Science, Technology and Medicine, London, 1997] free surface capturing volume of fluid (VoF) method and PLIC reconstruction to interpolate the volume fraction field during refinement and derefinement processes. The combination of high resolution adaptive hierarchical remeshing and CICSAM interface advection is shown to overcome the problems of interface smearing and high CPU intensivity inherent in most VoF schemes. The result is a combination of free surface tracking and free surface capturing in that the interface is effectively tracked by the adapting refinements in the quadtree grid. In this way, a sharp interface is achieved and the advantages of both free surface tracking and capturing are combined. The new method is applied to interface advection examples in translating, rotating and shearing flow fields, and the benefits of using adapting quadtree grids demonstrated.

title: Simulation of sharp gas–liquid interface using VOF method and adaptive grid local refinement around the interface
abstract: The volume of fluid (VOF) method is used to perform two-phase simulations (gas–liquid). The governing Navier–Stokes conservation equations of the flow field are numerically solved on two-dimensional axisymmetric or three-dimensional unstructured grids, using Cartesian velocity components, following the finite volume approximation and a pressure correction method. A new method of adaptive grid local refinement is developed in order to enhance the accuracy of the predictions, to capture the sharp gas–liquid interface and to speed up the calculations. Results are compared with experimental measurements in order to assess the efficiency of the method.

title: Spatially adaptive techniques for level set methods and incompressible flow
abstract: Since the seminal work of [Sussman, M, Smereka P, Osher S. A level set approach for computing solutions to incompressible two-phase flow.  on coupling the level set method of [Osher S, Sethian J. Fronts propagating with curvature-dependent speed: algorithms based on Hamilton–Jacobi formulations.  to the equations for two-phase incompressible flow, there has been a great deal of interest in this area. That work demonstrated the most powerful aspects of the level set method, i.e. automatic handling of topological changes such as merging and pinching, as well as robust geometric information such as normals and curvature. Interestingly, this work also demonstrated the largest weakness of the level set method, i.e. mass or information loss characteristic of most Eulerian capturing techniques. In fact, [Sussman M, Smereka P, Osher S. A level set approach for computing solutions to incompressible two-phase flow.  introduced a partial differential equation for battling this weakness, without which their work would not have been possible. In this paper, we discuss both historical and most recent works focused on improving the computational accuracy of the level set method focusing in part on applications related to incompressible flow due to both of its popularity and stringent accuracy requirements. Thus, we discuss higher order accurate numerical methods such as Hamilton–Jacobi WENO [Jiang G-S, Peng D. Weighted ENO schemes for Hamilton–Jacobi equations. , methods for maintaining a signed distance function, hybrid methods such as the particle level set method [Enright D, Fedkiw R, Ferziger J, Mitchell I. A hybrid particle level set method for improved interface capturing.  and the coupled level set volume of fluid method [Sussman M, Puckett EG. A coupled level set and volume-of-fluid method for computing 3d and axisymmetric incompressible two-phase flows. , and adaptive gridding techniques such as the octree approach to free surface flows proposed in [Losasso F, Gibou F, Fedkiw R. Simulating water and smoke with an octree data structure, ACM Trans Graph (SIGGRAPH Proc) 2004;23:457–62].

title: A second order accurate level set method on non-graded adaptive cartesian grids
abstract: We present a level set method on non-graded adaptive Cartesian grids, i.e. grids for which the ratio between adjacent cells is not constrained. We use quadtree and octree data structures to represent the grid and a simple algorithm to generate a mesh with the finest resolution at the interface. In particular, we present (1) a locally third order accurate reinitialization scheme that transforms an arbitrary level set function into a signed distance function, (2) a second order accurate semi-Lagrangian methods to evolve the linear level set advection equation under an externally generated velocity field, (3) a second order accurate upwind method to evolve the non-linear level set equation under a normal velocity as well as to extrapolate scalar quantities across an interface in the normal direction, and (4) a semi-implicit scheme to evolve the interface under mean curvature. Combined, we obtain a level set method on adaptive Cartesian grids with a negligible amount of mass loss. We propose numerical examples in two and three spatial dimensions to demonstrate the accuracy of the method.

title: Transient adaptivity applied to two-phase incompressible flows
abstract: An anisotropic adaptation process is applied to a three-dimensional incompressible two-phase flow solver. The solver uses a level set/finite element method on unstructured tetrahedral meshes. We show how the level set function can be used to build an anisotropic mesh with good properties. Some computations with a strong transient character and large densities ratios (1/1000) are presented. We show that the efficiency of the computations can be deeply enhanced by mesh adaptations.

title: Gerris: a tree-based adaptive solver for the incompressible Euler equations in complex geometries
abstract: An adaptive mesh projection method for the time-dependent incompressible Euler equations is presented. The domain is spatially discretised using quad/octrees and a multilevel Poisson solver is used to obtain the pressure. Complex solid boundaries are represented using a volume-of-fluid approach. Second-order convergence in space and time is demonstrated on regular, statically and dynamically refined grids. The quad/octree discretisation proves to be very flexible and allows accurate and efficient tracking of flow features. The source code of the method implementation is freely available.

title: Finite-volume-type VOF method on dynamically adaptive quadtree grids
abstract: This paper describes a finite-volume volume-of-fluid (VOF) method for simulating viscous free surface flows on dynamically adaptive quadtree grids. The scheme is computationally efficient in that it provides relatively fine grid resolution at the gas–liquid interface and coarse grid density in regions where flow variable gradients are small. Special interpolations are used to ensure volume flux conservation where differently sized neighbour cells occur. The numerical model is validated for advection of dyed fluid in unidirectional and rotating flows, and for two-dimensional viscous sloshing in a rectangular tank.

title: Simulating water and smoke with an octree data structure
abstract: We present a method for simulating water and smoke on an unrestricted octree data structure exploiting mesh refinement techniques to capture the small scale visual detail. We propose a new technique for discretizing the Poisson equation on this octree grid. The resulting linear system is symmetric positive definite enabling the use of fast solution methods such as preconditioned conjugate gradients, whereas the standard approximation to the Poisson equation on an octree grid results in a non-symmetric linear system which is more computationally challenging to invert. The semi-Lagrangian characteristic tracing technique is used to advect the velocity, smoke density, and even the level set making implementation on an octree straightforward. In the case of smoke, we have multiple refinement criteria including object boundaries, optical depth, and vorticity concentration. In the case of water, we refine near the interface as determined by the zero isocontour of the level set function.

title: Simulation of viscous water column collapse using adapting hierarchical grids
abstract: An adaptive hierarchical grid based method for predicting complex free surface flows is used to simulate collapse of a water column. Adapting quadtree grids are combined with a high-resolution interface capturing approach and pressure based coupling of the Navier Stokes equations. The Navier-Stokes flow solution scheme is verified for simulation of flow in a lid driven cavity at Re=1000. Two approaches to the coupling of the Navier-Stokes equations are investigated as are alternative face velocity and hanging node interpolations. Collapse of a water column as well as collapse of a water column and its subsequent interaction with an obstacle are simulated. The calculations are made on uniform and adapting quadtree grids, and the accuracy of the quadtree calculations is shown to be the same as those made on the equivalent uniform grids. Results are in excellent agreement with experimental and other numerical data. A sharp interface is maintained at the free surface. The new adapting quadtree-based method achieves a considerable saving in the size of the computational grid and CPU time in comparison with calculations made on equivalent uniform grids.

title: Three-dimensional numerical investigation of a droplet impinging normally onto a wall film
abstract: The paper presents a three-dimensional numerical investigation of a droplet impinging normally onto a wall film. The numerical method is based on the finite volume solution of the Navier-Stokes equations coupled with the volume of fluid method (VOF) and utilizing an adaptive local grid refinement technique for tracking more accurately the liquid-gas interface. The results are compared with available experimental data for integral quantities such as the lamella temporal development. Two mechanisms are identified leading to secondary droplet formation; in the initial and intermediate stages of splashing secondary droplet formation is according to Rayleigh instability while at later times surface tension effects contribute further to secondary atomization. Moreover, the influence of Weber number on the impingement process is investigated and correlations for the diameter and number of secondary droplets are proposed.

title: A balanced force refined level set grid method for two-phase flows on unstructured flow solver grids
abstract: A balanced force refined level set grid method for two-phase flows on structured and unstructured flow solver grids is presented. To accurately track the phase interface location, an auxiliary, high-resolution equidistant Cartesian grid is introduced. In conjunction with a dual-layer narrow band approach, this refined level set grid method allows for parallel, efficient grid convergence and error estimation studies of the interface tracking method. The Navier-Stokes equations are solved on an unstructured flow solver grid with a novel balanced force algorithm for level set methods based on the recently proposed method by Francois et al. [M.M. Francois, S.J. Cummins, E.D. Dendy, D.B. Kothe, J.M. Sicilian, M.W. Williams, A balanced-force algorithm for continuous and sharp interfacial surface tension models within a volume tracking framework,  (2006) 141-173] for volume of fluid methods on structured grids. To minimize spurious currents, a second order converging curvature evaluation technique for level set methods is presented. The results of several different test cases demonstrate the effectiveness of the proposed method, showing good mass conservation properties and second order converging spurious current magnitudes.

title: Adaptive vof with curvature-based refinement
abstract: Adaptive refinement is implemented in the context of the volume-of-fluid (VOF) methodology in order to study the efficacy of resolving interfaces adaptively based on the local value of curvature. The usual uniform mesh VOF implementation is modified slightly to ensure accurate advection of fluxes between cells at different resolutions. Normals and curvatures are calculated accurately via height functions. Results of a series of tests indicate that in most instances the use of adaptive refinement (when compared to uniform refinement with a similar number of cells) leads to more accurate VOF advection. The results also clearly show that curvature-based adaptive refinement leads to a distribution of errors along an interface that is nearly independent of curvature.

title: Analytical relations connecting linear interfaces and volume fractions in rectangular grids
abstract: The computational-geometric problems arising when a linear interface cuts a cube are considered. They are of interest in particular for the calculation of volume fractions or interface positions in three-dimensional interface calculations in the Volume of Fluid (VOF) methods. Typically, the normal vector is known. One then wants to compute the volume fraction knowing the interface position, or conversely the interface position knowing the volume fraction. Explicit expressions of general use are given, and the algorithms used to search for solutions are described in detail. Explicit formulas for cubic roots are found to be less than two thirds as time consuming as Newton–Raphson iterations.

title: Interface reconstruction with least‐square fit and split Eulerian–Lagrangian advection
abstract: Two new volume-of-fluid (VOF) reconstruction algorithms, which are based on a least-square fit technique, are presented. Their performance is tested for several standard shapes and is compared to a few other VOF/PLIC reconstruction techniques, showing in general a better convergence rate. The geometric nature of Lagrangian and Eulerian split advection algorithms is investigated in detail and a new mixed split Eulerian implicit–Lagrangian explicit (EI–LE) scheme is presented. This method conserves the mass to machine error, performs better than split Eulerian and Lagrangian algorithms, and it is only slightly worse than unsplit schemes. However, the combination of the interface reconstruction with the least-square fit and its advection with the EI–LE scheme appears superior to other existing approaches.

title: Interface reconstruction with least-squares fit and split advection in three-dimensional Cartesian geometry
abstract: In this paper we present and analyze different volume-of-fluid (VOF) reconstruction and advection algorithms that approximate the interface separating two immiscible fluids in the three-dimensional space. The paper describes the improvement of the reconstruction when a least-square fit algorithm, which minimizes a distance functional, is applied. Its performance is tested for several smooth surfaces against other simple reconstruction methods. Then Eulerian, Lagrangian and mixed split advection schemes are presented and analyzed. In particular, one advection method is discussed that conserves mass exactly for a divergence-free velocity field, thus allowing computations to machine precision.

title: High-order surface tension VOF-model for 3D bubble flows with high density ratio
abstract: An improved Volume of Fluid (VOF) method is presented which is applicable to high density ratio 3D flows for a large range of bubble Reynolds number (Re). The method is based on the Navier-Stokes equations for incompressible multi-phase flows which are discretized on a Cartesian staggered grid. The multi-grid technique together with the pressure-velocity coupling scheme for multi-phase flows have resulted in an efficient solver which nearly exponentially converge with the number of iterations. The convergence speed also shows negligible dependence on density ratio, viscosity ratio and Re. A second-order accurate, non-diffusive, mass conservative phase transport model is presented which does not suffer from unphysical over- or under-shoots of the phase variable. The high accurate normal, curvature and surface tension force model in combination with the high-order defect-correction scheme for multi-phase flows shows second-order global accuracy when applied to the transient bubble rise where the viscosity ratio is equal to one. In contrast, the commonly used viscosity model for VOF introduces a first order error for the same problem. The VOF method has been tested for different types of bubble flows at low Re and for path-oscillating and wobbling air bubbles (in water) with a diameter range of 1.82 6 mm. The numerical results agree quantitatively with the available experimental data. The investigations show that the proposed high accurate surface tension model can be used successfully for wobbling flows with bubble deformation while maintaining the mass of the phases. The error in mass conservation is directly proportional to the residual in solving the discrete problem.

title: Fully multidimensional flux-corrected transport algorithms for fluids
abstract: The theory of flux-corrected transport (FCT) developed by Boris and Book [ (1973) 38; 18 (1975) 248; 20 (1976) 397] is placed in a simple, generalized format, and a new algorithm for implementing the critical flux limiting stage' in multidimensions without resort to time splitting is presented. The new flux limiting algorithm allows the use of FCT techniques in multidimensional fluid problems for which time splitting would produce unacceptable numerical results, such as those involving incompressible or nearly incompressible flow fields. The “clipping” problem associated with the original one dimensional flux limiter is also eliminated or alleviated. Test results and applications to a two dimensional fluid plasma problem are presented.

title: A Method for Capturing Sharp Fluid Interfaces on Arbitrary Meshes
abstract: The paper describes a high resolution method (CICSAM) for the accurate capturing of fluid interfaces on meshes of arbitrary topology. It is based on the finite-volume technique and is fully conservative. The motion of the interface is tracked by the solution of a scalar transport equation for a phase-indicator field that is discontinuous at the interface and uniform elsewhere; no explicit interface reconstruction, which is perceived to be difficult to implement on unstructured meshes, is needed. The novelty of the method lies in the adaptive combination of high resolution discretisation schemes which ensure the preservation of the sharpness and shape of the interface while retaining boundedness of the field. The special implicit implementation presented herein makes it applicable to unstructured meshes and an extension to such grids is presented. The method is capable of handling interface rupture and coalescence. The paper outlines the methodology of CICSAM and its validation against academic test cases used to verify its accuracy.

title: A continuum method for modeling surface tension
abstract: A new method for modeling surface tension effects on fluid motion has been developed. Interfaces between fluids of different properties, or “colors,” are represented as transition regions of finite thickness, across which the color variable varies continuously. At each point in the transition region, a force density is defined which is proportional to the curvature of the surface of constant color at that point. It is normalized so that the conventional description of surface tension on an interface is recovered when the ratio of local transition region thickness to local radius of curvature approaches zero. The continuum method eliminates the need for interface reconstruction, simplifies the calculation of surface tension, enables accurate modeling of two- and three-dimensional fluid flows driven by surface forces, and does not impose any modeling restrictions on the number, complexity, or dynamic evolution of fluid interfaces having surface tension. Computational results for two-dimensional flows are given to illustrate the properties of the method.

title: A Coupled Level Set and Volume-of-Fluid Method for Computing 3D and Axisymmetric Incompressible Two-Phase Flows
abstract: We present a coupled level set/volume-of-fluid (CLSVOF) method for computing 3D and axisymmetric incompressible two-phase flows. This method combines some of the advantages of the volume-of-fluid method with the level set method to obtain a method which is generally superior to either method alone. We present direct comparisons between computations made with the CLSVOF method and computations made with the level set method, the volume-of-fluid method, and the boundary integral method. We also compare our computations to the exact solution for an oscillating ellipse due to Lamb and experimental results obtained for a rising gas bubble in liquid obtained by Hnat and Buckmaster. Our computational examples focus on flows in which surface tension forces and changes in topology are dominant features of the flow.

title: An adaptive coupled level-set/volume-of-fluid interface capturing method for unstructured triangular grids
abstract: We present an adaptive coupled level-set/volume-of-fluid (ACLSVOF) method for interfacial flow simulations on unstructured triangular grids. At each time step, we evolve both the level set function and the volume fraction. The level set function is evolved by solving the level set advection equation using a discontinuous Galerkin finite element method. The volume fraction advection is performed using a Lagrangian-Eulerian method. The interface is reconstructed based on both the level set and the volume fraction information. In particular, the interface normal vector is calculated from the level set function while the line constant is determined by enforcing mass conservation based on the volume fraction. Different from previous works, we have developed an analytic method for finding the line constant on triangular grids, which makes interface reconstruction efficient and conserves volume of fluid exactly. The level set function is finally reinitialized to the signed distance to the reconstructed interface. Since the level set function is continuous, the normal vector calculation is easy and accurate compared to a classic volume-of-fluid method, while tracking the volume fraction is essential for enforcing mass conservation. The method is also coupled to a finite element based Stokes flow solver. The code validation shows that our method is second order and mass is conserved very accurately. In addition, owing to the adaptive grid algorithm we can resolve complex interface changes and interfaces of high curvature efficiently and accurately.

title: Comparison of volume-of-fluid methods for surface tension-dominant two-phase flows
abstract: The capabilities of three volume-of-fluid methods for the calculation of surface tension-dominant two-phase flows are tested. The accurate calculation of the interface remains a problem for the volume-of-fluid method if the surface tension force plays an important role and the density ratios of the fluids in different phases are high. The result can be an artificial velocity field at the interface ("parasitic currents"), which can destabilize the interface significantly. The three different algorithms compared can be distinguished by their methods to compute the surface tension force, namely, the method using a kernel function for smoothing the discontinuity at the interface, a combined level-set and volume-of-fluid approach and a parabolic reconstruction of surface tension. The test cases consist of an equilibrium rod, a capillary wave and the Rayleigh-Taylor instability. The analytical solutions for each problem serve to examine the accuracy and the convergence behavior of each approach. Finally, the slow formation of a gas bubble at an underwater orifice was computed with the combined level-set and volume-of-fluid method and the results are compared with an analytical solution based on the Young-Laplace equation.

title: Motion of two superposed viscous fluids
abstract: The initial‐value problem associated with the small‐amplitude waves on the interface between two superposed viscous fluids is studied in the stable and in the unstable case (Rayleigh–Taylor instability). The results are illustrated with examples obtained by numerical and analytic inversion of the Laplace‐transformed solution for the wave amplitude. The diffusion of vorticity from the free surface into the fluids is also shown.

title: A Cascade of Structure in a Drop Falling from a Faucet
abstract: A drop falling from a faucet is a common example of a mass fissioning into two or more pieces. The shape of the liquid in this situation has been investigated by both experiment and computer simulation. As the viscosity of the liquid is varied, the shape of the drop changes dramatically. Near the point of breakup, viscous drops develop long necks that then spawn a series of smaller necks with ever thinner diameters. Simulations indicate that this repeated formation of necks can proceed ad infinitum whenever a small but finite amount of noise is present in the experiment. In this situation, the dynamical singularity occurring when a drop fissions is characterized by a rough interface.

title: Breakdown of scaling in droplet fission at high Reynolds number
abstract: In this paper we address the shape of a low-viscosity fluid interface near the breaking point. Experiments show that the shape varies dramatically as a function of fluid viscosity. At low viscosities, the interface develops a region with an extremely sharp slope, with the steepness of the slope diverging with vanishing viscosity. Numerical simulations demonstrate that this tip forms as a result of a convective instability in the fluid; in the absence of viscosity this instability results in a finite time singularity of the interface far before rupture (in which the interfacial curvature diverges). The dynamics before the instability roughly follow the scaling laws consistent with predictions based on dimensional analysis, though these scaling laws are violated at the instability. Since the dynamics after rupture is completely determined by the shape at the breaking point, the time dependences of recoiling do not follow a simple scaling law. In the process of demonstrating these results, we present detailed comparisons between numerical simulations and experimental drop shapes with excellent agreement.

title: Nonlinear dynamics and breakup of free-surface flows
abstract: Surface-tension-driven flows and, in particular, their tendency to decay spontaneously into drops have long fascinated naturalists, the earliest systematic experiments dating back to the beginning of the 19th century. Linear stability theory governs the onset of breakup and was developed by Rayleigh, Plateau, and Maxwell. However, only recently has attention turned to the nonlinear behavior in the vicinity of the singular point where a drop separates. The increased attention is due to a number of recent and increasingly refined experiments, as well as to a host of technological applications, ranging from printing to mixing and fiber spinning. The description of drop separation becomes possible because jet motion turns out to be effectively governed by one-dimensional equations, which still contain most of the richness of the original dynamics. In addition, an attraction for physicists lies in the fact that the separation singularity is governed by universal scaling laws, which constitute an asymptotic solution of the Navier-Stokes equation before and after breakup. The Navier-Stokes equation is thus continued uniquely through the singularity. At high viscosities, a series of noise-driven instabilities has been observed, which are a nested superposition of singularities of the same universal form. At low viscosities, there is rich scalingmore » behavior in addition to aesthetically pleasing breakup patterns driven by capillary waves. The author reviews the theoretical development of this field alongside recent experimental work, and outlines unsolved problems. « less

title: Universal pinching of 3D axisymmetric free-surface flow.
abstract: We consider the viscous motion of an axisymmetric column of fluid with a free surface. The Navier-Stokes equation forms a singularity as the height of the fluid neck goes to zero. Close to pinchoff, the solutions have a scaling form characterized by a set of universal exponents. The shape of the neck and its velocity field is described by scaling functions, which we predict without adjustable parameters.

title: The bifurcation of liquid bridges
abstract: Details of the shape of the liquid bridge joining a nascent water drop to its parent body are presented for times before, after and at the instant of bifurcation when the drop is created and also when the secondary droplet is formed. After the instant of bifurcation there is ‘unbalanced’ surface tension which gives an impulse to the rest of the fluid causing strong surface deformations. The major point of this work is to draw attention to the strong up–down asymmetry at each bifurcation point. The geometric similarity at each bifurcation instant supports the conjecture that the flow converges to just one similarity solution of the type described by Keller & Miksis (1983) in which only surface tension and inertia are important. Features of the flow before and after bifurcation are discussed.

title: On the separation of droplets from a liquid jet
abstract: The droplet separation from a liquid jet was investigated experimentally. Details of the shape of the thin liquid neck joining the droplet to its parent body were studied in terms of the uid viscosity and the jet diameter. As the viscosity increased, the neck rapidly elongated creating a long thread. Its nal diameter before rupture was approximately one micrometer and seems to be constant within wide range of parameters varied. One or multiple breakups of the micro-thread were observed, which produced micro-satellites, i.e. droplets in a micrometer range. The experimental results only partly connrmed the predictions of Eggers' (1993) similarity solution. The predicted shape of the pinch-oo region well overlaps the long thread observed for very viscous liquids. However, the nal jet diameter, retraction velocity of the thread and presence of multiple breakups diierentiate the experimental evidence from the model expectations.