Abstract and Applied Analysis

Volume 2015, Article ID 706034, 12 pages

http://dx.doi.org/10.1155/2015/706034

## Lattice Boltzmann Simulation of Multiple Bubbles Motion under Gravity

^{1}Institute of Fluid Mechanics, China Jiliang University, Hangzhou 310018, China^{2}Institute of Refrigeration and Cryogenics, Zhejiang University, Hangzhou 310027, China

Received 24 October 2014; Accepted 10 December 2014

Academic Editor: Shuyu Sun

Copyright © 2015 Deming Nie et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The motion of multiple bubbles under gravity in two dimensions is numerically studied through the lattice Boltzmann method for the Eotvos number ranging from 1 to 12. Two kinds of initial arrangement are taken into account: vertical and horizontal arrangement. In both cases the effects of Eotvos number on the bubble coalescence and rising velocity are investigated. For the vertical arrangement, it has been found that the coalescence pattern is similar. The first coalescence always takes place between the two uppermost bubbles. And the last coalescence always takes place between the coalesced bubble and the bottommost bubble. For four bubbles in a horizontal arrangement, the outermost bubbles travel into the wake of the middle bubbles in all cases, which allows the bubbles to coalesce. The coalescence pattern is more complex for the case of eight bubbles, which strongly depends on the Eotvos number.

#### 1. Introduction

Bubble motion is one of the most common gas-liquid flow phenomena and plays an important role in many industrial applications, such as cavitation in fluid machinery, nucleate boiling in reactors, and condenser/evaporator. In medical ultrasound imaging, small encapsulated bubbles called contrast agent are used to enhance the contrast. In thermal inkjet printing, vapor bubbles are used as actuators. They are occasionally used in the microfluidics applications as actuators [1]. The motion of multiple bubbles under gravity is complex due to bubble deformation, coalescence, and breakup. Understanding the dynamic interaction between bubbles is an important aspect of the design and operation of many industrial applications.

A number of investigations of the bubble motion in liquid have been conducted in the past [2–5] due to its scientific and engineering importance. Many numerical methods have been established to simulate the motion of gas bubbles in liquid, such as VOF (volume of fluid) [6] and level set method [7]. In recent decades, the lattice Boltzmann method (LBM) has proved to be a powerful numerical scheme for the simulation of multiphase flow which is based on mesoscopic particle dynamics. Its kinetic nature can provide many of the advantages of molecular dynamics, such as the introduction of repulsive interaction between particles without any boundary conditions for the interface, which makes it more useful than other more conventional methods for numerical modeling of multiphase fluid systems [8]. Several kinds of lattice Boltzmann model for simulating multiphase fluid have been established and applied to the simulations of gas bubbles under gravity successfully. These models include potential method proposed by Shan and Chen [9], color method proposed by Rothman and Keller [10], and free energy method proposed by Swift et al. [11] and improved by Zheng et al. [12]. Recently, Huang et al. [13] evaluated these models’ performance in modeling two-dimensional immiscible two-phase flow in porous media on the pore scale. The free energy method [12] was proved to be a good tool for the study of two-phase flows with high viscosity ratios and high density ratio. Takada et al. [8] used the LBM which was based on the free energy model [11] to simulate single bubble and two bubbles rising under gravity. Their results [8] of single bubble are consistent with those determined by VOF method. They [8] also simulated two bubbles rising in a circular tube full of stationary fluid. It has been found [8] that the two bubbles approach due to the wake formation and then coalesce into a single bubble eventually. Sankaranarayanan et al. [14] simulated the rise behavior of a single bubble in a periodic box by using an implicit formulation of LBM under the conditions of 1 < Eo (Eotvos number) < 10 and (Reynolds number) > 100. The effect of volume fraction on the rise characteristics was analyzed by changing the size of the box relative to that of the bubble. Recently, Amaya-Bower and Lee [15] presented a comprehensive study of the dynamics of a single bubble rising under gravitational force by using lattice Boltzmann method based on the Cahn-Hilliard diffuse interface approach. They [15] analyzed the dependence of terminal bubble shape and velocity on density ratio, viscosity ratio, surface tension, interface thickness, and domain size. In particular, their results [15] showed that the influence of viscosity ratio on the terminal velocity and shape of the bubble is significant, while the density ratio has little effect. Gupta and Kumar [16] used the lattice Boltzmann method based on potential method [9] to study the motion of multiple bubbles under gravity. Two kinds of initial arrangement were considered which are in-line arrangement and staggered arrangement in a vertical orientation, respectively. Their results [16] showed that as the Eotvos number increases, the uppermost bubble deforms the most and the bubbles coalesce eventually due to the wake behind the leading bubble. They [16] also showed the process of bubble coalescence for staggered bubbles in channels, in which lift forces come into play due to the presence of walls. More recently, Yu et al. [17] presented an adaptive lattice Boltzmann method to simulate the interaction between a pair of bubbles with spherical or ellipsoidal shapes under gravity. They observed both attractive and repulsive interactions between bubbles which depend on the relative position and the Reynolds number. They [17] also simulated a group of 14 bubbles and analyzed the effects of the bubble shape and Reynolds number on the spatial distribution of the bubbles.

As shown by Gupta and Kumar [16] and Yu et al. [17], in comparison with the case of single bubble, the motion of multiple bubbles under gravity could be much more complex due to hydrodynamic interactions between bubbles. Besides, the bubble coalescence and break-up could take place occasionally which have great effect on the motion of multiple bubbles. Furthermore, as far as we know, a detailed numerical study of the influence of bubble collision and coalescence on the rising velocity has not been undertaken. Thus, it is important to focus on the fundamental understanding of the bubble collision and coalescence when multiple bubbles are rising under gravity. To fulfil this task, the LBM proposed by Zheng et al. [12] is adopted in this work to study the rise behavior of multiple bubbles which are initially placed in in-line arrangement. The effects of Eotvos number which ranges from 1 to 12 are taken into account. This study will evaluate the coalescence pattern and rising velocity of the bubbles which not only depend on the computational parameters but also depend on the initial arrangements. In addition, the terminal velocity is compared with the result of single bubble under the same conditions to illustrate the influence of multiple bubbles motion.

#### 2. Numerical Method

The discrete lattice Boltzmann equations under external forces for the continuity and momentum equations are given by where is the density distribution function at the th microscopic velocity , is the equilibrium distribution function, is the external force (which is gravity in this work), is the time step, is the relaxation time, is the speed of sound, and are the weights related to the D2Q9 model. According to Zheng et al. [12], is the order parameter that is responsible for the gas-liquid interface. is the chemical potential, which is defined in the following. The macroscopic variables and are determined by the distribution function as follows:The variables and are defined aswhere and are the densities of liquid and gas phase, respectively. The lattice model of D2Q9 is used here and the equilibrium distribution functions are defined asAccording to Zheng et al. [12], the coefficients are defined as The discrete lattice Boltzmann equations for the interface capture equation are given by [12]where is a constant coefficient, which is determined byThe macroscopic variable, that is, the order parameter , is given byAccording to Zheng et al. [12], the chemical potential is computed usingwhere and and are the densities of the liquid phase and the gas phase, respectively. and are parameters related to the thickness of the interface layer and the surface tension coefficient , which are expressed asAccording to Zheng et al. [12] the lattice model of D2Q5 is used for and the equilibrium distribution functions are [12]For simplicity the D2Q5 is used here and the coefficients are taken asThe parameter is used to control the mobility. By Chapman-Enskog analysis, the Navier-Stokes equations and an interface evolution equation can be obtained by (1) and (6)where is the mobility, given by . The viscosity is .

#### 3. Validation

Three cases are carried out to validate the present computational code. The first is the Laplace law, which is given (for the two-dimensional case) bywhere is the pressure jump across the interface and is the bubble radius. In this work, the Laplace law is validated by calculating the pressure jump while varying the bubble radius from to . Other parameters are fixed at , , , and . The interface layer thickness is set to be for and for . The relaxation times are and . The results are shown in Figure 1(a), which shows good agreement between the numerical results and the analytical solution of (14).