Elsevier

Computers & Fluids

Volume 48, Issue 1, September 2011, Pages 98-112
Computers & Fluids

Isotropic color gradient for simulating very high-density ratios with a two-phase flow lattice Boltzmann model

https://doi.org/10.1016/j.compfluid.2011.04.001Get rights and content

Abstract

This study presents the integration of isotropic color gradient discretization into a lattice Boltzmann Rothman–Keller (RK) model designed for two-phase flow simulation. The proposed model removes one limitation of the RK model, which concerns the handling of O(1000) large density ratios between the fluids for a wide range of parameters. Taylor’s series expansions are used to characterize the difference between an isotropic gradient discretization and the commonly used anisotropic gradient. The proposed color gradient discretization can reduce, by one order of magnitude, the spurious current problem that affects the interface between the phases. A set of numerical tests is conducted to show that a rotationally invariant discretization enables widening of the parameter range for the surface tension. Surface tensions from O(10−2) to O(10−8), depending on the density ratio, are accurately simulated. An extreme density ratio of O(10,000) is successfully tested for a steady bubble with an error of 0.5% for Laplace’s law across a sharp interface, with a thickness of about 5–6 lattice units.

Introduction

For two decades now, a new method for examining fluids in the molecular state, rather than at the classical macroscopic level, has been developed and used as a tool for numerical flow simulation. This alternative, known as the lattice Boltzmann method (LBM), uses continuous distribution functions [1], and was originally built as an extension of the lattice gas automata [2]. Later, it was discovered that the method can also be regarded as a systematic truncation of the Boltzmann equation in the velocity space [3]. The advantages of the LBM are the way in which it addresses computational issues in complex geometries, and that it is simple to apply, because the streaming operator is linear.

One of the features of the LBM is its great flexibility in dealing with additional complex physics, which is difficult for the traditional CFD methods, based on the Navier–Stokes equations, to handle. Multiphase flow is a typical example of such complex physics. In this case, even some existing LBMs still face problems because of their limited parameter window with respect to surface tension and density ratio, which hampers its application in practical engineering problems. Before setting out the specific goal of this study, we briefly recall the techniques commonly used to predict multiphase flows within the lattice Boltzmann framework.

The lattice Boltzmann (LB) models for immiscible two-phase flows can generally be classified in five groups: Rothman–Keller (RK) [4], [5], [6], [7], [8], [9], [10], [11], [12], [13], Shan–Chen (SC) [14], Free Energy (FE) [15], [16], mean-field (MF) [17], [18], and field mediator (FM) [19].

The original RK model considers two types of fluids, red and blue, each of them with particle distributions following its own LB equation. Because there are two fluids involved, the collision step, intrinsic to any LB formulation, takes into account interactions among particles of the same color, as well as cross interactions between particles of different colors. These latter are related to the surface tension between the two phases, which requires the gradient of the color function. The phase separation is achieved with an additional recoloring step. For the original RK models, this step is seen as cumbersome to program. Recently, however, the ideas of Latva-Kokko and Rothman [7] have made it very easy and straightforward to implement, along with the additional capability of interface thickness adjustment. Work to introduce thermodynamics into the RK model has also been carried out, by Kono et al. [20].

The model developed by Shan and Chen [14] uses an interaction force between the particles to mimic microscopic interactions and automatically separate the concentrated and diluted phases. This spontaneous phase segregation is a feature that has attracted LB practitioners and contributed to the popularity of the SC model.

The free energy method of Swift et al. [15] describes an approach which, in equilibrium, leads to a steady state that can be associated with a free energy. The method uses collision rules, which ensure that the system evolves towards the minimum of an input free energy functional. This functional includes both the pure fluid part and the interface part. The first accounts for the equilibrium between two phases, while the second concerns surface tension. Like the SC model, this formulation is also capable of achieving automatic phase separation.

The MF methods simulate interparticle attraction in the same way as the Coulomb interaction is treated in the Vlasov equation [17]. These models can simulate thermodynamically consistent liquid–vapor systems, where the SC models fail [21].

The FM methods use null-mass particles, the only role of which is to invert the momentum of lattice particles in the transition layer to segregate fluids of different colors [19]. They have the advantage of being able to incorporate binary diffusivity, and this method can therefore be adapted to simulate miscible fluids.

A major problem common to these five methods is that they all appear to suffer from the presence of spurious currents at the interface between the fluids. These currents affect numerical stability, limit the density ratio, and reduce accuracy.

Spurious currents at a two-fluid interface in the RK method have been studied by Ginzbourg and Adler [6], who noted that an appropriate choice of the eigenvalues of the collision matrix leads to the elimination of these currents under certain circumstances. Unfortunately, the authors conclude that these unrealistic currents cannot be eliminated or reduced solely with the choice of eigenvalues. The spurious current phenomenon near an interface has also been mentioned by Halliday et al. in [5], where it was found that these currents are strongest at the interface between two fluids. To reduce them, Dupin et al. [10] proposed adjusting the surface tension between the fluids by disturbing the distribution functions in accordance with the directions of the lattice. This idea seemed promising, but has not been further explored.

In this paper, an improvement to the RK-type model, described in Refs. [4], [22], is proposed which reduces spurious currents. It should be noted that the methodology presented here can be applied to any RK-type LB model.

A characteristic that distinguishes the RK model from other models is the presence of a color gradient. For two-phase LB flows, the color gradient was introduced by Rothman and Keller [23] with their model based on cellular automata. It is by using the color gradient that the interface between the fluids is maintained, and evaluating the color gradient is the key point that we address in this study. Specifically, the goal of this work is to show that an isotropic discretization of the color gradient, that is, with a leading discretization error that is rotationally invariant [24], increases the robustness of the method. This is of primary importance, because it ensures less accumulation of directionally biased discretization errors, and one consequence of using it is a reduction in spurious currents. This ultimately means that it will be possible to solve problems with a higher density ratio and/or with lower surface tension than previously achieved with other RK-type models.

Another challenge facing LBMs designed to simulate two-phase flows is handling the high-density ratio between the phases. Ratios of O(160) were simulated by Kuzmin and Mohamad [25] using a multi-range, multi-relaxation time LB scheme. However, this improved SC method is still unable to simulate a O(1000) air-to-water density ratio. With the projection method, an FE scheme devised by Inamuro et al. [26], [27] succeeded in simulating density ratios of up to O(1000). Unfortunately, this method requires the solution of a costly Poisson equation at every time step. Other methods, such as the one proposed by Lee and Lin [28], can treat high-density ratios, like O(1000), using a special discretization of the LB equation. Although it seems quite complex, this discretization works well, according to the authors. Another promising approach has been developed by Zheng et al. [29], which can also achieve O(1000) density ratios and does not require the solution of a Poisson equation or the implementation of a complex treatment of derivatives. More recently, Becker et al. [30] succeeded in simulating density ratios as high as O(1000) using the LBM and the level set method. However, coupling the LBM and the level set method is not straightforward, because doing so requires knowledge from two different research fields.

The RK model has been shown to be capable of simulating complex two-phase flow problems, but with density ratios reaching only O(80) [4]. Our study will show that this limit can be substantially increased. The method is based on a modified recoloring operator presented by Leclaire et al. [4], which, when combined with the improvement presented in this work, can handle high-density ratios, up to O(10,000), with great accuracy for Laplace’s law. Specifically, the model is based on that of Reis and Phillips [22], with a variant of the recoloring operator from Latva-Kokko and Rothman [7]. This improvement is expected to allow high-density ratios to be addressed for all RK models.

When applying the RK method, anisotropic color gradients are mostly used [5], [11], [10], [9], [7], [6], [31], [13]. In fact, an isotropic color gradient is seldom used, and can only be found in [12], [32]. Unfortunately, in these studies, neither the advantage of using such a gradient discretization, nor the reason for doing so, is explained.

In contrast, for the three LBMs for multiphase flows, FE, MF, and SC, isotropic gradients have been used. For the FE methods, Tiribocchi et al. [33] show that there is a reduction in the spurious current by one order of magnitude when the gradient in the source term of their model is defined using a rotationally invariant discretization. Spurious currents were also reduced in the work of Pooley and Furtado [34] when an isotropic gradient is used, along with the addition of a special forcing term.

More recently, Chiappini et al. [35] proposed an MF model which completely eliminates spurious currents, but, unfortunately, the model does not respect the principle of mass conservation.

Similarly, isotropic discretizations have been introduced for the SC methods [36] for the gradient of the “effective mass function”, represented by the interparticle interactions. Again, it is shown that the more isotropy there is in the gradient, the lower the spurious currents. A detailed description of isotropic gradient stencils and weights can be found in the work of Sbragaglia et al. [37] for 2D and 3D spaces.

In this paper, we present the development of a second order color gradient discretization, while for higher order gradients, the stencils and weights developed by Sbragaglia et al. [37] will be used. This choice will make it possible to improve the accuracy of the RK model presented in [4]. To our knowledge, there is no study showing that the implementation of an isotropic color gradient in RK models provides a significant reduction in spurious currents, widens the range of application of the surface tension, or greatly increases the maximum density ratio handled between the two phases.

This paper is organized in two major sections. In the first, the LB model is given, along with theoretical details to illustrate the difference between isotropic and anisotropic gradients. In the second, a large number (hundreds) of numerical simulations are performed for a steady bubble, combining the various parameters of the model. These tests allow us to analyze the influence of gradient discretization on the solution and to demonstrate the need for an isotropic gradient, and its superiority, versus the standard anisotropic approximation. The primary goal of our paper is to demonstrate the usefulness of an isotropic gradient discretization for the RK model.

Section snippets

Lattice Boltzmann immiscible two-phase model

The current LB approach follows the model of Reis and Phillips [22], along with the improvement presented by Leclaire et al. [4] for the recoloring operator. For the sake of clarity, we recall the model described in Refs. [4], [22]. For this two-dimensional LB model, there are two sets of distribution functions, one for each fluid, moving on a D2Q9 grid with the velocity vectors ci. With θi=π4(4-i), the velocity vectors are defined as (Δx = Δy = 1):ci=(0,0),i=1[sin(θi),cos(θi)],i=2,4,6,8[sin(θi),

Numerical results

The focus of the numerical simulations is on the circular steady bubble, because it is the best test case for evaluating the isotropy of the numerical scheme. In this section, we concentrate on a steady (red) bubble of radius R immersed in a different (blue) fluid.

Only calculations aimed at stationary solutions are performed here, so it is convenient to introduce the stopping criterion common to all simulations:maxall distribution functions|(Nik)(n)-(Nik)(n-1)|ϵwith ϵ = 10−9 and n denoting the

Conclusion

In this study, we showed that the isotropy of the color gradient is crucial to RK-type lattice Boltzmann schemes for the solution of two-phase flow problems. The simulations revealed that, when the isotropy of the color gradient is poor, the surface tension result is not reliable. For a steady bubble problem, we showed that isotropic color gradients up to second order in space are not sufficient to treat low surface tension (bubble deformation), with an average density variation of O(10).

Acknowledgement

The authors thank the anonymous reviewers for their insightful and constructive comments. This work was supported by a grant from FQRNT “Le Fonds Québécois de la Recherche sur la Nature et les Technologies”.

References (41)

  • J. Tölke et al.

    Comput Fluids

    (2006)
  • L. Wu et al.

    Int J Multiphase Flow

    (2008)
  • L. Hao et al.

    J Power Sources

    (2009)
  • X. He et al.

    J Comput Phys

    (1999)
  • K. Kono et al.

    Comput Phys Commun

    (2000)
  • A. Kumar

    J Comput Phys

    (2004)
  • A. Kuzmin et al.

    Comput Math Appl

    (2010)
  • T. Inamuro et al.

    Future Gen Comput Syst

    (2004)
  • T. Inamuro et al.

    J Comput Phys

    (2004)
  • T. Lee et al.

    J Comput Phys

    (2005)
  • H.W. Zheng et al.

    J Comput Phys

    (2006)
  • J. Becker et al.

    Comput Math Appl

    (2009)
  • S. Hou et al.

    J Comput Phys

    (1997)
  • R.R. Nourgaliev et al.

    Int J Multiphase Flow

    (2003)
  • G. McNamara et al.

    Phys Rev Lett

    (1988)
  • U. Frisch et al.

    Phys Rev Lett

    (1986)
  • X. He et al.

    Phys Rev E

    (1997)
  • Leclaire S, Reggio M, Trepanier J. Numerical evaluation of two recoloring operators for an immiscible two-phase flow...
  • I. Halliday et al.

    Phys Rev E

    (1998)
  • I. Ginzbourg et al.

    Trans Porous Media

    (1995)
  • Cited by (0)

    View full text