- Email: [email protected]
Axisymmetric simulation of the interaction of a rising bubble with a rigid surface in viscous flow
International Journal of Multiphase Flow 52 (2013) 60–70
Contents lists available at SciVerse ScienceDirect
International Journal of Multiphase Flow j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / i j m u l fl o w
Axisymmetric simulation of the interaction of a rising bubble with a rigid surface in viscous flow Tong Qin a,⇑, Saad Ragab a, Pengtao Yue b a b
Department of Engineering Science and Mechanics, Virginia Tech, Blacksburg, VA 24061, USA Department of Mathematics, Virginia Tech, Blacksburg, VA 24061, USA
a r t i c l e
i n f o
Article history: Received 22 June 2012 Received in revised form 4 January 2013 Accepted 5 January 2013 Available online 12 January 2013 Keywords: Bubble Rigid wall Film drainage ALE method Moving mesh
a b s t r a c t The interaction between a rising deformable gas bubble and a solid wall in viscous liquids is investigated by direct numerical simulation via an arbitrary-Lagrangian–Eulerian (ALE) approach. The flow field is assumed to be axisymmetric. The bubble is driven by gravity only and the motion of the gas inside the bubble is neglected. Deformation of the bubble is tracked by a moving triangular mesh and the liquid motion is obtained by solving the Navier–Stokes equations in a finite element framework. To understand the mechanisms of bubble deformation as it interacts with the wall, the interaction process is studied as a function of two dimensionless parameters, namely, the Morton number (Mo) and Bond number (Bo). We study the range of Bo and Mo from (2, 6.5 106) to (16, 0.1). The film drainage process is also considered in this study. It is shown that the deformation of a bubble interacting with a solid wall can be classified into three modes depending on the values of Mo and Bo. Ó 2013 Elsevier Ltd. All rights reserved.
1. Introduction The interaction of a rising bubble with a rigid surface plays an important role in many technological processes (e.g. mineral flotation, fluidized beds and particle sedimentation). The mechanism of bubble-wall collision has not been fully understood. The significance of studying this problem lies not only the bubble deformation when the buoyancy-driven bubble rises, but also in the film drainage process that takes place when the rising bubble impacts a solid wall. For the motion of a bubble in a liquid, Moore (1965) calculated the terminal rising velocity of distorted gas bubbles in a liquid of small viscosity. Grace (1973) and Grace et al. (1976) developed a well-known graphical correlation to explain the shapes and velocities of single gas bubbles rising in infinite liquids. Hnat and Buckmaster (1976) also experimentally studied the shapes and terminal velocities of bubbles rising in liquids and their data are used as benchmarks for comparison with numerical studies (Gueyffier et al., 1999; Ohta et al., 2010). In addition, Bhaga and Weber (1981) observed the flow field around a rising bubble by using the hydrogen bubble tracer technique. They reported that the wake of bubble was closed and laminar for Mo > 4 103 and Reynolds number (Re) < 110 but was open and unsteady for Re > 110. Their work has been compared with numerical simulations by Gaudlitz and Adams (2009) who investigated the wake and shape of a rising
⇑ Corresponding author. E-mail address: [email protected] (T. Qin). 0301-9322/$ - see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijmultiphaseflow.2013.01.001
bubble on a zigzag path. Duineveld (1995) experimentally determined the velocity and shape of rising bubbles in pure water with an extremely low Mo (1011). They found the bubbles can approximately be considered as oblate spheroids with Re 102 and Weber number (We) 10. Many numerical simulations have been carried out for rising bubble, including groups of bubbles. Kushch et al. (2002) computed the motion of bubbles with a finite Weber number through a nearly inviscid liquid. Frank et al. (2006) and Yu and Fan (2008) both reported their research on numerical simulations of bubbles in a viscous liquid using lattice Boltzmann method and level set method, respectively. A review on direction numerical simulations of bubbly flows can be found in Tryggvason et al. (2006). In practice, researchers have paid much attention to the rising gas bubble interacting with surfaces, especially in water. In that case and depending on the bubble size, viscous effects are weak, and much of the kinetic energy can be transformed into surface energy of the bubble, a later release of which may result in a bubble rebound. Using a high-speed video camera, Tsao and Koch (1997) observed a bubble bouncing several times from a horizontal wall before viscosity dissipated the kinetic energy. Klaseboer et al. (2001) experimentally studied the rebound of a drop impinging on a rigid wall and presented a bouncing trajectory mode which included the film drainage effect. Then, Legendre et al. (2005) and Zenit and Legendre (2009) carried out experimental studies on bubble and drops bouncing on a solid wall in viscous liquids. Recently, Zawala et al. (2007) investigated the collision and bouncing of bubbles with various interfaces (water–gas and water–solid interfaces). The effect of the solid mate-
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
rial and surfactant on the bubble-wall interaction has been reported by Fujasova-Zednikova et al. (2010). The dynamics of a bubble or a drop approaching a surface have been considered in a number of theoretical studies. Yiantsios and Davis (1990) analyzed the shape of a drop driven by buoyancy towards a flat rigid surface and a deformable interface. Power (1997) studied the interaction of a deformable bubble with a rigid wall at small Re by solving the Stokes equation. Related numerical work has been reported by Shopov et al. (1990). They investigated the hydrodynamic behavior of a bubble near a rigid surface and the wall of a spherical container, and successfully simulated the dimple formation and the concavity in the rear part of the bubble. However, the lubrication layer formed between the wall and the gas–liquid interface is not thin enough to account for the whole process of film drainage. Mukundakrishnan et al. (2007) numerically studied the wall effects on a rising gas bubble in a liquid-filled finite cylinder. They described the bubble rising in liquids by Eötvös number (Eo) and Mo and investigated the bubble shape regimes by changing the cylinder radius. Many other numerical methods have been used for solving this multiphase problem, such as the volume-of-fluid (VOF) (Rabha and Buwa, 2010), front tracking (Tryggvason et al., 2001), level set (Minev et al., 2003), moving mesh (Wilkes et al., 1999), and phase-field (Yue et al., 2006). Nevertheless, a deformed bubble colliding with a rigid surface, especially in high-Mo liquid, has not been fully studied. To provide more details about a deformable bubble interacting with a rigid surface, the present work uses an arbitrary Lagrangian–Eulerian (ALE) formulation to simulate the bubble deformation and film drainage process. The ALE method can solve the movingboundary problems by using a moving mesh with grid points on the bubble surface (Hu et al., 2001). A parametric investigation which has been conducted revealed the existence of three modes of bubble-wall interactions. The prevailing mode of interaction depends on two parameters, namely the Morton number and the Bond number.
61
the outer boundary of the computational domain, no-slip, slip, or prescribed stress conditions can be applied. The governing equations are solved in axisymmetric geometry on a triangular mesh using a finite element method. The mesh moves and deforms according to an arbitrary Lagrangian–Eulerian (ALE) scheme as described by Hu (1996) and Hu et al. (2001). Further details on the numerical methods can be found in Yue et al. (2007). 3. Problem description and mesh refinement 3.1. Problem description In the present bubble-rigid surface interaction problem, the bubble is rising due to buoyancy only and the rigid surface is horizontal, so the overall motion of bubble happens only in the normal direction to the surface. Thus, this problem can be considered as an axisymmetric problem with the y axis being the axis of symmetry as shown in Fig. 1. The top surface of the computational domain is a rigid surface with no-slip boundary condition. For the side wall of the computational domain, we apply a slip boundary condition. At the bottom surface, we apply stress boundary condition, and specify the ambient pressure to be P0. The pressure in the bubble can be initially obtained by Pb = P0 qfgh + 2r/R. Initially, a spherical air bubble is released from rest at a height h = 5R above the bottom of computational domain, where R is the radius of the undeformed bubble. Then the bubble rises due to buoyancy and begins to deform. Since the rising distance is short, the change of bubble volume during the rising process is neglected. Thus, the bubble will obtain its terminal velocity U0 and reach its final deformed shape after the relaxation time. In this study, we are interested in the interaction between the rigid surface and the fully deformed bubble. Moreover, if the separation between the two sides (top and bottom) of the bubble reaches zero, the bubble will become toroidal in our simulation. For convenience, in the following discussions, we refer to this topological change as ‘bubble breakup’ or ‘bubble rupture’.
2. Governing equations and numerical methods In bubble-rigid surface interaction problems, it is necessary to pay attention to both fluid motion and bubble deformation. For the fluid, the governing equations are the continuity equation
$ u ¼ 0;
ð1Þ
and the momentum equation
qf
Du ¼ $p þ $ s þ qf f ; Dt
ð2Þ
where qf, u, p, s and f are fluid density, velocity, pressure, viscous stress tensor, and gravity. For the gas bubble, we do not consider gas diffusion which means that no bubble gains mass from surrounding fluids. The gas in the bubble is at thermodynamic and mechanical equilibrium and follows the ideal gas law:
e Pb V b ¼ mb RT;
ð3Þ
e the specific where Pb is the gas pressure, Vb is the bubble volume, R gas constant and T the absolute temperature. In this problem, we only consider the isothermal case, so T is a constant. On the bubble surface, we assume that there is no surfactant so the shear stress is zero. The normal stress is obtained from the Young–Laplace equation:
n ðpI þ sÞ ¼ ðPb þ K rÞn;
ð4Þ
where r is the surface tension, n is the unit normal to bubble surface pointing into the bubble, and K is the surface curvature. For
Fig. 1. The sketch of interaction between fully deformed bubble and rigid surface, y is axis of symmetry.
62
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
(a)
(b)
Fig. 2. Triangular mesh around a near-wall bubble, with refined elements between the bubble and rigid wall.
250
ues of Mo, such as pure water, researchers usually use Reynolds number (Re) and Weber number (We) to study the bubble energy (Tsao and Koch, 1997), bubble bouncing (Zenit and Legendre, 2009) because inertial effects are very important. For the high values of Mo, which correspond to highly viscous fluids, the Mo and Bo are used to examine the bubble shape (Grace, 1973; Mukundakrishnan et al., 2007) in that the viscous effects become more important. In this study, the focus is on the bubble shape evolution during the interaction after the bubble is fully developed, the process can be determined by two dimensionless groups: Morton number (Mo) and Bond number (Bo). The Morton number is defined as
U t (mm s-1)
200
150
100 Experiments Previous simulation Present work
50
0
2
4
6
8
10
d (mm) Fig. 3. Comparison of our simulations for the terminal velocity vs. bubble diameter in a Newtonian liquid for Mo = 3.769 104 with experiments reported by Maxworthy et al. (1996) and numerical prediction by Tsamopoulos et al. (2008).
3.2. Parameter determination To investigate the bubble-rigid wall interaction problem, two dimensionless parameters need to be determined. For the low val-
Mo ¼
g l4f Dq
q2f r3
ð5Þ
;
where g is the acceleration of gravity, lf and qf are the viscosity and density of the surrounding fluid, Dq is the difference in density of the air and fluid, and r is the surface tension between the air bubble and the fluid. In our problem, air density qair in the bubble is neglected, so Mo can be rewritten as
Mo ¼
g l4f
qf r3
ð6Þ
:
The Bond number is defined as
Bo ¼
qf gR2 : r
(a)
ð7Þ
(b)
Fig. 4. (a) Comparison between coarse grids (left) and fine grids (right). (b) Bubble shapes: solid line for the coarse grid and dashed line for the fine grid.
63
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0
-0.05
(a)
(b) v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
(c)
(d) v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
(e)
(f) v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
(g)
(h)
Fig. 5. Vertical velocity component contours of interaction between a rising bubble and a rigid wall at Mo = 6.5 102, Bo = 6: (a) t = 0; (b) t = 5; (c) t = 10; (d) t = 15; (e) t = 20; (f) t = 30; (g) t = 50 and (h) t = 100.
It is important to note that the parameters Re, We, Mo, Bo are not independent of each other. Grace et al. (1976) gave a famous diagram to show the relationship among Mo, Bo and Re. From this diagram, each of the three parameters can be uniquely determined from the other two. We number can be also obtained by
We ¼ Re2
rffiffiffiffiffiffiffi Mo : Bo
ð8Þ
In the bubble-rigid wall interaction problem, we focus on the bubble shape evolution as it approaches the rigid wall. The simulations assume highly viscous liquids (Mo > 6.5 106). Therefore, the bubble rebound is unlikely.
64
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
20
RigidWall
y
yt and yb
19.5
Top surfaceyt Bottom surfaceyb
yt
yb
19
18.5
18 10
x
15
20
25
30
Time Fig. 6. A sketch of the top (yt) and bottom (yb) surfaces of the bubble. Fig. 8. The positions of the two sides of the bubble along y-axis for Mo = 6.5 102 and Bo = 6.
Surface distance (y t-yb)
2
Oh. As Oh increases, the viscous dissipation starts to dominate and the bubble may not bounce back. In our study, if we consider the value of CAm as 0.5 and deq as 2R, for qair qf, the Ca/St⁄ can be written as
1.5
rffiffiffiffiffiffiffi Ca Mo 9l2 ¼ ¼ 9 : St Bo Rqf r 1
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi The minimum value of 9 Mo=Bo, in our simulations, is larger than the rebound condition Ca/St⁄ < 102. Thus for the parameter range in our simulations, the bubble will not rebound. Although the bubble does not rebound as a whole, we do observe that, in some cases, some part of the bubble surface gains a downward velocity after interacting with the wall.
0.5
0
0
10
20
30
3.3. Mesh refinement
40
Time Fig. 7. The distance between two sides of the bubble for Mo = 6.5 10
2
and Bo = 6.
Zenit and Legendre (2009) carefully studied the coefficient of restitution for air bubbles colliding against solid walls. They indicated that the coefficient of restitution can be described by a function
rffiffiffiffiffiffiffi! Ca ¼ exp 30 ; St
ð9Þ
where Ca = U0l/r is the capillary number. The modified Stokes number St⁄ is defined as
St ¼
ð11Þ
ðqair þ C AM qf Þdeq U 0 ; 9l
ð10Þ
where CAM(v) is the added mass coefficient which depends on the bubble aspect ratio v (for a sphere, v = 1 and CAm = 0.5), deq is the pffiffiffiffiffiffiffiffiffiffiffiffiffiffi bubble equivalent diameter. It should be noted that Ca=St (the Ohnesorge number (Oh)) relates the viscous forces to inertial and surface tension forces. In their study, most rebounds are observed for Ca/St⁄ < 102. It means that the rebound only occurs at a small
There are two important processes in the problem of interaction between a fully developed bubble and a rigid surface. The first process is the rising of a bubble in highly viscous liquids. Since the interest is in the bubble shape evolution, it is important to resolve the bubble shape and obtain accurate terminal velocity. The detailed discussion will be presented in the next section. The second process is the collision between the deformed bubble and a rigid surface. It requires the computational mesh to be fine enough between the bubble surface and the rigid wall for detailed analysis of film drainage. In Fig. 2a, a case of fully deformed bubble colliding with a rigid surface is shown for Mo = 0.1 and Bo = 8.5. As can be seen more clearly in Fig. 2b, the mesh in the thin film between the two surfaces is adaptively refined. Even the thinnest part of the film contains at leat three layers of interior nodes to simulate the film drainage process.
4. Results and discussion First we will demonstrate that our numerical method produces accurate bubble shape, flow field structure and bubble terminal velocities at moderate Mo through comparison with known solutions. Grid independence is also shown by using coarse and fine grids. Then, we fix a certain Mo to discuss the bubble shape evolution by changing Bo. Three modes for bubble-rigid wall interaction
65
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
clearly that a dimple appears on the near-wall side of the bubble as shown in Fig. 5e. The formation of the dimple makes the smallest film thickness appears at the bubble rim rather than at its center. When the film between the rim of the bubble and rigid surface is thin enough to rupture, the bubble will attach to the rigid surface. However, due to limitations of the present numerical method, we cannot capture the film rupture. The simulation stops when the thickness of the film is smaller than 1/100 of the bubble radius. To investigate the bubble deformation and the interaction between the bubble and the rigid wall, we denote the instantaneous position of the intersections of bubble surface with the axis of symmetry as yb and yt for the bottom and top surface, respectively (see Fig. 6). Fig. 7 shows the temporal evolution of the vertical distance between the top and bottom sides of the bubble (yt–yb) along the axis of symmetry. The figure shows several phases of deformation of the bubble. Initially, because the bubble is released from rest and then begins to deform, the corresponding vertical distance deceases sharply from twice the bubble radius until the bubble is fully deformed after t 10. Then, since the bubble has been fully developed and is still far away from the rigid wall, the vertical distance changes very little. Afterwards, with the fully deformed bubble approaching the rigid wall, the vertical distance deceases sharply again and then rises a little. Finally, the vertical distance deceases very slowly until it reaches a constant value, which is caused by the film drainage process between the bubble and rigid wall. At this point, it is important to note that there exists a little increase in the vertical distance as the bubble comes near the rigid wall. To clarify the cause for this increase, we show the positions of both sides of the bubble versus time in Fig. 8. When the bubble interacts with the rigid wall, it can be clearly seen that the position of top surface almost becomes a constant whereas the bottom surface recedes away from the rigid wall. This leads to the little increase of separation between the top and bottom surfaces. The cause for the overshoot of the height of the bubble bottom surface around t = 20, as shown in Fig. 8, is the storage and release of the bubble surface energy. When the bubble interacts with the rigid wall, the surrounding liquid below the bubble still has an upward velocity. Then, a part of the kinetic energy of this liquid is transformed to the surface energy of the bubble. As the upward velocity vanishes, the bottom surface tends to restore to its original shape. This leads to the bottom surface receding away from the rigid wall. To illustrate this phenomenon, Fig. 9 presents the flow field with velocity contours for Mo = 6.5 102 and Bo = 6 at
have been distinguished. Finally, the boundaries between the three modes are presented in a (Mo, Bo)-plane. 4.1. Numerical method validation First, we compare our results with the experimental data by Maxworthy et al. (1996) who measured bubble terminal velocities as a function of bubble diameter in the mixture of distilled water with glycerin (Mo = 3.769 104) and with the numerical predictions by Tsamopoulos et al. (2008) who conducted their numerical simulation under the same condition by using a mixed finite element method in Fig. 3. We observe very good agreements between our simulations and the two previous studies. The code used in this study has been extensively validated by Yue et al. (2007). The bubble shape, terminal velocity, and the streamlines in the wake of a rising bubble match very well with the experiments by Hnat and Buckmaster (1976). The grid and time convergence studies were also presented in Yue et al. (2007). Beside those, we compare the coarse grid which is used in the following simulations with the fine one in Fig. 4. As shown in Fig. 4b, the bubble shapes predicted by the two grids are indistinguishable. 4.2. Effects of Bond number In our simulations, we release the bubble far below a solid surface at a distance of 15R. The bubble rises toward the solid wall by buoyancy only, and reaches a terminal velocity before interacting with the solid wall. Through numerical experiments, it is found that if we fix Mo and change Bo, there will be three modes of interaction. To explain the three modes, we set Mo = 6.5 102 and change Bo from 6 to 14. All parameters are non-dimensionalized by bubble radius, gravitational acceleration and fluid viscosity. 4.2.1. Mode one For small Bo, which means the ratio of buoyancy force to surface tension is small, the deformability of the bubble is weak. When the bubble rises towards a solid surface, a liquid film forms between the bubble and the solid wall. The high pressure inside the film causes a slight dimple on the bubble’s upper surface. Fig. 5 shows a time sequence of a gas bubble approaching a solid wall for Mo = 6.5 102 and Bo = 6. The bubble is released at t = 0 with zero velocity (see Fig. 5a). It then rises due to buoyancy force and reaches its terminal velocity (Vt = 0.85) and deformation (see Fig. 5c). As the bubble approaches the rigid surface, it can be seen
-2
-1
0
1
2
v
v
0.55 0.51 0.47 0.43 0.39 0.35 0.31 0.27 0.23 0.19 0.15 0.11 0.07 0.03 -0.01 -0.05
0.26 0.24 0.22 0.2 0.18 0.16 0.14 0.12 0.1 0.08 0.06 0.04 0.02 0 -0.02 -0.04
-2
-1
0
x
x
(a)
(b)
1
Fig. 9. Vertical velocity component contours for Mo = 6.5 102 and Bo = 6 at: (a) t = 18 and (b) t = 21.
66
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(a)
(b) v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(c)
(d) v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(e)
(f) v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(g)
(h)
Fig. 10. Vertical velocity component contours of interaction between a rising bubble and a rigid wall at Mo = 102, Bo = 10: (a) t = 0; (b) t = 2.5; (c) t = 5; (d) t = 7.5; (e) t = 10; (f) t = 15; (g) t = 20 and (h) t = 30.
t = 18 and t = 21. It is obvious that the liquid below the bubble has an upward velocity at t = 18 and a downward velocity at t = 21. 4.2.2. Mode two As Bo increases, the inertia of the liquid becomes more important and the bubble is more easily deformed. When Bo reaches a critical value, the top and bottom sides of the bubble collide with
each other and the bubble may form a torus. Fig. 10 shows the deformation phases of a rising bubble and its interaction with a rigid surface for Mo = 6.5 102 and Bo = 10. In the rising stage, the fully developed bubble shape is a little flatter than that in mode one (Fig. 10c) and the terminal velocity increases to 0.88. In Fig. 10h, the final stage of interaction is depicted where a bubble with a large dimple is formed and the bubble will break up instantaneously.
67
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
2
Surface distance yt-yb
v 0.24 0.22 0.2 0.18 0.16 0.14 0.12 0.1 0.08 0.06 0.04 0.02 0 -0.02
1.5
1
0.5 -2
0
2
x 0
0
5
10
15
20
25
30
35
Time Fig. 11. The distance between two sides of the bubble for Mo = 6.5 102 and Bo = 10.
20
yt and yb
19.5
19 Top surface yt Bottom surface yb 18.5
18 10
15
20
25
30
Time Fig. 12. The positions of the two sides of the bubble along y-axis for Mo = 6.5 102 and Bo = 10.
To describe the process leading to the rupture of the bubble, a temporal evolution of the vertical distance between the top and bottom surfaces of the bubble along y-axis is displayed in Fig. 11. Before the bubble approaches the rigid wall, the vertical distance has a similar behavior as that in mode one. As the bubble interacts with the rigid wall, the vertical distance first increases (around t = 20), which is still the same as that in mode one. However, the distance decreases to zero right after this increase. This is caused by the deformation of the top surface of the bubble. As shown in Fig. 12, the change in the position of the bottom surface is similar to that in mode one. However, the top surface bounces back quickly after the bubble collides with the rigid wall. The reason is that the bubble is easier to deform and therefore a larger dimple is formed on the upper surface. If the upper surface reaches the lower surface, the bubble ruptures. The above analysis is consistent with the flow field shown in Fig. 13, in which the surrounding liquid near the upper surface of the bubble has a downward velocity while the lower surface is almost stationary.
Fig. 13. Vertical velocity component contours for Mo = 6.5 102 and Bo = 10 at t = 26.
4.2.3. Mode three As Bo is further increased, we expect a situation in which the bubble breaks up before it collides with the rigid wall due to increased bubble deformability. We denote this situation as mode three. As an example, a typical process of bubble-rigid surface interaction with Mo = 6.5 102 and Bo = 14 is presented in Fig. 14. For the high Bo in this mode, the ratio of the buoyancy force to surface tension is much higher than in the previous two modes, which leads to a significant bubble deformation. The fully developed bubble is flatter and the edge of the bubble becomes sharper (Fig. 14c) and the terminal velocity increases to 0.9. It should be noted that, the bubble will also break up in this mode (Fig. 14h), however, it is different from the rupture in mode two. To explain the reason, Fig. 15 shows a temporal evolution of the vertical distance between the top and bottom sides of the bubble along the axis of symmetry. Similar to the other modes, the bubble gets fully developed and reaches its terminal velocity very soon (t 7). However, different from the modes one and two, the vertical distance sharply decreases to zero without a temporary increase as the bubble interacts with the rigid wall. In Fig. 16, it can be clearly seen that when the bubble comes close to the rigid wall, the rising of the top surface slows down but the bottom surface continues to rise. Hence, after a while, the bottom surface will catch up with the top one and then the bubble breaks up. Fig. 17 shows the vertical velocity contours of the flow field and the surface of bubble. It can be observed that the surrounding liquid, above the top surface of the bubble is almost immobile but the liquid under the bottom surface still has an upward velocity. This will cause the two sides of the bubble to touch and the bubble ruptures. It is interesting to note that, as shown in Fig. 13, a dimple forms at the top surface of the bubble. At the edge of the film, the squeezing effect maximizes due to the smallest film thickness. As a consequence, a portion of the fluid there is squeezed inward into the film, which causes flow reversal inside the film. In Fig. 17, as Bo gets higher, the bubble has a higher terminal velocity and the flow inertia in the wake plays a more important role. Consequently, the bubble deforms more and ruptures before the dimple on the upper surface forms. Thus, there is no chance to observe the reversal flow in the film. It should also be noted that the film thickness in Fig. 17 is much thicker than the minimum thickness in Fig. 13 due to the earlier rupture of bubble. 4.3. Effects of Morton number In the last section, a moderate Morton number (Mo) case has been investigated and three modes have been identified by numer-
68
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(b)
(a)
v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(c)
(d) v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(f)
(e)
v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(g)
(h)
Fig. 14. Vertical velocity component contours of interaction between a rising bubble and a rigid wall at Mo = 6.5 102, Bo = 14: (a) t = 0; (b) t = 2; (c) t = 4; (d) t = 6; (e) t = 8; (f) t = 10; (g) t = 14 and (h) t = 18.
ical simulations. Also, we find that the three modes can be distinguished from each other by a transition Bound number (Bo). In practice, due to flow instability and the roughness of the rigid wall, a precise transition Bo may be difficult to find. However, for convenience, it is reasonable to assume the existence of a transition Bo between two neighboring modes of interaction. Such a transition Bo will be helpful in understanding the physical mechanisms and
achieving a criterion for mode transition at different Mo. Here we should point out that the transition between different modes is regular and there is no sudden change. For the cases listed in Table 1, we have carefully examined the variation of transition Bo for different Mo. The transition Bo is established by incrementing Bo by 0.1 each step to investigate the problem and reduce the number of the time consuming simu-
69
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
2
Surface distance yt-yb
v 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1.5
1
0.5
-4 0
-3
-2
-1
0
1
2
3
4
x 0
5
10
15
20
Time
Fig. 17. Vertical velocity component contours for Mo = 6.5 102 and Bo = 14 at t = 17.5.
Fig. 15. The distance between two sides of the bubble for Mo = 6.5 102 and Bo = 14. Table 1 Transition Bond number for different Morton numbers.
20 Top surfacey t Bottom surfaceyb
yt and yb
19.5
19
Mo
Transition Bo between modes one and two
Transition Bo between modes two and three
1.0 101 6.5 102 2.0 102 8.4 103 4.0 103 9.8 104 1.0 104 3.3 105 6.5 106
8.7 8.3 7.7 7.4 6.9 5.8 4.1 3.4 2.7
15.3 13.9 11.0 9.3 8.1 6.4 4.6 3.7 2.8
16
18.5
Mode 2 12
14
16
18
Time Fig. 16. The positions of the two sides of the bubble along y-axis for Mo = 6.5 102 and Bo = 14.
lations. In our results, Mo varies from 1.0 101 to 6.5 106 and the corresponding transition Bo between mode one and mode two decreases from 8.7 to 2.7, and similarly the one between mode two and mode three drops from 15.3 to 2.8. Furthermore, the difference between two transition Bond numbers at the same Mo is also lessened from 6.6 to 0.1 with the decreasing Mo. That means that the region of mode two shrinks as the Mo becomes small. In other words, the existance of the three modes is more distinct at a high Mo. To clearly present the results, we describe the boundaries between different modes on the Mo-Bo chart. In Fig. 18, we represent Mo on x-direction with log scale and display Bo in y-direction. The two curves shown in this figure give the boundaries between different modes of interaction. The top one distinguishes mode two from mode three and the bottom one separates mode one from mode two. In this figure, for large Mo (>102), it is easy to distinguish the three modes from each other. When Mo is large, the viscosity of the fluid plays an important role. The surface tension,
Mode 3
12
20
Bo number
18 10
8
4
0 -6 10
Mode 1
-5
10
-4
10
-3
10
-2
10
-1
10
0
10
Mo number Fig. 18. Three modes diagram in the (Mo, Bo)-plane.
which determines Bo, needs to vary in a wide range to obtain different modes of bubble deformation. However, for small Mo number (<104), the surface tension dominates the bubble deformation. The gas bubble either breaks up at a high Bo or remains unbroken at a low Bo, which causes the region of mode two to decrease. Especially for the Mo < 105, the two curves merge into one and the region of mode two disappears. In addition, for the case of
70
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
Mo < 105 and in the range of Bo given in Fig. 18, the shape of the bubble remains spherical, see Grace et al. (1976) and Tsamopoulos et al. (2008). There may exist more modes, for example, the bubble rebound and vibration as observed experimentally (Tsao and Koch, 1997; Zenit and Legendre, 2009; Hendrix et al., 2012). However, the bubble deformation in these studies is small and the corresponding Mo is beyond the range considered in the present study. It should be noted that Mo depends only on the material properties of the viscous liquid. However, Bo depends on both the properties of the viscous liquid and the bubble size. Thus, in practice, we can predict the critical bubble size in a given viscous liquid using the above diagram, according to Mo and Bo. For example, if a solid container is filled with cyclohexanol, the properties of which are q = 962 kg/m3, l = 0.0575 Pa s and r = 0.0329 N/m (Haynes (2011-2012)), we obtain Morton number Mo = 3 103. From the above diagram, we can determine the transition Bo between mode one and mode two to be Bo = 6.6. Thus, we can predict that the radius of the bubble which can attach to the upper wall of the solid container without break up should be smaller than 5 mm at its undeformed state. 5. Conclusion In this article, direct numerical simulations have been performed by using an arbitrary-Lagrangian–Eulerian method in order to study the interaction between a deformable bubble and a rigid wall and to clarify the bubble deformation as it approaches the rigid wall. It is found that three modes of bubble-rigid wall interaction exist as Bo changes at a moderate Mo (6.5 106 < Mo < 0.1). The first mode prevails at small Bo where the bubble deformation is small. For this mode, the bubble is hard to break up and will bounce back and eventually attach to the rigid wall. In the second mode, the bubble may break up after it collides with the rigid wall, which is determined by the film drainage. In the third mode, which prevails at high Bo, the bubble breaks up due to the bottom surface catches up the top surface during the interaction with a rigid wall. We also found that the three modes are easily distinguishable when Mo is large. When Mo gets small, the region of the second mode diminishes and finally disappears. We have to point out that the current work assumes axisymmetric flow. For the bubble rising problem, according to the experiments by Bhaga and Weber (1981), the flow should be axisymmetric in the Reynolds number range (Re < 50) in this work. However, for the more complicated bubble-rigid surface interaction problem, the flow axisymmetry and actual existence of the three modes should be established by future 3D simulations or experiments. Acknowledgments The authors thank the anonymous reviewers for helpful comments. P. Yue acknowledges the support by NSF-DMS 0907788. References Bhaga, D., Weber, M.E., 1981. Bubbles in viscous liquids: shapes, wakes and velocities. J. Fluid Mech. 105, 61–85. Duineveld, P.C., 1995. The rise velocity and shape of bubbles in pure water at high reynolds number. J. Fluid Mech. 292, 325–332. Frank, X., Funfschilling, D., Midoux, N., Li, H.Z., 2006. Bubbles in a viscous liquid: lattice boltzmann simulation and experimental validation. J. Fluid Mech. 546, 113–122.
Fujasova-Zednikova, M., Vobecka, L., Vejrazka, J., 2010. Effect of solid material and surfactant presence on interactions of bubbles with horizontal solid surface. Can. J. Chem. Eng. 88, 473–481. Gaudlitz, D., Adams, N.A., 2009. Numerical investigation of rising bubble wake and shape variations. Phys. Fluids 21, 122102. Grace, J.R., 1973. Shapes and velocities of bubbles rising in infinite liquids. Trans. Instn. Chem. Eng. 51, 116–120. Grace, J.R., Wairegi, T., Nguyen, T.H., 1976. Shapes and velocities of single drops and bubbles moving freely through immiscible liquids. Trans. Instn. Chem. Eng. 54, 167–173. Gueyffier, D., Li, J., Nadim, A., Scardovelli, R., Zaleski, S., 1999. Volume-of-fluid interface tracking with smoothed surface stress methods for three-dimensional flows. J. Comput. Phys. 152, 423–456. Haynes, W.M., 2011-2012. CRC Handbook of Chemistry and Physics, 92nd ed. CRC Press. Hendrix, M.H.W., Manica, R., Klaseboer, E., Chan, D.Y.C., Ohl, C., 2012. Spatiotemporal evolution of thin liquid films during impact of water bubbles on glass on a micrometer to nanometer scale. Phys. Rev. Lett. 108, 247803. Hnat, J.G., Buckmaster, J.D., 1976. Spherical cap bubbles and skirt formation. Phys. Fluids 19, 182–194. Hu, H.H., 1996. Direct simulation of flows of solid–liquid mixtures. Int. J. Multiphase Flow 22, 335–352. Hu, H.H., Patankar, N.A., Zhu, M.Y., 2001. Direct numerical simulations of fluid–solid systems using the arbitrary lagrangian-eulerian technique. J. Comput. Phys. 169, 427–462. Klaseboer, E., Chevaillier, J., Mate, A., Masbernat, O., Gourdon, C., 2001. Model and experiments of a drop impinging on an immersed wall. Phys. Fluids 13, 45–57. Kushch, V.I., Sangani, A.S., Spelt, P.D.M., Koch, D.L., 2002. Finite-weber-number motion of bubbles through a nearly inviscid liquid. J. Fluid Mech. 460, 241–280. Legendre, D., Daniel, C., Guiraud, P., 2005. Experimental study of a drop bouncing on a wall in a liquid. Phys. Fluids 17, 097105. Maxworthy, T., Gnann, C., Kurten, M., Durst, F., 1996. Experiments on the rise of air bubbles in clean viscous liquids. J. Fluid Mech. 321, 421–441. Minev, P.D., Chen, T., Nandakumar, K., 2003. A finite element technique for multifluid imcompressible flow using Eulerian grids. J. Comput. Phys. 187, 255– 273. Moore, D.W., 1965. The velocity of rise of distorted gas bubbles in a liquid of small viscosity. J. Fluid Mech. 23, 749–766. Mukundakrishnan, K., Quan, S., Eckmann, D.M., Ayyaswamy, P.S., 2007. Numerical study of wall effects on buoyant gas-bubble rise in a liquid-filled finite cylinder. Phys. Rev. E 76, 036308. Ohta, M., Imura, T., Yoshida, Y., Sussman, M., 2010. A computational study of the effect of initial bubble conditions on the motion of a gas bubble rising in viscous liquids. Int. J. Multiphase Flow 31, 223–237. Power, H., 1997. The interaction of a deformable bubble with a rigid wall at small reynolds number: a general approach via integral equations. Eng. Anal. Bound. Elem. 19, 291–297. Rabha, S.S., Buwa, V.V., 2010. Volume-of-fluid (vof) simulations of rise of single/ multiple bubbles in sheared liquids. Chem. Eng. Sci. 65, 527–537. Shopov, P.J., Minev, P.D., Bazhlekov, I.B., Zapryanov, Z.D., 1990. Interaction of a deformable bubble with a rigid wall at moderate reynolds numbers. J. Fluid Mech. 219, 241–271. Tryggvason, G., Bunner, B., Esmaeeli, A., Juric, D., Al-Rawahi, N., Tauber, W., Han, J., Nas, S., Jan, Y.J., 2001. A front-tracking method for the computations of multiphase flow. J. Comput. Phys. 169, 708–759. Tryggvason, G., Esmaeeli, A., Lu, J., Biswas, S., 2006. Direct numerical simulations of gas/liquid multiphase flows. Fluid Dyn. Res. 38, 660–681. Tsamopoulos, J., Dimakopoulos, Y., Chatzidai, N., Karapetsas, G., Pavlidis, M., 2008. Steady bubble rise and deformation in newtonian and viscoplastic fluids and conditions for bubble entrapment. J. Fluid Mech. 601, 123–164. Tsao, H., Koch, D.L., 1997. Observations of high reynolds number bubbles interacting with a rigid wall. Phys. Fluids 9, 44–56. Wilkes, E.D., Phillips, S.D., Basaran, O.A., 1999. Computational and experimental analysis of dynamics of drop formation. Phys. Fluids 11, 3577–3598. Yiantsios, S.G., Davis, R.H., 1990. On the buoyancy-driven motion of a drop towards a rigid surface or a deformable interface. J. Fluid Mech. 217, 547–573. Yu, Z., Fan, L., 2008. Direct simulations of the buoyant rise of bubbles in infinite liquid using level set method. Can. J. Chem. Eng. 86, 267–275. Yue, P., Feng, J.J., Bertelo, C.A., Hu, H.H., 2007. An arbitrary Lagrangian–Eulerian method for simulating bubble growth in polymer foaming. J. Comput. Phys. 226, 2229–2249. Yue, P., Zhou, C., Feng, J.J., Ollivier-Gooch, C.F., Hu, H.H., 2006. Phase-field simulations of interfacial dynamics in viscoelastic fluids using finite elements with adaptive meshing. J. Comput. Phys. 219, 47–67. Zawala, J., Krasowska, M., Dabros, T., Malysa, K., 2007. Influence of bubble kinetic energy on its bouncing during collisions with various interfaces. Can. J. Chem. Eng. 85, 669–678. Zenit, R., Legendre, D., 2009. The coefficient of restitution for air bubbles colliding against solid walls in viscous liquids. Phys. Fluids 21, 083306.
Contents lists available at SciVerse ScienceDirect
International Journal of Multiphase Flow j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / i j m u l fl o w
Axisymmetric simulation of the interaction of a rising bubble with a rigid surface in viscous flow Tong Qin a,⇑, Saad Ragab a, Pengtao Yue b a b
Department of Engineering Science and Mechanics, Virginia Tech, Blacksburg, VA 24061, USA Department of Mathematics, Virginia Tech, Blacksburg, VA 24061, USA
a r t i c l e
i n f o
Article history: Received 22 June 2012 Received in revised form 4 January 2013 Accepted 5 January 2013 Available online 12 January 2013 Keywords: Bubble Rigid wall Film drainage ALE method Moving mesh
a b s t r a c t The interaction between a rising deformable gas bubble and a solid wall in viscous liquids is investigated by direct numerical simulation via an arbitrary-Lagrangian–Eulerian (ALE) approach. The flow field is assumed to be axisymmetric. The bubble is driven by gravity only and the motion of the gas inside the bubble is neglected. Deformation of the bubble is tracked by a moving triangular mesh and the liquid motion is obtained by solving the Navier–Stokes equations in a finite element framework. To understand the mechanisms of bubble deformation as it interacts with the wall, the interaction process is studied as a function of two dimensionless parameters, namely, the Morton number (Mo) and Bond number (Bo). We study the range of Bo and Mo from (2, 6.5 106) to (16, 0.1). The film drainage process is also considered in this study. It is shown that the deformation of a bubble interacting with a solid wall can be classified into three modes depending on the values of Mo and Bo. Ó 2013 Elsevier Ltd. All rights reserved.
1. Introduction The interaction of a rising bubble with a rigid surface plays an important role in many technological processes (e.g. mineral flotation, fluidized beds and particle sedimentation). The mechanism of bubble-wall collision has not been fully understood. The significance of studying this problem lies not only the bubble deformation when the buoyancy-driven bubble rises, but also in the film drainage process that takes place when the rising bubble impacts a solid wall. For the motion of a bubble in a liquid, Moore (1965) calculated the terminal rising velocity of distorted gas bubbles in a liquid of small viscosity. Grace (1973) and Grace et al. (1976) developed a well-known graphical correlation to explain the shapes and velocities of single gas bubbles rising in infinite liquids. Hnat and Buckmaster (1976) also experimentally studied the shapes and terminal velocities of bubbles rising in liquids and their data are used as benchmarks for comparison with numerical studies (Gueyffier et al., 1999; Ohta et al., 2010). In addition, Bhaga and Weber (1981) observed the flow field around a rising bubble by using the hydrogen bubble tracer technique. They reported that the wake of bubble was closed and laminar for Mo > 4 103 and Reynolds number (Re) < 110 but was open and unsteady for Re > 110. Their work has been compared with numerical simulations by Gaudlitz and Adams (2009) who investigated the wake and shape of a rising
⇑ Corresponding author. E-mail address: [email protected] (T. Qin). 0301-9322/$ - see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijmultiphaseflow.2013.01.001
bubble on a zigzag path. Duineveld (1995) experimentally determined the velocity and shape of rising bubbles in pure water with an extremely low Mo (1011). They found the bubbles can approximately be considered as oblate spheroids with Re 102 and Weber number (We) 10. Many numerical simulations have been carried out for rising bubble, including groups of bubbles. Kushch et al. (2002) computed the motion of bubbles with a finite Weber number through a nearly inviscid liquid. Frank et al. (2006) and Yu and Fan (2008) both reported their research on numerical simulations of bubbles in a viscous liquid using lattice Boltzmann method and level set method, respectively. A review on direction numerical simulations of bubbly flows can be found in Tryggvason et al. (2006). In practice, researchers have paid much attention to the rising gas bubble interacting with surfaces, especially in water. In that case and depending on the bubble size, viscous effects are weak, and much of the kinetic energy can be transformed into surface energy of the bubble, a later release of which may result in a bubble rebound. Using a high-speed video camera, Tsao and Koch (1997) observed a bubble bouncing several times from a horizontal wall before viscosity dissipated the kinetic energy. Klaseboer et al. (2001) experimentally studied the rebound of a drop impinging on a rigid wall and presented a bouncing trajectory mode which included the film drainage effect. Then, Legendre et al. (2005) and Zenit and Legendre (2009) carried out experimental studies on bubble and drops bouncing on a solid wall in viscous liquids. Recently, Zawala et al. (2007) investigated the collision and bouncing of bubbles with various interfaces (water–gas and water–solid interfaces). The effect of the solid mate-
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
rial and surfactant on the bubble-wall interaction has been reported by Fujasova-Zednikova et al. (2010). The dynamics of a bubble or a drop approaching a surface have been considered in a number of theoretical studies. Yiantsios and Davis (1990) analyzed the shape of a drop driven by buoyancy towards a flat rigid surface and a deformable interface. Power (1997) studied the interaction of a deformable bubble with a rigid wall at small Re by solving the Stokes equation. Related numerical work has been reported by Shopov et al. (1990). They investigated the hydrodynamic behavior of a bubble near a rigid surface and the wall of a spherical container, and successfully simulated the dimple formation and the concavity in the rear part of the bubble. However, the lubrication layer formed between the wall and the gas–liquid interface is not thin enough to account for the whole process of film drainage. Mukundakrishnan et al. (2007) numerically studied the wall effects on a rising gas bubble in a liquid-filled finite cylinder. They described the bubble rising in liquids by Eötvös number (Eo) and Mo and investigated the bubble shape regimes by changing the cylinder radius. Many other numerical methods have been used for solving this multiphase problem, such as the volume-of-fluid (VOF) (Rabha and Buwa, 2010), front tracking (Tryggvason et al., 2001), level set (Minev et al., 2003), moving mesh (Wilkes et al., 1999), and phase-field (Yue et al., 2006). Nevertheless, a deformed bubble colliding with a rigid surface, especially in high-Mo liquid, has not been fully studied. To provide more details about a deformable bubble interacting with a rigid surface, the present work uses an arbitrary Lagrangian–Eulerian (ALE) formulation to simulate the bubble deformation and film drainage process. The ALE method can solve the movingboundary problems by using a moving mesh with grid points on the bubble surface (Hu et al., 2001). A parametric investigation which has been conducted revealed the existence of three modes of bubble-wall interactions. The prevailing mode of interaction depends on two parameters, namely the Morton number and the Bond number.
61
the outer boundary of the computational domain, no-slip, slip, or prescribed stress conditions can be applied. The governing equations are solved in axisymmetric geometry on a triangular mesh using a finite element method. The mesh moves and deforms according to an arbitrary Lagrangian–Eulerian (ALE) scheme as described by Hu (1996) and Hu et al. (2001). Further details on the numerical methods can be found in Yue et al. (2007). 3. Problem description and mesh refinement 3.1. Problem description In the present bubble-rigid surface interaction problem, the bubble is rising due to buoyancy only and the rigid surface is horizontal, so the overall motion of bubble happens only in the normal direction to the surface. Thus, this problem can be considered as an axisymmetric problem with the y axis being the axis of symmetry as shown in Fig. 1. The top surface of the computational domain is a rigid surface with no-slip boundary condition. For the side wall of the computational domain, we apply a slip boundary condition. At the bottom surface, we apply stress boundary condition, and specify the ambient pressure to be P0. The pressure in the bubble can be initially obtained by Pb = P0 qfgh + 2r/R. Initially, a spherical air bubble is released from rest at a height h = 5R above the bottom of computational domain, where R is the radius of the undeformed bubble. Then the bubble rises due to buoyancy and begins to deform. Since the rising distance is short, the change of bubble volume during the rising process is neglected. Thus, the bubble will obtain its terminal velocity U0 and reach its final deformed shape after the relaxation time. In this study, we are interested in the interaction between the rigid surface and the fully deformed bubble. Moreover, if the separation between the two sides (top and bottom) of the bubble reaches zero, the bubble will become toroidal in our simulation. For convenience, in the following discussions, we refer to this topological change as ‘bubble breakup’ or ‘bubble rupture’.
2. Governing equations and numerical methods In bubble-rigid surface interaction problems, it is necessary to pay attention to both fluid motion and bubble deformation. For the fluid, the governing equations are the continuity equation
$ u ¼ 0;
ð1Þ
and the momentum equation
qf
Du ¼ $p þ $ s þ qf f ; Dt
ð2Þ
where qf, u, p, s and f are fluid density, velocity, pressure, viscous stress tensor, and gravity. For the gas bubble, we do not consider gas diffusion which means that no bubble gains mass from surrounding fluids. The gas in the bubble is at thermodynamic and mechanical equilibrium and follows the ideal gas law:
e Pb V b ¼ mb RT;
ð3Þ
e the specific where Pb is the gas pressure, Vb is the bubble volume, R gas constant and T the absolute temperature. In this problem, we only consider the isothermal case, so T is a constant. On the bubble surface, we assume that there is no surfactant so the shear stress is zero. The normal stress is obtained from the Young–Laplace equation:
n ðpI þ sÞ ¼ ðPb þ K rÞn;
ð4Þ
where r is the surface tension, n is the unit normal to bubble surface pointing into the bubble, and K is the surface curvature. For
Fig. 1. The sketch of interaction between fully deformed bubble and rigid surface, y is axis of symmetry.
62
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
(a)
(b)
Fig. 2. Triangular mesh around a near-wall bubble, with refined elements between the bubble and rigid wall.
250
ues of Mo, such as pure water, researchers usually use Reynolds number (Re) and Weber number (We) to study the bubble energy (Tsao and Koch, 1997), bubble bouncing (Zenit and Legendre, 2009) because inertial effects are very important. For the high values of Mo, which correspond to highly viscous fluids, the Mo and Bo are used to examine the bubble shape (Grace, 1973; Mukundakrishnan et al., 2007) in that the viscous effects become more important. In this study, the focus is on the bubble shape evolution during the interaction after the bubble is fully developed, the process can be determined by two dimensionless groups: Morton number (Mo) and Bond number (Bo). The Morton number is defined as
U t (mm s-1)
200
150
100 Experiments Previous simulation Present work
50
0
2
4
6
8
10
d (mm) Fig. 3. Comparison of our simulations for the terminal velocity vs. bubble diameter in a Newtonian liquid for Mo = 3.769 104 with experiments reported by Maxworthy et al. (1996) and numerical prediction by Tsamopoulos et al. (2008).
3.2. Parameter determination To investigate the bubble-rigid wall interaction problem, two dimensionless parameters need to be determined. For the low val-
Mo ¼
g l4f Dq
q2f r3
ð5Þ
;
where g is the acceleration of gravity, lf and qf are the viscosity and density of the surrounding fluid, Dq is the difference in density of the air and fluid, and r is the surface tension between the air bubble and the fluid. In our problem, air density qair in the bubble is neglected, so Mo can be rewritten as
Mo ¼
g l4f
qf r3
ð6Þ
:
The Bond number is defined as
Bo ¼
qf gR2 : r
(a)
ð7Þ
(b)
Fig. 4. (a) Comparison between coarse grids (left) and fine grids (right). (b) Bubble shapes: solid line for the coarse grid and dashed line for the fine grid.
63
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0
-0.05
(a)
(b) v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
(c)
(d) v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
(e)
(f) v
v
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05
(g)
(h)
Fig. 5. Vertical velocity component contours of interaction between a rising bubble and a rigid wall at Mo = 6.5 102, Bo = 6: (a) t = 0; (b) t = 5; (c) t = 10; (d) t = 15; (e) t = 20; (f) t = 30; (g) t = 50 and (h) t = 100.
It is important to note that the parameters Re, We, Mo, Bo are not independent of each other. Grace et al. (1976) gave a famous diagram to show the relationship among Mo, Bo and Re. From this diagram, each of the three parameters can be uniquely determined from the other two. We number can be also obtained by
We ¼ Re2
rffiffiffiffiffiffiffi Mo : Bo
ð8Þ
In the bubble-rigid wall interaction problem, we focus on the bubble shape evolution as it approaches the rigid wall. The simulations assume highly viscous liquids (Mo > 6.5 106). Therefore, the bubble rebound is unlikely.
64
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
20
RigidWall
y
yt and yb
19.5
Top surfaceyt Bottom surfaceyb
yt
yb
19
18.5
18 10
x
15
20
25
30
Time Fig. 6. A sketch of the top (yt) and bottom (yb) surfaces of the bubble. Fig. 8. The positions of the two sides of the bubble along y-axis for Mo = 6.5 102 and Bo = 6.
Surface distance (y t-yb)
2
Oh. As Oh increases, the viscous dissipation starts to dominate and the bubble may not bounce back. In our study, if we consider the value of CAm as 0.5 and deq as 2R, for qair qf, the Ca/St⁄ can be written as
1.5
rffiffiffiffiffiffiffi Ca Mo 9l2 ¼ ¼ 9 : St Bo Rqf r 1
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi The minimum value of 9 Mo=Bo, in our simulations, is larger than the rebound condition Ca/St⁄ < 102. Thus for the parameter range in our simulations, the bubble will not rebound. Although the bubble does not rebound as a whole, we do observe that, in some cases, some part of the bubble surface gains a downward velocity after interacting with the wall.
0.5
0
0
10
20
30
3.3. Mesh refinement
40
Time Fig. 7. The distance between two sides of the bubble for Mo = 6.5 10
2
and Bo = 6.
Zenit and Legendre (2009) carefully studied the coefficient of restitution for air bubbles colliding against solid walls. They indicated that the coefficient of restitution can be described by a function
rffiffiffiffiffiffiffi! Ca ¼ exp 30 ; St
ð9Þ
where Ca = U0l/r is the capillary number. The modified Stokes number St⁄ is defined as
St ¼
ð11Þ
ðqair þ C AM qf Þdeq U 0 ; 9l
ð10Þ
where CAM(v) is the added mass coefficient which depends on the bubble aspect ratio v (for a sphere, v = 1 and CAm = 0.5), deq is the pffiffiffiffiffiffiffiffiffiffiffiffiffiffi bubble equivalent diameter. It should be noted that Ca=St (the Ohnesorge number (Oh)) relates the viscous forces to inertial and surface tension forces. In their study, most rebounds are observed for Ca/St⁄ < 102. It means that the rebound only occurs at a small
There are two important processes in the problem of interaction between a fully developed bubble and a rigid surface. The first process is the rising of a bubble in highly viscous liquids. Since the interest is in the bubble shape evolution, it is important to resolve the bubble shape and obtain accurate terminal velocity. The detailed discussion will be presented in the next section. The second process is the collision between the deformed bubble and a rigid surface. It requires the computational mesh to be fine enough between the bubble surface and the rigid wall for detailed analysis of film drainage. In Fig. 2a, a case of fully deformed bubble colliding with a rigid surface is shown for Mo = 0.1 and Bo = 8.5. As can be seen more clearly in Fig. 2b, the mesh in the thin film between the two surfaces is adaptively refined. Even the thinnest part of the film contains at leat three layers of interior nodes to simulate the film drainage process.
4. Results and discussion First we will demonstrate that our numerical method produces accurate bubble shape, flow field structure and bubble terminal velocities at moderate Mo through comparison with known solutions. Grid independence is also shown by using coarse and fine grids. Then, we fix a certain Mo to discuss the bubble shape evolution by changing Bo. Three modes for bubble-rigid wall interaction
65
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
clearly that a dimple appears on the near-wall side of the bubble as shown in Fig. 5e. The formation of the dimple makes the smallest film thickness appears at the bubble rim rather than at its center. When the film between the rim of the bubble and rigid surface is thin enough to rupture, the bubble will attach to the rigid surface. However, due to limitations of the present numerical method, we cannot capture the film rupture. The simulation stops when the thickness of the film is smaller than 1/100 of the bubble radius. To investigate the bubble deformation and the interaction between the bubble and the rigid wall, we denote the instantaneous position of the intersections of bubble surface with the axis of symmetry as yb and yt for the bottom and top surface, respectively (see Fig. 6). Fig. 7 shows the temporal evolution of the vertical distance between the top and bottom sides of the bubble (yt–yb) along the axis of symmetry. The figure shows several phases of deformation of the bubble. Initially, because the bubble is released from rest and then begins to deform, the corresponding vertical distance deceases sharply from twice the bubble radius until the bubble is fully deformed after t 10. Then, since the bubble has been fully developed and is still far away from the rigid wall, the vertical distance changes very little. Afterwards, with the fully deformed bubble approaching the rigid wall, the vertical distance deceases sharply again and then rises a little. Finally, the vertical distance deceases very slowly until it reaches a constant value, which is caused by the film drainage process between the bubble and rigid wall. At this point, it is important to note that there exists a little increase in the vertical distance as the bubble comes near the rigid wall. To clarify the cause for this increase, we show the positions of both sides of the bubble versus time in Fig. 8. When the bubble interacts with the rigid wall, it can be clearly seen that the position of top surface almost becomes a constant whereas the bottom surface recedes away from the rigid wall. This leads to the little increase of separation between the top and bottom surfaces. The cause for the overshoot of the height of the bubble bottom surface around t = 20, as shown in Fig. 8, is the storage and release of the bubble surface energy. When the bubble interacts with the rigid wall, the surrounding liquid below the bubble still has an upward velocity. Then, a part of the kinetic energy of this liquid is transformed to the surface energy of the bubble. As the upward velocity vanishes, the bottom surface tends to restore to its original shape. This leads to the bottom surface receding away from the rigid wall. To illustrate this phenomenon, Fig. 9 presents the flow field with velocity contours for Mo = 6.5 102 and Bo = 6 at
have been distinguished. Finally, the boundaries between the three modes are presented in a (Mo, Bo)-plane. 4.1. Numerical method validation First, we compare our results with the experimental data by Maxworthy et al. (1996) who measured bubble terminal velocities as a function of bubble diameter in the mixture of distilled water with glycerin (Mo = 3.769 104) and with the numerical predictions by Tsamopoulos et al. (2008) who conducted their numerical simulation under the same condition by using a mixed finite element method in Fig. 3. We observe very good agreements between our simulations and the two previous studies. The code used in this study has been extensively validated by Yue et al. (2007). The bubble shape, terminal velocity, and the streamlines in the wake of a rising bubble match very well with the experiments by Hnat and Buckmaster (1976). The grid and time convergence studies were also presented in Yue et al. (2007). Beside those, we compare the coarse grid which is used in the following simulations with the fine one in Fig. 4. As shown in Fig. 4b, the bubble shapes predicted by the two grids are indistinguishable. 4.2. Effects of Bond number In our simulations, we release the bubble far below a solid surface at a distance of 15R. The bubble rises toward the solid wall by buoyancy only, and reaches a terminal velocity before interacting with the solid wall. Through numerical experiments, it is found that if we fix Mo and change Bo, there will be three modes of interaction. To explain the three modes, we set Mo = 6.5 102 and change Bo from 6 to 14. All parameters are non-dimensionalized by bubble radius, gravitational acceleration and fluid viscosity. 4.2.1. Mode one For small Bo, which means the ratio of buoyancy force to surface tension is small, the deformability of the bubble is weak. When the bubble rises towards a solid surface, a liquid film forms between the bubble and the solid wall. The high pressure inside the film causes a slight dimple on the bubble’s upper surface. Fig. 5 shows a time sequence of a gas bubble approaching a solid wall for Mo = 6.5 102 and Bo = 6. The bubble is released at t = 0 with zero velocity (see Fig. 5a). It then rises due to buoyancy force and reaches its terminal velocity (Vt = 0.85) and deformation (see Fig. 5c). As the bubble approaches the rigid surface, it can be seen
-2
-1
0
1
2
v
v
0.55 0.51 0.47 0.43 0.39 0.35 0.31 0.27 0.23 0.19 0.15 0.11 0.07 0.03 -0.01 -0.05
0.26 0.24 0.22 0.2 0.18 0.16 0.14 0.12 0.1 0.08 0.06 0.04 0.02 0 -0.02 -0.04
-2
-1
0
x
x
(a)
(b)
1
Fig. 9. Vertical velocity component contours for Mo = 6.5 102 and Bo = 6 at: (a) t = 18 and (b) t = 21.
66
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(a)
(b) v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(c)
(d) v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(e)
(f) v
v
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1 0.95 0.9 0.85 0.8 0.75 0.7 0.65 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
(g)
(h)
Fig. 10. Vertical velocity component contours of interaction between a rising bubble and a rigid wall at Mo = 102, Bo = 10: (a) t = 0; (b) t = 2.5; (c) t = 5; (d) t = 7.5; (e) t = 10; (f) t = 15; (g) t = 20 and (h) t = 30.
t = 18 and t = 21. It is obvious that the liquid below the bubble has an upward velocity at t = 18 and a downward velocity at t = 21. 4.2.2. Mode two As Bo increases, the inertia of the liquid becomes more important and the bubble is more easily deformed. When Bo reaches a critical value, the top and bottom sides of the bubble collide with
each other and the bubble may form a torus. Fig. 10 shows the deformation phases of a rising bubble and its interaction with a rigid surface for Mo = 6.5 102 and Bo = 10. In the rising stage, the fully developed bubble shape is a little flatter than that in mode one (Fig. 10c) and the terminal velocity increases to 0.88. In Fig. 10h, the final stage of interaction is depicted where a bubble with a large dimple is formed and the bubble will break up instantaneously.
67
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
2
Surface distance yt-yb
v 0.24 0.22 0.2 0.18 0.16 0.14 0.12 0.1 0.08 0.06 0.04 0.02 0 -0.02
1.5
1
0.5 -2
0
2
x 0
0
5
10
15
20
25
30
35
Time Fig. 11. The distance between two sides of the bubble for Mo = 6.5 102 and Bo = 10.
20
yt and yb
19.5
19 Top surface yt Bottom surface yb 18.5
18 10
15
20
25
30
Time Fig. 12. The positions of the two sides of the bubble along y-axis for Mo = 6.5 102 and Bo = 10.
To describe the process leading to the rupture of the bubble, a temporal evolution of the vertical distance between the top and bottom surfaces of the bubble along y-axis is displayed in Fig. 11. Before the bubble approaches the rigid wall, the vertical distance has a similar behavior as that in mode one. As the bubble interacts with the rigid wall, the vertical distance first increases (around t = 20), which is still the same as that in mode one. However, the distance decreases to zero right after this increase. This is caused by the deformation of the top surface of the bubble. As shown in Fig. 12, the change in the position of the bottom surface is similar to that in mode one. However, the top surface bounces back quickly after the bubble collides with the rigid wall. The reason is that the bubble is easier to deform and therefore a larger dimple is formed on the upper surface. If the upper surface reaches the lower surface, the bubble ruptures. The above analysis is consistent with the flow field shown in Fig. 13, in which the surrounding liquid near the upper surface of the bubble has a downward velocity while the lower surface is almost stationary.
Fig. 13. Vertical velocity component contours for Mo = 6.5 102 and Bo = 10 at t = 26.
4.2.3. Mode three As Bo is further increased, we expect a situation in which the bubble breaks up before it collides with the rigid wall due to increased bubble deformability. We denote this situation as mode three. As an example, a typical process of bubble-rigid surface interaction with Mo = 6.5 102 and Bo = 14 is presented in Fig. 14. For the high Bo in this mode, the ratio of the buoyancy force to surface tension is much higher than in the previous two modes, which leads to a significant bubble deformation. The fully developed bubble is flatter and the edge of the bubble becomes sharper (Fig. 14c) and the terminal velocity increases to 0.9. It should be noted that, the bubble will also break up in this mode (Fig. 14h), however, it is different from the rupture in mode two. To explain the reason, Fig. 15 shows a temporal evolution of the vertical distance between the top and bottom sides of the bubble along the axis of symmetry. Similar to the other modes, the bubble gets fully developed and reaches its terminal velocity very soon (t 7). However, different from the modes one and two, the vertical distance sharply decreases to zero without a temporary increase as the bubble interacts with the rigid wall. In Fig. 16, it can be clearly seen that when the bubble comes close to the rigid wall, the rising of the top surface slows down but the bottom surface continues to rise. Hence, after a while, the bottom surface will catch up with the top one and then the bubble breaks up. Fig. 17 shows the vertical velocity contours of the flow field and the surface of bubble. It can be observed that the surrounding liquid, above the top surface of the bubble is almost immobile but the liquid under the bottom surface still has an upward velocity. This will cause the two sides of the bubble to touch and the bubble ruptures. It is interesting to note that, as shown in Fig. 13, a dimple forms at the top surface of the bubble. At the edge of the film, the squeezing effect maximizes due to the smallest film thickness. As a consequence, a portion of the fluid there is squeezed inward into the film, which causes flow reversal inside the film. In Fig. 17, as Bo gets higher, the bubble has a higher terminal velocity and the flow inertia in the wake plays a more important role. Consequently, the bubble deforms more and ruptures before the dimple on the upper surface forms. Thus, there is no chance to observe the reversal flow in the film. It should also be noted that the film thickness in Fig. 17 is much thicker than the minimum thickness in Fig. 13 due to the earlier rupture of bubble. 4.3. Effects of Morton number In the last section, a moderate Morton number (Mo) case has been investigated and three modes have been identified by numer-
68
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(b)
(a)
v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(c)
(d) v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(f)
(e)
v
v
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1
(g)
(h)
Fig. 14. Vertical velocity component contours of interaction between a rising bubble and a rigid wall at Mo = 6.5 102, Bo = 14: (a) t = 0; (b) t = 2; (c) t = 4; (d) t = 6; (e) t = 8; (f) t = 10; (g) t = 14 and (h) t = 18.
ical simulations. Also, we find that the three modes can be distinguished from each other by a transition Bound number (Bo). In practice, due to flow instability and the roughness of the rigid wall, a precise transition Bo may be difficult to find. However, for convenience, it is reasonable to assume the existence of a transition Bo between two neighboring modes of interaction. Such a transition Bo will be helpful in understanding the physical mechanisms and
achieving a criterion for mode transition at different Mo. Here we should point out that the transition between different modes is regular and there is no sudden change. For the cases listed in Table 1, we have carefully examined the variation of transition Bo for different Mo. The transition Bo is established by incrementing Bo by 0.1 each step to investigate the problem and reduce the number of the time consuming simu-
69
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
2
Surface distance yt-yb
v 0.6 0.55 0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1
1.5
1
0.5
-4 0
-3
-2
-1
0
1
2
3
4
x 0
5
10
15
20
Time
Fig. 17. Vertical velocity component contours for Mo = 6.5 102 and Bo = 14 at t = 17.5.
Fig. 15. The distance between two sides of the bubble for Mo = 6.5 102 and Bo = 14. Table 1 Transition Bond number for different Morton numbers.
20 Top surfacey t Bottom surfaceyb
yt and yb
19.5
19
Mo
Transition Bo between modes one and two
Transition Bo between modes two and three
1.0 101 6.5 102 2.0 102 8.4 103 4.0 103 9.8 104 1.0 104 3.3 105 6.5 106
8.7 8.3 7.7 7.4 6.9 5.8 4.1 3.4 2.7
15.3 13.9 11.0 9.3 8.1 6.4 4.6 3.7 2.8
16
18.5
Mode 2 12
14
16
18
Time Fig. 16. The positions of the two sides of the bubble along y-axis for Mo = 6.5 102 and Bo = 14.
lations. In our results, Mo varies from 1.0 101 to 6.5 106 and the corresponding transition Bo between mode one and mode two decreases from 8.7 to 2.7, and similarly the one between mode two and mode three drops from 15.3 to 2.8. Furthermore, the difference between two transition Bond numbers at the same Mo is also lessened from 6.6 to 0.1 with the decreasing Mo. That means that the region of mode two shrinks as the Mo becomes small. In other words, the existance of the three modes is more distinct at a high Mo. To clearly present the results, we describe the boundaries between different modes on the Mo-Bo chart. In Fig. 18, we represent Mo on x-direction with log scale and display Bo in y-direction. The two curves shown in this figure give the boundaries between different modes of interaction. The top one distinguishes mode two from mode three and the bottom one separates mode one from mode two. In this figure, for large Mo (>102), it is easy to distinguish the three modes from each other. When Mo is large, the viscosity of the fluid plays an important role. The surface tension,
Mode 3
12
20
Bo number
18 10
8
4
0 -6 10
Mode 1
-5
10
-4
10
-3
10
-2
10
-1
10
0
10
Mo number Fig. 18. Three modes diagram in the (Mo, Bo)-plane.
which determines Bo, needs to vary in a wide range to obtain different modes of bubble deformation. However, for small Mo number (<104), the surface tension dominates the bubble deformation. The gas bubble either breaks up at a high Bo or remains unbroken at a low Bo, which causes the region of mode two to decrease. Especially for the Mo < 105, the two curves merge into one and the region of mode two disappears. In addition, for the case of
70
T. Qin et al. / International Journal of Multiphase Flow 52 (2013) 60–70
Mo < 105 and in the range of Bo given in Fig. 18, the shape of the bubble remains spherical, see Grace et al. (1976) and Tsamopoulos et al. (2008). There may exist more modes, for example, the bubble rebound and vibration as observed experimentally (Tsao and Koch, 1997; Zenit and Legendre, 2009; Hendrix et al., 2012). However, the bubble deformation in these studies is small and the corresponding Mo is beyond the range considered in the present study. It should be noted that Mo depends only on the material properties of the viscous liquid. However, Bo depends on both the properties of the viscous liquid and the bubble size. Thus, in practice, we can predict the critical bubble size in a given viscous liquid using the above diagram, according to Mo and Bo. For example, if a solid container is filled with cyclohexanol, the properties of which are q = 962 kg/m3, l = 0.0575 Pa s and r = 0.0329 N/m (Haynes (2011-2012)), we obtain Morton number Mo = 3 103. From the above diagram, we can determine the transition Bo between mode one and mode two to be Bo = 6.6. Thus, we can predict that the radius of the bubble which can attach to the upper wall of the solid container without break up should be smaller than 5 mm at its undeformed state. 5. Conclusion In this article, direct numerical simulations have been performed by using an arbitrary-Lagrangian–Eulerian method in order to study the interaction between a deformable bubble and a rigid wall and to clarify the bubble deformation as it approaches the rigid wall. It is found that three modes of bubble-rigid wall interaction exist as Bo changes at a moderate Mo (6.5 106 < Mo < 0.1). The first mode prevails at small Bo where the bubble deformation is small. For this mode, the bubble is hard to break up and will bounce back and eventually attach to the rigid wall. In the second mode, the bubble may break up after it collides with the rigid wall, which is determined by the film drainage. In the third mode, which prevails at high Bo, the bubble breaks up due to the bottom surface catches up the top surface during the interaction with a rigid wall. We also found that the three modes are easily distinguishable when Mo is large. When Mo gets small, the region of the second mode diminishes and finally disappears. We have to point out that the current work assumes axisymmetric flow. For the bubble rising problem, according to the experiments by Bhaga and Weber (1981), the flow should be axisymmetric in the Reynolds number range (Re < 50) in this work. However, for the more complicated bubble-rigid surface interaction problem, the flow axisymmetry and actual existence of the three modes should be established by future 3D simulations or experiments. Acknowledgments The authors thank the anonymous reviewers for helpful comments. P. Yue acknowledges the support by NSF-DMS 0907788. References Bhaga, D., Weber, M.E., 1981. Bubbles in viscous liquids: shapes, wakes and velocities. J. Fluid Mech. 105, 61–85. Duineveld, P.C., 1995. The rise velocity and shape of bubbles in pure water at high reynolds number. J. Fluid Mech. 292, 325–332. Frank, X., Funfschilling, D., Midoux, N., Li, H.Z., 2006. Bubbles in a viscous liquid: lattice boltzmann simulation and experimental validation. J. Fluid Mech. 546, 113–122.
Fujasova-Zednikova, M., Vobecka, L., Vejrazka, J., 2010. Effect of solid material and surfactant presence on interactions of bubbles with horizontal solid surface. Can. J. Chem. Eng. 88, 473–481. Gaudlitz, D., Adams, N.A., 2009. Numerical investigation of rising bubble wake and shape variations. Phys. Fluids 21, 122102. Grace, J.R., 1973. Shapes and velocities of bubbles rising in infinite liquids. Trans. Instn. Chem. Eng. 51, 116–120. Grace, J.R., Wairegi, T., Nguyen, T.H., 1976. Shapes and velocities of single drops and bubbles moving freely through immiscible liquids. Trans. Instn. Chem. Eng. 54, 167–173. Gueyffier, D., Li, J., Nadim, A., Scardovelli, R., Zaleski, S., 1999. Volume-of-fluid interface tracking with smoothed surface stress methods for three-dimensional flows. J. Comput. Phys. 152, 423–456. Haynes, W.M., 2011-2012. CRC Handbook of Chemistry and Physics, 92nd ed. CRC Press. Hendrix, M.H.W., Manica, R., Klaseboer, E., Chan, D.Y.C., Ohl, C., 2012. Spatiotemporal evolution of thin liquid films during impact of water bubbles on glass on a micrometer to nanometer scale. Phys. Rev. Lett. 108, 247803. Hnat, J.G., Buckmaster, J.D., 1976. Spherical cap bubbles and skirt formation. Phys. Fluids 19, 182–194. Hu, H.H., 1996. Direct simulation of flows of solid–liquid mixtures. Int. J. Multiphase Flow 22, 335–352. Hu, H.H., Patankar, N.A., Zhu, M.Y., 2001. Direct numerical simulations of fluid–solid systems using the arbitrary lagrangian-eulerian technique. J. Comput. Phys. 169, 427–462. Klaseboer, E., Chevaillier, J., Mate, A., Masbernat, O., Gourdon, C., 2001. Model and experiments of a drop impinging on an immersed wall. Phys. Fluids 13, 45–57. Kushch, V.I., Sangani, A.S., Spelt, P.D.M., Koch, D.L., 2002. Finite-weber-number motion of bubbles through a nearly inviscid liquid. J. Fluid Mech. 460, 241–280. Legendre, D., Daniel, C., Guiraud, P., 2005. Experimental study of a drop bouncing on a wall in a liquid. Phys. Fluids 17, 097105. Maxworthy, T., Gnann, C., Kurten, M., Durst, F., 1996. Experiments on the rise of air bubbles in clean viscous liquids. J. Fluid Mech. 321, 421–441. Minev, P.D., Chen, T., Nandakumar, K., 2003. A finite element technique for multifluid imcompressible flow using Eulerian grids. J. Comput. Phys. 187, 255– 273. Moore, D.W., 1965. The velocity of rise of distorted gas bubbles in a liquid of small viscosity. J. Fluid Mech. 23, 749–766. Mukundakrishnan, K., Quan, S., Eckmann, D.M., Ayyaswamy, P.S., 2007. Numerical study of wall effects on buoyant gas-bubble rise in a liquid-filled finite cylinder. Phys. Rev. E 76, 036308. Ohta, M., Imura, T., Yoshida, Y., Sussman, M., 2010. A computational study of the effect of initial bubble conditions on the motion of a gas bubble rising in viscous liquids. Int. J. Multiphase Flow 31, 223–237. Power, H., 1997. The interaction of a deformable bubble with a rigid wall at small reynolds number: a general approach via integral equations. Eng. Anal. Bound. Elem. 19, 291–297. Rabha, S.S., Buwa, V.V., 2010. Volume-of-fluid (vof) simulations of rise of single/ multiple bubbles in sheared liquids. Chem. Eng. Sci. 65, 527–537. Shopov, P.J., Minev, P.D., Bazhlekov, I.B., Zapryanov, Z.D., 1990. Interaction of a deformable bubble with a rigid wall at moderate reynolds numbers. J. Fluid Mech. 219, 241–271. Tryggvason, G., Bunner, B., Esmaeeli, A., Juric, D., Al-Rawahi, N., Tauber, W., Han, J., Nas, S., Jan, Y.J., 2001. A front-tracking method for the computations of multiphase flow. J. Comput. Phys. 169, 708–759. Tryggvason, G., Esmaeeli, A., Lu, J., Biswas, S., 2006. Direct numerical simulations of gas/liquid multiphase flows. Fluid Dyn. Res. 38, 660–681. Tsamopoulos, J., Dimakopoulos, Y., Chatzidai, N., Karapetsas, G., Pavlidis, M., 2008. Steady bubble rise and deformation in newtonian and viscoplastic fluids and conditions for bubble entrapment. J. Fluid Mech. 601, 123–164. Tsao, H., Koch, D.L., 1997. Observations of high reynolds number bubbles interacting with a rigid wall. Phys. Fluids 9, 44–56. Wilkes, E.D., Phillips, S.D., Basaran, O.A., 1999. Computational and experimental analysis of dynamics of drop formation. Phys. Fluids 11, 3577–3598. Yiantsios, S.G., Davis, R.H., 1990. On the buoyancy-driven motion of a drop towards a rigid surface or a deformable interface. J. Fluid Mech. 217, 547–573. Yu, Z., Fan, L., 2008. Direct simulations of the buoyant rise of bubbles in infinite liquid using level set method. Can. J. Chem. Eng. 86, 267–275. Yue, P., Feng, J.J., Bertelo, C.A., Hu, H.H., 2007. An arbitrary Lagrangian–Eulerian method for simulating bubble growth in polymer foaming. J. Comput. Phys. 226, 2229–2249. Yue, P., Zhou, C., Feng, J.J., Ollivier-Gooch, C.F., Hu, H.H., 2006. Phase-field simulations of interfacial dynamics in viscoelastic fluids using finite elements with adaptive meshing. J. Comput. Phys. 219, 47–67. Zawala, J., Krasowska, M., Dabros, T., Malysa, K., 2007. Influence of bubble kinetic energy on its bouncing during collisions with various interfaces. Can. J. Chem. Eng. 85, 669–678. Zenit, R., Legendre, D., 2009. The coefficient of restitution for air bubbles colliding against solid walls in viscous liquids. Phys. Fluids 21, 083306.