Hostname: page-component-7c8c6479df-r7xzm Total loading time: 0 Render date: 2024-03-26T14:25:15.011Z Has data issue: false hasContentIssue false

The shape and motion of gas bubbles in a liquid flowing through a thin annulus

Published online by Cambridge University Press:  21 September 2018

Q. Lei*
Affiliation:
Department of Earth Science and Engineering, Imperial College London, SW7 2AZ, UK
Z. Xie
Affiliation:
Department of Earth Science and Engineering, Imperial College London, SW7 2AZ, UK Department of Chemical Engineering, Imperial College London, SW7 2AZ, UK School of Engineering, Cardiff University, CF24 3AA, UK
D. Pavlidis
Affiliation:
Department of Earth Science and Engineering, Imperial College London, SW7 2AZ, UK
P. Salinas
Affiliation:
Department of Earth Science and Engineering, Imperial College London, SW7 2AZ, UK
J. Veltin
Affiliation:
TNO, 2628 CK, The Netherlands
O. K. Matar
Affiliation:
Department of Chemical Engineering, Imperial College London, SW7 2AZ, UK
C. C. Pain
Affiliation:
Department of Earth Science and Engineering, Imperial College London, SW7 2AZ, UK
A. H. Muggeridge
Affiliation:
Department of Earth Science and Engineering, Imperial College London, SW7 2AZ, UK
A. J. Gyllensten
Affiliation:
Statoil ASA, 3936, Norway
M. D. Jackson
Affiliation:
Department of Earth Science and Engineering, Imperial College London, SW7 2AZ, UK
*
Email address for correspondence: q.lei12@imperial.ac.uk

Abstract

We study the shape and motion of gas bubbles in a liquid flowing through a horizontal or slightly inclined thin annulus. Experimental data show that in the horizontal annulus, bubbles develop a unique ‘tadpole-like’ shape with a semi-circular cap and a highly stretched tail. As the annulus is inclined, the bubble tail tends to vanish, resulting in a significant decrease of bubble length. To model the bubble evolution, the thin annulus is conceptualised as a ‘Hele-Shaw’ cell in a curvilinear space. The three-dimensional flow within the cell is represented by a gap-averaged, two-dimensional model, which achieved a close match to the experimental data. The numerical model is further used to investigate the effects of gap thickness and pipe diameter on the bubble behaviour. The mechanism for the semi-circular cap formation is interpreted based on an analogous irrotational flow field around a circular cylinder, based on which a theoretical solution to the bubble velocity is derived. The bubble motion and cap geometry is mainly controlled by the gravitational component perpendicular to the flow direction. The bubble elongation in the horizontal annulus is caused by the buoyancy that moves the bubble to the top of the annulus. However, as the annulus is inclined, the gravitational component parallel to the flow direction becomes important, causing bubble separation at the tail and reduction in bubble length.

Type
JFM Papers
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution, and reproduction in any medium, provided the original work is properly cited.
Copyright
© 2018 Cambridge University Press

1 Introduction

Simultaneous flow of gas and liquid through a conduit, over a wide range of flow rates, exhibits a ‘slug flow’ pattern, which consists of the pseudo-periodic appearance of large bubbles (often called Taylor bubbles) separated by liquid slugs (Fabre & Liné Reference Fabre and Liné1992). The shape and motion of Taylor bubbles in tubes have been extensively investigated over the past decades (Davies & Taylor Reference Davies and Taylor1950; White & Beardmore Reference White and Beardmore1962; Zukoski Reference Zukoski1966; Collins Reference Collins1967b ; Collins et al. Reference Collins, de Moraes, Davidson and Harrison1978; Reinelt Reference Reinelt1987b ; Vanden-Broeck Reference Vanden-Broeck1984, Reference Vanden-Broeck1992; Viana et al. Reference Viana, Pardo, Yánez, Trallero and Joseph2003; Figueroa-Espinoza & Fabre Reference Figueroa-Espinoza and Fabre2011; Fabre Reference Fabre2016). In contrast, much less effort has been devoted to studying Taylor bubbles in annuli, which are, however, very relevant for many industrial applications such as energy dissipation facilities for chemical/nuclear reactors, casing-tubing annuli for oil production and double-pipe heat exchangers for geothermal systems (Kelessidis & Dukler Reference Kelessidis and Dukler1990; Das et al. Reference Das, Das, Purohit and Mitra1998; Ekberg et al. Reference Ekberg, Ghiaasiaan, Abdel-Khalik, Yoda and Jeter1999).

In an annulus, flow is confined between the inner surface of an outer pipe and the outer surface of an inner pipe. If the thickness of the gap between the two surfaces is much smaller than the perimeter of the pipe, the annulus is then similar to a Hele-Shaw cell. The problem of two-phase flow in a Hele-Shaw cell between two parallel flat plates has been the focus of many previous theoretical, experimental and numerical studies (e.g. Taylor & Saffman Reference Taylor and Saffman1959; Collins Reference Collins1965a , Reference Collins1967a ; Maxworthy Reference Maxworthy1986; Kopf-Sill & Homsy Reference Kopf-Sill and Homsy1988; Roig et al. Reference Roig, Roudet, Risso and Billet2012; Cueto-Felgueroso & Juanes Reference Cueto-Felgueroso and Juanes2014). Large bubbles in a Hele-Shaw cell (the bubble characteristic length, $d$ , is much greater than the gap thickness between the plates, $h$ ) behave very differently from bubbles in an unbounded liquid due to the presence of wall confinement and wettability effects (Thompson Reference Thompson1968; Eck & Siekmann Reference Eck and Siekmann1978; Saffman & Tanveer Reference Saffman and Tanveer1989; Cueto-Felgueroso & Juanes Reference Cueto-Felgueroso and Juanes2014). The bubble dynamics in a Hele-Shaw cell can be characterised by the bubble Reynolds number $Re=\unicode[STIX]{x1D70C}_{l}u_{\infty }d/\unicode[STIX]{x1D707}_{l}$ ( $u_{\infty }$ is the bubble velocity relative to the ambient liquid having a density $\unicode[STIX]{x1D70C}_{l}$ and a dynamic viscosity $\unicode[STIX]{x1D707}_{l}$ ). In addition, the gap Reynolds number $Re(h/d)^{2}$ is often used to compare the in-plane inertial effect and the cross-gap viscous effect (Thompson Reference Thompson1968; Eck & Siekmann Reference Eck and Siekmann1978; Bush & Eames Reference Bush and Eames1998; Roig et al. Reference Roig, Roudet, Risso and Billet2012). When $Re(h/d)^{2}\ll 1$ , the flow corresponds to the classical Hele-Shaw flow regime where the inertia is negligible (Saffman & Taylor Reference Saffman and Taylor1958; Taylor & Saffman Reference Taylor and Saffman1959; Eck & Siekmann Reference Eck and Siekmann1978; Maxworthy Reference Maxworthy1986; Tanveer Reference Tanveer1986, Reference Tanveer1987; Tanveer & Saffman Reference Tanveer and Saffman1987; Kopf-Sill & Homsy Reference Kopf-Sill and Homsy1988; Meiburg Reference Meiburg1989; Saffman & Tanveer Reference Saffman and Tanveer1989; Maruvada & Park Reference Maruvada and Park1996; Cueto-Felgueroso & Juanes Reference Cueto-Felgueroso and Juanes2014). When $Re(h/d)^{2}$ is close to or larger than 1, the inertia becomes important, as has been observed for bubbles in inclined or vertical Hele-Shaw cells (Collins Reference Collins1965a ,Reference Collins b ; Grace & Harrison Reference Grace and Harrison1967; Lazarek & Littman Reference Lazarek and Littman1974; Maneri & Zuber Reference Maneri and Zuber1974; Hills Reference Hills1975; Vanden-Broeck Reference Vanden-Broeck1984, Reference Vanden-Broeck1992; Couet & Strumolo Reference Couet and Strumolo1987; Bush Reference Bush1997; Kelley & Wu Reference Kelley and Wu1997; Bush & Eames Reference Bush and Eames1998; Huisman, Ern & Roig Reference Huisman, Ern and Roig2012; Roig et al. Reference Roig, Roudet, Risso and Billet2012; Wang et al. Reference Wang, Klaasen, Degrève, Blanpain and Verhaeghe2014, Reference Wang, Klaasen, Degrève, Mahulkar, Heynderickx, Reyniers, Blanpain and Verhaeghe2016; Filella, Ern & Roig Reference Filella, Ern and Roig2015; Piedra, Ramos & Herrera Reference Piedra, Ramos and Herrera2015).

In this paper, we study bubble behaviour in another type of Hele-Shaw cell, i.e. thin annuli. We provide an insight, based on experimental observations, numerical simulations and theoretical analysis, into simultaneous gas–liquid flow through such systems. We focus on the regime where inertia becomes important, i.e. $Re(h/d)^{2}$ is close to or larger than unity. Laboratory experiments are conducted to measure the bubble distribution and evolution in a horizontal or slightly inclined thin annulus (§ 2). A gap-averaged, two-dimensional (2-D) formulation of the Navier–Stokes equations is then derived to represent the three-dimensional (3-D) flow dynamics (§ 3). The numerical model predictions are compared with the experimental data for different inclinations. Further numerical simulations are then performed to elucidate the roles of gap thickness and pipe diameter in bubble formation (§ 4). In addition, a theoretical analysis is presented to interpret the flow mechanisms that govern bubble dynamics in thin annuli (§ 5). Finally, a brief discussion is given and conclusions are drawn (§ 6).

2 Experimental apparatus and procedure

The experimental apparatus consisted of air and water supply systems, an annulus formed between two concentric pipes (comprising three connected subsections and initially filled with water), and a collection tank (figure 1). The parameters of the experimental apparatus and working fluids are as given in table 1. During the experiments, air and dyed water were continuously pumped into a mixer before entering the test section. Visual observation showed the two fluids were well segregated by the time they entered the annulus. Three high-speed cameras were mounted above the annulus subsections for visualisation and recording purposes. The outer pipe was made of Plexiglas to facilitate visualisation, while the inner stainless steel tubing was spray coated to create a white surface which maximised the colour contrast between the air and dyed water. The videos were taken at a frame rate of 17 frames per second with a resolution of $1280\times 720$ pixels, and a colour depth of 24 bits per pixel. The data were imported into MATLAB for post-processing to extract bubble size and velocity data.

Figure 1. Schematic illustration of the experimental apparatus.

3 Mathematical formulation and numerical model

3.1 Governing equations

Figure 2. Schematic of (a) the thin annulus and (b) the equivalent Hele-Shaw cell. (c) An illustration showing how the direction of gravitational acceleration changes along the equivalent Hele-Shaw cell to ensure the flow replicates that seen in the actual annulus.

Table 1. Parameters of the experimental apparatus and fluids.

As shown in figure 2(a), we define a 3-D curvilinear coordinate system ( $x,y,z$ ) based on the mid-plane of the thin annulus (i.e. the distance of this mid-plane from both the outer and inner walls is $h/2$ ), where $x$ is defined to be along the pipe length ( $0\leqslant x\leqslant L$ ), $y$ is around the perimeter ( $-W/2\leqslant y\leqslant W/2$ with $y=0$ and $y=\pm W/2$ corresponding to the top and bottom of the pipe, respectively), and $z$ is the direction across the gap ( $-h/2\leqslant z\leqslant h/2$ with $z=h/2$ and $z=-h/2$ defining the outer and inner channel walls, respectively). The system is therefore equivalent to a Hele-Shaw cell (viz. figure 2 b) under a distorted gravity field (viz. figure 2 c).

We use the approximation that the flow velocity in the $z$ direction is negligible ( $u_{z}=0$ ) given the very high aspect ratio (i.e. $h/W\ll 1$ ). By eliminating the terms involving $u_{z}$ , the governing equations of the incompressible two-phase, air–water flow are obtained, including the mass and momentum conservation equations given respectively by

(3.1) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}=0, & \displaystyle\end{eqnarray}$$
(3.2) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70C}\left(\frac{\unicode[STIX]{x2202}\boldsymbol{u}}{\unicode[STIX]{x2202}t}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}\right)=-\unicode[STIX]{x1D735}p+\unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D707}[\unicode[STIX]{x1D735}\boldsymbol{u}+(\unicode[STIX]{x1D735}\boldsymbol{u})^{\text{T}}]+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}z}\left(\unicode[STIX]{x1D707}\frac{\unicode[STIX]{x2202}\boldsymbol{u}}{\unicode[STIX]{x2202}z}\right)+\unicode[STIX]{x1D70C}\boldsymbol{g}+\unicode[STIX]{x1D70E}\unicode[STIX]{x1D705}\unicode[STIX]{x1D6FF}\boldsymbol{n}, & \displaystyle\end{eqnarray}$$

and an equation for the volume fraction of air given as

(3.3) $$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x2202}t}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D6FC}=0. & & \displaystyle\end{eqnarray}$$

Here, $\boldsymbol{u}=[u_{x},u_{y}]$ is the velocity vector, $\unicode[STIX]{x1D735}=[\unicode[STIX]{x2202}/\unicode[STIX]{x2202}x,\unicode[STIX]{x2202}/\unicode[STIX]{x2202}y]$ is the 2-D gradient, $t$ is the time, $\unicode[STIX]{x1D6FC}$ is the volume fraction of air, the bulk density is $\unicode[STIX]{x1D70C}=\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D70C}_{a}+(1-\unicode[STIX]{x1D6FC})\unicode[STIX]{x1D70C}_{w}$ , $p$ is the pressure, the bulk dynamic viscosity is $\unicode[STIX]{x1D707}=\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D707}_{a}+(1-\unicode[STIX]{x1D6FC})\unicode[STIX]{x1D707}_{w}$ , the gravitational acceleration vector is $\boldsymbol{g}=[g\sin \unicode[STIX]{x1D6FD},g\cos \unicode[STIX]{x1D6FD}\sin (2y/d)]$ , $\unicode[STIX]{x1D70E}$ is the surface tension, $\unicode[STIX]{x1D705}$ is the interface curvature, $\unicode[STIX]{x1D6FF}$ is the Dirac delta function and $\boldsymbol{n}$ is the interface outward-pointing unit normal. Surface tension force is treated as a continuous surface force (Brackbill, Kothe & Zemach Reference Brackbill, Kothe and Zemach1992). Note (3.1)–(3.3) are written in a form following Bush (Reference Bush1997) that eases the derivation of gap-averaged, 2-D equations in the following § 3.2.

3.2 Gap-averaging process

We assume the velocity profile in the annulus to be parabolic in the $z$ -direction (Gondret & Rabaud Reference Gondret and Rabaud1997):

(3.4) $$\begin{eqnarray}\displaystyle \boldsymbol{u}=\frac{3}{2}\left[1-\left(\frac{2z}{h}\right)^{2}\right]\bar{\boldsymbol{u}}, & & \displaystyle\end{eqnarray}$$

where $\bar{\boldsymbol{u}}$ is the gap-averaged velocity. This is consistent with our calculation $\unicode[STIX]{x1D70C}_{w}v_{i}(2h)/\unicode[STIX]{x1D707}_{w}=1134$ , which reveals the flow to be laminar (Beavers, Sparrow & Magnuson Reference Beavers, Sparrow and Magnuson1970).

We derive a coupled set of gap-averaged equations by integrating equations (3.1)–(3.3) along the $z$ direction (Roig et al. Reference Roig, Roudet, Risso and Billet2012):

(3.5) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\bar{\boldsymbol{u}}=0, & \displaystyle\end{eqnarray}$$
(3.6) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70C}\left(\frac{\unicode[STIX]{x2202}\bar{\boldsymbol{u}}}{\unicode[STIX]{x2202}t}+\frac{6}{5}\bar{\boldsymbol{u}}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\bar{\boldsymbol{u}}\right)=-\unicode[STIX]{x1D735}p+\unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D707}[\unicode[STIX]{x1D735}\bar{\boldsymbol{u}}+(\unicode[STIX]{x1D735}\bar{\boldsymbol{u}})^{\text{T}}]-\frac{12\unicode[STIX]{x1D707}}{h^{2}}\bar{\boldsymbol{u}}+\unicode[STIX]{x1D70C}\boldsymbol{g}+\unicode[STIX]{x1D70E}\unicode[STIX]{x1D705}\unicode[STIX]{x1D6FF}\boldsymbol{n}, & \displaystyle\end{eqnarray}$$
(3.7) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x2202}t}+\bar{\boldsymbol{u}}\boldsymbol{\cdot }\unicode[STIX]{x1D6FB}\unicode[STIX]{x1D6FC}=0. & \displaystyle\end{eqnarray}$$

We acknowledge the fact that (3.5)–(3.7) may not represent fully all relevant dynamics, e.g. 3-D effects local to the air–water interface (Oliveira & Meiburg Reference Oliveira and Meiburg2011). It is, nonetheless, of interest to determine the extent to which our approximate model can capture as many of the phenomena observed experimentally as possible.

3.3 Curvature and contact angle

The shape of the air–water interface is affected by the wettability of the channel walls. To incorporate this effect, the curvature in (3.6) may be replaced by (Thompson, Juel & Hazel Reference Thompson, Juel and Hazel2014)

(3.8) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D705}=\unicode[STIX]{x1D705}_{l}+\unicode[STIX]{x1D705}_{t}\approx \unicode[STIX]{x1D705}_{l}+\frac{\cos \unicode[STIX]{x1D703}_{o}+\cos \unicode[STIX]{x1D703}_{i}}{h}, & & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D705}_{l}$ is the lateral curvature (on the $x$ $y$ plane) which is estimated using the diffuse-interface approach based on the volume fraction field (Xie et al. Reference Xie, Pavlidis, Salinas, Percival, Pain and Matar2016), $\unicode[STIX]{x1D705}_{t}$ is the transverse curvature (across the gap) which is estimated from the meniscus profile (figure 3 a), $\unicode[STIX]{x1D703}_{o}$ is the local contact angle at the outer wall and $\unicode[STIX]{x1D703}_{i}$ is that at the inner wall. The local contact angle along the contact line can vary significantly from the advancing value (at the bubble tail) to the receding value (at the bubble front). The experimental measurements by Antonini et al. (Reference Antonini, Carmona, Pierce, Marengo and Amirfazli2009) showed that the angle can be interpolated using a third-order polynomial function:

(3.9) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D703}=2\frac{\unicode[STIX]{x1D703}_{min}-\unicode[STIX]{x1D703}_{max}}{\unicode[STIX]{x03C0}^{3}}\unicode[STIX]{x1D719}^{3}-3\frac{\unicode[STIX]{x1D703}_{min}-\unicode[STIX]{x1D703}_{max}}{\unicode[STIX]{x03C0}^{2}}\unicode[STIX]{x1D719}^{2}+\unicode[STIX]{x1D703}_{min}, & & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D703}$ corresponds to the local contact angle at either the outer or inner wall (i.e. $\unicode[STIX]{x1D703}_{o}$ or $\unicode[STIX]{x1D703}_{i}$ , respectively), $\unicode[STIX]{x1D703}_{max}$ and $\unicode[STIX]{x1D703}_{min}$ denote the advancing and receding contact angles of the corresponding surface, respectively, $\unicode[STIX]{x1D719}$ is the angle between the local velocity vector $\bar{\boldsymbol{u}}$ and the interface normal $\boldsymbol{n}$ (figure 3 b). Figure 4 shows the variation of the local contact angle $\unicode[STIX]{x1D703}$ at the outer and inner walls with respect to $\unicode[STIX]{x1D719}$ . Contact angle hysteresis is promoted by roughness of the channel wall surface (Dussan Reference Dussan1979); we also discuss in § 6 the potential effect of roughness on the presence/absence of thin water films on the channel walls.

Figure 3. Schematic of (a) the cross-gap meniscus profile governed by the local contact angles, $\unicode[STIX]{x1D703}_{o}$ and $\unicode[STIX]{x1D703}_{i}$ , at the outer and inner walls, respectively, and (b) the moving contact line in the $x$ $y$ plane with the local advancing/receding state determined by the intersection angle $\unicode[STIX]{x1D719}$ between the velocity vector $\bar{\boldsymbol{u}}$ and the interface normal $\boldsymbol{n}$ .

Figure 4. Interpolation of contact angles around a contact line at the outer or inner wall of the annulus using a third-order polynomial fitting (Antonini et al. Reference Antonini, Carmona, Pierce, Marengo and Amirfazli2009);  $\unicode[STIX]{x1D703}$ denotes the local contact angle and $\unicode[STIX]{x1D719}$ denotes the intersection angle between the local velocity vector and the interface normal as shown in figure 3.

3.4 Model set-up and numerical methods

A numerical experiment was designed as shown in figure 5. The length and width of the 2-D domain were the pipe length $L$ and perimeter $W$ , respectively. This domain was originally filled with water that had an initial velocity $v_{i}$ along the $x$ direction. The inlet condition was defined by the inflow velocity $v_{i}$ uniformly imposed at the left boundary, such that the air flowed into the domain through the middle part ( $W/5$ wide), corresponding to the inflow rate condition of $q_{a}/(q_{a}+q_{w})=1/5$ (table 1), while the water was injected from the remaining portion of the left boundary. This is consistent with the experimental observation that air and water were well segregated in the mixer before entering the annulus. A free-slip condition was used for the two longitudinal sides of the domain, considering that the air phase is away from the pipe bottom and the flow across the pipe bottom is negligible. A hydrostatic pressure condition was imposed along the outlet to account for the gravity-induced pressure gradient across the channel.

Figure 5. The 2-D model set-up for numerical simulation.

The numerical model was solved based on a mixed control-volume finite-element method, which has been validated against a series of benchmark cases of single rising bubbles, coalescence of two bubbles and droplet impacts (Xie et al. Reference Xie, Pavlidis, Percival, Gomes, Pain and Matar2014, Reference Xie, Pavlidis, Salinas, Percival, Pain and Matar2016, Reference Xie, Hewitt, Pavlidis, Salinas, Pain and Matar2017). The computational domain was discretised into an unstructured grid of triangular $\text{P}_{1}\text{DG}\text{-}\text{P}_{2}$ elements (linear discontinuous velocity between elements and quadratic continuous pressure between elements) (Cotter et al. Reference Cotter, Ham, Pain and Reich2009). A finite volume discretisation of the mass conservation equation and a linear discontinuous Galerkin discretisation of the momentum equation were used with an adaptive implicit/explicit time stepping scheme (Xie et al. Reference Xie, Pavlidis, Salinas, Percival, Pain and Matar2016). Within each time step, the equations were iterated upon using a pressure projection method until all equations are simultaneously balanced (Pavlidis et al. Reference Pavlidis, Xie, Percival, Gomes, Pain and Matar2014; Xie et al. Reference Xie, Pavlidis, Percival, Gomes, Pain and Matar2014). Interface dynamics were captured through a compressive advection-based volume-of-fluid approach, which used a novel and mathematically rigorous nonlinear Petrov–Galerkin method for maintaining sharp interfaces (Pavlidis et al. Reference Pavlidis, Gomes, Xie, Percival, Pain and Matar2016). Surface tension was modelled using a diffused-interface formulation that can accurately estimate the lateral curvature based on the volume fraction field (Xie et al. Reference Xie, Pavlidis, Salinas, Percival, Pain and Matar2016). Anisotropic mesh adaptation (Pain et al. Reference Pain, Umpleby, de Oliveira and Goddard2001) was used to place a finer mesh around the interface during its evolution as required to capture the dynamics (Xie et al. Reference Xie, Pavlidis, Percival, Gomes, Pain and Matar2014, Reference Xie, Pavlidis, Salinas, Percival, Pain and Matar2016). The criteria for mesh refinement, based on the volume fraction field, included an interpolation error bound of 0.1, a minimum element size of $W/60$ and maximum element size of $L/6$ and $W/20$ in the $x$ and $y$ directions, respectively. The 2-D simulation outputs were projected onto the original 3-D Cartesian coordinate system to facilitate visual comparison with experimental data.

4 Results

4.1 Experimental observation and numerical simulation

Figures 69 show the formation of air bubbles in the water flowing through the annulus at different inclinations. Despite the limitations of the gap-averaged 2-D numerical model in representing some 3-D effects such as the detailed transverse shape of the air–water interface and the cross-gap velocity/gravity component, very similar flow behaviour in the lateral directions (i.e. $x$ and $y$ directions) was observed in the laboratory experiments and numerical simulations.

In the horizontal case (figure 6), air was conveyed by long bubbles having a tadpole-like shape with a semi-circular cap and a highly stretched tail connected via an apparent neck. These bubbles, separated by liquid slugs, were forced to flow at the top of the annulus due to buoyancy. The tail of a bubble was initially connected to the inlet (figure 6 a) and was stretched as the bubble advanced until it finally disconnected from the inlet (figure 6 b). The bubble then contracted as a result of surface tension (figure 6 c) and continued to move with an essentially steady shape through the co-current flowing liquid towards the outlet (figure 6 d). When the bubble was separated from the inlet, another bubble began to form, resulting in the appearance of a sequence of similarly shaped/sized bubbles moving along the top of the annulus.

Figure 6. The formation and evolution of bubbles in the horizontal annulus (top view): (a) generation of the bubble near the inlet (in annulus subsection 1), (b) disconnection of the bubble (in annulus subsection 1), (c) stabilisation of the bubble after slight contraction (in annulus subsection 2) and (d) translation of the steady bubble to the outlet (in annulus subsection 3). Each figure panel includes the experimental observation, the simulation result and the adaptive unstructured mesh used in the simulation.

In the $1.9^{\circ }$ inclined case (figure 7), bubbles still exhibited a tadpole-like shape but with a much shorter tail, whilst the cap shape remained similar to that in the horizontal case. The bubble detached quickly from the inlet after its formation (figure 7 a). Due to the constant inflow of air, another bubble formed, disconnected and migrated following the leading one. In this way, a stream of distinct bubbles developed as in the horizontal case. However, if the spacing between two bubbles was sufficiently small, the trailing bubble might catch up (figure 7 b) and merge with its upstream neighbour (figure 7 c). The larger bubble formed might further coalesce with another bubble upstream or downstream. Finally, the bubbles that formed relatively stable shapes were transported to the outlet by the flowing water (figure 7 d).

Figure 7. The formation and evolution of bubbles in the $1.9^{\circ }$ inclined annulus (top view): (a) generation of bubbles near the inlet (in annulus subsection 1), (b) interaction of bubbles, with the trailing bubble catching up with the leading bubble (in annulus subsection 2), (c) coalescence of the trailing and leading bubbles (in annulus subsection 2) and (d) translation of steady bubbles to the outlet (in annulus subsection 3). Each figure panel includes the experimental observation, the simulation result and the adaptive unstructured mesh used in the simulation.

There appears to be three stages in the coalescence process. (i) Interacting stage (0.0–0.4 s in figure 8). The trailing bubble with a larger velocity caught up with the leading bubble. As the two bubbles became closer, the trailing bubble experienced a marked change in shape, with its cap becoming more pointed with an increased local curvature at the front. The geometry of the leading bubble only changed slightly with the tail becoming shorter and flatter. The changing shape of the two bubbles was caused by the increased pressure gradient between the two bubble interfaces as they got closer together, during which time the liquid between the two bubbles was gradually pushed out. (ii) Merging stage (0.5–0.6 s in figure 8). The cap of the trailing bubble invaded the rear of the leading bubble. The liquid film between the two bubbles was ruptured and drained as the two bubbles united. The merged bubble exhibited a complex, tortuous shape with the cap and tail connected by a narrow neck and a wider body. (iii) Stabilising stage (0.7–1.0 s in figure 8). The merged larger bubble underwent deformation and adjusted itself from a complex, unsteady shape to the stable tadpole-like form.

Figure 8. Simulation results showing the coalescence process for two bubbles in the $1.9^{\circ }$ inclined annulus.

Figure 9. The formation and evolution of bubbles in the $4.6^{\circ }$ inclined annulus (top view): (a) generation and interaction of bubbles close to the inlet (in annulus subsection 1), (b) coalescence of the trailing, middle and leading bubbles (in annulus subsection 1), (c) stabilisation of merged larger bubbles (in annulus subsection 2) and (d) translation of steady bubbles to the outlet (in annulus subsection 3). Each figure panel includes the experimental observation, the simulation result and the adaptive unstructured mesh used in the simulation.

In the $4.6^{\circ }$ inclined case (figure 9), bubbles became more compact and almost lost the tail structure although the cap shape remained similar to that in the horizontal case. Bubbles separated from the inlet quickly and again formed a train (figure 9 a). The bubbles could merge at an earlier stage by virtue of their smaller spacing and simultaneous coalescence of three bubbles could occur (figure 9 a,b). During the coalescence of three bubbles, the leading bubble remained a semi-circular cap, whereas its rear was flattened; the middle bubble was between the leading and trailing bubbles and showed a sharpened front and a flattened rear; the trailing bubble displayed a narrowed front invading into the back of middle bubble. As merging occurred, the liquid films between bubbles broke and dissipated, such that a united larger bubble was formed with a highly tortuous, unsteady shape. The coalesced larger bubble gradually stabilised itself (figure 9 c) and then migrated to the exit (figure 9 d).

The statistics of the normalised length, cap diameter and terminal velocity of stabilised bubbles in the cases of different inclinations are compared in figure 10. Data from ten fully developed bubbles were collected in the region 2–5 m downstream of the inlet, where the bubbles were not influenced by either inlet or outlet effects. As expected from the discussion above, the bubble length reduces significantly as the inclination of the annulus increases (figure 10 a), due to the deterioration of tail structure. However, the cap diameter, measured based on the curvature of the cap front, only shows a slight increase with the inclination of the annulus (figure 10 b). The terminal velocity of the bubbles also increases slightly with inclination (figure 10 c). An excellent quantitative match with the experimental measurements was achieved by the numerical model. The pseudo-periodic nature of the slug flow and variation of bubble properties as observed in the experiment was also well captured in the simulation (figure 10). Although the inlet condition was stationary, the slug flow pattern varied constantly.

Figure 10. Variation of bubble properties with the inclination angle $\unicode[STIX]{x1D6FD}$ : (a) normalised bubble length $l/D$ , (b) normalised bubble cap diameter $d/D$ and (c) normalised bubble terminal velocity $u/v_{i}$ .

4.2 Effects of gap thickness and pipe diameter

We used the numerical model to further elucidate the effects of gap thickness and pipe diameter on bubbles in thin annuli. First, we set the pipe diameter to be $D=0.137$  m, the air-to-water inflow rate ratio to be $q_{a}/q_{w}=1/4$ and the superficial inflow velocity to be $v_{i}=0.162~\text{m}~\text{s}^{-1}$ . We varied the gap thickness $h$ from 0.001 to 0.005 m, for which the flow remained to be laminar (Beavers et al. Reference Beavers, Sparrow and Magnuson1970) and cross-gap gravity effect was not significant. As shown in figure 11(a), when $h=0.001$  m (i.e. $h/D=0.0073$ ), the bubbles in the horizontal and inclined annuli are very similar, having a slender oval shape with a semi-circular front and a gradually narrowed rear (no discernible neck existed). In the horizontal case, as $h$ becomes larger, the tadpole-like shape becomes more pronounced with a much longer tail and a more distinct neck where the tail joins the cap. If the annulus is inclined, as $h$ increases, the tail tends to vanish and the bubble length decreases significantly. In contrast to the marked variation of bubble length with $h$ , the cap diameter and terminal velocity are much less sensitive to $h$ (figure 11 b,c).

Figure 11. Variation of bubble properties with the gap thickness $h$ while the pipe diameter $D$ is fixed to be $0.137$ m: (a) normalised bubble length $l/D$ , (b) normalised bubble cap diameter $d/D$ and (c) normalised bubble terminal velocity $u/v_{i}$ .

Next, we explored the effect of changing pipe diameter $D=0.05$ –0.4 m (the pipe length $L$ also varied with $D$ according to a fixed ratio $L/D=43.8$ ), keeping $h=0.0035$  m and the same inflow conditions as before. Note that, in figure 12, the bubble length and cap diameter in each case are normalised by the corresponding pipe diameter. As shown in figure 12(a), the ratio of bubble length to the pipe diameter generally does not vary significantly with the change of $D$ . However, in the horizontal annulus with $D=0.4$  m (i.e. $h/D=0.0087$ ), the bubble has a very long, extremely stretched, wavy tail, which does not disconnect from the inlet; the bubble length, therefore, corresponds to the whole pipe length. With the decrease of $D$ (i.e. increase of $h/D$ ), the cap diameter relative to the pipe diameter increases while the bubble velocity decreases (figure 12 b,c).

Figure 12. Variation of bubble properties with the pipe diameter $D$ while the gap thickness $h$ is fixed to be 0.0035 m: (a) normalised bubble length $l/D$ , (b) normalised bubble cap diameter $d/D$ and (c) normalised bubble terminal velocity $u/v_{i}$ .

5 Flow mechanisms

We use the numerical model to investigate the pressure and velocity fields around the bubbles, which determine the bubble motion, shape and dynamics, and can thus provide clues for interpreting the mechanisms that control these behaviours. Figure 13(a) shows the typical steady bubbles for three different inclinations, and figures 13(b) and 13(c) give their pressure and velocity fields, respectively. To facilitate comparison, the absolute pressure at the foremost point of each bubble is chosen as the reference pressure in each case, so that the pressure contour gives the differential pressure between the absolute pressure and the reference pressure. It can be seen that inside the bubble, the pressure is generally constant due to the very low air density. The pressure in the liquid behind the bubble is higher than that in front of the bubble, driving the bubble to flow towards the outlet. The velocity of air inside the bubble is significantly higher than that of the surrounding liquid.

Figure 13. Simulation results showing (a) typical steady bubbles in annuli with different inclinations, and their surrounding (b) pressure and (c) velocity fields. To facilitate comparison, the absolute pressure at the frontal point of each bubble is chosen as the reference pressure for that case, so that the pressure contour gives the differential pressure between the absolute pressure and the reference pressure.

Figure 14. (a) The unrolled geometry of a bubble in the $x$ $y$ plane, and the views of (b) a transversal cross-section in the $y$ $z$ plane and (c) a longitudinal cross-section in the $x$ $z$ plane.

We explore the mechanisms controlling the bubble shape and motion through a simple analysis. We define two polar coordinate systems: one based on the unrolled 2-D bubble cap geometry in the $x$ $y$ plane (figure 14 a), and the other based on the transverse cross-section of the pipe in the $y$ $z$ plane (figure 14 b). By applying Bernoulli’s equation to the steady flow around the bubble cap, assumed to be irrotational with the surface tension being neglected, the following relationship is obtained:

(5.1) $$\begin{eqnarray}\displaystyle \frac{u_{\unicode[STIX]{x1D713}}^{2}}{2}=g\left[\frac{D}{2}(1-\cos \unicode[STIX]{x1D6FE})\cos \unicode[STIX]{x1D6FD}+\frac{d}{2}(1-\cos \unicode[STIX]{x1D713})\sin \unicode[STIX]{x1D6FD}\right], & & \displaystyle\end{eqnarray}$$

where $u_{\unicode[STIX]{x1D713}}$ is the relative velocity of the surrounding liquid to the bubble at a point with a polar angle of $\unicode[STIX]{x1D713}$ , and $\unicode[STIX]{x1D6FE}$ is the polar angle in the coordinate system of the transverse cross-section. We assume that flow around the bubble cap is similar to that seen around a circular cylinder to a first approximation (Collins Reference Collins1967a ). This flow field can be solved using the stream function and velocity potential, with the following relationship obtained (Wilkes Reference Wilkes2006):

(5.2) $$\begin{eqnarray}\displaystyle \frac{u_{\unicode[STIX]{x1D713}}^{2}}{u_{\infty }^{2}\sin ^{2}\unicode[STIX]{x1D713}}=4, & & \displaystyle\end{eqnarray}$$

where $u_{\infty }$ is the far-field relative velocity of the flow to the bubble. By combining (5.1), (5.2) and the relation $(d/2)\sin \unicode[STIX]{x1D713}=(D/2)\unicode[STIX]{x1D6FE}$ (see figure 14 a,b), we obtain

(5.3) $$\begin{eqnarray}\displaystyle u_{\infty }^{2}=\left(\frac{1-\cos \unicode[STIX]{x1D6FE}}{\unicode[STIX]{x1D6FE}^{2}}\frac{d}{D}\cos \unicode[STIX]{x1D6FD}+\frac{1-\cos \unicode[STIX]{x1D713}}{\sin ^{2}\unicode[STIX]{x1D713}}\sin \unicode[STIX]{x1D6FD}\right)\frac{gd}{4}. & & \displaystyle\end{eqnarray}$$

In the limit of small $\unicode[STIX]{x1D6FE}$ and $\unicode[STIX]{x1D713}$ , equation (5.3) becomes

(5.4) $$\begin{eqnarray}\displaystyle u_{\infty }\approx \frac{1}{2}\sqrt{\left(\frac{d}{D}\cos \unicode[STIX]{x1D6FD}+\sin \unicode[STIX]{x1D6FD}\right)\frac{gd}{2}}. & & \displaystyle\end{eqnarray}$$

If $r\rightarrow \infty$ , $\unicode[STIX]{x1D6FD}=90^{\circ }$ , then solution (5.4) converges to that for plane bubbles rising through a quiescent liquid between vertical, parallel flat plates, i.e. $0.5(gd/2)^{1/2}$ (Davies & Taylor Reference Davies and Taylor1950; Collins Reference Collins1965a ).

For the simultaneous gas–liquid flow here, the bubble terminal velocity $u$ is calculated as (Griffith & Wallis Reference Griffith and Wallis1961)

(5.5) $$\begin{eqnarray}\displaystyle u=u_{\infty }+v_{i}\approx \frac{1}{2}\sqrt{\left(\frac{d}{D}\cos \unicode[STIX]{x1D6FD}+\sin \unicode[STIX]{x1D6FD}\right)\frac{gd}{2}}+v_{i}, & & \displaystyle\end{eqnarray}$$

where $v_{i}$ is the superficial inflow velocity (table 1). We substitute the measured diameters of bubble cap from the experimental/numerical results into (5.5) to derive the predicted bubble terminal velocity $u_{pred}$ . The predicted values are compared to the actual velocity $u_{meas}$ measured in the experiment and simulation. It can be seen that our simplified analysis gives a good estimation of the bubble velocity in thin annuli for different inclination angles, pipe diameters and gap thicknesses (figure 15).

Figure 15. Comparison of the predicted terminal velocity $u_{pred}$ from (5.5) and the measured terminal velocity $u_{meas}$ from experiments and simulations for different configurations of pipe diameter $D$ , gap thickness $h$ and inclination angle $\unicode[STIX]{x1D6FD}$ . The velocity values are normalised by the superficial inflow velocity $v_{i}$ .

We then identify the following dimensionless groups of the annulus system:

(5.6a-c ) $$\begin{eqnarray}\displaystyle \hat{Fr}_{\Vert }^{-1}=\frac{\sqrt{gD\sin \unicode[STIX]{x1D6FD}}}{v_{i}},\quad \hat{Fr}_{\bot }^{-1}=\frac{\sqrt{gD\cos \unicode[STIX]{x1D6FD}}}{v_{i}},\quad \hat{\unicode[STIX]{x1D716}}=\displaystyle \frac{h}{D}, & & \displaystyle\end{eqnarray}$$

where the two inverse Froude numbers, $\hat{Fr}_{\Vert }^{-1}$ and $\hat{Fr}_{\bot }^{-1}$ , are related to the effects of gravity parallel and perpendicular to the flow direction, respectively, and $\hat{\unicode[STIX]{x1D716}}$ characterises the effects of confinement associated with the annulus. We also characterise the bubble shape and motion based on the bubble aspect ratio $l/d$ , the bubble Reynolds number $Re=\unicode[STIX]{x1D70C}_{w}u_{\infty }d/\unicode[STIX]{x1D707}_{w}$ and the gap Reynolds number $Re(h/d)^{2}$ . Here, $u_{\infty }$ is used because the bubble behaviour is governed by the relative motion between the two fluids.

As shown in figure 16(a), as $\hat{Fr}_{\Vert }^{-1}$ increases, the bubble shape becomes less elongated, i.e. $l/d$ approaches unity, since the gravitational component parallel to the flow direction causes an increased pressure at the bubble rear, which tends to separate the bubble and reduce the bubble length. As $\hat{Fr}_{\bot }^{-1}$ increases, $Re$ increases (figure 16 b), indicating inertia becomes more dominant. This is because, for the annuli with small inclinations to the horizontal, the gravitational component perpendicular to the flow direction mainly controls the flow field around the bubble cap, as can be expected from (5.1). As the annulus confinement decreases (i.e. $\hat{\unicode[STIX]{x1D716}}$ increases), the gap Reynolds number $Re(h/d)^{2}$ increases (figure 16 c), due to the reduced viscous effect across the gap, reaching a plateau with increasing $\hat{\unicode[STIX]{x1D716}}$ .

Figure 16. (a) Variation of the bubble aspect ratio $l/d$ with the inverse Froude number $\hat{Fr}_{\Vert }^{-1}$ . (b) Variation of the bubble Reynolds number $Re$ with the inverse Froude number $\hat{Fr}_{\bot }^{-1}$ . (c) Variation of the gap Reynolds number $Re(h/d)^{2}$ with the confinement ratio $\hat{\unicode[STIX]{x1D716}}$ . Refer to figure 15 for the legend.

6 Discussion and conclusions

In this paper, we presented an experimental and numerical investigation into the shape and motion of gas bubbles in a liquid flowing through a horizontal or slightly inclined thin annulus. We have focused on the regime where inertia is important (the gap Reynolds number is close to or larger than unity) and the bubble behaviour is governed by the complex interplay of inertia, gravity, viscosity and surface tension. Bubbles in the horizontal annulus developed a unique ‘tadpole-like’ shape featuring a semi-circular cap and a highly stretched tail. As the annulus was inclined with respect to the horizontal, the length of the bubble decreased. We developed a gap-averaged, 2-D numerical model to represent the 3-D flow dynamics, which achieved a close match to the experimental data for different small inclinations. The numerical model was used to further elucidate the effects of gap thickness and pipe diameter on the bubble evolution in thin annuli. We found that the bubble velocity is strongly correlated to the cap structure, but is independent of the bubble length, as has also been reported for bubbles in tubes (Davies & Taylor Reference Davies and Taylor1950; Zukoski Reference Zukoski1966; Fabre & Liné Reference Fabre and Liné1992; Viana et al. Reference Viana, Pardo, Yánez, Trallero and Joseph2003) and between parallel flat plates (Grace & Harrison Reference Grace and Harrison1967; Maneri & Zuber Reference Maneri and Zuber1974; Hills Reference Hills1975; Couet & Strumolo Reference Couet and Strumolo1987). We reported that the elongated bubble shape in horizontal annuli is due to the buoyancy which causes the bubbles to spread along the top of the annulus. The gravitational component along the flow direction, which increases as the annulus is inclined, impinges the liquid slug and causes a reduction of the bubble tail. The gravitational component perpendicular to the flow direction controls the bubble motion and the cap structure. These mechanisms produce the unique tadpole-like shape with a sharp tail tip because of the cross-sectional curvature of the annulus channel, in contrast to the bubble shape with a tongue-like rear seen in flows between parallel flat plates (Maneri & Zuber Reference Maneri and Zuber1974; Couet & Strumolo Reference Couet and Strumolo1987).

It is remarkable that the 2-D numerical model well captured the bubble evolution and interaction behaviour as observed in the experimental annulus. It is still worth mentioning that some complex 3-D effects were not fully represented in the gap-averaged formulation, such as the detailed transverse shape of the bubble cap and tail within the gap as well as the influence of the cross-gap velocity/gravity component (Oliveira & Meiburg Reference Oliveira and Meiburg2011), which require direct three-dimensional but computationally very expensive numerical simulations. Such limitations of the 2-D model may explain the discrepancy of the bubble geometry, e.g. in the region near the bubble neck (figure 6), between our experimental and numerical results. The gap-averaged 2-D model assumed an absence of liquid films between the bubble and channel walls. We examined different mathematical formulations that account for the presence of liquid films (Pitts Reference Pitts1980; McLean & Saffman Reference McLean and Saffman1981; Park & Homsy Reference Park and Homsy1984; Reinelt Reference Reinelt1987a ; Kopf-Sill & Homsy Reference Kopf-Sill and Homsy1988), which however resulted in a poor prediction of the experimentally observed bubble patterns (the simulated bubbles were much smaller, probably due to an overestimation of the curvature of the air–water interface across the gap). Saffman & Tanveer (Reference Saffman and Tanveer1989) also reported that the thin film hypothesis, compared to the contact angle assumption, gave a less consistent prediction of bubble shapes in Hele-Shaw cells. Furthermore, we did not see any evidence of such films in the experiments. The thickness of thin films (if present under no-gravity conditions) is estimated to be $\unicode[STIX]{x1D6FF}/h\approx 0.67(\unicode[STIX]{x1D707}_{w}u_{\infty }/\unicode[STIX]{x1D70E})^{2/3}$ (Park & Homsy Reference Park and Homsy1984; Reinelt Reference Reinelt1987a ; Klaseboer, Gupta & Manica Reference Klaseboer, Gupta and Manica2014), which may however be reduced further in the top region of the annulus under gravity-induced drainage downwards along the peripheral direction (Tso & Sugawara Reference Tso and Sugawara1990; Paras & Karabelas Reference Paras and Karabelas1991). Thus, we expect that $\unicode[STIX]{x1D6FF}$ may be close to the scale of the surface roughness of the Plexiglas and stainless steel walls, i.e. around 0.01 mm (Li et al. Reference Li, Su, Kuhn, Meusel, Ammann, Shao, Pöschl and Cheng2018). We therefore do not expect the films to persist in a stable way. The observed contact angle hysteresis in our study may also be attributed to the presence of pronounced surface roughness effects (Dussan Reference Dussan1979). The phenomenon that air wets the channel wall has also been found between the air bubble and the upper wall of a horizontal/slightly inclined tube (Fabre & Liné Reference Fabre and Liné1992) or Hele-Shaw cell (Maneri & Zuber Reference Maneri and Zuber1974). We suggest that a good approximation of the gap-averaged 2-D model to 3-D flow holds if the cross-gap effects are minor.

Acknowledgements

We would like to acknowledge Statoil ASA for funding Q.L., M.D.J., C.C.P., O.K.M. and A.H.M., and for granting permission to publish this work. We thank K. Årland from Statoil ASA for project management. O.K.M. and C.C.P. acknowledge the funding provided by the Engineering and Physical Sciences Research Council (EPSRC) through the Programme Grant MEMPHIS (grant EP/K003976/1). D.P. is grateful to the EPSRC (grant EP/M012794/1). P.S. thanks ExxonMobil and the EPSRC (grant EP/R005761/1) for funding. No data were generated in the course of this work. For further information, please contact the corresponding author Q.L. (email: q.lei12@imperial.ac.uk or qinghua.lei@erdw.ethz.ch), the AMCG Group (http://www.imperial.ac.uk/earth-science/research/research-groups/amcg/), the NORMS group (http://www.imperial.ac.uk/earth-science/research/research-groups/norms) or the Matar Fluids Group (http://www.imperial.ac.uk/matar-fluids-group).

References

Antonini, C., Carmona, F. J., Pierce, E., Marengo, M. & Amirfazli, A. 2009 General methodology for evaluating the adhesion force of drops and bubbles on solid surfaces. Langmuir 25 (11), 61436154.Google Scholar
Beavers, G. S., Sparrow, E. M. & Magnuson, R. A. 1970 Experiments on the breakdown of laminar flow in a parallel-plate channel. Intl J. Heat Mass Transfer 13, 809815.Google Scholar
Brackbill, J. U., Kothe, D. B. & Zemach, C. 1992 A continuum method for modeling surface tension. J. Comput. Phys. 100 (2), 335354.Google Scholar
Bush, J. W. M. 1997 The anomalous wake accompanying bubbles rising in a thin gap: a mechanically forced Marangoni flow. J. Fluid Mech. 352, 283303.Google Scholar
Bush, J. W. M. & Eames, I. 1998 Fluid displacement by high Reynolds number bubble motion in a thin gap. Intl J. Multiphase Flow 24 (3), 411430.Google Scholar
Collins, R. 1965a A simple model of the plane gas bubble in a finite liquid. J. Fluid Mech. 22, 763771.Google Scholar
Collins, R. 1965b Structure and behaviour of wakes behind two-dimensional air bubbles in water. Chem. Engng Sci. 20, 851853.Google Scholar
Collins, R. 1967a The cycloidal-cap bubble: a neglected solution in the theory of large plane gas bubbles in liquids. Chem. Engng Sci. 22, 8997.Google Scholar
Collins, R. 1967b The effect of a containing cylindrical boundary on the velocity of a large gas bubble in a liquid. J. Fluid Mech. 28 (1), 97112.Google Scholar
Collins, R., de Moraes, F. F., Davidson, J. F. & Harrison, D. 1978 The motion of a large gas bubble rising through liquid flowing in a tube. J. Fluid Mech. 89, 497514.Google Scholar
Cotter, C. J., Ham, D. A., Pain, C. C. & Reich, S. 2009 LBB stability of a mixed Galerkin finite element pair for fluid flow simulations. J. Comput. Phys. 228 (2), 336348.Google Scholar
Couet, B. & Strumolo, G. S. 1987 The effects of surface tension and tube inclination on a two-dimensional rising bubble. J. Fluid Mech. 184, 114.Google Scholar
Cueto-Felgueroso, L. & Juanes, R. 2014 A phase-field model of two-phase Hele-Shaw flow. J. Fluid Mech. 758, 522552.Google Scholar
Das, G., Das, P. K., Purohit, N. K. & Mitra, A. K. 1998 Rise velocity of a Taylor bubble through concentric annulus. Chem. Engng Sci. 53 (5), 977993.Google Scholar
Davies, R. M. & Taylor, G. I. 1950 The mechanics of large bubbles rising through extended liquids and through liquids in tubes. Proc. R. Soc. Lond. A 200 (1062), 375390.Google Scholar
Dussan, E. B. 1979 On the spreading of liquids on solid surfaces: static and dynamic contact lines. Ann. Rev. Fluid Mech. 11, 371400.Google Scholar
Eck, W. & Siekmann, J. 1978 On bubble motion in a Hele-Shaw cell, a possibility to study two-phase flows under reduced gravity. Ing.-Arch. 47, 153168.Google Scholar
Ekberg, N. P., Ghiaasiaan, S. M., Abdel-Khalik, S. I., Yoda, M. & Jeter, S. M. 1999 Gas–liquid two-phase flow in narrow horizontal annuli. Nucl. Engng Des. 192, 5980.Google Scholar
Fabre, J. 2016 A long bubble rising in still liquid in a vertical channel: a plane inviscid solution. J. Fluid Mech. 797, R4.Google Scholar
Fabre, J. & Liné, A. 1992 Modeling of two-phase slug flow. Annu. Rev. Fluid Mech. 24, 2146.Google Scholar
Figueroa-Espinoza, B. & Fabre, J. 2011 Taylor bubble moving in a flowing liquid in vertical channel: transition from symmetric to asymmetric shape. J. Fluid Mech. 679, 432454.Google Scholar
Filella, A., Ern, P. & Roig, V. 2015 Oscillatory motion and wake of a bubble rising in a thin-gap cell. J. Fluid Mech. 778, 6088.Google Scholar
Gondret, P. & Rabaud, M. 1997 Shear instability of two-fluid parallel flow in a Hele-Shaw cell. Phys. Fluids 9 (11), 32673274.Google Scholar
Grace, J. R. & Harrison, D. 1967 The influence of bubble shape on the rising velocities of large bubbles. Chem. Engng Sci. 22, 13371347.Google Scholar
Griffith, P. & Wallis, G. B. 1961 Two-phase slug flow. Trans. ASME J. Heat Transfer 83, 307318.Google Scholar
Hills, J. H. 1975 The two-dimensional elliptical cap bubble. J. Fluid Mech. 68, 503512.Google Scholar
Huisman, S. G., Ern, P. & Roig, V. 2012 Interactions and coalescence of large bubbles rising in a thin gap. Phys. Rev. E 85, 027302.Google Scholar
Kelessidis, V. C. & Dukler, A. E. 1990 Motion of large gas bubbles through liquids in vertical concentric and eccentric annuli. Intl J. Multiphase Flow 16 (3), 375390.Google Scholar
Kelley, E. & Wu, M. 1997 Path instability of rising air bubbles in a Hele-Shaw cell. Phys. Rev. Lett. 79 (7), 12651268.Google Scholar
Klaseboer, E., Gupta, R. & Manica, R. 2014 An extended Bretherton model for long Taylor bubbles at moderate capillary numbers. Phys. Fluids 26, 032107.Google Scholar
Kopf-Sill, A. R. & Homsy, G. M. 1988 Bubble motion in a Hele-Shaw cell. Phys. Fluids 31, 1826.Google Scholar
Lazarek, G. M. & Littman, H. 1974 The pressure field due to a large circular capped air bubble rising in water. J. Fluid Mech. 66, 673687.Google Scholar
Li, G., Su, H., Kuhn, U., Meusel, H., Ammann, M., Shao, M., Pöschl, U. & Cheng, Y. 2018 Technical note: influence of surface roughness and local turbulence on coated-wall flow tube experiments for gas uptake and kinetic studies. Atmos. Chem. Phys. 18 (4), 26692686.Google Scholar
Maneri, C. C. & Zuber, N. 1974 An experimental study of plane bubbles rising at inclination. Intl J. Multiphase Flow 1, 623645.Google Scholar
Maruvada, S. R. K. & Park, C. W. 1996 Retarded motion of bubbles in Hele-Shaw cells. Phys. Fluids 8 (12), 32293233.Google Scholar
Maxworthy, T. 1986 Bubble formation, motion and interaction in a Hele-Shaw cell. J. Fluid Mech. 173, 95114.Google Scholar
McLean, J. W. & Saffman, P. G. 1981 The effect of surface tension on the shape of fingers in a Hele Shaw cell. J. Fluid Mech. 102, 455469.Google Scholar
Meiburg, E. 1989 Bubbles in a Hele-Shaw cell: numerical simulation of three-dimensional effects. Phys. Fluids A 1 (6), 938946.Google Scholar
Oliveira, R. M. & Meiburg, E. 2011 Miscible displacements in Hele-Shaw cells: three-dimensional Navier–Stokes simulations. J. Fluid Mech. 687, 431460.Google Scholar
Pain, C. C., Umpleby, A. P., de Oliveira, C. R. E. & Goddard, A. J. H. 2001 Tetrahedral mesh optimisation and adaptivity for steady-state and transient finite element calculations. Comput. Meth. Appl. Mech. Engng 190 (29–30), 37713796.Google Scholar
Paras, S. V. & Karabelas, A. J. 1991 Properties of the liquid layer in horizontal annular flow. Intl J. Multiphase Flow 17 (4), 439454.Google Scholar
Park, C.-W. & Homsy, G. M. 1984 Two-phase displacement in Hele Shaw cells: theory. J. Fluid Mech. 139, 291308.Google Scholar
Pavlidis, D., Gomes, J. L. M. A., Xie, Z., Percival, J. R., Pain, C. C. & Matar, O. K. 2016 Compressive advection and multi-component methods for interface-capturing. Intl J. Numer. Meth. Fluids 80 (4), 256282.Google Scholar
Pavlidis, D., Xie, Z., Percival, J. R., Gomes, J. L. M. A., Pain, C. C. & Matar, O. K. 2014 Two- and three-phase horizontal slug flow simulations using an interface-capturing compositional approach. Intl J. Multiphase Flow 67, 8591.Google Scholar
Piedra, S., Ramos, E. & Herrera, J. R. 2015 Dynamics of two-dimensional bubbles. Phys. Rev. E 91 (6), 063013.Google Scholar
Pitts, E. 1980 Penetration of fluid into a Hele-Shaw cell: the Saffman–Taylor experiment. J. Fluid Mech. 97, 5364.Google Scholar
Reinelt, D. A. 1987a Interface conditions for two-phase displacement in Hele-Shaw cells. J. Fluid Mech. 183, 219234.Google Scholar
Reinelt, D. A. 1987b The rate at which a long bubble rises in a vertical tube. J. Fluid Mech. 175, 557565.Google Scholar
Roig, V., Roudet, M., Risso, F. & Billet, A.-M. 2012 Dynamics of a high-Reynolds-number bubble rising within a thin gap. J. Fluid Mech. 707, 444466.Google Scholar
Saffman, P. G. & Tanveer, S. 1989 Prediction of bubble velocity in a Hele-Shaw cell: thin film and contact angle effects. Phys. Fluids A 1 (2), 219223.Google Scholar
Saffman, P. G. & Taylor, G. I. 1958 The penetration of a fluid into a porous medium or Hele-Shaw cell containing a more viscous liquid. Proc. R. Soc. Lond. A 245 (1242), 312329.Google Scholar
Tanveer, S. 1986 The effect of surface tension on the shape of a Hele-Shaw cell bubble. Phys. Fluids 29 (11), 35373548.Google Scholar
Tanveer, S. 1987 New solutions for steady bubbles in a Hele-Shaw cell. Phys. Fluids 30 (3), 651658.Google Scholar
Tanveer, S. & Saffman, P. G. 1987 Stability of bubbles in a Hele-Shaw cell. Phys. Fluids 30 (9), 26242635.Google Scholar
Taylor, G. I. & Saffman, P. G. 1959 A note on the motion of bubbles in a Hele-Shaw cell and porous medium. Q. J. Mech. Appl. Math. 12, 265279.Google Scholar
Thompson, A. B., Juel, A. & Hazel, A. L. 2014 Multiple finger propagation modes in Hele-Shaw channels of variable depth. J. Fluid Mech. 746, 123164.Google Scholar
Thompson, B. W. 1968 Secondary flow in a Hele-Shaw cell. J. Fluid Mech. 372, 2544.Google Scholar
Tso, C. P. & Sugawara, S. 1990 Film thickness prediction in a horizontal annular two-phase flow. Intl J. Multiphase Flow 16 (5), 867884.Google Scholar
Vanden-Broeck, J.-M. 1984 Rising bubbles in a two-dimensional tube with surface tension. Phys. Fluids 27 (11), 26042607.Google Scholar
Vanden-Broeck, J.-M. 1992 Rising bubbles in a two-dimensional tube: asymptotic behavior for small values of the surface tension. Phys. Fluids A 4 (11), 23322334.Google Scholar
Viana, F., Pardo, R., Yánez, R., Trallero, J. L. & Joseph, D. D. 2003 Universal correlation for the rise velocity of long gas bubbles in round pipes. J. Fluid Mech. 494, 379398.Google Scholar
Wang, X., Klaasen, B., Degrève, J., Blanpain, B. & Verhaeghe, F. 2014 Experimental and numerical study of buoyancy-driven single bubble dynamics in a vertical Hele-Shaw cell. Phys. Fluids 26 (12), 123303.Google Scholar
Wang, X., Klaasen, B., Degrève, J., Mahulkar, A., Heynderickx, G., Reyniers, M.-F., Blanpain, B. & Verhaeghe, F. 2016 Volume-of-fluid simulations of bubble dynamics in a vertical Hele-Shaw cell. Phys. Fluids 28, 053304.Google Scholar
White, E. T. & Beardmore, R. H. 1962 The velocity of rise of single cylindrical air bubbles through liquids contained in vertical tubes. Chem. Engng Sci. 17 (5), 351361.Google Scholar
Wilkes, J. O. 2006 Fluid Mechanics for Chemical Engineers with Microfluidics and CFD, 2nd edn. Prentice Hall.Google Scholar
Xie, Z., Hewitt, G. F., Pavlidis, D., Salinas, P., Pain, C. C. & Matar, O. K. 2017 Numerical study of three-dimensional droplet impact on a flowing liquid film in annular two-phase flow. Chem. Engng Sci. 166, 303312.Google Scholar
Xie, Z., Pavlidis, D., Percival, J. R., Gomes, J. L. M. A., Pain, C. C. & Matar, O. K. 2014 Adaptive unstructured mesh modelling of multiphase flows. Intl J. Multiphase Flow 67, 104110.Google Scholar
Xie, Z., Pavlidis, D., Salinas, P., Percival, J. R., Pain, C. C. & Matar, O. K. 2016 A balanced-force control volume finite element method for interfacial flows with surface tension using adaptive anisotropic unstructured meshes. Comput. Fluids 138, 3855.Google Scholar
Zukoski, E. E. 1966 Influence of viscosity, surface tension, and inclination angle on motion of long bubbles in closed tubes. J. Fluid Mech. 25, 821837.Google Scholar
Figure 0

Figure 1. Schematic illustration of the experimental apparatus.

Figure 1

Figure 2. Schematic of (a) the thin annulus and (b) the equivalent Hele-Shaw cell. (c) An illustration showing how the direction of gravitational acceleration changes along the equivalent Hele-Shaw cell to ensure the flow replicates that seen in the actual annulus.

Figure 2

Table 1. Parameters of the experimental apparatus and fluids.

Figure 3

Figure 3. Schematic of (a) the cross-gap meniscus profile governed by the local contact angles, $\unicode[STIX]{x1D703}_{o}$ and $\unicode[STIX]{x1D703}_{i}$, at the outer and inner walls, respectively, and (b) the moving contact line in the $x$$y$ plane with the local advancing/receding state determined by the intersection angle $\unicode[STIX]{x1D719}$ between the velocity vector $\bar{\boldsymbol{u}}$ and the interface normal $\boldsymbol{n}$.

Figure 4

Figure 4. Interpolation of contact angles around a contact line at the outer or inner wall of the annulus using a third-order polynomial fitting (Antonini et al.2009); $\unicode[STIX]{x1D703}$ denotes the local contact angle and $\unicode[STIX]{x1D719}$ denotes the intersection angle between the local velocity vector and the interface normal as shown in figure 3.

Figure 5

Figure 5. The 2-D model set-up for numerical simulation.

Figure 6

Figure 6. The formation and evolution of bubbles in the horizontal annulus (top view): (a) generation of the bubble near the inlet (in annulus subsection 1), (b) disconnection of the bubble (in annulus subsection 1), (c) stabilisation of the bubble after slight contraction (in annulus subsection 2) and (d) translation of the steady bubble to the outlet (in annulus subsection 3). Each figure panel includes the experimental observation, the simulation result and the adaptive unstructured mesh used in the simulation.

Figure 7

Figure 7. The formation and evolution of bubbles in the $1.9^{\circ }$ inclined annulus (top view): (a) generation of bubbles near the inlet (in annulus subsection 1), (b) interaction of bubbles, with the trailing bubble catching up with the leading bubble (in annulus subsection 2), (c) coalescence of the trailing and leading bubbles (in annulus subsection 2) and (d) translation of steady bubbles to the outlet (in annulus subsection 3). Each figure panel includes the experimental observation, the simulation result and the adaptive unstructured mesh used in the simulation.

Figure 8

Figure 8. Simulation results showing the coalescence process for two bubbles in the $1.9^{\circ }$ inclined annulus.

Figure 9

Figure 9. The formation and evolution of bubbles in the $4.6^{\circ }$ inclined annulus (top view): (a) generation and interaction of bubbles close to the inlet (in annulus subsection 1), (b) coalescence of the trailing, middle and leading bubbles (in annulus subsection 1), (c) stabilisation of merged larger bubbles (in annulus subsection 2) and (d) translation of steady bubbles to the outlet (in annulus subsection 3). Each figure panel includes the experimental observation, the simulation result and the adaptive unstructured mesh used in the simulation.

Figure 10

Figure 10. Variation of bubble properties with the inclination angle $\unicode[STIX]{x1D6FD}$: (a) normalised bubble length $l/D$, (b) normalised bubble cap diameter $d/D$ and (c) normalised bubble terminal velocity $u/v_{i}$.

Figure 11

Figure 11. Variation of bubble properties with the gap thickness $h$ while the pipe diameter $D$ is fixed to be $0.137$ m: (a) normalised bubble length $l/D$, (b) normalised bubble cap diameter $d/D$ and (c) normalised bubble terminal velocity $u/v_{i}$.

Figure 12

Figure 12. Variation of bubble properties with the pipe diameter $D$ while the gap thickness $h$ is fixed to be 0.0035 m: (a) normalised bubble length $l/D$, (b) normalised bubble cap diameter $d/D$ and (c) normalised bubble terminal velocity $u/v_{i}$.

Figure 13

Figure 13. Simulation results showing (a) typical steady bubbles in annuli with different inclinations, and their surrounding (b) pressure and (c) velocity fields. To facilitate comparison, the absolute pressure at the frontal point of each bubble is chosen as the reference pressure for that case, so that the pressure contour gives the differential pressure between the absolute pressure and the reference pressure.

Figure 14

Figure 14. (a) The unrolled geometry of a bubble in the $x$$y$ plane, and the views of (b) a transversal cross-section in the $y$$z$ plane and (c) a longitudinal cross-section in the $x$$z$ plane.

Figure 15

Figure 15. Comparison of the predicted terminal velocity $u_{pred}$ from (5.5) and the measured terminal velocity $u_{meas}$ from experiments and simulations for different configurations of pipe diameter $D$, gap thickness $h$ and inclination angle $\unicode[STIX]{x1D6FD}$. The velocity values are normalised by the superficial inflow velocity $v_{i}$.

Figure 16

Figure 16. (a) Variation of the bubble aspect ratio $l/d$ with the inverse Froude number $\hat{Fr}_{\Vert }^{-1}$. (b) Variation of the bubble Reynolds number $Re$ with the inverse Froude number $\hat{Fr}_{\bot }^{-1}$. (c) Variation of the gap Reynolds number $Re(h/d)^{2}$ with the confinement ratio $\hat{\unicode[STIX]{x1D716}}$. Refer to figure 15 for the legend.