The Effect of Liquid Viscosity on the Rise Velocity of Taylor Bubbles in Small Diameter Bubble Column

The rise velocity of Taylor bubbles in small diameter bubble column was measured via cross-correlation between two planes of time-averaged void fraction data obtained from the electrical capacitance tomography (ECT). This was subsequently compared with the rise velocity obtained from the high-speed camera, manual time series analysis and likewise empirical models. The inertia, viscous and gravitational forces were identified as forces, which could influence the rise velocity. Fluid flow analysis was carried out using slug Reynolds number, Froude number and inverse dimensionless viscosity, which are important dimensionless parameters influencing the rise velocity of Taylor bubbles in different liquid viscosities, with the parameters being functions of the fluid properties and column diameter. It was found that the Froude number decreases with an increase in viscosity with a variation in flow as superficial gas velocity increases with reduction in rise velocity. A dominant effect of viscous and gravitational forces over inertia forces was obtained, which showed an agreement with Stokes law, where drag force is directly proportional to viscosity. Hence, the drag force increases as viscosity increases (5 < 100 < 1000 < 5000 mPa s), leading to a decrease in the rise velocity of Taylor bubbles. It was concluded that the rise velocity of Taylor bubbles decreases with an increase in liquid viscosity and, on the other hand, increases with an increase in superficial gas velocity.


Introduction
Slug flow is characterized by Taylor bubbles, which has large pockets of bullet shaped bubbles occupying almost the entire cross-section of the column. The Taylor bubble is surrounded by a thin film of liquid, and below, it is the liquid slugs, which are agglomerate of small bubbles. Zukoski [1], Tomiyama et al. [2] and Mandal and Das [3] described the Taylor bubble length to exceed 1.5 times of the tube diameter or its diameter is greater than 60% of the tube diameter. The rise velocity of a single isolated Taylor bubble is dependent on inertia and drag forces [4]. A number of parameters affect the rise velocity of Taylor bubbles through a stagnant liquid; such parameters include density of liquid, surface tension of liquid, liquid viscosity, acceleration due to gravity, diameter of bubbles etc. [5]. Mao and Dukler [6] explained that in a situation whereby the liquid is flowing, the rise velocity of a Taylor bubble must depend on the velocity of the liquid flowing upstream as well as the rise due to buoyancy. A typical example of bubble rising through a stagnant liquid as taken from a high-speed video camera (from the current study) is shown in Figure 1.

Background
The rise velocity of Taylor bubbles is otherwise known as structure velocity. It can also be defined as the velocity of periodic structures in the slug [7,8]. Some researchers have carried out studies on the rise velocity of Taylor bubbles through stagnant liquid. Mao and Dukler [9] from their experimental results explained that for a wide range of viscosity and surface tension, the rise velocity can be expressed in terms of a constant Froude number.
Hence, the rise velocity of Taylor bubble is given as: where g is the acceleration due to gravity and D is the diameter of the tube. The rise velocity of Taylor bubbles in stagnant liquids was first studied by Dumitrescu [10] and Davies and Taylor [11] in which observed bubbles were of characteristics shape referred to as Dumitrescu or Taylor bubbles. Griffith and Wallis [12] eventually proposed the name as Taylor bubbles. Dumitrescu [10] carried out a study on the rise velocity of bubbles using water in a vertical tube and it was established both theoretically and experimentally that the bubble velocity was: This is synonymous to the proposition of Mao and Dukler [9], where Froude's number is given as 0.35.
Dumitrescu [10] assumed that the bubble would have a spherical nose, solving simultaneously the flow around the bubble and the asymptotic film which eventually led to the bubble velocity [13] given as Eq. (3). Davies and Taylor [11] gave the bubble velocity as: after solving the problem using different assumptions. Nicklin et al. [7] later postulated that the Davies and Taylor [11] solution was not unique but should tend to the limiting value given as: Eqs. (3) and (4) proposed by Dumitrescu [10] and Davies and Taylor [11], respectively, assume that the Taylor bubble was obtained from a gas of zero density [13]. Neal [14] proposed that if the bubble density is significant, the bubble velocity is given as: where c is approximately 0.35, Δρ ¼ ρ L À ρ G , ρ L and ρ G are the liquid and gas densities respectively.
Brown [15] from his experimental studies found that the solutions of Dumitrescu [10] and Davies and Taylor [11] were not suitable for high viscosity liquids, that they only describe the behaviour of gas bubbles in low viscosity liquids [13]. So, Brown [15] gave the bubble velocity as: (8) and (9) where N LB is the liquid viscosity number and μ L is the liquid viscosity. Zukoski [1] proposed an expression for velocity of large bubbles in a closed horizontal pipe with large diameter (neglecting surface tension effects) given as: where the Froude number is 0.54. A correlation for the bubble rise velocity was proposed by Griffith and Wallis [12] based on the studies on vertical slug flow given as: where U SG and U SL are the superficial gas and liquid velocities respectively, and K 1 = 0.35.
They investigated the effect of different velocity profiles in the liquid slug by varying K 2 [16].
Nicklin et al. [7] from their vertical slug experiments proposed the rise velocity of a Taylor bubble in the liquid in a vertical tube as: (12) where U o is the translational velocity in a stagnant liquid or velocity of bubble propagating into stagnant liquid, given as: where C o is the distribution coefficient which is close to 1.2 for fully developed turbulent flow (low viscous liquid e.g. water) and close to 2 for laminar flow (high viscous liquid e.g. 1000 mPa s silicone oil) [7,[17][18][19].
From the experimental work reported by Sylvester [20], the rise velocity of Taylor bubble was presented as: where D is the pipe diameter. They proposed C o and C 1 to be 1.2 and 0.35, respectively.
From the experiment carried out by Bendiksen [21] in a vertical tube with flowing liquid, the distribution coefficient, C o was obtained to be 1.2 for Reynold's number in the range 5000-110,000, i.e. low viscous liquid.
Nicklin et al. [7] interpreted their equation as: Rise velocity of Taylor bubble (structure velocity) is equal to the velocity of the liquid at the tip of the bubble nose plus the rise velocity of bubble in a stagnant liquid (translational velocity).
Mao and Dukler [6] called this effective upstream velocity (rise velocity of Taylor bubble), the centreline velocity of the liquid.
White and Beadmore [22] carried out an experimental investigation on the rise velocity of Taylor bubbles through liquids in a vertical tube using three dimensionless parameters: Froude number, Fr, Eotvos number, E o and Morton number, M o . A recent review on vertical gas-liquid slug flow which highlights previous works on the rise velocity of Taylor bubbles is provided by Morgado et al. [23]. They discussed experimental, theoretical and numerical methods of investigating the rise velocity of Taylor bubbles, where the so-called numerical methods involve the use of empirical correlations. The limitations of these studies that have been addressed in the current study include (1) low range of liquid viscosities and more emphasis on low viscosities rather than high viscosities, (2) consideration of column or pipe diameter greater than 50 mm, (3) limited exploration of the effect of forces such as surface tension, inertia, gravitational and viscous forces acting on Taylor bubble, (4) limited exploration of the relationship between fluid dimensionless parameters and the Taylor bubble rise velocity, and (5) detailed comparison between different methods for obtaining the Taylor bubble rise velocity.
The rising Taylor bubble in a stagnant liquid as observed from the high-speed camera and ECT instrument 3D image (from current study) can be seen in

Fluid flow studies using dimensionless numbers
A significant number of dimensionless parameters have been identified to be of relevance in fluid flow studies. Examples of such include bond number, capillary number, drag coefficient, Froude number, inverse dimensionless viscosity, Reynolds number and Weber number to mention a few. In this study, the Reynolds number, inverse dimensionless viscosity and Froude number were used to analyse the experimental results to clearly explain the effect of liquid viscosity on structure velocity.

Reynolds number, Re
The Reynolds number gives a measure of the ratio of inertia forces to viscous forces. Hence, it can be used to depict the competitive interplay between the effect of inertia forces and viscous forces [24].

Reynolds number, Re ¼
Inertia forces Viscous forces (15) Reynolds number can also be used to characterize flow regimes into laminar or turbulent flow. The occurrence of laminar flow is at low Reynolds number in which viscous forces dominate. This is characterized by smooth, constant fluid motion. Turbulent flow on the other hand is at high Reynolds numbers which is associated with chaotic eddies, vortices and other flow instabilities [25].
The slug Reynolds number which is the Reynolds number of the rising slug in the gas-liquid mixture [26][27][28] is expressed as: where μ = dynamic viscosity of the fluid, ρ = density of the fluid, D = diameter of column and U M is the mixture velocity [26,29].

Inverse dimensionless viscosity, N f
According to Lu and Prosperetti [30], the inverse dimensionless viscosity is directly proportional to the fourth root of the Eotvos number raised to a power of three and inversely proportional to the fourth root of the Morton number. Morton number, M o is used alongside with Eotvos number, E o to characterize the shape of bubbles or drops moving in a surrounding fluid or continuous phase.
The inverse dimensionless viscosity is given as: where Morton number is given as: Eotvos number is given as:

Froude number
Apart from the effect of viscous force, gravitational force also affects the rise velocity of Taylor bubble through the liquid. Froude number is a dimensionless parameter which gives a relationship between inertia and gravitational forces. It describes different flow regimes of open channel flow as in the case of the bubble column in this study.

Experimental arrangements
The bubble column experimental set-up consists of a 50 mm internal diameter and 1.6 m long perspex column in a vertical orientation. At the bottom of the column is a single nozzle gas distributor through which gas is introduced into the column. A phantom high-speed camera was used to obtain the video of the gasliquid flow in the column. A frame rate of 1000 pictures per second (pps) and exposure time of 100 μs was used. The geometry specified from the high-speed camera setting gives the image width by image height as 512 by 512 pixel.
Fitted midway to the column is the twin-planes electrical capacitance tomography (ECT) sensor with an interplanar spacing of 30 mm. The 8-electrode system consists of measurement and driven guard electrodes, which is connected to the electrical capacitance tomography processor box, TFLR 5000-20. This sensor electronics gives 28 measurements which are relayed to the computer where image reconstruction occurs, and the data are acquired and processed to obtain the liquid holdup (which is the fraction of liquid in the gas-liquid mixture). The void fraction otherwise known as gas holdup is hence obtained from this.
The ECT is located about 0.7 m above the nozzle, while the liquid level is located about 0.095 m above the ECT sensor. On injecting the gas, the gas flows into the bubble column through the single nozzle gas injector with an orifice diameter of 6.8 mm. A range of silicone oil with viscosities 5, 100, 1000 and 5000 mPa s was used.
The liquid holdup obtained from the ECT was used to obtain the structure velocity (rise velocity of Taylor bubbles) via cross-correlation between two planes -plane 1 and plane 2 putting into consideration the distance between the two planes, 30 mm. The ECT sensor used is shown in Figure 3.
Further details of the experimental arrangements are given in Kajero et al. [32,33].

Rise velocity of Taylor bubbles
In this study, the rise velocity of Taylor bubbles was obtained from ECT (via cross-correlation between signals from planes 1 and 2 as shown in Figure 3), manual time series analysis and the high-speed camera. This rise velocity of Taylor bubbles is also referred to as structure velocity which is from the Taylor bubble periodic structures velocity at real-time measurements.

Cross-correlation
Correlation is the measure of the degree of linear relationship between two variables. Cross-correlation is a statistical method of estimating the degree to which two variables (in this case, time series data sets) are correlated.
The structure velocity was computed from the cross-sectional time averaged void fraction data measured by the ECT for both planes 1 and 2. The cross-correlation between the signals obtained from the two planes gave the structure velocity.
Given two functions x t ð Þ and y t ð Þ, the cross-correlation function, R xy τ ð Þ between them is given as: The correlation coefficient function is expressed as: where τ is the time delay, T is the record length (period), C xy τ ð Þ is the crosscovariance function, C xx 0 ð Þ and C yy 0 ð Þ are auto-covariance functions for x and y, respectively when time delay is zero, μ x and μ y are mean of the corresponding series, and R xx 0 ð Þ and R yy 0 ð Þ are the auto-correlation functions at a time delay of zero [34].
In this experimental work, the two functions x t ð Þ and y t ð Þare time series data of planes 1 and 2, respectively.
The important parameters required for the computation of structure velocity using cross-correlation include: a. Void fraction data for planes 1 and 2.
b. Number of data points.
c. Sampling frequency of data. d. Distance between two planes (planes 1 and 2).
The time taken for the bubbles to travel between the two planes is calculated which then leads to the calculation of the structure velocity. This is done via an Excel Visual Basic Macro program used for the analysis of the time series data [35].
The time series, upstream and downstream with the corresponding correlation are shown in Figure 4. The time delay which is the time taken for the signal to travel between the two planes 1 and 2 is in the interval À 1≤ τ ≤1, where τ is the time delay. Figure 4 indicates periodic structures of short slugs defined as advanced form of spherical cap bubbles gradually developing into clearly distinct slugs. These periodic structures are identified to be void waves in Taylor bubbles and liquid slugs in slug flow. Cross correlation of the time series data from the two axial locations can give  From Figure 5, the y-axis on the cross-correlation plot is the correlation coefficient. This is in the range of À0.3 to +0.9 (though generally falls between À1 and +1). There could be perfect positive correlation (correlation coefficient of +1) or perfect negative correlation (correlation coefficient of À1). A positive correlation indicates that if a signal moves either up or down, the other signal will move in the same direction, while for a negative correlation, if a signal moves either up or down, the other signal will move by an equal amount in the opposite direction. When the correlation is 0, the movement of the signals gives no correlation and is completely random.
The mean velocity for the slug structure can hence be defined as distance between centres of measurement electrodes for two planes divided by the time delay [36][37][38].
The schematics of the ECT Sensor signals generation for measurement of velocity in two-phase flows is shown in Figure 6.

Manual time series analysis
The structure velocity can be obtained manually by estimating the time displacement of output signals from both planes 1 and 2 via the peaks of the time series plots for both planes. Since the distance between the centres of the measurement electrodes for both planes are known, the structure velocity is hence computed using Eq. (22).
The time displacement from the time series data of both planes 1 and 2 used in the calculation of the structure velocity alongside with distance between centres of   Figure 7. This is compared with that obtained via cross correlation later in this paper.

High-speed camera estimate
The camera control software of the Phantom High-Speed Camera can be used to obtain an estimate of the structure velocity. This is done via a playback in which the distance covered by the bubble at a given time is computed which gives a corresponding value for the velocity. A typical window of the camera control software is shown in Figure 8.

Effect of viscosity on the rise velocity (structure velocity)
The effect of liquid viscosity on the rise velocity (structure velocity) has been studied by making a comparison between the respective structure velocities obtained from the ECT for the range of viscosities considered (i.e. 5, 100, 1000 and 5000 mPa s). The physical properties of the liquids used are given in Table 1.
A plot of structure velocity versus superficial gas velocity for all the viscosities considered is given in Figure 9 which shows that structure velocity increases with an increase in superficial gas velocity which is in agreement with the observations of Abdulkaldir et al. [39] and decreases with increase in viscosity as shown in Figure 10 (obtained from ECT Plot3d Image reconstruction software). The structure velocity of 5 and 100 mPa s is found to be approximately the same due to similar void fraction data values. The variation from small to bigger spherical cap and developing slug in 5 and 100 mPa s, and the slug flow in 1000 and 5000 mPa s (as shown in Figure 9) has been discussed by Kajero et al. [33].
This can be explained using the slug Reynolds number, a dimensionless parameter (Eq. 16).
A plot of slug Reynolds number versus superficial gas velocity is made at various viscosities as shown in Figure 11, with an indication of laminar flow. Figure 11 reveals that as viscosity increases, slug Reynolds number decreases tending towards zero. This can be explained as follow: i. Occurrence and prevalence of laminar flow as viscosity increases: According to Bendiksen [21], Reynolds number in the range 5000-110,000 (for low viscous fluids) give turbulent flow. As the Reynolds numbers of the viscosities considered are less than 5000, laminar flow prevails. For large slug Reynolds number, viscous effect will be negligible, while for small slug Reynolds number, viscous effect will be dominant [40]. So, since the slug Reynolds number for all the viscosities are small, viscous effect will be dominant.
ii. Dominating effect of viscous forces over inertia forces: Inertia forces are forces acting due to motion of bubbles through the liquid. It opposes any  force that could resist motion. Viscous forces are forces acting due to the viscous nature of the liquid. Hence, from the plot, it can be inferred that viscous forces have a domineering effect over inertia forces [41].
The forces acting on the Taylor bubble are shown in Figure 12. These forces have an influence on its rise velocity. Surface tension force helps to hold the bubbles together due to the cohesive force existing between them. This accounts for why 'surface tension force' was indicated at the centre of the bubble in Figure 12. Based on the proposition of White and Beardmore [22], the effect of surface tension force can be neglected when Eotvos number is greater than 70. So, since for all the viscosities considered, Eotvos number is greater than 70, its effect on the rise velocity of Taylor bubbles can be neglected.
The dominating effect of viscous forces over inertia forces can be further confirmed using the Inverse dimensionless viscosity according to White and Beardmore [22].
From the various viscosities considered, the dimensionless property numbers are given as follows ( Table 2).
It will be observed that as viscosity increases, Morton number, M o increases while the Eotvos number, E o decreases, which culminates in the decrease of the dimensionless inverse viscosity, N f . This confirms the proposition of Fabre and Line [42]. As the dimensionless inverse viscosity decreases, viscous effect dominates [31].
White and Beardmore [22] proposed that viscous effects come into play when the square of dimensionless inverse viscosity, N f is less than 3 Â 10 5 . i.e.

Drag force, F D
From Stokes law, drag force is given as: (24) based on the fact that the flow is laminar. Drag force is the force due to the resistance provided by the fluid to the motion of a body through it.
From Eq. (24), Drag force is directly proportional to viscosity of the fluid. Drag coefficient, C D is given as: On substituting Eq. (25) into (24) given that Re ¼ ρ L Vd=μ L , Eq. (26) is used to compute drag force, F D which is plotted against superficial gas velocity for all viscosities considered as shown in Figure 13(a). This can be interpreted as follows: i. Drag force is directly proportional to viscosity, which confirms Eq. (26) according to Stokes law. Drag force is the force due to the resistance provided by the viscous nature of the silicone oil fluids to the motion of the large bubble through them. So, as viscosity increases, there is a linear increase in drag force.
ii. The increase in drag force shows an opposition to the motion of the bubbles through the liquid acting parallel to direction of relative motion, hence causing the bubbles to rise at much slower rate. This obviously leads to a decrease in the rise velocity [42].
To further confirm the effect of drag force, drag force was plotted against void fraction, which is the volume fraction of gas in the gas-liquid mixture as shown in Figure 13(b). As void fraction increases, the drag force increases. This is probably due to the increase in the volume fraction of gas (void fraction); hence, an increase in the bubble size relative to the column diameter. This increase in bubble size is due to coalescence and it subsequently requires an increase in drag force to be effective enough to oppose the motion of the bubbles through the liquid. Hence, it can be said that the drag force is exponentially proportional to the void fraction. The drag force was found to be low at low viscosities (5 and 100 mPa s) and increases as viscosity increases. This is due to viscous effect as explained earlier on. As viscosity increases, void fraction increases; hence, an increase in drag force subsequently leads to hindered rise velocity.
The relationship between drag force and drag coefficient is shown in Figure 14. As viscosity decreases, the curves of the respective viscosities are tending towards zero. As drag coefficient increases, the drag force decreases tending towards zero. This gives an inverse relationship, which fits well into a power law expression, with the inverse proportionality constant increasing with increase in viscosity as shown in Table 3.
A plot of drag coefficient against Reynolds number on a log-log plot agrees with Stokes law, which gives an inverse relationship between the drag coefficient and Reynolds number as shown in Figure 15, with purple, green, blue and red for 5, 100, 1000 and 5000 mPa s respectively. This further confirms the fact that the drag force is directly proportional to the liquid viscosity. At low Reynolds number, drag coefficient is high, while at high Reynolds number, drag coefficient is low. Hence, the structure velocity of the latter is greater than the former.  A plot of Froude number versus the superficial gas velocity is shown in Figure 16. The following can be inferred from the plot:   gas velocities, Fr< 1, this implies subcritical flow. The viscous effect makes the flow very slow, with low rise velocity of Taylor bubbles.
2. As viscosity increases, Froude number decreases which indicates the dominant effect of gravitational force over inertia force. This also has a retarding effect on the rise velocity of Taylor bubbles. Hence, this shows that the combined effect of the viscous and gravitational force causes a decrease in the rise velocity of bubbles as viscosity increases. For the high viscous liquids, 1000 and 5000 mPa s, as superficial gas velocity increases, the Froude number tends to be the same. This can be seen at 0.242 and 0.361 m/s. This implies that the dominant effect of gravitational force over inertia force tends to be the same.

Comparison of structure velocity computation methods
The structure velocity obtained from ECT was compared with the manual estimate from the time series [34], high-speed camera and empirical models such as modified Viana et al. [43], modified De Cachard and Delhaye [44], model and Nicklin et al. [7] model. These are shown in Figure 17.
The modified form of De Cachard and Delhaye [44] model is given as: where Γ is given as: Condition for m: The modified form of Viana et al. [43] model is given as: where Fr is obtained from Eqs. (27) to (32).
G ¼ m=C (37) and the parameters a, b, … … … , l ð Þ are: a ¼ 0: 34 On the average, a good agreement exists between the structure velocity from the ECT and that measured from the high-speed camera. A reasonably fair agreement exists between the former and that estimated from the time series. The modified Viana et al. [43], modified De Cachard and Delhaye [44], and Nicklin et al. [7] model gave roughly similar pattern with structure velocity from cross-correlation, manual time series analysis and high-speed camera. Modified Viana et al. [43] and modified De Cachard and Delhaye [44] models showed good agreement with the velocity obtained from cross-correlation, manual time series analysis and high-speed camera at low superficial gas velocity for 1000 and 5000 mPa s, while Nicklin et al. [7] over predicts it.
The variations observed in the agreement could be due to the viscosities of the liquids used. Viana et al. [43] used silicone oil of viscosity range 1-3900 mPa s. De Cachard and Delhaye [44] and Nicklin et al. [7] used water. Also, from the structure velocity plots using cross-correlation, the distribution coefficient is in the range 1.07 to 1.6, while that of modified Viana et al. [43] and De Cachard and Delhaye [44] is approximately 2.25. A distribution coefficient of 2.0 was used for Nicklin et al. [7] as proposed for laminar flow.
The video technique of determining the rise velocity of bubbles gave errors of 4.3, 4.6, 7.3 and 11.5% for 5, 100, 1000 and 5000 mPa s, respectively when compared with cross-correlation technique.

Conclusions
From the foregoing, the following conclusions can be drawn: 1. The forces acting on a Taylor bubble in a 50 mm column diameter include inertia force, surface tension force, viscous force and gravitational force. Surface tension force can be neglected based on Eotvos number greater than 70. The remaining forces have an influence on the rise velocity of Taylor bubble. Furthermore, viscous and gravitational forces were observed to have dominant effect over inertia forces, hence causing the rise velocity of Taylor bubbles to decrease as viscosity increases.
2. The dimensionless parameters: Froude number, Reynolds number and inverse dimensionless viscosity all played vital roles in affecting the rise velocity of Taylor bubbles in various viscosities for a 50 mm diameter column. The dimensionless parameters being functions of the fluid properties and column diameter. Froude number helped to categorize the flow in the four viscosities considered into subcritical (slow and tranquil flow, Fr < 1), with lower rise velocity of Taylor bubbles, critical (Fr = 1) and supercritical flow (fast rapid flow, Fr > 1), with higher rise velocity of Taylor bubbles.
3. The rise of large bubbles through the liquid in the column agrees with Stokes law where drag force is directly proportional to viscosity and an inverse relationship exists between drag coefficient and Reynolds number, as superficial gas velocity increases. The drag force is also exponentially proportional to the void fraction and it retards the motion of the Taylor bubbles through the liquid. An inverse relationship in the form of Power law expression also exists between the drag force and drag coefficient. The drag coefficient was high at low Reynolds number but low at high Reynolds number, which contributed to the rise velocity of Taylor bubbles decreasing from low to higher viscosity liquids.
4. The rise velocity of Taylor bubbles increases with an increase in superficial gas velocity for each viscosity considered (i.e. 5, 100, 1000, 5000 mPa s).
5. The comparison between the rise velocity of Taylor bubbles obtained from the ECT, high-speed camera, cross correlation, manual time series, Nickel et al. [7] model, modified models of De Cachard and Delhaye [44] and Viana et al. [43] gave reasonably fair agreement.

Dedication
This publication is in loving memory of Late Prof. Barry Azzopardi.