Experimental Study of Single Taylor Bubble Rising in Stagnant and Downward Flowing Non-Newtonian Fluids in Inclined Pipes

: An experimental investigation of single Taylor bubbles rising in stagnant and downward ﬂowing non-Newtonian ﬂuids was carried out in an 80 ft long inclined pipe (4 ◦ , 15 ◦ , 30 ◦ , 45 ◦ from vertical) of 6 in. inner diameter. Water and four concentrations of bentonite–water mixtures were applied as the liquid phase, with Reynolds numbers in the range 118 < Re < 105,227 in countercurrent ﬂow conditions. The velocity and length of Taylor bubbles were determined by differential pressure measurements. The experimental results indicate that for all ﬂuids tested, the bubble velocity increases as the inclination angle increases, and decreases as liquid viscosity increases. The length of Taylor bubbles decreases as the downward ﬂow liquid velocity and viscosity increase. The bubble velocity was found to be independent of the bubble length. A new drift velocity correlation that incorporates inclination angle and apparent viscosity was developed, which is applicable for non-Newtonian ﬂuids with the Eötvös numbers ( E 0 ) ranging from 3212 to 3405 and apparent viscosity ( µ app ) ranging from 0.001 Pa · s to 129 Pa · s. The proposed correlation exhibits good performance for predicting drift velocity from both the present study (mean absolute relative difference is 0.0702) and a database of previous investigator’s results (mean absolute relative difference is 0.09614).


Introduction
Gas bubble entry into a wellbore during oil and gas drilling operations can occur under a variety of conditions. Their removal is usually required before drilling can continue. In most operations where sedimentary formations are being drilled with wellbore fluid pressure higher than that of the reservoir (i.e., overbalanced), gas bubbles are removed from the wellbore at surface. Attempting to push the gas back into the reservoir via bullheading will inevitably result in fracturing the rock being drilled-thereby jeopardizing the ability to drill further. However, when drilling through highly fractured and vugular carbonate reservoirs, wellbore fluid pressure is balanced with reservoir pressure. This is because the highly productive fracture/vug network intersects the wellbore with width dimensions too large to be plugged by drilled cuttings or drilling mud filter cake. Such conditions make it relatively easy for gas bubbles to enter a well bore and migrate upward. However, under these same conditions it is also possible to bullhead gas bubbles down a wellbore and back into the reservoir without hindering continued drilling operations. The fluid injection rate necessary to push a gas bubble downward depends on multiple factors-a key one being fluid rheology. In this work, Bingham plastic fluids of varying plastic viscosities, µ p , and yield points, τ y , were tested to determine their impact on gas bubble migration rates under static and downward-moving (i.e., countercurrent) fluid conditions.
When gas bubbles enter a wellbore, the flow pattern can be characterized as gasliquid slug flow. In a vertical conduit, a typical slug unit contains a bullet-shaped nose with a cylindrical body, also known as Taylor bubble, followed by liquid slug. In an inclined conduit, the Taylor bubble is non-concentric, hugging the top side of the flow path. To make sure of the success of a bullheading operation, modeling two-phase slug flow in wellbores is essential. Numerical models, including drift-flux models and twofluid models [1], are commonly used to represent two and three-phase flow in pipes and wellbores, which requires closure relationships to robustly solve the partial differential equations. One of the closure relationships required by slug flow is the translational velocity of Taylor bubbles. Although the translational velocity determined by the experiments of a single Taylor bubble rising in co-current flow fluids is extensively investigated, few of them were carried out in countercurrent flow conditions [2,3]. Additionally, to calculate translational velocity, the drift velocity of a Taylor bubble measured in stagnant liquid should be determined. However, most of the drift velocity correlations were developed from experiments with Newtonian fluids. These correlations may not be directly applicable to the non-Newtonian behavior of drilling fluids in wellbore flow.
Although not considered in this paper, the shape and rise velocity of Taylor bubbles through annuli have been determined experimentally and numerically by several researchers. Kelessidis and Dukler [4] studied the motion of Taylor bubble in vertical concentric and eccentric annuli. They assumed the bubble in an annulus was a 2D bubble of uniform thickness, and indicated that the asymmetric bubbles in an annulus take an elliptic shape which results in higher rise velocities. Based on the work of Kelessidis and Dukler, Das et al. [5] considered the 3D shape of the bubbles and took account of the film thickness. Then, the rise velocity was proposed as a function of the annulus dimensions. Agarwal et al. [6] extended the work of Das to evaluate the shape of Taylor bubbles in annuli with extremely small inner diameter. The rise velocity for such situations showed a better agreement with the prediction model proposed by Das. Lei et al. [7] investigated the shape and motion of a Taylor bubble in a liquid flowing through a thin annulus; a 2D gap-averaged numerical model was developed to evaluate the effects of gap thickness and pipe diameter on the bubble motion in thin annuli. The bubble velocity was found to be highly dependent on the cap structure, whereas it was independent of the bubble length. Although Taylor bubble rising in the annulus is dominant in the drilling process, the motion of Taylor bubbles in pipes is still significant when there is no drillstring in the wellbore (i.e., when tripping out). Moreover, experimental study of Taylor bubbles rising in downward-flowing non-Newtonian fluids in inclined pipes is limited, it is necessary to investigate the movement of Taylor bubbles in the pipe to lay a foundation for the subsequent research on the movement in the annulus.
In the present study, experimental work focused on single Taylor bubbles rising in stagnant and downward flowing water and Bingham plastic fluids in near-vertical and inclined pipes. Four concentrations of bentonite-water mixtures were considered. The velocity and length of Taylor bubbles were determined by measuring differential pressures along the flow path. The effect of liquid rheology and inclination angles on bubble velocity and length are discussed in both stagnant and countercurrent flow conditions. A closure relationship for the drift velocity in terms of apparent viscosity and inclination angle is proposed. The performance of the proposed drift correlation is then compared to that of other existing correlations using our experimental results. Finally, the proposed drift correlation is applied to the experimental data from other investigators to further gauge its performance.

Drift Velocity Correlation
The upward translational velocity of a Taylor bubble (U T ) for all inclinations is modeled by Nicklin et al. [8]: where C 0 is the distribution coefficient which represents the impact of the velocity and concentration profile. The value of constant C 0 is approximately 1.2 for fully developed Energies 2021, 14, 578 3 of 28 turbulent flow and 2.0 for laminar flow. The research of several authors [9][10][11][12] shows that the distribution coefficient significantly depends on liquid viscosity and void fraction. Equation (1) and its associated C 0 values listed above are normally associated with cocurrent flow. In this study, however, we investigate countercurrent flow, which will result in C 0 values much different from those stated above. To represent countercurrent conditions, we can recast Equation (1) as: where U L is the mean downward liquid phase velocity, and U d is the upward drift velocity measured in the stagnant liquid. The movement of Taylor bubbles has been studied numerous times by researchers in different experiment conditions, and many correlations have been proposed to predict the drift velocity. When the gas viscosity is negligible, Dumitrescu [13] and Davies and Taylor [14] formulated the bubble velocity inside a vertical (i.e., θ = 0 • ) tube as: where D is the pipe diameter, g is the gravitational acceleration, and Fr θ is the dimensionless Froude number accounting for the liquid inertia. Other notable dimensionless groups such as the E .. otv ..
os number (E 0 ) and Morton number (M 0 ) have been applied to describe the ratio of gravitational to interfacial forces and viscous to interfacial forces, respectively: where ρ L , σ GL , and µ L are the liquid phase density, gas-liquid interface tension, and viscosity, respectively. As suggested by the work of Viana [15] and Morgado [16], the gas-liquid interface tension (σ GL ) effect on the motion of Taylor bubbles can be negligible when E 0 > 40 and M 0 < 10 -5 The inverse viscosity number N f can be obtained by the Eötvös number (E 0 ) and Morton number (M 0 ): Bendiksen [17] studied the elongated bubbles rising in water in inclined pipes and proposed a simplified correlation for drift velocity by accounting for vertical and horizontal components: where θ is the well bore inclination relative to vertical. Hasan and Kabir [18] correlated the drift velocity through the experimental data gathered at various pipe inclination angles, defined as: Gokcal et al. [19] extended the analysis of Benjamin [20] to evaluate the drift velocity for a horizontal vase, and considered the effect of viscosity on vertical drift velocity: Note that Equation (15) is not dimensionally consistent, because taking the logarithm of a dimensional quantity is not permitted. Here, µ L is the magnitude of its dimensional quantity.
Jeyachandra et al. [21] investigated the effects of high oil viscosities and pipe diameter on the drift velocity for horizontal and upward inclined pipes by extending the work of Gokcal et al. [19] for different pipe diameters and viscosity range. The proposed correlation was formulated in terms of the Eötvös number and inverse viscosity number: Based on the work in refs 19,21, Moreiras et al. [22] proposed a dimensionless closure relationship for drift velocity applicable for a board range of liquid viscosity, inclination angles, and pipe diameters: The drift velocity correlation proposed by Lizarraga-Garcia et al. [23] was extracted from a large database generated by numerical simulations covering a wide range of fluid properties and pipe inclination angles: Fr v = 0. 34 The above correlations were built based on experimental or numerical data. Although the applications of those correlations covered a large range of liquid properties, inclination angles, and pipe diameters, none of them were built for non-Newtonian fluids; research about the movement of a single Taylor bubble in non-Newtonian fluid is still necessary. Carew et al. [24] reported the rise velocity of slug bubbles in Newtonian and non-Newtonian fluids, considering the effect of power-law rheology and inclination angles, and a semi-theoretical correlation applicable for a high Reynolds number (Re) was derived: . .
γ e f f n−1 (27) log 10 Fr θ = −0.2 1.08 − log 10 Re 0.8 Majumdar and Das [25] explored the dynamics of Taylor bubbles rising through powerlaw fluids using CFD (Computational Fluid Dynamics) and a semi-analytical technique; the expression for the drift velocity of Taylor bubble rising through a power-law fluid is given as: From Equations (26) and (27), we know that the drift velocity correlations developed for non-Newtonian fluid are applied for a particular rheology model, and its applicability vanishes for fluids having a yield point.

Theoretical Part
Shosho and Ryan [26] studied the velocity of long bubbles in terms of the Froude, E .. otv .. os, and Morton numbers in Newtonian and non-Newtonian fluids in inclined tubes. Kamışlı [27] proposed one-dimensional approximate equations for long bubbles in capillary vertical and inclined tube filled with power-law fluid; he pointed out that the dimensionless bubble radius decreases with decreasing the power-law index. Sousa et al. [28,29] described the flow field of single Taylor bubble rising in carboxymethyl cellulose (CMC) and polyacrylamide (PAA) polymer by the method of particle image velocimetry (PIV) and shadowgraphy. Rabenjafimanantsoa et al. [30] performed experiments on the dynamics of Taylor bubbles rising in polyanionic cellulose (PAC) in water by using PIV and differential pressure techniques; a more stable and flat trailing edge of the Taylor bubble was found in high viscosity fluids. Araújo et al. [31] focused on the CFD study about the rise of individual Taylor bubbles through inelastic non-Newtonian fluids in stagnant, co-current, and countercurrent conditions; the influence of shear-thinning and shear thickness rheology on Taylor bubble was discussed.
The studies of Taylor bubbles rising in downward flowing non-Newtonian fluids are more focused on the bubble shape and stability. Griffith and Wallis. [32] noted that the stability of the bubble changes and becomes eccentric as the downward flowing liquid velocity increased to a certain point. Lu and Prosperetti [33] studied the instability of Taylor bubbles rising against an incoming liquid stream owing to the flattening of the bubble nose as the liquid flowed downward. Fabre and Figueroa-Espinoza [2] pointed Energies 2021, 14, 578 6 of 28 out that above some critical liquid velocity, the bubble symmetry was broken, and the bubble nose moved close to the tube wall, resulting in a step-function increase in the Taylor bubble upward drift velocity due to its more hydrodynamically efficient shape. Similar to the work of Fabre and Figueroa-Espinoza [2], the relationship between the translational velocity and Taylor bubble shape was investigated in a downward liquid pipe flow by Fershtman et al. [3]. They pointed out that, at downward liquid velocity exceeding a critical value, three different modes of bubble motion were observed (symmetric, asymmetric, and transition between symmetric and asymmetric).
The literature review shows that single Taylor bubbles rising in stagnant Newtonian fluids have been vastly studied. However, the influence of inclination angles, rheological properties, and downward flowing velocities of non-Newtonian fluids, especially having a yield point, on the characteristics of a single Taylor bubble's translational velocity, drift velocity, and length has not been so deeply investigated.

Experimental Facility
Experiments were conducted on the University of Tulsa Low Pressure-Ambient Temperature flow loop (LPAT). Figure 1 shows a general view and schematic view of the experimental facility.
The studies of Taylor bubbles rising in downward flowing non-Newtonian fluids ar more focused on the bubble shape and stability. Griffith and Wallis. [32] noted that th stability of the bubble changes and becomes eccentric as the downward flowing liquid velocity increased to a certain point. Lu and Prosperetti [33] studied the instability of Tay lor bubbles rising against an incoming liquid stream owing to the flattening of the bubbl nose as the liquid flowed downward. Fabre and Figueroa-Espinoza [2] pointed out tha above some critical liquid velocity, the bubble symmetry was broken, and the bubble nos moved close to the tube wall, resulting in a step-function increase in the Taylor bubbl upward drift velocity due to its more hydrodynamically efficient shape. Similar to th work of Fabre and Figueroa-Espinoza [2], the relationship between the translational ve locity and Taylor bubble shape was investigated in a downward liquid pipe flow by Fershtman et al. [3]. They pointed out that, at downward liquid velocity exceeding a crit ical value, three different modes of bubble motion were observed (symmetric, asymmet ric, and transition between symmetric and asymmetric).
The literature review shows that single Taylor bubbles rising in stagnant Newtonian fluids have been vastly studied. However, the influence of inclination angles, rheologica properties, and downward flowing velocities of non-Newtonian fluids, especially having a yield point, on the characteristics of a single Taylor bubble's translational velocity, drif velocity, and length has not been so deeply investigated.

Experimental Facility
Experiments were conducted on the University of Tulsa Low Pressure-Ambien Temperature flow loop (LPAT). Figure 1 shows a general view and schematic view of th experimental facility. The test section was approximately 80 ft long and consisted of a 6 in. inner diamete (ID) transparent acrylic pipe. One end of the flow loop was attached to a vertically mova ble platform, while the other was connected to a trolley, which enabled the experimenta operator to obtain any desired inclination angle between 4° to 90° from vertical [34][35][36][37] Air as the gas phase was provided from an in-house air compressor (working capacity 0 125 psi, maximum 825 scfm capacity) to a buffer accumulator, through a regulator valve which filled the buffer accumulator with air at a fixed pressure. There were four air valve corresponding to four injection ports at different locations ( Figure 2). Port A1 was at th bottom of the flow loop, which was used to release the air in stagnant liquid, we called static condition. Ports A2 and A4 were installed approximately 20 ft from each end of th acrylic pipe, and Port A3 was located in the middle of the other two, as shown in th schematic view of the facility (Figure 1). Those three ports were used to inject air whil The test section was approximately 80 ft long and consisted of a 6 in. inner diameter (ID) transparent acrylic pipe. One end of the flow loop was attached to a vertically movable platform, while the other was connected to a trolley, which enabled the experimental operator to obtain any desired inclination angle between 4 • to 90 • from vertical [34][35][36][37]. Air as the gas phase was provided from an in-house air compressor (working capacity 0-125 psi, maximum 825 scfm capacity) to a buffer accumulator, through a regulator valve, which filled the buffer accumulator with air at a fixed pressure. There were four air valves corresponding to four injection ports at different locations ( Figure 2). Port A1 was at the bottom of the flow loop, which was used to release the air in stagnant liquid, we called static condition. Ports A2 and A4 were installed approximately 20 ft from each end of the acrylic pipe, and Port A3 was located in the middle of the other two, as shown in the schematic view of the facility (Figure 1). Those three ports were used to inject air while the fluid was flowing downward, which was the so-called countercurrent flow condition. If necessary, air trapped at the top of the test section could be discharged directly to the atmosphere through a vent hose. the fluid was flowing downward, which was the so-called countercurrent flow condition. If necessary, air trapped at the top of the test section could be discharged directly to the atmosphere through a vent hose.  (Figure 3b) was used to control pressure in the system. Maintaining adequately high pressures was necessary to prevent the formation of vacuum bubble in the test apparatus-especially at inclinations near vertical. For higher-rheology fluids, a triplex pump supplied by an oil-field service company was used. Liquid flow rate was measured using a Micro-Motion TM Emerson (St. Louis, MO, USA) mass flow meter. Figure 4 shows the flow rig built for this study. A 75 hp centrifugal mud pump (maximum capacity of 650 gpm) was used to supply water and low-rheology fluids, which were injected from the top of the flow loop. Figure 3a,b present the pump connection (inlet) and mud return (outlet), respectively. Fluids were pumped into the riser pipe ( Figure 3a) to the top of the text section to create downward flow. At the bottom of the test section, a gate valve (Figure 3b) was used to control pressure in the system. Maintaining adequately high pressures was necessary to prevent the formation of vacuum bubble in the test apparatus-especially at inclinations near vertical.
Energies 2021, 14, x 7 of 29 the fluid was flowing downward, which was the so-called countercurrent flow condition. If necessary, air trapped at the top of the test section could be discharged directly to the atmosphere through a vent hose.  For higher-rheology fluids, a triplex pump supplied by an oil-field service company was used. Liquid flow rate was measured using a Micro-Motion TM Emerson (St. Louis, MO, USA) mass flow meter. Figure 4 shows the flow rig built for this study. For higher-rheology fluids, a triplex pump supplied by an oil-field service company was used. Liquid flow rate was measured using a Micro-Motion TM Emerson (St. Louis, MO, USA) mass flow meter. Figure 4 shows the flow rig built for this study.  The applied instruments include: 1. Ten fiber optic pressure transducers, sensitive enough to detect small pressure variations with an accuracy of ±0.1%, were installed approximately 5 ft apart for differential pressure and velocity measurements. The pressure and time resolutions of the sensors were 0.001 psi and 0.1 sec, respectively. 2. Three non-contact nuclear densitometers (X96S: Valencia, CA, USA) using a gamma radiation technique, which covered the entire pipe cross section area, for holdup and velocity measurements, were located approximately 20 ft from each end of the test section, and in its middle (spaced about 21 ft apart). The resolutions in density and time were 0.001 ppg and 0.1 sec, respectively. All densitometers were pre-calibrated by the full liquid/gas phase based on two-phase flow systems. 3. Eight high-speed Amcrest security video cameras (Houston, TX, USA) were used for visualization. The non-Newtonian fluid applied in the experiment was opaque; therefore, cameras were only used in the water test. Figure 5 shows a part of the instruments positioned along the test section.   The applied instruments include: 1.
Ten fiber optic pressure transducers, sensitive enough to detect small pressure variations with an accuracy of ±0.1%, were installed approximately 5 ft apart for differential pressure and velocity measurements. The pressure and time resolutions of the sensors were 0.001 psi and 0.1 sec, respectively.

2.
Three non-contact nuclear densitometers (X96S: Valencia, CA, USA) using a gamma radiation technique, which covered the entire pipe cross section area, for holdup and velocity measurements, were located approximately 20 ft from each end of the test section, and in its middle (spaced about 21 ft apart). The resolutions in density and time were 0.001 ppg and 0.1 sec, respectively. All densitometers were pre-calibrated by the full liquid/gas phase based on two-phase flow systems.

3.
Eight high-speed Amcrest security video cameras (Houston, TX, USA) were used for visualization. The non-Newtonian fluid applied in the experiment was opaque; therefore, cameras were only used in the water test.  The applied instruments include: 1. Ten fiber optic pressure transducers, sensitive enough to detect small pressure variations with an accuracy of ±0.1%, were installed approximately 5 ft apart for differential pressure and velocity measurements. The pressure and time resolutions of the sensors were 0.001 psi and 0.1 sec, respectively. 2. Three non-contact nuclear densitometers (X96S: Valencia, CA, USA) using a gamma radiation technique, which covered the entire pipe cross section area, for holdup and velocity measurements, were located approximately 20 ft from each end of the test section, and in its middle (spaced about 21 ft apart). The resolutions in density and time were 0.001 ppg and 0.1 sec, respectively. All densitometers were pre-calibrated by the full liquid/gas phase based on two-phase flow systems. 3. Eight high-speed Amcrest security video cameras (Houston, TX, USA) were used for visualization. The non-Newtonian fluid applied in the experiment was opaque; therefore, cameras were only used in the water test. Figure 5 shows a part of the instruments positioned along the test section.  Data acquisition was performed with LabView TM (Austin, TX, USA) software which automatically recorded data, including mud flow rate, mixture density, inclination angle, pressure and differential pressure along the test section, and pressure and temperature at the inlet and outlet of the test section.

Experimental Procedure and Test Matrix
In this test, bentonite-water mixtures were used as the liquid phase. The bentonitewater mixtures were firstly prepared in a 100-barrel mud tank, then circulated into the system for at least 30 min to fully mix, waiting for at least 24 h for its full hydration. Rheology properties were checked before every new daily test.
In the static test, the fluid was firstly circulated at a high flow rate to make sure that the flow condition was homogeneous and to flush out the air in the test section. Then, the gate valve at the bottom was gradually closed while maintaining sufficient pressure on the system with the pump to prevent the formation of vacuum bubbles. At this point, the whole system was shut-in to avoid any effect of gas expansion on the bubble migration. Using the injection port A1, a finite burst of air was released. The pressure of the air accumulation chamber (Figure 2) was set anywhere between 20 to 50 psi above the pressure at port A1. Air was released for a 3 sec interval. When the bubble reached the top of the loop, the test was complete and data recording was stopped. The bubble was then purged from the system in preparation for another test.
In the countercurrent flow test, for each inclination angle, fluids were downwardly flowing at a given velocity; the Taylor bubble was released at port A2 and its rise velocity was measured. The bubble was then purged from the system in preparation for another test at a different countercurrent rate.
Water and four concentrations of bentonite-water mixtures were applied in the tests. The physical properties and test matrix are described in Table 1. All the experiments were carried out in inclined pipes at 4 • , 15 • , 30 • , and 45 • from vertical. In the present study, the gas-liquid interface tension of various concentrations of bentonite in water was applied from the experimental results of Dadashev et al. [38]-the range of the Eötvös number was 3212 < E 0 < 3405; therefore, the effect of σ GL on the movement of the Taylor bubble was not discussed in this work [15]. 1 Plastic viscosity (µ p ) is defined as µ p = θ 600 − θ 300 ; 2 Yield point was obtained by τ y = θ 300 − µ p , where θ 600 and θ 300 are the dial readings at 600 rpm and 300 rpm, respectively, using a Fann 35 rotational viscometer with an R1B1 bob (using F1 and F2 springs as appropriate); 3 Reynolds number (R e ) is defined as R e = ρ L U L D µp [39].
In Table 1, "10 lb/bbl BW" means 10 pounds per barrel bentonite-water mixture; "S" and "C" represent the tests conducted in static or countercurrent flow conditions, respectively; and µ p and τ y are the plastic viscosity and yield point in the Bingham plastic model, respectively. Rheological behavior is presented in Appendix B. As can be seen, countercurrent flow tests carried out in water were in a turbulent regime (R e > 2100) and a laminar flow regime (R e < 2100) in non-Newtonian fluids.

Determination of Bubble Velocity and Length
The main objective of this study was to investigate the migration of a single Taylor bubble in stagnant and downward non-Newtonian fluid flow in inclined pipes. The main hydrodynamic behaviors including translational velocity (U T ), drift velocity (U d ), and the length of the Taylor bubble (L TB ) in the countercurrent flow condition were determined.
Velocities of Taylor bubbles are determined by the time for the Taylor bubble to travel between a fixed distance. Then, with the known of velocity, the length of Taylor bubble can be obtained by the time required for the bubble to pass through its own length ( Figure 6).

Determination of Bubble Velocity and Length
The main objective of this study was to investigate the migration of a single Taylo bubble in stagnant and downward non-Newtonian fluid flow in inclined pipes. The ma hydrodynamic behaviors including translational velocity ( ), drift velocity ( ), and th length of the Taylor bubble ( ) in the countercurrent flow condition were determined Velocities of Taylor bubbles are determined by the time for the Taylor bubble to trav between a fixed distance. Then, with the known of velocity, the length of Taylor bubb can be obtained by the time required for the bubble to pass through its own length (Figu 6). In the present study, differential pressure (∆P) was applied to determine the veloci and length of Taylor bubbles. Figure 7 shows the differential pressure measured between two pressure transdu ers. ∆Pn measures the pressure difference between pressure transducer ∆Pn and ∆Pn while ∆Pn+1 measures the pressure difference between ∆Pn+1 and ∆Pn+2. Considering a Ta lor bubble moving inside the pipe, when its nose passes the pressure transducer PTn, th hydrostatic pressure at PTn decreases, the differential pressure ∆Pn starts decreasing sy chronously. ∆Pn+1 remains stable until the bubble moves passed pressure transducer ∆Pn Then, the same trend happens on ∆Pn+1. In the present study, differential pressure (∆P) was applied to determine the velocity and length of Taylor bubbles. Figure 7 shows the differential pressure measured between two pressure transducers. ∆P n measures the pressure difference between pressure transducer ∆P n and ∆P n+1 , while ∆P n+1 measures the pressure difference between ∆P n+1 and ∆P n+2 . Considering a Taylor bubble moving inside the pipe, when its nose passes the pressure transducer PT n , the hydrostatic pressure at PT n decreases, the differential pressure ∆P n starts decreasing synchronously. ∆P n+1 remains stable until the bubble moves passed pressure transducer ∆P n+1 . Then, the same trend happens on ∆P n+1 .

Determination of Bubble Velocity and Length
The main objective of this study was to investigate the migration of a single Taylor bubble in stagnant and downward non-Newtonian fluid flow in inclined pipes. The main hydrodynamic behaviors including translational velocity ( ), drift velocity ( ), and the length of the Taylor bubble ( ) in the countercurrent flow condition were determined. Velocities of Taylor bubbles are determined by the time for the Taylor bubble to travel between a fixed distance. Then, with the known of velocity, the length of Taylor bubble can be obtained by the time required for the bubble to pass through its own length ( Figure  6). In the present study, differential pressure (∆P) was applied to determine the velocity and length of Taylor bubbles. Figure 7 shows the differential pressure measured between two pressure transducers. ∆Pn measures the pressure difference between pressure transducer ∆Pn and ∆Pn+1, while ∆Pn+1 measures the pressure difference between ∆Pn+1 and ∆Pn+2. Considering a Taylor bubble moving inside the pipe, when its nose passes the pressure transducer PTn, the hydrostatic pressure at PTn decreases, the differential pressure ∆Pn starts decreasing synchronously. ∆Pn+1 remains stable until the bubble moves passed pressure transducer ∆Pn+1. Then, the same trend happens on ∆Pn+1.  Thus, the time between two sudden drops on the differential pressure curves is the duration of the Taylor bubble traveling between a fixed distance. The velocity of the Taylor bubble is obtained by averaging the velocities measured across the various possible pairings of pressure transducers:

Length of Taylor Bubbles
To determine the length of a bubble, two different scenarios are discussed. The first scenario is when the length of Taylor bubble is shorter than the distance between two pressure transducers, as can be seen in Figure 8a. As the bubble reaches PT m , differential pressure ∆P m starts decreasing reaching a minimum when the bubble tail reaches PT m . As the tail of the bubble moves past PT m , there is a slight increase in ∆P m (shown between regions (1) and (2)). This is an artifact of the higher-pressure region behind the bubble being sensed by PT m . The fluid moving past the bubble creates a pressure reduction that is recovered in the wake of the upwardly migrating bubble. Region (2) is the period where the entire bubble length resides between PT m and PT m+1 . The transition between regions (2) and (3) indicates when the bubble's nose reaches PT m+1 .

Length of Taylor Bubbles
To determine the length of a bubble, two different scenarios are discussed. The first scenario is when the length of Taylor bubble is shorter than the distance between two pressure transducers, as can be seen in Figure 8a. As the bubble reaches PTm, differential pressure ∆Pm starts decreasing reaching a minimum when the bubble tail reaches PTm. As the tail of the bubble moves past PTm, there is a slight increase in ∆Pm (shown between regions (1) and (2)). This is an artifact of the higher-pressure region behind the bubble being sensed by PTm. The fluid moving past the bubble creates a pressure reduction that is recovered in the wake of the upwardly migrating bubble. Region (2) is the period where the entire bubble length resides between PTm and PTm+1. The transition between regions (2) and (3) indicates when the bubble's nose reaches PTm+1.
Then, a second scenario is considered where a bubble's length is larger than the distance between two pressure transducers (Figure 8b). We start with the bubble reaching PTm. Region (1) signifies when the bubble's nose is moving between PTm and PTm+1. Region (2) signifies when the bubble spans both PTm and PTm+1. Additionally, the transition between regions (2) and (3) signifies when the tail of the bubble passes PTm. In each scenario, Δtm defines the time required for a Taylor bubble to pass PTm. Thus, with the known bubble velocity and transit time across the sensor, the average length of Taylor bubble can be obtained: (1) Time, t (s) Time, t (s) Figure 8. Differential pressure to determine the length of Taylor bubble. (a) The length of the Taylor bubble is shorter than the distance of two pressure transducers; (b) the length of the Taylor bubble is larger than the distance between two pressure transducers.
Then, a second scenario is considered where a bubble's length is larger than the distance between two pressure transducers (Figure 8b). We start with the bubble reaching PT m . Region (1) signifies when the bubble's nose is moving between PT m and PT m+1 . Region (2) signifies when the bubble spans both PT m and PT m+1 . Additionally, the transition between regions (2) and (3) signifies when the tail of the bubble passes PT m .
In each scenario, ∆t m defines the time required for a Taylor bubble to pass PT m . Thus, with the known bubble velocity and transit time across the sensor, the average length of Taylor bubble can be obtained:

Results and Discussion
The experimental data collected from both static tests and countercurrent flow tests were processed to obtain the velocity and length of the Taylor bubble. The propagation of a Taylor bubble through various stagnant liquids is presented in Figure 9, where the drift velocity is plotted as a function of inclination angles. Each line corresponds to a particular fluid. As can be seen, the dependence of drift velocity on inclination angle for non-Newtonian fluids is observed to be similar to Newtonian fluids. In each case, the drift velocity increases as the inclination angle increases. This result is consistent with Zukoski [40], who verified that the drift velocity would reach the maximum between 40 • and 60 • , and then decrease as the inclination approaches horizontal. Additionally, increasing the concentration of bentonite in our fluid, which increases plastic viscosity and the yield point, decreased drift velocity. This reduction was more significant at a lower inclination angle. It is noteworthy that an unexpected case was detected in 48.5 lb/bbl BW at 4 • , where the Taylor bubble stopped moving soon after it entered the pipe. We will discuss this behavior in more detail in the subsequent section.

Results and Discussion
The experimental data collected from both static tests and countercurrent flow tests were processed to obtain the velocity and length of the Taylor bubble.

Influence of and on
The propagation of a Taylor bubble through various stagnant liquids is presented in Figure 9, where the drift velocity is plotted as a function of inclination angles. Each line corresponds to a particular fluid. As can be seen, the dependence of drift velocity on inclination angle for non-Newtonian fluids is observed to be similar to Newtonian fluids. In each case, the drift velocity increases as the inclination angle increases. This result is consistent with Zukoski [40], who verified that the drift velocity would reach the maximum between 40° and 60°, and then decrease as the inclination approaches horizontal. Additionally, increasing the concentration of bentonite in our fluid, which increases plastic viscosity and the yield point, decreased drift velocity. This reduction was more significant at a lower inclination angle. It is noteworthy that an unexpected case was detected in 48.5 lb/bbl BW at 4°, where the Taylor bubble stopped moving soon after it entered the pipe. We will discuss this behavior in more detail in the subsequent section. The lower drift velocities in a more viscous fluid can be explained by the higher drag forces imposed on the bubble's flanks and by changes in the shape of the bubble's nose. Figures 10a-d illustrate the air-liquid mixture density (ρm) change with time as the bubble migrates past the densitometers at 45° from vertical. As the Taylor bubble nose reached the densitometer, the mixture density dropped because of the lower density of air. It also shows that there was a flow development between densitometers D1 and D2 after the gas bubble was released. The mixture density can be obtained by the relationship: The lower drift velocities in a more viscous fluid can be explained by the higher drag forces imposed on the bubble's flanks and by changes in the shape of the bubble's nose. Figure 10a-d illustrate the air-liquid mixture density (ρ m ) change with time as the bubble migrates past the densitometers at 45 • from vertical. As the Taylor bubble nose reached the densitometer, the mixture density dropped because of the lower density of air. It also shows that there was a flow development between densitometers D1 and D2 after the gas bubble was released. The mixture density can be obtained by the relationship: where α g is the void fraction defined as the fraction of the cross-sectional area of the channel that is occupied by the gas phase. In fluids with lower concentrations of bentonite (Figure 10a,b), the mixture density gradually decreased during the passage of the bubble's nose. However, for higher concentrations of bentonite (Figure 10c,d), a steeper decrease in the mixture density occurred, indicating a steeper increase in α g . This indicates that the curvature of the Taylor bubble's nose is blunted with a significant increase in viscosity and yield point. This is also suggested by the experimental work conducted by Carew et al. [24], who reported that an increase in the liquid viscosity caused the nose region of the bubbles to become blunter and led to a decrease in the bubble rise velocities. (Figures 10a,b), the mixture density gradually decreased during the passage of the bubble's nose. However, for higher concentrations of bentonite (Figure 10c,d), a steeper decrease in the mixture density occurred, indicating a steeper increase in . This indicates that the curvature of the Taylor bubble's nose is blunted with a significant increase in viscosity and yield point. This is also suggested by the experimental work conducted by Carew et al. [24], who reported that an increase in the liquid viscosity caused the nose region of the bubbles to become blunter and led to a decrease in the bubble rise velocities. In this study, we used apparent viscosity as a singular measure of a non-Newtonian fluid's rheological state, which was useful in the formulation of our generalized bubble drift correlation. The bentonite-water mixtures used in this study exhibited non-Newtonian rheological behavior having a yield point, which can be represented by the Herschel-Bulkley model: where is the shear stress, is the yield point, is the consistency index, is the shear rate, and is the flow index. For a single Taylor bubble rising in stagnant liquid, the characteristic shear rate around the bubble can be scaled as [18]: In this study, we used apparent viscosity as a singular measure of a non-Newtonian fluid's rheological state, which was useful in the formulation of our generalized bubble drift correlation. The bentonite-water mixtures used in this study exhibited non-Newtonian rheological behavior having a yield point, which can be represented by the Herschel-Bulkley model: where τ is the shear stress, τ y is the yield point, k is the consistency index, . γ is the shear rate, and n is the flow index. For a single Taylor bubble rising in stagnant liquid, the characteristic shear rate around the bubble can be scaled as [18]: .
Thus, the apparent viscosity is given by: In this study, bentonite-water mixtures were represented by a Bingham plastic model, having a yield point. Thus, the k in Equation (36) was the plastic viscosity (µ p ), and n became unity. Therefore, the apparent viscosity for a Bingham plastic fluid can be written as: From Equation (39), we know that when the pipe diameter is fixed, µ app is then governed by both fluid properties and drift velocity.
Based on Equation (39), the quantitative influence of apparent viscosity on drift velocity can be obtained-which is illustrated in Figure 11. For all fluids studied, the drift velocity decreased dramatically with the increase in apparent viscosity, and this reduction became even more significant when the inclination angle decreased. The apparent viscosity for Bingham plastic fluids in terms of plastic viscosity and yield point became significant when the plastic viscosity and yield point increased with the concentration of bentonite, making it harder for the liquid film to flow around the bubble.
Thus, the apparent viscosity is given by: In this study, bentonite-water mixtures were represented by a Bingham plastic model, having a yield point. Thus, the in Equation (36) was the plastic viscosity ( ), and became unity. Therefore, the apparent viscosity for a Bingham plastic fluid can be written as: From Equation (39), we know that when the pipe diameter is fixed, is then governed by both fluid properties and drift velocity.
Based on Equation (39), the quantitative influence of apparent viscosity on drift velocity can be obtained-which is illustrated in Figure 11. For all fluids studied, the drift velocity decreased dramatically with the increase in apparent viscosity, and this reduction became even more significant when the inclination angle decreased. The apparent viscosity for Bingham plastic fluids in terms of plastic viscosity and yield point became significant when the plastic viscosity and yield point increased with the concentration of bentonite, making it harder for the liquid film to flow around the bubble.

Proposed Drift Velocity Correlation
Based on the variation of drift velocity with the inclination angles and apparent viscosity ranging from 0.001 • < < 129 • , the proposed correlation (3212 < < 3405) combines Bendiksen et al. [17] in terms of and a yield-power-law correlation in terms of [41,42]. The mathematical form of the proposed drift velocity correlation is given as: This proposed correlation is similar to Bendiksen et al. [17], and the proportional value of the horizontal drift velocity component is modified as 0.45, which was inspired Figure 11. Influence of apparent viscosity on drift velocity (data set in Appendix A).

Proposed Drift Velocity Correlation
Based on the variation of drift velocity with the inclination angles and apparent viscosity ranging from 0.001 Pa·s < µ app < 129 Pa·s , the proposed correlation (3212 < E 0 < 3405) combines Bendiksen et al. [17] in terms of θ and a yield-power-law correlation in terms of µ app [41,42]. The mathematical form of the proposed drift velocity correlation is given as: This proposed correlation is similar to Bendiksen et al. [17], and the proportional value of the horizontal drift velocity component is modified as 0.45, which was inspired by Gokcol et al. [19] and suggested by Bhagwat and Ghajar [10]. Fitted values of these parameters are given as a = 0.99416, b = −0.18458, c = 0.21337, d = 0.48996, and e = 0.02661, respectively.
Performance of Proposed Correlation Figure 12 presents the performance of the proposed and comparative correlations against measured drift velocity (drift velocity measured in 48.5 lb/bbl at 4 • , which is zero, is excluded). The correlations developed by Carew et al. [24] and Majumdar and Das. [25] for power-law fluids were not applied for performance-comparison. Those correlations utilize n and k and are, therefore, not aligned with our Bingham plastic model. The dashed lines are the ±15%, ±20%, and 50% error bands showing in Figure 12. It is obvious that the proposed correlation shows good performance, slightly under-predicting the drift velocity. The major differences among correlations occurred in the relatively low value of Fr corresponding to non-Newtonian fluids. The Gokcal et al. [19], Moreiras et al. [22] and Jeyachandra et al. [21] correlations show better prediction when Fr was high and tend to overestimate the drift velocity when Fr < 0.4, which is thought to be caused by over-predicting the influence of viscosity on drift velocity. Bendiksen [17] and Hasan and Kabir [13] correlations were developed by water, which could have a relatively better performance for low-viscosity fluids.
against measured drift velocity (drift velocity measured in 48.5 lb/bbl at 4°, whic is excluded). The correlations developed by Carew et al. [24] and Majumdar and for power-law fluids were not applied for performance-comparison. Those cor utilize n and k and are, therefore, not aligned with our Bingham plastic model. Th lines are the ±15%, ±20%, and 50% error bands showing in Figure 12. It is obviou proposed correlation shows good performance, slightly under-predicting the dr ity. The major differences among correlations occurred in the relatively low va corresponding to non-Newtonian fluids. The Gokcal et al. [19], Moreiras et al. Jeyachandra et al. [21] correlations show better prediction when was high an overestimate the drift velocity when < 0.4, which is thought to be caused by dicting the influence of viscosity on drift velocity. Bendiksen [17] and Hasan a [13] correlations were developed by water, which could have a relatively bette mance for low-viscosity fluids.
Mean absolute relative difference, ε 2 : Fr θi, cal − Fr θi,mea Fr θi,mea (42) Standard deviation, ε 3 : As shown in Table 2, the proposed correlation predicted 90% of data points within ±15% error bands and 85% within ±10%. The proposed correlation had a mean absolute relative difference value of 0.0702. The Gokcal et al. [19] correlation gave the best performance among the comparative correlations within ±15%, but fewer data points were within the ±15% error band compared with the Jeyachandra et al. [21] correlation. The poor performance of the Hasan and Kabir [18] correlation is thought to be because of the lowviscosity fluid which they used to develop the correlation. Based on the performance of different correlations for predicting drift velocity, we know that although the behaviors of drift velocity against inclination angles between Newtonian and non-Newtonian fluids are similar, those correlations developed for Newtonian fluids may not be applicable for non-Newtonian fluids.

Validation of the Proposed Correlation
An additional database collected from Carew et al. [24], Sousa et al. [28], Eghorieta et al. [43], and Livinus et al. [44] was applied to validate the applicability of the proposed correlation. The database is summarized in Table 3 (inclination angles are measured from vertical). The non-Newtonian fluids (power-law) experimental data reported by Carew et al. [24] and Sousa et al. [28] are for different concentrations of CMC and Carbopol, covering a range of pipe diameters and liquid viscosities. The Newtonian fluids used by Eghorieta et al. [43], and Livinus et al. [44] were water and oil, respectively. According to the rheological model for Newtonian and power-law fluids, the apparent viscosity for a Newtonian fluid is the same as liquid viscosity; the yield point in a power-law fluid is zero. Thus, Equation (38) can be written as: Power-law : The performance of the proposed correlation against the new database is shown in Figure 13. It is evident that the proposed correlation exhibits good experimental performance for both Newtonian and non-Newtonian fluid data points, but slightly overpredicted by Carew et al. [24], who obtained data from small diameter pipes. The proposed correlation shows a small value of mean absolute relative difference of 0.09614. This is attributed to the apparent viscosity used in this study, which can be applied for both Newtonian and non-Newtonian fluids.

Newtonian:
= ( Power-law: = ( The performance of the proposed correlation against the new database is shown Figure 13. It is evident that the proposed correlation exhibits good experimental perf mance for both Newtonian and non-Newtonian fluid data points, but slightly over-p dicted by Carew et al. [24], who obtained data from small diameter pipes. The propos correlation shows a small value of mean absolute relative difference of 0.09614. This attributed to the apparent viscosity used in this study, which can be applied for both Ne tonian and non-Newtonian fluids.

Influence of on
In Figure 14, the translational velocity is plotted against inclination angles (cur with = 0 m/s for reference). As can be seen, the translational velocity increases with increase in inclination angle; this behavior is similar to the results obtained in static c dition because of an increase in the drift velocity.

Influence of θ on U T
In Figure 14, the translational velocity is plotted against inclination angles (curves with U L = 0 m/s for reference). As can be seen, the translational velocity increases with the increase in inclination angle; this behavior is similar to the results obtained in static condition because of an increase in the drift velocity. Figure 15a,b indicate the influence of liquid viscosity induced by different concentrations of bentonite-water solutions on bubble velocity under fixed downward liquid velocities. As can be seen, in agreement with the observations in stagnant liquid test, an increase in the liquid viscosity results in a reduction in the bubble velocity for a given countercurrent flow liquid velocity. This reduction in the translational velocity is qualitatively the same as the inclination angle increases.

Distribution coefficient, C
Once the measured in countercurrent flow and measured from static test are obtained, the distribution coefficient can be calculated by Equation (2) in the form of = ( − )/ . Figures 16a,b plot against the inclination angles for both water and non-Newtonian fluids.
tends to increase with the increase in inclination angle.

Distribution coefficient, C 0
Once the U T measured in countercurrent flow and U d measured from static test are obtained, the distribution coefficient C o can be calculated by Equation (2) in the form of Figure 16a,b plot C 0 against the inclination angles for both water and non-Newtonian fluids. C 0 tends to increase with the increase in inclination angle.

Distribution coefficient, C
Once the measured in countercurrent flow and measured from static tes obtained, the distribution coefficient can be calculated by Equation (2) in the for = ( − )/ . Figures 16a,b plot against the inclination angles for both w and non-Newtonian fluids.
tends to increase with the increase in inclination angl The general increasing trend of with wellbore inclination is the result of creasing more rapidly than as increases. is measured under static fluid co tions; in a near-vertical wellbore, a Taylor bubble's shape is relatively concentric, minimizing (because it is hydrodynamically inefficient). As wellbore inclination parts from vertical, buoyancy forces lateral to the bubble increasingly deform its shap make the bubble is more hydrodynamically efficient, resulting in significantly highe values.
, on the other hand, is the translational velocity of a Taylor bubble that, bec of countercurrent flow, is always in a non-concentric (and, thus, a hydrodynamically cient) shape, regardless of the flow conduit's inclination for all our experiments. Hen progressively increasing drives greater increases in measured under static co tions (see Figure 9) than measured under constant conditions (see Figure 14). The general increasing trend of C o with wellbore inclination is the result of U d increasing more rapidly than U T as θ increases. U d is measured under static fluid conditions; in a near-vertical wellbore, a Taylor bubble's shape is relatively concentric, thus minimizing U d (because it is hydrodynamically inefficient). As wellbore inclination departs from vertical, buoyancy forces lateral to the bubble increasingly deform its shape to make the bubble is more hydrodynamically efficient, resulting in significantly higher U d values. U T , on the other hand, is the translational velocity of a Taylor bubble that, because of countercurrent flow, is always in a non-concentric (and, thus, a hydrodynamically efficient) shape, regardless of the flow conduit's inclination for all our experiments. Hence, a progressively increasing θ drives greater increases in U d measured under static conditions (see Figure 9) than U T measured under constant U L conditions (see Figure 14).
According to Zuber and Findley [45], the value of C 0 could be less than unity when the ratio of the volumetric gas concentration evaluated at the wall, α w , and at the centerline, α c , of the duct is greater than one, α w > α c . This occurs when a Taylor bubble flows inside an inclined pipe, where a loss of symmetry is observed due to the component of buoyancy force leading to the liquid film thickness between the top wall of the pipe and Taylor bubble decreasing when inclination angle is increased from vertical. Although Zuber and Findley's insight comes from an analysis based on a co-current flow perspective, their finding correctly applies to our countercurrent test results, regardless of inclination. This is because our imposed countercurrent flow conditions always result in a non-concentric Taylor bubble that tracks along the flow conduit wall.
Additionally, under co-current flow conditions, C 0 values obtained for non-Newtonian fluid is generally greater than that in Newtonian fluid (e.g., water), which is due to the different flow profiles between Newtonian fluids and non-Newtonian fluids ( Figure 17). For Newtonian fluid (Poiseuille flow), the flow profile is parabolic, while for shear-thinning non-Newtonian fluid, such as a Bingham plastic, the profile of the flow is flatter in the middle, called the plug region, and decays faster towards the wall. For a Taylor bubble rising in co-current flow liquid in a vertical orientation, the leading nose will follow the highest velocity, which is at the center of the pipe. In contrast, the bubble will tend to follow the lowest liquid velocity which is at the wall in countercurrent flow liquid. As shown in Figure 16, under countercurrent conditions, there is little difference in the values of C 0 between Newtonian and non-Newtonian fluids. Again, this is because our imposed countercurrent flow conditions always result in a non-concentric Taylor bubble that tracks along the conduit wall. A break from this trend only occurs at higher bentonite concentrations where viscous forces are so high that plug flow occurs (i.e., when C 0 = 1 and U d = 0 m/s).
bubble rising in co-current flow liquid in a vertical orientation, the leading low the highest velocity, which is at the center of the pipe. In contrast, the bu to follow the lowest liquid velocity which is at the wall in countercurrent f shown in Figure 16, under countercurrent conditions, there is little differenc of between Newtonian and non-Newtonian fluids. Again, this is becaus countercurrent flow conditions always result in a non-concentric Taylor bub along the conduit wall. A break from this trend only occurs at higher ben trations where viscous forces are so high that plug flow occurs (i.e., when = 0 m/s).

Length of Taylor Bubbles
The influence of downward flowing velocity and inclination angle on Taylor bubble (the volume of the Taylor bubble is the same in all cases) in w ized in Figure 18a-d. The reason why all pipes appear to be vertically orient the camera was mounted right above the pipe. As shown in Figure 18, the b asymmetric with a fluctuating nose under countercurrent flow condition. locities of water from left to right for each figure are 0.173 m/s, 0.346 m/s, respectively. In higher downward flowing velocities, bubble surface waves nose and are swept down to the tail. It is obvious that for each inclination an of the Taylor bubbles decrease with the increase in downward liquid flo cause of increasing liquid film thickness/holdup and increasing entrainm the wake region. However, no direct relationships between the Taylor bub the inclination angle were observed.

Length of Taylor Bubbles
The influence of downward flowing velocity and inclination angle on the length of Taylor bubble (the volume of the Taylor bubble is the same in all cases) in water is visualized in Figure 18a-d. The reason why all pipes appear to be vertically orientated is because the camera was mounted right above the pipe. As shown in Figure 18

Influence of Liquid Viscosity on L TB
The influence of liquid viscosity induced by higher concentration of bentonite on L TB under fixed downward liquid velocities is presented in Figure 20a (U L = 0.218 m/s) and Figure 20b (U L = 0.5 m/s). As can be seen, an increase in the liquid viscosity results in a reduction in the bubble length for a given countercurrent flow velocity. As discussed before, the reduction in the translational velocity as the liquid viscosity increases is qualitatively the same as the inclination angle increases. However, the decrease in the bubble length due to viscosity increase is different as the inclination angle increases, which indicates that the translational velocity is independent of bubble length. As can be seen, an increase in the liquid viscosity results in a reduction in the bubble length for a given countercurrent flow velocity. As discussed before, the reduction in the translational velocity as the liquid viscosity increases is qualitatively the same as the inclination angle increases. However, the decrease in the bubble length due to viscosity increase is different as the inclination angle increases, which indicates that the translational velocity is independent of bubble length.

Unexpected Case: The Stagnant Taylor Bubble
In the static test of 48.5lb/bbl BW at 4°, the Taylor bubble stopped moving soon after it entered the pipe (Figure 21). This was because of the thixotropic properties of bentonite-

Unexpected Case: The Stagnant Taylor Bubble
In the static test of 48.5lb/bbl BW at 4 • , the Taylor bubble stopped moving soon after it entered the pipe (Figure 21). This was because of the thixotropic properties of bentonitewater mixtures. When a thixotropic material is sheared, the buildup and breakdown of gel structure processes compete, and a dynamic equilibrium eventually results. For high concentrations of bentonite, a strong gelled structure is formed over little time when not subject to shearing. This provides resistance to the bubble movement, eventually resulting in a stagnant bubble.

Unexpected Case: The Stagnant Taylor Bubble
In the static test of 48.5lb/bbl BW at 4°, the Taylor bubble stopped moving soon after it entered the pipe (Figure 21). This was because of the thixotropic properties of bentonitewater mixtures. When a thixotropic material is sheared, the buildup and breakdown of gel structure processes compete, and a dynamic equilibrium eventually results. For high concentrations of bentonite, a strong gelled structure is formed over little time when not subject to shearing. This provides resistance to the bubble movement, eventually resulting in a stagnant bubble.  The effect of gel strength can be directly obtained by comparing the drift velocities at different positions along the pipe. Figure 22a-d show the drift velocity measured at different pressure transducers against the distance it travels, each line representing a particular inclination angle. As can be seen, in the lower concentrations, the drift velocity remains relatively constant during the bubble movement in all inclination angles, while a slight decrease in drift velocity can be observed in 45 lb/bbl bentonite at 4 • as the bubble migrates (Figure 22c), then a great reduction in drift velocity is shown in 48.5 lb/bbl bentonite at even higher inclination angles (Figure 22d), which is owing to the lower drift velocity at low inclination angles corresponding to a longer stationary time, eventually resulting in a stronger gel structure.
Gel structure breakdown can be induced by shearing. Figure 23 shows the propagation of a Taylor bubble against a range of downward liquid velocities (0.073 m/s < U L < 0.5 m/s) at 4 • in 48.5 lb/bbl bentonite-water mixtures after it has been circulated. Notice that at U L = 0 m/s the bubble migration stops. This is due to the fluid gelling tendency. However, under a shearing condition (U L > 0 m/s) the gel structure is broken, and the bubble migration becomes dependent on U L . The relationship between the U L and U T under countercurrent flow condition can be fitted in a linear function (red dash line), as suggested in Equation (2). The second constant (0.172 m/s) is normally interpreted as the drift velocity. In this case, however, the strong gelling nature of the fluid invalidates this interpretation.  (Figure 22c), then a great reduction in drift velocity is shown in 48.5 lb/bbl bentonite at even higher inclination angles (Figure 22d), which is owing to the lower drift velocity at low inclination angles corresponding to a longer stationary time, eventually resulting in a stronger gel structure. Gel structure breakdown can be induced by shearing. Figure 23 shows the propagation of a Taylor bubble against a range of downward liquid velocities (0.073 m/s < < 0.5 m/s) at 4° in 48.5 lb/bbl bentonite-water mixtures after it has been circulated. Notice that at = 0 m/s the bubble migration stops. This is due to the fluid gelling tendency. However, under a shearing condition ( > 0 m/s) the gel structure is broken, and the bubble migration becomes dependent on . The relationship between the and under countercurrent flow condition can be fitted in a linear function (red dash line), as suggested in Equation (2). The second constant (0.172 m/s) is normally interpreted as the drift velocity. In this case, however, the strong gelling nature of the fluid invalidates this interpretation.

Conclusions
The movement of single Taylor bubble rising in stagnant and downward flowing non-Newtonian fluids in inclined pipes has been experimentally studied. The following conclusions can be drawn:

Conclusions
The movement of single Taylor bubble rising in stagnant and downward flowing non-Newtonian fluids in inclined pipes has been experimentally studied. The following conclusions can be drawn:

1.
A flow rig has been built for the purpose of investigating the movement of a single Taylor bubble rising in stagnant and downward flowing bentonite-water mixture in inclined pipes (4 • < θ < 45 • from vertical).

2.
The measured drift velocity (U d ) data presented in this study can contribute to improve the knowledge of a single Taylor bubble rising in stagnant non-Newtonian fluids having a yield point.

3.
The experimental results indicate that for all fluids tested, the bubble velocity increased as the inclination angles (θ) were increased, while velocity decreased with an increase in plastic viscosity (µ L ) and yield point (τ y ). The length of the Taylor bubble (L TB ) decreased as the downward flowing liquid velocity (U L ) and τ y /µ p were increased. The bubble velocity was found to be independent of the L TB .

4.
A reduction in drift velocity along its migration was detected, because increasing the concentration of bentonite results in a strong gel strength buildup. Eventually, a stagnant Taylor bubble can occur.

5.
The drift velocity is formulated by combining a function of θ and a yield-powerlaw correlation in terms of apparent viscosity (µ app ), making it applicable for non-Newtonian fluids with the Eötvös number (E 0 ) ranging from 3212 to 3405 and apparent viscosity (µ app ) ranging from 0.001 Pa·s to 129 Pa·s. The proposed correlation shows good performance for predicting drift velocity from both the present study (mean absolute relative difference is 0.0702) and a new database (mean absolute relative difference is 0.09614).
Author Contributions: Y.L. performed the experiments, analyzed the data, and wrote the paper; E.R.U. conceived, designed, and performed the experiments, and edited the manuscript; E.M.O. designed the experiments, led the project, and advised on the whole process of the manuscript preparation. Each author has contributed to the present paper. All authors have read and agreed to the published version of the manuscript.
Funding: This research was supported by the Chevron Corporation.

Institutional Review Board Statement: Not applicable
Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.

Data Availability Statement:
The data presented in this study are available in Appendix A.

Conflicts of Interest:
The authors declare no conflict of interest. Table A1 summarizes the experimental data for single Taylor bubbles rising in stagnant water and Bingham plastic fluids in inclined pipe, to provide traceability from the experimental measurements to the final conclusions. Appendix B Figure A1 displays the measured shear-stress/shear-rate relationship for the fluids used in this study. All data were measured using a Fann 35 rotational viscometer at rotational speeds of 3, 6, 100, 200, 300 and 600 rpm. The dotted line y-axis intercept for each fluid represents its theoretical yield point, τ y , for a Bingham plastic fluid. The curves for 30, 45 and 48.5 lb/bbl BW mixtures, at initial inspection, suggest that τ y does not actually represent a physical characteristic of these fluids. This intuitive conclusion, however, would not be correct. Under the dynamic conditions of a rotational viscometer, the low-shear-rate behavior shown below is common for a bentonite-water mixture. However, under static conditions, bentonite particles quickly build a gel structure that results in an actual τ y reasonably approximated by the theoretical Bingham plastic τ y . Thus, the yield point of a fluid is certainly a meaningful physical characteristic for understanding Taylor bubble behavior under static fluid conditions, and the shear stresses observed at low shear rates are material to dynamic fluid conditions. Considering that both static and dynamic fluid conditions play a role in this analysis, our method for estimating apparent viscosity based on a Bingham plastic model is a reasonable tool for use in our correlation. Furthermore, our preference for using the theoretical Bingham Plastic τ y instead of a measured τ y in our analysis stems from the difficulty of obtaining consistent ultra-low-shear measurements for high-viscosity fluids. The theoretical Bingham Plastic τ y is an easier-to-use and more consistent proxy for direct τ y laboratory measurements.

Appendix A
fluid conditions play a role in this analysis, our method for estimating apparent viscosit based on a Bingham plastic model is a reasonable tool for use in our correlation. Furthe more, our preference for using the theoretical Bingham Plastic y instead of a measured in our analysis stems from the difficulty of obtaining consistent ultra-low-shear measure ments for high-viscosity fluids. The theoretical Bingham Plastic y is an easier-to-use an more consistent proxy for direct y laboratory measurements.