Numerical Study of the Lift Force, Velocities and Pressure Distribution of a Single Air Bubble and Two Interacting Air Bubbles Rising in Quiescent Liquid

doi:10.4236/ojfd.2020.101003

Open Journal of Fluid Dynamics
Vol.10 No.01(2020), Article ID:98684,21 pages
10.4236/ojfd.2020.101003

Chao Guan^1*, Shinichiro Yanase¹^#, Koji Matsuura², Toshinori Kouchi¹, Yasunori Nagata¹

●How to Cite this Article

¹Graduate School of Natural Science and Technology, Okayama University, Okayama, Japan

²Department of Biomedical Engineering, Faculty of Engineering, Okayama University of Science, Okayama, Japan

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

Received: January 19, 2020; Accepted: March 2, 2020; Published: March 5, 2020

ABSTRACT

The study of multiphase flow consisting of liquid and air bubbles has been attracting the interest of many researchers. Numerical methods for such a system are, however, facing difficulty in numerical accuracy and a heavy computational load. In this paper, we made corrections to the modified force-coupling method in our previous papers and applied it to the numerical studies of a single air bubble rising near a vertical wall and two interacting air bubbles rising in line in quiescent liquid. Corrections were made to the effective ranges of the force-coupling method. The calculation results showed that the lift force acting on an air bubble obtained by the experimental data was more accurately reproduced than those by our previous method. We accurately calculated the time evolution of the velocities of interacting two air bubbles rising in line obtained in the previous experiments and resolved the physical mechanism of the relative movement of two bubbles. We also found the present method is much quicker and needs much smaller memory capacity than other methods, such as the volume of fluid method.

Keywords:

Force-Coupling Method, Multiphase Flow, Bubble Motion, Two Bubble Interaction

1. Introduction

The study of multiphase flow of liquid (water) and air bubbles has been attracting the interest of many researchers. Experimental investigations of the flow including microbubbles are now well developed regarding engineering applications, such as drag reduction of pipe flow [1], cleaning of polymer ink [2] and separation of carbon fibers [3] using microbubbles. However, numerical investigations of the liquid flow containing air bubbles are facing difficulty in numerical accuracy and a heavy computational load. Multiphase flow has been computationally studied by the volume of fluid method (VOF) [4] [5], the discrete element method (DEM) [6] [7], the immersed boundary method (IBM) [8], the moving particle semi-implicit method (MPS) [9], and the theoretical model equations taking into account of inertial and added-mass body acceleration forces [10] [11]. Although VOF can be applied to many complex configurations, it needs much CPU time and heavy memory capacity. IBM proposed by Peskin [8] has been largely applied to solve the liquid motion including solid particles. But it cannot be used to solve particle-liquid flow if the particle density is less than 0.3 of the liquid density [12]. MPS also has a problem in a heavy calculation load.

Takemura et al. performed experimental studies of a rising air bubble [13] and a solid-like air bubble [14] near a vertical wall. There are some works which studied multiphase flow including plural bubbles. At low Reynolds numbers, Katz and Meneveau [15] conducted experiments of the interacting air bubbles rising in a quiescent distilled water, and observed the velocities of the two bubbles rising in line. Watanabe and Sanada [16] and Kusuno et al. [17] performed experiments of two rising bubbles in a quiescent silicon oil. Kusuno et al. [17] observed that a trailing bubble collided with a leading bubble, but two bubbles did not coalesce. Gumulya et al. [5] numerically studied two rising bubbles which collided and coalesced using VOF.

The force-coupling method (FCM) originally proposed by Maxey and Patel [18], which we call the original FCM (OFCM) in this paper, requires no moving boundaries of particles, but represents the center position of a particle by the point-body force acting from the particle to the liquid [19] [20] [21] [22]. OFCM has been developed to calculate multiphase flow consisting of liquid and solid particles. However, the applicability of OFCM to air bubbles was notclear [19]. Very recently, Guan et al. [23] [24] proposed the modification of OFCM in order to simulate the multiphase flow composed of liquid and air bubbles, where they called it the modified FCM (MFCM), and obtained fairly good agreement between the numerical results and the experimental data by Takemura et al. [13] [14]. However, inconsistency remained in the definition of the effective ranges for the force terms of MFCM in the previous papers [23] [24], where the values valid for flow containing solid particles were used for that containing air bubbles. In the present paper, we proposed the calculation method with a more appropriate definition of the effective ranges. Specifically, the effective range of the force monopole term and that of the force dipole term of OFCM were changed to be applied to air bubbles, which had not been made in the previous papers [23] [24]. As a result, we reproduced more accurately the results of the experimental data of Takemura et al. [13] than those in the previous study [23] for a single air bubble rising near a vertical wall, and the calculated velocities of rising interacting two air bubbles in line were very close to those in the experimental study by Katz and Meneveau [15].

2. Formulation of the Modified FCM with Corrections

In OFCM and MFCM, the incompressible Navier-Stokes equations with a body-force term, $f_{b, i}$ , representing the interaction of the liquid with bubbles,

$\frac{\partial u_{i}}{\partial t} = - u_{j} \frac{\partial u_{i}}{\partial x_{j}} - \frac{1}{ρ_{f}} \frac{\partial p}{\partial x_{i}} + ν \frac{\partial^{2} u_{i}}{\partial x_{j}^{2}} + g_{i} + f_{b, i}, (i = 1, 2, 3)$ (1)

where $i = 1, 2, 3$ indicate the x, y, z-direction components, respectively, and the equation of continuity,

$\frac{\partial u_{i}}{\partial x_{i}} = 0,$ (2)

are the basic equations, where the Einstein’s summation convention is used. $u_{i}$ is the velocity vector of the liquid, t the time, $ρ_{f}$ the density of the liquid, p the pressure, $ν$ the kinematic viscosity of the liquid, and $g_{i}$ the gravitational acceleration vector. $f_{b, i}$ denotes the body-force acting on the liquid from N bubbles, which is given by

$f_{b, i} = \sum_{n = 1}^{N} {F_{i}^{(n)} Δ_{M}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{M}^{(n)}) + G_{i j}^{(n)} \frac{\partial}{\partial x_{j}} Δ_{D}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{D}^{(n)})} .$ (3)

$f_{b, i}$ consists of two parts, the force monopole term (FMT, the first term on the right-hand side of Equation (3)) and the force dipole term (FDT, the second term on the right-hand side of Equation (3)). Here, $x_{i}$ is the position vector in the liquid, and $Y_{i}^{(n)}$ is the position vector of the center of the n-th bubble. $σ_{M}^{(n)}$ and $σ_{D}^{(n)}$ are the effective ranges of FMT and FDT around the center of the particle, respectively, of the n-th bubble, and $Δ_{M}^{(n)}$ and $Δ_{D}^{(n)}$ are the force distribution of FMT and FDT, respectively, given by Maxey and Patel [18] and Lomholt et al. [20] as

$Δ_{M}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{M}^{(n)}) = \frac{1}{{[2 π {(σ_{M}^{(n)})}^{2}]}^{3 / 2}} \exp [- \frac{{(x_{i} - Y_{i}^{(n)})}^{2}}{2 {(σ_{M}^{(n)})}^{2}}],$ (4)

$Δ_{D}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{D}^{(n)}) = \frac{1}{{[2 π {(σ_{D}^{(n)})}^{2}]}^{3 / 2}} \exp [- \frac{{(x_{i} - Y_{i}^{(n)})}^{2}}{2 {(σ_{D}^{(n)})}^{2}}] .$ (5)

$σ_{M}^{(n)}$ and $σ_{D}^{(n)}$ are related to the radius of the n-th particle $R^{(n)}$ as $σ_{M}^{(n)} = R^{(n)} / \sqrt{π}$ and $σ_{D}^{(n)} = R^{(n)} / {(6 \sqrt{π})}^{1 / 3}$ in OFCM [18] [19] [20] [21] [22]. However, because these values of $σ_{M}^{(n)}$ and $σ_{D}^{(n)}$ were obtained by considering the theoretical solution of the Stokes approximation of the flow around a solid particle, they are not guaranteed to be applied to an air bubble.

Therefore, we proposed in this paper that $σ_{M}^{(n)} = R^{(n)} / \sqrt{1.88 π}$ and $σ_{D}^{(n)} = R^{(n)} / {(3 \sqrt{π})}^{1 / 3}$ , where the former was obtained in consideration of the theoretical solution of the Stokes approximation of the flow around an air bubble, while the latter of the theoretical solution of dissipation rate for an air bubble by Batchelor [25] (§4.11). Please refer Appendix for the evaluation of $σ_{M}^{(n)}$ and $σ_{D}^{(n)}$ .These are the corrections newly introduced in the present paper to OFCM [18] [19] [20] [21] [22] and MFCM [23] [24]. In this paper, we assumed that the radius of each bubble, $R^{(n)}$ , was the same, so that we put $R^{(n)} = R$ .

In Equation (3), FMT represents the reaction to the drag force acting on the bubble for its translational motion. $F_{i}^{(n)}$ is the force acting on the liquid from the n-th bubble, which is given by

$F_{i}^{(n)} = \frac{4}{3} π {(R^{(n)})}^{3} (ρ_{f} - ρ_{b}) (g_{i} - \frac{d U_{i}^{(n)}}{d t}),$ (6)

where $ρ_{b}$ is the density of the bubble, $U_{i}^{(n)}$ the velocity of the n-th bubble. $U_{i}^{(n)}$ was given by the following equation in terms of the Gaussian distribution function.

$U_{i}^{(n)} = ∭ u_{i} (x_{i}, t) Δ_{M}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{M}^{(n)}) d x_{1} d x_{2} d x_{3},$ (7)

in OFCM [15]. Equation (7) means that the particle (bubble) moves with the liquid velocity averaged over the region, where the body force of FMT is effective. In the present study, we used the similar method as Guan et al. [23], and derived the following equation for $U_{i}^{(n)}$ according to Mei et al. [26] :

$U_{i}^{(n)} = ∭ [u_{i} (x_{i}, t) Δ_{M}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{M}^{(n)}) \times \frac{1}{1 + {[\frac{8}{R e_{\infty}^{(n)}} + 0.5 (1 + \frac{3.315}{{(R e_{\infty}^{(n)})}^{0.5}})]}^{- 1}}] d x_{1} d x_{2} d x_{3},$ (8)

in order to apply MFCM over a wide range of the bubble Reynolds number. Here, $R e_{\infty}^{(n)}$ is the bubble Reynolds number determined by the terminal velocity of a bubblerising in quiescent liquid, $U_{\infty}^{(n)}$ , given by

$R e_{\infty}^{(n)} = \frac{2 R^{(n)} U_{\infty}^{(n)}}{ν} = \frac{2}{3} \frac{g_{i} {(R^{(n)})}^{3}}{ν^{2}} \cdot \frac{1}{1 + {[\frac{8}{R e_{\infty}^{(n)}} + 0.5 (1 + \frac{3.315}{{(R e_{\infty}^{(n)})}^{0.5}})]}^{- 1}} .$ (9)

Note that in Guan et al. [23], the equation corresponding to Equation (8) (Equation (12) in Guan et al. [23] ) was developed for the flow including a solid particle. This is a new correction introduced in the present paper for the air bubble case.

In Equation (3), FDT represents the influence of the velocity gradient around the bubble. $G_{i j}^{(n)}$ in Equation (3) is $G_{i j}^{(n)} = T_{i j}^{(n)} + S_{i j}^{(n)}$ , where $T_{i j}^{(n)}$ is an anti-symmetric and $S_{i j}^{(n)}$ is symmetric parts, respectively. $T_{i j}^{(n)}$ is related to the torque to the fluid from the bubble and is given by the equation,

$T_{i j}^{(n)} = \frac{1}{2} ϵ_{i j k} T_{k}^{(n) (e x t)},$ (10)

where $ϵ_{i j k}$ is the Eddington’s alternating sign tensor, and $T_{k}^{(n) (e x t)}$ is calculated by use of $ρ_{f}$ , $ρ_{b}$ , and $R^{(n)}$ [23]. $S_{i j}^{(n)}$ represents a stresslet acting on the fluid and has been already given by the following equation for an air bubble [23].

$S_{i j}^{(n)} = \frac{8}{3} π ν ρ_{f} {(R^{(n)})}^{3} E_{i j}^{(n)},$ (11)

where $E_{i j}^{(n)}$ is the rate of strain, and determined to satisfy the constraint, for each particle using an iteration method in order to keep the bubble shape a sphere.

$E_{i j}^{(n)} = \frac{1}{2} \int^{} (\frac{\partial u_{i}}{\partial x_{j}} + \frac{\partial u_{j}}{\partial x_{i}}) Δ_{D}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{D}^{(n)}) d x_{1} d x_{2} d x_{3} = 0,$ (12)

The calculation domain of the present study was a cube with 0.05-m-length of each side, as shown in Figure 1. The size of the computational domain was the same to the experimental ones [13]. The gravitational acceleration, $g = | g | = 9.81 m / s^{2}$ , is directed to the negative x-direction. In the numerical calculation, the Fourier transform was used in the x-direction with the collocation method, where the number of Fourier modes was $n_{x} = 128$ . The finite-difference method was used with equally spaced staggered grids, where the number of grids $n_{y} = 128$ in the y-direction and $n_{z} = 128$ in the z-direction. The kinematic viscosity, $ν$ , and the bubble radius, R, were taken to be equal to those of the experimental conditions.

The calculation procedure of the MFCM is summarized as follows:

1) Give the initial conditions.

2) Get $F_{i}^{(n)}$ from the initial conditions (Equation (6)).

3) Obtain the FMT from $Δ_{M}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{M}^{(n)})$ , and $F_{i}^{(n)}$ (Equations (4), and (6)).

Figure 1. Computational domain and coordinate system.

4) Get $G_{i j}^{(n)}$ from $Δ_{D}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{D}^{(n)})$ , $T_{i j}^{(n)}$ , $S_{i j}^{(n)}$ , and $E_{i j}^{(n)}$ (Equations (5), (10), (11), and (12)).

5) Obtain the FDT from $Δ_{D}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{D}^{(n)})$ , and $G_{i j}^{(n)}$ (Equations (5) and Step 4).

6) Acquire $f_{b, i}$ from the FMT and the FDT (Equation (3)).

7) Solve the Navier-Stokes equations including $f_{b, i}$ and obtain the velocity field $u_{i} (x_{i}, t)$ (Equations (1), and (2)).

8) Obtain the bubble velocity $U_{i}^{(n)}$ from $u_{i} (x_{i}, t)$ and $Δ_{M}^{(n)} (x_{i} - Y_{i}^{(n)}, σ_{M}^{(n)})$ (Equations (8), (9) and Step 3).

9) Solve the next-time step bubble position from $U_{i}^{(n)}$ and $Δ t$ .

10) Return to Step 2.

3. Comparison of the Numerical Results with the Experimental Data for a Single Air Bubble Rising near a Vertical Wall

To show better accuracy of the proposed method applied to the bubble motion in quiescent liquid than our previous method [23], comparison was conducted between the results of the numerical calculation using the present method and the experimental data by Takemura et al. [13] for the motion of a single air bubble rising near a vertical wall. Since the number of a concerning bubble was unity, $N = 1$ was set. Figure 2 shows the test sections in the y-z plane and the position of the bubble in the experiment [13] and the present numerical calculation, where $L (t)$ denotes the distance from the bubble center to the nearest vertical wall. Note that L is a function of t. The initial position of the bubble center was put at the point $(0, 0, 0.025 - L_{0})$ , where $L_{0}$ is $L (0) = 2 R$ .

Takemura et al. [13] presented the values related to the lift force by observing the motion of the bubble in their paper. We also calculated those values from the motion of the bubble obtained by our simulations for comparison. They introduced the wall distance Reynolds number, $R e_{L}^{(1)}$ , defined in terms of $L (t)$ as

$R e_{L}^{(1)} = \frac{L (t) U_{\infty}^{(1)}}{ν} .$ (13)

Figure 2. Experimental test section and computational domain.

Note that $R e_{L}^{(1)}$ is also a function of t. Takemura et al. [13] also evaluated the inverse Fourier transform of the lift force on a bubble, $I_{L}$ , using the following equation:

$I_{L} = \frac{2}{R e_{\infty}^{(1)}} \cdot \frac{U_{3}^{(1)}}{U_{1}^{(1)}} \cdot {1 + {[\frac{8}{R e_{\infty}^{(1)}} + \frac{1}{2} (1 + \frac{3.315}{{(R e_{\infty}^{(1)})}^{0.5}})]}^{- 1} + \frac{R}{L (t)} R e_{L}^{(1)} I_{D}},$ (14)

where $I_{D}$ is the inverse Fourier transform of the horizontal component of the drag force. We evaluated $I_{L}$ in the same manner of Takemura et al. [13] and compared it with the experimental ones. Note $I_{D}$ was theoretically evaluated under the assumption of the Oseen approximation [13] [27].

Figure 3(a) and Figure 3(b) show the values of $I_{L}$ as a function of $R e_{L}^{(1)}$ for $R e_{\infty}^{(1)} = 2.1$ and 15, respectively. The black symbols indicate the experimental data of Takemura et al. [13]. The black dashed lines show $I_{L}$ by Guan et al. [23] in which the numerical simulations of this case were conducted without the correction of the effective ranges, while the black solid lines show $I_{L}$ of the present numerical simulations using MFCM with the correction of the effective ranges. It is found that the present results are more accurate than those in the

(a)(b)

Figure 3. The values of $I_{L}$ as a function of $R e_{L}^{(1)}$ for $R e_{\infty}^{(1)} = 2.1$ and 15. (a) $R e_{\infty}^{(1)} = 2.1$ ; (b) $R e_{\infty}^{(1)} = 15$ .

previous study [23]. Note that the maximum discrepancy between the experimental data [13] and the present numerical results for $R e_{\infty}^{(1)} = 2.1$ was 5%and that for $R e_{\infty}^{(1)} = 15$ was 28%. The present numerical predictions were much better than the previous ones [23]. Thus, the present MFCM proved to be applicable to the calculation of multiphase flow consisting of water and an air bubble.

4. Comparison of the Numerical Results with the Experimental Data for Two Interacting Air Bubbles Rising in Quiescent Liquid

Then, we applied the present MFCM with the corrections of the effective ranges to the case of two interacting bubbles rising in quiescent liquid and compared the results with the experimental data of Katz and Meneveau [15]. For this case, $N = 2$ was set. Figure 4 shows the initial positions of two bubbles on the vertical line, $y = z = 0$ , in the computational domain of the present study. The initial position of the lower bubble center was set at the point $(0, 0, 0)$ , and that of the upper bubble center at the point $(0, 0, S_{0})$ , where $S_{0}$ was the initial value of S, the distance between two bubbles, and was taken to be equal to the experimental value of Katz and Meneveau [15]. Note that two bubbles moved only vertically.

Figures 5(a)-(c) show the dimensionless bubble-rising velocities, $U_{1}^{(1)} (t) / U_{\infty}^{(1)}$ and $U_{1}^{(2)} (t) / U_{\infty}^{(2)}$ , as a function of a dimensionless distance between two bubbles, S/R, for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ , 19 and 35, respectively. Solid and blank circle symbols indicate the dimensionless bubble velocities of the experimental value [15], and red long-dashed and dash-dotted lines are those obtained by the present calculations.

Figure 5(a) shows that the results by the present MFCM are in good agreement with the experimental values for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ . It is interesting that the velocities of the two bubbles became faster than the terminal velocity of the one bubble case when two bubbles were present. In the experiment [15] two bubbles

Figure 4. Initial position of two bubbles in the computational domain.

Figure 5. Change of the velocities of two bubbles as a function of the nondimensional distance between two bubbles for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ , 19, and 35. (a) $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ ; (b) $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 19$ ; (c) $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 35$ .

finally coalesced, while in the present numerical calculations they were stuck and moved together. In Figure 6, the change of the distance between two bubbles is plotted as a function of time.

For $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 19$ and 35, the present numerical simulations could reproduce the experimental data [15] very well if two bubbles were not very close as shown in Figure 5(b) and Figure 5(c). However, when two bubbles were very close, discrepancy between the present calculation and the experimental data

Figure 6. Change of the distance between two bubbles as a function of time for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ .

happened for these bubble Reynolds numbers, because if two bubbles approached, deformation of the bubble surface and coalescence of bubbles might occur, which was not considered in the present calculation. Note also that the experiments of Katz and Meneveau [15] were on many (more than two) rising air bubbles aligned like a chain, which might cause a little difference between the present simulations and their experiments.

To understand the physical mechanism of the acceleration of the bubble velocity and the approach of two bubbles, we plot the pressure contours (the unit is Pascal [Pa]), the pressure distribution along the line $y = z = 0$ and the vertical velocity distribution for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ in Figure 7. When two bubbles were apart from each other at t = 0.05 s as shown in Figure 7(a1) and Figure 7(a2), the pressure distribution and vertical velocity distribution were similar for two bubbles. If two bubbles came closer at t = 2.3 s as shown in Figure 7(b1) and Figure 7(b2), the high pressure region of the upper side of the lower bubble shrunk because the lower bubble entered the wake region of the upper bubble, which caused the acceleration of the lower bubble velocity. When two bubbles were very close at t = 2.9 s as shown in Figure 7(c1) and Figure 7(c2), pressure variation between two bubbles became smaller, which implies that two bubbles tended to move together. Remark that the maximum rising velocity is $2 U_{\infty} (= 2 U_{\infty}^{(1)} = 2 U_{\infty}^{(2)})$ , which is realized if two bubbles coalesce to form a single spherical bubble whose volume is twice of each bubble. In Figure 7(a1) and Figure 7(a2) to Figure 7(c1) and Figure 7(c2), we plotted two circles on the abscissa which show the radius, R, of each bubble. It is reasonable that the pressure peak occurred at the head surface of each bubble and the valley at the tail. Therefore, we found that the present MFCM could approximately exhibit the boundary of the bubble because the effective ranges of FMT and FDT, both of which are proportional to R, were properly set in the present method.

(a1) (a2) (b1) (b2) (c1) (c2)

Figure 7. Pressure contours and velocity distribution in the plane $y = 0$ , and pressure distribution along the line $y = z = 0$ for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ . Circles plotted on the abscissa show the surface of each bubbles. (a1) Pressure contours and vertical velocity distribution around two bubbles in the plane $y = 0$ at $S / R = 15$ and $t = 0.05 s$ for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ . (a2) Pressure distribution along the line $y = z = 0$ , at $S / R = 15$ and $t = 0.05 s$ for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ . (b1) Pressure contours and vertical velocity distribution around two bubbles in the plane $y = 0$ at $S / R = 8$ $t = 2.3 s$ for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ . (b2) Pressure distribution along the line $y = z = 0$ , at $S / R = 8$ and $t = 2.3 s$ for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ . (c1) Pressure contours and vertical velocity distribution around two bubbles in the plane $y = 0$ at $S / R = 3$ $t = 2.9 s$ for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ . (c2) Pressure distribution along the line $y = z = 0$ , at $S / R = 3$ and $t = 2.9 s$ for $R e_{\infty}^{(1)} = R e_{\infty}^{(2)} = 3$ .

5. Comparison of CPU Time and Memory Capacity between MFCM and VOF

We conducted a calculation of mutually interacting two air bubbles rising in linetreated in § 4 by using MFCM as well as VOF [28] to show the quickness of the present MFCM for the simulation of multiphase flow. For both VOF and MFCM calculations, we used the same computer as shown in Table 1. In the present case, we did not use the parallel computing system.

Figure 8 shows the typical results for both calculations. For this case, the bubble Reynolds number of both the bubbles was set to be 3. Figure 8 shows the distance between the two bubbles once increased and gradually decreased. Both calculations agreed with each other. We found that VOF required 28 grids (using regular grid) in one dimension within the diameter of a bubble for an accurate calculation, while MFCM required only 6 grids to capture the accurate motion of the bubble. An increase in the number of the grid for VOF required much larger computer memory and CPU time for calculations.

Required CPU time and memory for both calculations are also listed in Table 2. We found that VOF required over 50 times CPU time of MFCM and over 100 times memory of MFCM. Table 2 clearly shows the advantage of MFCM to study the multiphase flow including many air bubbles numerically with keeping the ability of the accurate prediction of bubble orbits. Remark that the study of rising two bubbles in line were dealt with numerically by VOF [4] [5], MPS [9] and the use of the theoretical model equations [10] [11], but we applied the present MFCM to this case, and found that the present method needed a very small calculation load, and resulted in very accurate solutions and a good physical insight into the interaction of two bubbles.

Table 1. The specifications of the computer used in this study.

Table 2. Required CPU time and memory for two bubbles calculation using MFCM and VOF.

Figure 8. Change of the distance between two bubbles as a function of time for.

6. Conclusion

In the present paper, we changed the effective ranges of the force monopole and force dipole terms in the modified force-coupling method (MFCM) proposed by Guan et al. [23] to calculate multiphase flow composed of liquid (water) and air bubbles more accurately. First, we performed numerical simulations of a single air bubble rising near a vertical wall. We found that the present method reproduced the experimental data more accurately than our previous method. Second, we conducted numerical simulations of mutually interacting two air bubbles rising in line. We showed that the present method accurately reproduced the experimental data over a wide range of the bubble Reynolds number. The good agreement proved that our new method is proper. Third, we compared the CPU time and memory needed for the present method with those for VOF. We found that the present method was much quicker and needed much smaller memory capacity. This indicates that our new method has the potential to provide accurate physical insights into the interaction of many bubbles with a much smaller computational load.

Acknowledgements

This study was supported by a grant-in-aid for Scientific Research for Fundamental Sciences (No. 18K12069 to K. M.) from The Ministry of Education, Culture, Sports, Science and Technology (MEXT) of Japan.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Cite this paper

Guan, C., Yanase, S., Matsuura, K., Kouchi, T. and Nagata, Y. (2020) Numerical Study of the Lift Force, Velocities and Pressure Distribution of a Single Air Bubble and Two Interacting Air Bubbles Rising in Quiescent Liquid. Open Journal of Fluid Dynamics, 10, 31-51. https://doi.org/10.4236/ojfd.2020.101003

References

1. Shatat, M.M.E., Yanase, S., Takami, T. and Hyakutake, T. (2009) Drag Reduction Effects of Micro-Bubbles in Straight and Helical Pipes. Journal of Fluid Science and Technology, 4, 156-167. https://doi.org/10.1299/jfst.4.156

2. Matsuura, K., Ogawa, S., Kasaki, S., Koyama, K., Kodama, M. and Yanase S. (2015) Cleaning Polymer Ink from a Glass Substrate Using Microbubbles Generated by Hydrogen Bubble Method. Separation and Purification Technology, 142, 242-250.https://doi.org/10.1016/j.seppur.2015.01.009

3. Matsuura, K., Uchida, T., Guan, C. and Yanase, S. (2018) Separation of Carbon Fibers in Water Using Microbubbles Generated by Hydrogen Bubble Method. Separation and Purification Technology, 190, 190-194.https://doi.org/10.1016/j.seppur.2017.08.065

4. Hirt, C.W. and Nichols, B.D. (1981) Volume of Fluid (VOF) Method for the Dynamics of Free Boundaries. Journal of Computational Physics, 39, 201-225.https://doi.org/10.1016/0021-9991(81)90145-5

5. Gumulya, M., Utikar, R.P., Evans, G.M., Joshi, J.B. and Pareek, V. (2017) Interaction of Bubbles Rising Inline in Quiescent Liquid. Chemical Engineering Science, 166, 1-10. https://doi.org/10.1016/j.ces.2017.03.013

6. Cundall, P.A. and Strack, O.D.L. (1979) A Discrete Numerical Model for Granular Assemblies. Geotechnique, 29, 47-65. https://doi.org/10.1680/geot.1979.29.1.47

7. Fuchimoto, T., Yanase, S., Mizushima, J. and Senda, J. (2009) Dynamics of Vortex Rings in the Spray from a Swirl Injector. Fluid Dynamics Research, 41, Article ID: 045503. https://doi.org/10.1088/0169-5983/41/4/045503

8. Peskin, C.S. (1972) Flow Patterns around Heart Valves: A Digital Computer Method for Solving the Equations of Motion. PhD Thesis, Albert Einstein College of Medicine, University Microfilms, 72-30.

9. Chen, R.H., Tian, W.X., Su, G.H., Qiu, S.Z., Ishiwatari, Y. and Oka, Y. (2011) Numerical Investigation on Coalescence of Bubble Pairs Rising in a Stagnant Liquid. Chemical Engineering Science, 66, 5055-5063.https://doi.org/10.1016/j.ces.2011.06.058

10. Zhang, J. and Fan, L-S. (2003) On the Rise Velocity of an Interactive Bubble in Liquids. Chemical Engineering Journal, 92, 169-176.https://doi.org/10.1016/S1385-8947(02)00189-4

11. Ramírez-Munoz, J., Baz-Rodríguez, S., Salinas-Rodríguez, E., Castellanos-Sahagún, E. and Puebla, H. (2013) Forces on Aligned Rising Spherical Bubbles at Low-to-Moderate Reynolds Number. Physics of Fluids, 25, Article ID: 093303.https://doi.org/10.1063/1.4822183

12. Kempe, T. and Frohlich, J. (2012) An Improved Immersed Boundary Method with Direct Forcing for the Simulation of Particle Laden Flows. Journal of Computational Physics, 231, 3663-3684. https://doi.org/10.1016/j.jcp.2012.01.021

13. Takemura, F., Takagi, S., Magnaudet, J. and Matsumoto, Y. (2002) Drag and Lift Forces on a Bubble Rising Near a Vertical Wall in a Viscous Liquid. Journal of Fluid Mechanics, 461, 277-300. https://doi.org/10.1017/S0022112002008388

14. Takemura, F., Takagi, S. and Matsumoto, Y. (2002b) Lift Force on a Bubble Rising Near a Wall in Water below Reynolds Number of 20. Transactions of the Japan Society of Mechanical Engineers, Series B, 68, 1684-1690 (In Japanese).https://doi.org/10.1299/kikaib.68.1684

15. Katz, J. and Meneveau, C. (1996) Wake-Induced Relative Motion of Bubbles Rising in Line. International Journal of Multiphase Flow, 22, 239-258.https://doi.org/10.1016/0301-9322(95)00081-X

16. Watanabe, M. and Sanada, T. (2006) In-Line Motion of a Pair of Bubbles in a Viscous Liquid. JSME International Journal Series B, 49, 410-418.https://doi.org/10.1299/jsmeb.49.410

17. Kusuno, H., Yamamoto, H. and Sanada, T. (2019) Lift Force Acting on a Pair of Clean Bubbles Rising in-Line. Physics of Fluids, 31, Article ID: 072105.https://doi.org/10.1063/1.5100183

18. Maxey, M.R. and Patel, B.K. (2001) Localized Force Representations for Particles Sedimenting in Stokes Flow. International Journal of Multiphase Flow, 27, 1603-1626. https://doi.org/10.1016/S0301-9322(01)00014-3

19. Xu, J., Maxey, M.R. and Karniadakis, G.E. (2002) Numerical Simulation of Turbulent Drag Reduction Using Micro-Bubbles. Journal of Fluid Mechanics, 468, 271-281. https://doi.org/10.1017/S0022112002001659

20. Lomholt, S., Stenum, B. and Maxey, M.R. (2002) Experimental Verification of the Force Coupling Method for Particulate Flows. International Journal of Multiphase Flow, 28, 225-246. https://doi.org/10.1016/S0301-9322(01)00045-3

21. Lomholt, S. and Maxey, M.R. (2004) Force-Coupling Method for Particulate Two-Phase Flow: Stokes Flow. Journal of Computational Physics, 184, 381-405.https://doi.org/10.1016/S0021-9991(02)00021-9

22. Yeo, K. and Maxey, M.R. (2010) Simulation of Concentrated Suspensions Using the Force-Coupling Method. Journal of Computational Physics, 229, 2401-2421.https://doi.org/10.1016/j.jcp.2009.11.041

23. Guan, C., Yanase, S., Matsuura, K., Kouchi, T. and Nagata, Y. (2017) Application of the Modified Force-Coupling Method of Tracing the Trajectories of Spherical Bubbles with Solid-Like and Slip Surfaces. Open Journal of Fluid Dynamics, 7, 657-672.https://doi.org/10.4236/ojfd.2017.74043

24. Guan, C., Yanase, S., Matsuura, K., Kouchi, T. and Nagata, Y. (2018) Numerical Study on Motion of a Single Spherical Bubble Near Vertical Wall in the Quiescent Fluid Using the Modified Force-Coupling Method. Nagare, 37, 281-289 (In Japanese).

25. Batchelor, G.K. (1967) An Introduction to Fluid Dynamics. Chapter 4, Cambridge University Press, Cambridge, UK, 174-255.

26. Mei, R., Klausner, J.F. and Lawrence, C.J. (1994) A Note on the History Force on Spherical Bubble at Finite Reynolds Number. Physics of Fluids, 6, 418-420.https://doi.org/10.1063/1.868039

27. Vasseur, P. and Cox, R.G. (1977) The Lateral Migration of Spherical Particle Sedimenting in Stagnant Bounded Fluid. Journal of Fluid Mechanics, 80, 561-591.https://doi.org/10.1017/S0022112077001840

28. OpenFOAM. The Open Source CFD Toolbox. https://www.openfoam.com

Appendix

In this section, we describe the derivation of and for air bubble. that can be applied to air bubble is derived from the following equations. The Stokes equation of the incompressible flow (without considering the inertia term of the Navier-Stokes equation), which acts only on the force of the FMT, is given as the following equation when considering the steady flow,

(A1)

Fourier transform of the above equation gave the following equation,

(A2)

Substituting the Oseen operator and assuming holds, the following equation was obtained,

(A3)

Inverse Fourier transform of the above equation gives the following equation,

(A4)

where, and are given as follows [18],

(A5)

(A6)

where erf represents the error function. Considering the coordinate origin, , the following equation was obtained,

(A7)

(A8)

In the case of Maxey and Patel [18], the velocity of movement of the solid spherical particle under the action of only gravity is considered to be equivalent to the velocity of the spherical particle predicted by Stokes law, and substituting the Oseen operator into the Stokes law equation of air bubble gives,

(A9)

at the coordinate origin, assuming, the above equation is as follows,

(A10)

Substituting Equation (A4) into Equation (A10) gives,

(A11)

Using Equation (A7) and Equation (A8), the above equation is changed to,

(A12)

Solving the above equation gave the following equation,

(A13)

where, assuming that, the above equation is given as follows,

(A14)

therefore,

(A15)

From the above equation, the effective ranges of FMT for air bubble is, but considering other forces (e.g. the added mass force), we found that which was more approximate.

that can be applied to air bubble is derived from the following equations. The Stokes equation of the incompressible flow (without considering the inertia term of the Navier-Stokes equation) acting on the force of the FMT and the FDT is given as the following equation when considering the steady flow.