A mass-conserving volume-of-fluid method: Volume tracking and droplet surface-tension in incompressible isotropic turbulence
Introduction
The computational methods to perform fully-resolved simulations of multiphase flows fall into two main categories: 1. interface capturing methods (ICM), such as the volume of fluid (VoF) [1] and the level-set [2] methods; 2. interface tracking methods (ITM), such as the front tracking [3] and the immersed boundary [4] methods. Both ICM and ITM are designed to compute (capture or track) sharp interfaces, and can be used on a fixed Cartesian mesh. Mass conservation, the ability to compute strong deformations of the interface, and interface topology changes due to break-up and coalescence are highly desirable features in the simulations of multiphase flows, e.g. gas–liquid and liquid–liquid. Interface topology changes need ad hoc modeling in ITM [3], whereas they are directly captured by the ICM [5].
Among the ICM, the level-set is a widely used method where a distance function to the interface is advected with the flow throughout the computational domain in place of the interface itself, thus avoiding the need to advect a discontinuous function. The drawback of the level-set method is that it does not conserve mass even when coupled with a particle tracking method [6]. In the VoF method, instead, the advection equation of the volume fraction is directly solved for, thus, the VoF method could potentially conserve mass exactly. Limits to this potential may only come from the numerics adopted. A color function representing the volume fraction of the reference phase is advected geometrically and the interface is reconstructed from this function typically with a piecewise linear representation.
VoF advection schemes can be broadly classified as either direction split or unsplit schemes. Split advection schemes consist of a sequence of one-dimensional advection and reconstruction steps in each coordinate direction, thus they are algorithmically straightforward to implement compared to multidimensional (unsplit) schemes. On the other hand, unsplit methods have the advantage of only requiring one advection and reconstruction step per time step, however the advection step is often algorithmically complex. This is because unsplit methods require either the computation of a flux polyhedron for each cell face [7] or calculation of polyhedra volumes with non-planar surfaces requiring triangulation [8]. In our experience, these geometric calculations are the computational bottleneck of the VoF advection scheme. Thus, the computational savings of the unsplit versus the split method (if any) are likely to be small. Furthermore, unsplit methods do not necessarily conserve mass to zero machine precision [8]. For these reasons, we chose to adopt a split mass-conserving VoF advection approach.
This paper presents a numerical methodology to perform DNS of droplet-laden incompressible turbulent flows with uniform density and viscosity. This is a necessary step to verify that the numerical method is accurately simulating the motion of finite-size droplets, i.e., volume tracking and surface tension computations, before we can develop a coupled flow solver where density and viscosity variations occur between the droplet and the surrounding fluid (such a flow solver has been developed in Ref. [9]). We adopted the volume of fluid (VoF) method because of its potential to conserve mass with zero-machine accuracy. The VoF advection is performed through a spatially split approach, i.e., the Eulerian implicit - Eulerian algebraic - Lagrangian explicit (EI-EA-LE) algorithm which was originally proposed by Scardovelli et al. [10]. We chose EI-EA-LE over the EILE-3D algorithm [11] because EI-EA-LE requires half the number of advection and reconstruction steps and does not require the calculation of three two-dimensional divergence-free velocity fields. Thus, EI-EA-LE is at least two times computationally faster than EILE-3D. The original EI-EA-LE algorithm [10] is globally mass-conserving but generates wisps, and does not conserve the mass of the individual volumes tracked in the flow. We have improved this method with the addition of a redistribution and a wisp suppression algorithm. Our method is consistent (i.e., the volume of fluid function, C, satisfies the condition ) and wisps-free and, thus, conserves mass both globally and locally within each volume tracked. We also present and analyze the numerical treatment of the surface tension force within a projection method combined with our VoF method. The surface tension effects are treated with the continuum surface force (CSF) approach proposed by Brackbill et al. [12] and adapted to the VoF method by Francois et al. [13].
We present the governing equations for droplet-laden flows in Section 2, the projection method in Section 3, the VoF interface reconstruction and advection algorithms in Section 4, and the method to compute the interface curvature in Section 5. In Section 6, we present the numerical results for the volume-tracking test-cases (zero surface tension), and for the coupled droplet-flow cases (non-zero surface tension). In Section 7, we give the concluding remarks.
Section snippets
Governing equations
The non-dimensional governing equations for a droplet-laden incompressible flow are the continuity and momentum (Navier–Stokes) equations,where is the fluid velocity, is the pressure, is the density, is the viscosity, and is the force per unit volume due to the surface tension,In Eq. (3), is the surface tension coefficient, is the curvature of the interface between the two fluids (i.e., droplet and
Projection method
The Navier–Stokes equations, (6), (7), are solved numerically using the projection method. Time integration of (7) from time to is performed using the second-order Adams–Bashforth scheme,In Eq. 8, is a non-divergence-free fluid velocity approximating , andThe divergence-free condition on the updated fluid velocity is imposed by solving the Poisson equation for pressure,and by updating the fluid
Volume-of-fluid (VoF) method
This section describes the numerical method for tracking the droplets in space and time with the volume-of-fluid method. In the volume of fluid (VoF) method, the sharp interface between the two phases (e.g. liquid and gas) is determined using the VoF function, C, that represents the volume fraction of the reference phase in each computational cell, e.g. in the gas, and in the liquid (). The VoF method requires two steps, the VoF interface reconstruction and the VoF advection which
Curvature computation
This section describes the numerical method for computing the curvature of the interface between the two fluids that is needed for calculating the surface tension force, (Eqs. (7), (12)). In order to compute the interface curvature from the VoF function, C, the height function technique by Cummins et al. [19] has been adopted and investigated. A correction based on the local orientation of the interface normal is performed to minimize the error on a spherical interface [20]. The technique is
Results
In this section, we first report the numerical results on the accuracy of the VoF interface reconstruction (Section 6.1) and VoF advection (Section 6.2), described in Section 4, for several test-cases where volumes of initially spherical shapes are tracked without surface tension in specified velocity fields. This assesses the properties of the VoF method independently from the two-way coupling between the droplet and the surrounding fluid by setting the surface tension force to zero. Then, in
Concluding remarks
We have presented a coupled volume-of-fluid/projection method to simulate droplet-laden incompressible turbulent flows with uniform density and viscosity. We adopted a three-dimensional direction-split volume of fluid (VoF) method for tracking volumes accurately in incompressible velocity fields. Our VoF method relies on the MYC interface reconstruction [11] and the EI–EA–LE advection scheme [10] to which we added a redistribution and wisp suppression algorithms (Section 4.2). EI–EA–LE has the
Acknowledgments
This work was supported by the National Science Foundation CAREER Award, Grant No. OCI-1054591. The numerical simulations were performed in part on Hyak, high-performance computer cluster at the University of Washington, Seattle, and in part on the XSEDE computational resources provided by the National Institute for Computational Sciences (NICS) at the Oak Ridge National Laboratory, under XRAC Grant No. TG-CTS100024. We specifically acknowledge the assistance of the XSEDE Extended Collaborative
References (26)
- et al.
Modelling merging and fragmentation in multiphase flows with SURFER
J Comput Phys
(1994) - et al.
Fronts propagating with curvature-dependent speed: algorithms based on Hamilton–Jacobi formulations
J Comput Phys
(1988) - et al.
A front-tracking method for viscous, incompressible, multi-fluid flows
J Comput Phys
(1992) - et al.
Direct numerical simulation of evaporating droplets
J Comput Phys
(2008) - et al.
A hybrid particle level set method for improved interface capturing
J Comput Phys
(2002) - et al.
Combination of 3D unsplit forward and backward volume-of-fluid transport and coupling to the level set method
J Comput Phys
(2013) - et al.
Interface reconstruction with least-squares fit and split advection in three-dimensional Cartesian geometry
J Comput Phys
(2007) - et al.
A continuum method for modeling surface tension
J Comput Phys
(1992) - et al.
A balanced-force algorithm for continuous and sharp interfacial surface tension models within a volume tracking framework
J Comput Phys
(2006) - et al.
Analytical and geometrical tools for 3D volume of fluid methods in general grids
J Comput Phys
(2008)
A new volume of fluid advection algorithm: the stream scheme
J Comput Phys
Estimating curvature from volume fractions
Comput Struct
On reducing interface curvature computation errors in the height function technique
J Comput Phys
Cited by (50)
A novel steepness-adjustable harmonic volume-of-fluid method for interface capturing
2024, Journal of Computational PhysicsShear-thinning stimulative fluid breakup in 3D pore-throat
2024, Chemical Engineering Research and DesignTwo-phase free-surface flow interaction with moving bodies using a consistent, momentum preserving method
2023, Journal of Computational PhysicsComparison of turbulent drop breakup in an emulsification device and homogeneous isotropic turbulence: Insights from numerical experiments
2023, Colloids and Surfaces A: Physicochemical and Engineering AspectsA review of VOF methods for simulating bubble dynamics
2022, Progress in Nuclear Energy
- 1
Ph.D. student, William E. Boeing Department of Aeronautics & Astronautics, University of Washington, Seattle, USA.
- 2
Assistant Professor, William E. Boeing Department of Aeronautics & Astronautics, University of Washington, Seattle, USA.