Two-phase flow and morphology of the gas–liquid interface for bubbles or droplets in different microchannels

A_s

Generalized asymmetry factor

Ca

Capillary number

D

Diameter of the channel, m

K

Constant system parameter

L

Given length in the flow direction, m

L_b

Bubble length, m

L₀

The critical droplet length, m

P₁

Outlet pressures of the branch channels, Pa

P₂

Outlet pressures of the branch channels, Pa

Q

Flow rate in the main channel, m³ s⁻¹

R

Radius of the channel, m

Re

Reynolds number

t

Time, s

t_c

Critical time, s

$U ¯$

The mean velocity of droplet, m s⁻¹

U_b

The velocity of bubble, m s⁻¹

$U ¯ x$

The mean velocity of solute, m s⁻¹

V_total

Total droplet volume, m³

V_1fianl

Sub-droplet volumes in the branch channels, m³

V_2fianl

Sub-droplet volumes in the branch channels, m³

w

Width of the microchannel, m

w_b

Width of the branch channel, m

We

Weber number

3D

Three dimensional

Greek symbols

α

Channel aspect ratio

β

The eccentricity ratio

δ

Minimum bubble neck thickness, m

δ₀

Droplet neck thickness, m

δ_c

Critical bubble neck thickness, m

ε

Bubble liquid film thickness, m

ε_H

Bubble liquid film thickness on the side wall, m

ε_V

Bubble liquid film thickness on the bottom and top walls, m

μ

Liquid viscosity, Pa s

ρ

Liquid density, kg m⁻³

σ

Molecular diffusion coefficient

ω

Molecular diffusion coefficient

As a typical gas–liquid two-phase fluid, the bubble flow is used in oil and gas field development,^1–3 contaminated soil remediation,⁴ and medical fields.^5–7 In the oil and gas field development, a foam fluid, formed by the packing bubbles, is widely used as the fracturing fluid in the development of low-permeability oil and gas reservoirs. It has the characteristics of less formation damage, lower filtration loss, and stronger prop-carrying capacity. The foam fracturing and the prop-carrying foam fracturing fluids flow in porous medium of the stratum fracture. In the field of soil remediation, some pollutants (toxic pesticides, fertilizer particles, and heavy metals) are retained in the formation pores, due to the capillary force. Compared with the traditional surfactant flushing, the potential advantages of the foam flushing include less injected liquid, more uniform contact between the flushing agent and the pollutants, and stronger ability to adsorb the pollutants. Thus, the foam flushing improves the purification efficiency of the soil.⁸ In addition, the foam can also carry some nutrient particles to repair the soil.⁹ In the medical field, the study of confined bubble hydrodynamics can lay theoretical basis for the gas embolism therapy,^5,7 microcirculation and oxygen transport in blood vessels,^10,11 and targeted microbubble drug delivery.^12,13 A few studies have shown that the microbubbles in the gas embolism therapy keep in blood vessels with a diameter of less than or equal to 20 μm to block the flow of the red blood cells, thereby causing tumor cell necrosis.⁷ However, due to the complexities of the blood vessels and the structural evolution of bubbles, this technology has not been extensively promoted and effectively applied in the medical field. In addition, the fatal bubble embolism in the vascular system, as well as the bubble movement in the bronchus during the treatment of pulmonary infection, is also related to the bubble flow in complex microchannels.

In microchannels, the bubble flow involves film formation,^14–17 breakup,^17–20 and merging.^21–23 In a straight capillary, the flow pressure drop increases significantly, which causes different flow characteristics in the microchannel, when compared with the bubble flow in a conventional channel. For example, the occurrence of steam core countercurrent can lead to the increases of system resistance and inlet pressure.²⁴ Bubble transport in straight capillaries has been extensively studied, and the relationship between the thin film and the capillary number is generally considered. At present, T-junction,^25–27 Y-junction,^28,29 and their combination are usually used to generate specific bubbles. In general, bubbles (or droplets) break into the two branches at the bifurcation of the channel, forming a front cover at each branch. Under the influence of the viscous shear stress of the continuous phase, the two parts flow downstream and the bubbles (or droplets) are stretched. When the bubble (or droplet) reaches the maximum extension, it will reduce its surface energy by breaking into two sub-bubbles (sub-droplets) or retracting the whole bubble (or droplet) into one of the branch channels.^30,31

However, the research achievements of bubble dynamics at the bifurcation are still insufficient, and the unified conclusions are not built on the quantitative laws of bubble splitting. In the experimental study of bubble dynamics in a microfluidic T-junction with the square section, the bubble breakup and non-breakup behaviors were observed by Fu et al.²⁵ Their results showed that the bubble breakup would occur in the bifurcated channel when the bubble length reached a critical value. At the same time, the authors proposed a model to predict the bubble length. Some scholars have also studied the shape, size, and formation mechanism of the bubbles in the Y-shaped microchannels with different bifurcation angles.^32–34 In a bubble flow experiment²⁹ for studying continuous liquid flow in a Y-junction microchannel, the researchers observed that the bubble length increased with the increase in the gas–liquid velocity ratio, and there were two behaviors at the bifurcation: breakup and non-breakup. Furthermore, for the bubbles with a greater length than the channel diameter, they observed the breakup behavior when Re < 20. In practice, when the bubbles split at the bifurcation, the sub-bubbles entering the two branch channels may be same or different. The homogeneity of bubble breakup increases with the increase in the bubble driving pressure and decreases with increasing bifurcation angle.⁵ 

At present, many experimental studies are concerned with the symmetric bubble transport and breakup behavior in microchannels. However, the asymmetric behavior of bubbles in microchannel is also common in practical applications. Therefore, some scholars^35,36 have studied the asymmetric behavior of bubbles at Y-junctions. For example, when bubbles move in a symmetric Y-junction, the pressure gradients in the upper, the lower, and the middle parts of the inner wall are similar, while the bubble movement in an asymmetric model will result in a larger pressure gradient in the middle and the upper parts. Qamar et al.³⁷ also studied the bubble breakup in a Y-junction and observed the asymmetric behavior when the pressure at the each outlet of the channel was different.

Studying the evolution mechanism of a flowing Taylor bubble in a obstructed channel can be used for understanding the applications in chemistry, biology, and engineering fields, such as drug delivery,^38–40 cell focus,⁴¹ and separation of tumor cells.⁴² Nowadays, extensive research types have been conducted on the structural evolution of bubbles when encountering obstacles of different shapes, including triangles, squares, and circles. Tao et al.⁴³ found that the triangular barrier in the T-junction microchannel could induce to transfer the cavities generated by the liquid asynchrony to a stable bubble on the target wall. Furthermore, as for the bubble breakup, one can obtain different length ratios of sub-bubbles by adjusting the height and position of the obstacle when setting a square obstacle at the T-junction bifurcation.⁴⁴ When the square obstacle is located in the center of the straight pipe, the Taylor bubble will symmetrically break into two parts. However, when the obstacle leans to one side of the channel, the bubble will move to the other side in an non-breakup way.⁴⁵ In the porous medium channel of cylindrical obstacles, Baudouin et al.⁴ quantitatively expressed the relationship between the geometry of the medium and the bubble flow. They found that the intermittency of bubble flow was related to the evolution of the preferred flow path geometry, and the initial average size, as well as the bubble number, had a significant influence on the breakup efficiency.

Experimental studies²⁵ have shown that the evolution laws of other two-phase interface fluid, such as droplets, can be also applicable to bubbles, and to obtain a relatively complete research progress, we combine the related research results. Therefore, in this review, we also involve the related evolution behavior of the droplet flow. In Secs. II–IV, we will introduce the morphological evolution and flow mechanics of the two-phase interface fluid in the single capillary, the bifurcated capillary, and the obstructed channels. At last, we will summarize the status and provide suggestions for the potential development of the two-phase flow.

Taylor bubbles passing through a horizontal single capillary have been widely studied.^46–48 In the microchannel, bubble will adjust its shape to meet the restrictions.^49,50 The width of the lubricating film between the bubble and the capillary wall changes with the flow velocity, and the mechanism has been studied.^51,52 Concurrently, the mass transfer mechanism between the liquid and the bubble interface in the microchannel has also been studied.^53–56 Taylor^57,58 discussed the dispersion of soluble substances in a viscous liquid when flowing in a circular pipe under the condition of laminar flow. He found that the solute was dispersed with an apparent diffusion coefficient $R 2 U ¯ x 2 / 48 ω$⁠, where R, $U ¯ x$⁠, and $ω$ were the radius of the channel, the average velocity, and the molecular diffusion coefficient, respectively. However, this coefficient was valid when 4 L/R >  $U ¯ x R / ω$ > 6.9, where L was a given length in the flow direction. In 1956, Aris⁵⁹ found that the growth of variance for the solute distribution was proportional to the sum of the molecular diffusion coefficient and the Taylor diffusion coefficient, which removed the restriction. In 1972, Gupta et al.⁶⁰ discovered that the flow of the passive scalars in the circular pipes could enhance the axial diffusion. In 2005, Ng⁶¹ extended the Aris–Taylor dispersion theory to the situation that the solid wall absorbed the solute. During this period, although the relevant research was developed, there was still no systematic method that could fully summarize the mass transfer mechanism at the gas–liquid interface.^62,63 At present, Picchi et al.⁶⁴ extended the Aris–Taylor dispersion theory, revealing the influencing mechanism of the passive scalar diffusion around Taylor bubbles. In the future, we can further obtain the analytical expression of the relationship among the effective velocity, the diffusion, and the mass transfer coefficients in the advection–diffusion–mass transfer equation to accurately understand the dominant transport system of bubbles.

In the related study, the measurement of the liquid film thickness around the bubble is important for understanding the hydrodynamics of the bubble flow in microchannels. The film thickness plays a significant role in the analysis of the heat and mass transfer. The relevant formulas of the liquid film thickness are summarized in Table I. Fairbrother and Stubbs⁶⁵ first realized that the velocity of slender bubbles was different from the average velocity of the fluid in channels. They observed that the wetting film would be deposited on the microchannel wall when the wetting viscous liquid was replaced by bubbles. As shown in Table I, in their experiment, the liquid film thickness was proportional to Ca^1/2. Taylor⁶⁶ further confirmed the Fairbrother and Stubbs' conclusion by measuring the fraction of liquid remained on the channel wall when the air was discharged. At the same time, as Ca was close to two, the dimensionless liquid film thickness ε/R would approach 1/3, where ε was the bubble liquid film thickness and R was the radius of capillary. Marchessault et al.⁶⁷ obtained the thickness model in the range Ca = 7 × 10⁻⁶–2 × 10⁻⁴ by studying the two-phase flow of glycerol aqueous solution, water, and air. Their results are slightly different from the previous models.

Subsequently, Bretherton⁶⁸ assumed the liquid film peristalsis and calculated the thickness of the flat film of the spherical bubble by using the lubrication theory. It was found that the thickness of the middle part of bubble was constant, and the liquid film was proportional to Ca^2/3. However, this relation is only valid at a low Ca. When Ca < 5 × 10⁻³, the calculation accuracy is within 10%. If Ca is close to 2 × 10⁻⁴, there is a smaller error between the experimental thickness and the theoretical thickness. Aussillous and Quéré⁷⁰ proposed an empirical “repair” and obtained the generalized formula for calculating the thickness of liquid film in the middle part, which could be applicable to the situation for Ca of 10⁻³–1.4. Afterward, Evert Klaseboer et al.⁷⁴ carried out the theoretical derivation and optimized the “empirical model.” They obtained the relevant formula for measuring the liquid film thickness by analyzing the extended version of Bretherton's Taylor bubble model. This formula is similar with the correlation obtained by Aussillous and Quéré.⁷⁰ 

In order to measure the thickness of liquid film in vertical microchannels, Irandoust and Andersson⁶⁹ proposed an empirical equation, which was quite different from the relevant formulas. The thickness estimated by the new empirical formula was much larger than that proposed by the previous scholars.^65–68 An experimental study on the two-phase flow of water, FC-40, tetradecane, and dodecane was conducted by Grimes et al.⁷² They summarized the relevant formula for measuring the film thickness around the water droplets. They found that the Bretherton's theoretical correlation, in which Ca was defined by the carrier viscosity, underestimated the liquid film thickness, so that they proposed to define Ca by the liquid droplet viscosity. Han and Shikazono⁷³ studied the two-phase flow of air, FC-40, ethanol, and water in a wide range of Re based on a few experiments, and obtained two correlations for measuring the liquid film thickness. Their studies focused on the horizontal microchannels and showed that the film thickness was a function of Re, Ca, and We and remained constant for Re > 2000. Ni et al.⁷⁵ used the commercial software package, COMSOL Multiphysics, to derive a correlation formula that was same with that of Han and Shikazono.⁷³ Etminan et al.⁷⁶ also conducted a few studies on the thickness of liquid film. They showed that the thickness of liquid film was almost unchanged while the length of the flat film area increased with the flow of bubbles.

Rectangular microchannels are widely used in engineering fields, so the studies on the distribution characteristics of bubble liquid film in rectangular microchannels are also of significance. In Matin and Moghaddam's study on the liquid film thickness of bubbles in rectangular channels, they found that their liquid film distribution was different from that in circular channels.⁷⁷ As shown in Fig. 1, since circular cross sections are isotropically symmetric, the liquid forms a uniformly thick film around the bubble. In contrast, square and rectangular sections are not isotropically symmetric, resulting a non-uniform film of liquid around the bubbles. At the corners of the square and rectangular sections, the thickness of the liquid film is greater, while at the four sides and the center, the film thickness is less. This difference is mainly caused by the different pressures generated by the fluid flow in the different cross sections. In circular cross sections, the fluid flow and circumferential pressure distribution are uniform, and thus, the liquid film distribution is also uniform and consistent. In square and rectangular cross sections, as the fluid flow is slower in the corners, the liquid film accumulates more in these areas and forms a thicker liquid film, whereas in the sides and the center part, the fluid flow is faster and the liquid film is thinner.^78–80 The prediction of liquid film thickness in Table I is applicable to the circular channels, as shown in Fig. 1(a). In a square channel, the cross section of bubble is still symmetrical on the diagonal and the liquid film thickness is the same near each of the walls, despite the difficulty in determining the liquid film distribution at the corners, as shown in Fig. 1(b). However, in a rectangular channel, the liquid film thickness varies in the direction of the channel axis in Figs. 1(c) and 1(d), and it also depends on the degree of constraint (i.e., channel aspect ratio α). The confinement in the rectangular channel squeezes the bubble into a non-axisymmetric cross-sectional shape, while the surface tension restores the gas–liquid interface in an axisymmetric state. Comparing Fig. 1(c) with Fig. 1(d), we can get a result that increasing α leads to increasing ε_H, as well as decreasing ε_V, where ε_H is the liquid film thickness on the sidewall of channel and ε_V is the liquid film thickness on the bottom and the top walls. Due to the deposition of liquid in corner of the section, the central liquid film thickness in the rectangular channel is thinner than that in the circular channel,^81,82 when the bubble size is the same. Bartkus and Kuznetsov⁸³ used the LIF method to study the liquid film distribution of slender bubbles in the rectangular microchannels. Their study showed that the calculated liquid film thickness was in good agreement with Taylor's law when Ca > 0.05.

In practice, besides the liquid film thickness around the bubble, there is more information about the bubble in the single capillary, which is needed to be explored. For instance, the gas–liquid interface equation of the bubble is still lacking. Furthermore, smaller bubbles whose size is lower than the channel diameter show distorted morphologies due to the gravitational force and other constraints. Their detailed gas–liquid interface morphologies for the related two-phase flow are also expected to be further present.

Y-junction microchannel is a common channel structure in porous medium of the formation. Studying the flow and breakup of bubbles in a Y-junction microchannel will help to accurately control the bubble parameters and increase the production efficiency of crude oil in reservoirs. The flow and breakup of bubbles are closely related to the microchannel structures and the physical properties of the fluid. They are subject to the combined actions of viscous shear force, interfacial tension, and inertial force. At a low continuous-phase velocity, bubble breakup is mainly acted by a squeezing mechanism, but at a high continuous-phase velocity, it is a shearing mechanism.¹⁷ In addition, a few of scholars have conducted studies to get the breakup mechanisms of droplets (bubbles). Droplet breakup in microchannels can be generally divided into three states, including non-breakup, tunnel breakup, and obstruction breakup.^17,84–87 The difference between the tunnel breakup and the obstruction breakup is in a visible gap between the bubble (droplet) and the channel wall during the breakup process.^88,89 Some scholars have also classified the flow of bubbles in Y-junction microchannels into three types: homogeneous breakup, non-homogeneous breakup, and non-breakup.³⁷ Nagargoje and Gupta²⁹ conducted a detailed study on these three types, as shown in Fig. 2. In the non-breakup state, the bubbles expanded and deformed due to the inertia of the continuous phase, but did not split near the bifurcation, as shown in Fig. 2(a). In the breakup state, when the bubble fully reached the bifurcation, it was compressed and necked at the bifurcation due to the inertia of the liquid phase, causing bubble breakup, as shown in Figs. 2(b) and 2(c).

When the breakup of the bubble occurs at the Y-junction, the process can be divided into three stages: extension, extrusion, and breakup stages.⁸⁵ As shown in Fig. 3, the process of bubble breakup occurring in the Y-junction microchannel can be further divided into three stages: extension stage, extrusion stage, and breakup stage. The flow of bubbles from the main channel to the bifurcation is the expansion phase where the bubbles deform and block the channels. The deformation of the bubble at the bifurcation is due to the pressure of the upstream fluid. In the extrusion stage, the bubble head will divide into two parts and enter the bifurcation channel separately, while the curvature of the tail will become larger. Under the action of the upstream fluid, the bubble expands to the two bifurcation channels, respectively, and the bubble tail's curvature keeps increasing until breakup, which is the breakup stage. After breakup, the bubble tail will be rounded again under the action of interfacial tension. In addition, if the bubbles do not break, bubble stretching and squeezing will occur at the bifurcation. The extension stage is similar with that when the breakup occurs. The difference is mainly in the extrusion stage, where the volume of the bubble entering the two branch microchannels is different. Because of the interfacial tension, the bubble will be automatically pulled into the larger side.⁹⁰ Pan et al.⁸⁵ observed that bubble breakup at Y-junctions could be categorized as tunnel–tunnel breakup, obstruction–obstruction breakup, tunnel–obstruction breakup, and non-breakup based on the presence or absence of gaps between the sub-bubbles and the walls of the microchannels. As shown in Fig. 4(a), in the tunnel–tunnel breakup, there are gaps between the bubbles and the walls, which are mainly affected by the squeezing and the shear force of liquid. Since the liquid passes through the tunnel with a large velocity gradient, the shear force is more obvious in tunnels. As shown in Fig. 4(b), in the obstruction–obstruction breakup, there is not gap between the bubbles and the walls, which is mainly affected by the squeezing force of liquid. The shear force only appears at the ends of the contact interface between the bubbles and walls, so there is not a shear force in this case. As shown in Fig. 4(c), in the tunnel–obstruction breakup, one sub-bubble has a gap with the wall, while the other sub-bubble does not. The former is subjected to the squeezing force and the shear, while the latter is only subjected to the squeezing force. As shown in Fig. 4(d), in the non-breakup, the large bubble drags the small bubble out of the upper channel and into the lower channel under the action of surface tension. At this situation, the surface tension is greater than the shear force and the squeezing force.

The uniformity of bubble breakup is also a study theme. The length of two sub-bubbles formed by the bubble breakup increases with the increases of the gas flow rate and the bubble length, but decreases with the increases of the liquid flow rate and the liquid viscosity. The asymmetry of bubble breakup decreases with the increases of the liquid velocity and the liquid viscosity.⁹⁴ At the same time, the uniformity of breakup increases with the rise of driving pressure and decreases with the increase in the bifurcation angle of a Y-junction.⁵ Some scholars found that the phase breakup at the bifurcation depended on the flow pattern at the entrance.⁹⁵ For a plug flow, the gas phase preferentially flows into the branch. In an annular flow, the liquid preferentially flows into the branch. Xu et al.⁹⁶ experimentally observed the effects of the channel bifurcation angle (30°–150°), the capillary number, and the bubble length on the bubble breakup. They found that the curve of bubble neck was similar to a parabola rather than an arc in the extrusion process, and proposed a quadrilateral model, which could reflect the change of gas–liquid interface during the bubble breakup. However, in their study, the neck curve of the bubble did not follow the parabola model in the fast pinch stage.

At present, although the bubble breakup transition at the bifurcation of the symmetric Y-junction has been widely studied, the research about the breakup mechanism of the two-phase interface fluid in an asymmetric Y-junction is still insufficient, especially for the droplet breakup. For asymmetric bifurcated microchannels with different branch angles, some researchers have observed three breakup states during the experimental processes, and the distribution is approximatively a function of the capillary number and the length of droplet.⁹⁷ Many scholars have found that once the droplet length is greater than the circumference of the channel, the droplet will break at the T-junction bifurcation. Similar correlations have been also demonstrated between the critical capillary number and the droplet length in other bifurcation structures.^30,98,99 Wang et al.¹⁰⁰ further proved that the length and the velocity of droplets greatly affected the breakup, and the asymmetry of the droplet breakup was positively related to the length and the velocity. They also found that, at Ca < 0.01, the critical droplet length L₀/w ≈ 2.5, and the increase in Ca would also cause increasing critical length. From the above studies, it can be seen that the length and the velocity of droplet play a key role in the breakup. In addition, some scholars¹⁰¹ also found that the vortex of droplet would inhibit its deformation and breakup, and they obtained the change rules of the droplet neck thickness (δ₀/w) at different stages. As shown in Fig. 5, with the droplet moving downstream, the two emerging front caps penetrate into the branch, while the droplet neck keeps contracting. The droplet breakup process can be divided into three stages by two critical conditions. Initially, the droplet neck thickness is larger than the channel height, and the droplet neck contacts the top and bottom walls, with only a thin film between the droplet and the two walls. Thus, the droplet neck thickness decreases in a two-dimensional fashion, which is defined as the squeezing phase. The first critical state occurs when the droplet neck thickness equals the channel height (176 ms). After that, the droplet neck detaches from the top and bottom walls and starts to contract in a direction perpendicular to the top view. The second stage of the droplet neck contraction as a three-dimensional shrinkage is called the transition stage. When the rear cap of the droplet reverses its normal direction (220 ms), the last phase is entered. In the pinch-off phase, the back cap changes from convex to concave and the thread breaks quickly, with a duration of less than 8 ms in Fig. 5.

In the study of droplet breakup mechanism in a unilateral Y-junction,³⁴ the critical condition of droplet breakup could be expressed as L₀/w_b = aCa^−0.21, and the pre-factor a changes with the bifurcation angle and the viscosity ratio, where w_b is the width of the branches. When the bifurcation angle is reduced or the viscosity ratio is increased, the cutting effect of sharp corner or the vortex-related resistance will reduce and the droplets are easier to break.

Based on the application of bubble or droplet microfluidics in chemical engineering, medicine, and materials synthesis, it is possible to achieve the control of bubble size. Currently, the bubble size can be successfully tailored by a specific T-junction, which has become one of the most fundamental structures for customizing bubbles. It is no doubt that the behaviors of bubbles (or droplets) at T-junctions are highly non-linear.¹⁰² The bubble (or droplet) breakup at the bifurcation mainly depends on their sizes, the capillary number, and the viscosity.^88,89,103 We can use the breakup of bubbles (or droplets) in a symmetric T-junction to increase the productivity. An asymmetric T-junction can be used to manipulate the sizes of bubbles or droplets. Link et al.⁹¹ used a T-junction to break a droplet into two sub-droplets. They found that the flow of droplets could be adjusted to induce symmetrical or asymmetric breakup. In addition, an analytical model was proposed by Garstecki et al.¹⁰⁴ to explain the transition between the breakup and the non-breakup in terms of the stretching of droplet in an extensional flow, based on the classical Rayleigh–Plateau instability.

As mentioned in Sec. III A, it has been confirmed that the capillary number and the droplet length play a crucial role in the breakup of droplets in a bifurcated structure. When Ca is between 4 × 10⁻⁴ and 0.2, Jullien et al.¹⁰⁵ observed the tunnel breakup and the obstruction breakup at the T-junction bifurcation. They found that the critical length of droplet breakup had no relationship with the droplet flow, but changed with the two-phase viscosity ratio. When the droplet length exceeded the critical length, it would present the obstruction breakup. In addition, Ménétrier-Deremble et al.⁹⁷ observed two different flow patterns for droplet breakup: “direct breakup” and “retarded breakup,” and they concluded a critical length controlling the transition between the two flow patterns.

As for the bubble breakup at the T-junction, a few of scholars have conducted extensive studies on the dynamic structural evolution of bubble. It has been found that the dynamic characteristics of bubble breakup at T-junctions are related to the constraint of microchannels and the capillary number.^25,104 Samie et al.³⁰ obtained four flow patterns: non-breakup, non-tunnel breakup, tunnel breakup, and unstable breakup at a series of T-junction bifurcations with different ratios of the branch width. They investigated the quantitative relationship between the distribution of droplet volume and the branch width ratio to obtain a prediction formula for the droplet breakup in T-junction microchannels. Similarly, Fu et al.²⁵ used high-speed photography and microscopic particle image velocimetry techniques to study the breakup of nitrogen bubble in glycerol aqueous solution at the symmetric T-junction bifurcation when Ca was within the range of 0.001–0.1. Under the action of surfactants, by changing the flow rates of the gas and the liquid, they discovered four different patterns of bubble behaviors: symmetric breakup of bubbles type I (BSI), symmetric breakup of bubbles type II (BSII), symmetric breakup of bubbles type III (BSIII), and non-breakup of bubbles (NB). As shown in Fig. 6, the BSI bubble is mainly controlled by the pressure in the liquid phase, and the BSII bubble is controlled by the pressure, as well as the viscous force. In the BSIII bubble, there is a scaling law of the minimum bubble neck in the bubble breakup process. In the study of Fu et al.,²⁵ it was found that the bubbles eventually broke up when the bubble neck reached a critical thickness δ_c (δ_c/w =0.32). The existence of the critical width was found both experimentally and numerically for droplets breakup in T-junctions. The critical width, δ_c/w =0.3, was obtained experimentally¹⁰⁵ for the droplet breakup and δ_c/w = 0.5 obtained theoretically⁹² from a 2D approach, which suggested that the 2D theory could be approximatively used to explain the present experiments. The deviation between these values might stem from the difference in the two-phase fluid systems, as well as in the size of microchannels. The experimental results of Fu et al.²⁵ also proved the mechanism of pressurization to explain the 2D breakup dynamics at the T-junction. Furthermore, some scholars¹⁰⁶ have also studied the breakup mechanism of droplets at the T-junction in 3D, and their conclusion is consistent with that of the 2D approach.

In a subsequent study, Afkhami¹⁰⁷ found that a minimum film thickness, whose value was close to 1.08Ca^2/5, existed at the droplet breakup at a low capillary number. Lu et al.¹⁰⁸ found that the thinning of the bubble neck experienced a slow-breakup stage and a fast-breakup stage by examining the evolution of the minimum thickness δ of the bubble neck with the remaining time t before the final pinch-off for bubble breakup as shown in Fig. 7. The rate of the bubble breakup first decreases and then increases when the width of the bubble neck is less than a critical value of δ_c/w = 0.59, where δ_c is the critical neck thickness for the breaking bubble between the slow-breakup and the fast-breakup stages. When the width of the neck is larger than the critical value δ_c, the breakup is driven by the continuous phase and this stage is the slow-breakup stage, and when the width of the neck is less than the critical value δ_c, the breakup is driven by the surface tension and this stage is the rapid-breakup stage. As shown in Fig. 7(a), when the pump stops operating at some time during the bubble breaking process, the bubble neck would return to the original state during the slow-breakup stage or continue to break in the case of the fast-breakup stage. This shows that the slow-breakup is reversible and the fast-breakup is irreversible. As shown in Fig. 7(b), the minimum width of the bubble neck shows a power-law link with the time, and the index is 0.22 at the beginning and 0.5 when it is close to the final clip.¹⁰⁸ 

Meanwhile, some scholars used the multiphase Lattice Boltzmann model to explore the effects of pressure, shear, and viscosity on the mechanism of droplet (or bubble) breakup in T-junctions. For example, Chen et al.⁸⁹ found that the tunnel effect could reduce the deformation rate of droplets and even induce the non-breakup phenomenon. From this study, it can be obtained that the shear action in the tunnel has little contribution to the elongation of droplet in the tunnel flow pattern. The vortexes in droplets play an important role in the breakup. Especially, when the intensity of the vortexes exceeds the critical value, the droplet will not break. At the equal shear action of a continuous flow, a larger droplet viscosity can limit the generation of vortex, which is easier to induce the droplet to break. For the bubble breakup, by means of a numerical simulation, Zhang et al.¹⁰⁹ found that a lower gas velocity, smaller bubble size, and uneven pressure distribution in the branch channels were the main factors causing non-breakup behavior of the bubbles.

The sizes of the sub-bubbles can be controlled by adjusting the resistance of branch channels. The resistance includes the inherent one of channels and the one caused by bubbles. At the same time, bubbles are also affected by the bubble collision at the junction.^30,112–114 Therefore, the resistance not only depends on the liquid viscosity, the cross section, and the length of microchannels,¹¹⁵ but also is decided by the number and the size of bubbles. The flow resistance increases with the increase in the bubble number in microchannels. In addition, the quantitative additional resistance of a bubble also relies on its size. The resistance of a long bubble is related to the capillary number and can be expressed by the Bretherton equation⁶⁸ and the Ratulowski–Chang equation.¹¹⁶ The resistance of a short bubble can be expressed as a characteristic length, but this has not been completely documented in previous literatures, which is needed to be further studied. When using a T-junction to predict the sizes of bubbles, we are also supposed to consider the resistance caused by bubbles.¹¹³ 

At present, a lot of efforts have been invested in the research of bubble and droplet breakup mechanism in symmetric T-junctions, but the understanding in asymmetric T-junctions is still insufficient. The studies on the breakup mechanism and the flow state in asymmetric T-junctions are of significance, since droplets and bubbles mostly flow in asymmetric structures in practical applications. There are four different flow patterns in an asymmetric T-junction: permanent obstruction breakup, partly obstruction breakup, permanent tunnel breakup, and non-breakup. Compared with symmetric T-junctions, asymmetric T-junctions are more possible to cause bubble breakup.¹⁰³ The bubble breakup process in asymmetric T-junctions can be divided into three stages: the extrusion, the transition, and the pinching. Fluid surface velocity is one of the key factors affecting the thinning rate of bubble neck and decides the bubble breakup time. In addition, the increase in liquid viscosity will accelerate the thinning of bubble neck in extrusion and transition stages, but it is opposite in the pinching stage.¹¹⁷ Liu et al.¹¹⁸ found that there was a sharp decrease in the pressure drop at the gas–liquid interface in the extrusion stage, and the pressure drop at the front interface area was far greater than that at the depression area.

Some scholars tried to achieve a deeper understanding about the droplet breakup in a T-junction with asymmetric branches through the analysis and calculation.¹¹⁹ Bedram and Moosavi¹²⁰ studied the sub-droplets produced by T-junctions with different branch channel widths, and they found that a larger capillary number could cause wider range for the droplet size. Samie et al.³⁰ found that increasing ratio of the branch channel width, as well as decreasing capillary number, would reduce the length of the sub-droplets. Therefore, by applying a high width ratio and a small capillary number, we can produce the sub-droplets with the similar size, which improves the efficiency of droplet breakup. In the following study, Cheng et al.¹²¹ explored the influence of length ratio of the branch channels on droplet rupture from a 2D perspective, and they found that, due to the asymmetric channel structure, only through 3D simulation could one more accurately predict the influence of the length ratio of the branch channels. The flow characteristics and the breakup mechanism of bubbles in asymmetric structures have not been perfectly revealed, and a few of mechanism studies are needed to fill the content.¹¹³ 

Many scholars have studied the flow characteristics of bubbles in obstructed microchannels. Shu et al.¹²² found that when a square obstacle was located at the center of a microchannel, bubbles would split into two same sub-bubbles. In practice, as the obstacle is deflected toward one side of the channel wall, bubbles may not break into two parts. In addition, the bubble will not break at a lower flow rate. On the contrary, if the flow rate is higher, the bubble flowing through the obstacle will be squeezed by the channel wall and shear force of the liquid phase, which induces bubble to break. At the same time, a stagnation point is formed behind the narrowest part of channel, causing vortices in the liquid and promoting bubble breakup.^123,124

Not coincidentally, the breakup of liquid droplets in the obstructed microchannels has also been studied. Li et al.¹²⁵ studied the influence of a cylinder obstacle and obtained the same conclusion with Shu et al.,¹²² as shown in Fig. 8. As shown in Fig. 8(a), the droplet initially deforms and spreads radially shaped like a curved moon face suspended from the cylinder. Afterward, the curved moon face is pulled downward by gravity and fills the two gaps between the cylinder and the wall. The liquid film layer between the rear and the tip of the droplet thins while surface tension enhances its contraction to the sides. As in Fig. 8(b), the droplet breaks up into two asymmetric sub-droplets and the gap between the cylinder and the right wall decreases, resulting in greater pressure developing on the right side and making more liquid inclined to invade the left gap. A further increase in β to 2/7 showed that the left gap anterior meniscus moved faster and the right gap small meniscus pulled backward, as shown in Fig. 8(c). As a result, the droplet does not break and the unbroken droplet falls faster as it passes through the cylinder. As shown in the Fig. 9, when the obstacle lies in the center of channel, the β is 0 and the droplet breaks up into two symmetric sub-droplets. As the obstacle gets closer to the right wall, β increases, and the volume ratio of the left and right sub-droplets also increases, until the eccentricity ratio reaches a critical value, the droplet will not break. In addition, the viscosity ratio also has a great influence on the shape, breakup position, and liquid film thickness of the sub-droplets. When the viscosity is relatively large, the liquid film thickness of the deposition decreases quickly. At the same time, the wettability also plays an important role in the dynamic behavior of droplets. As shown in Fig. 8(d), when the droplet is wetting the cylinder, there are two bubbles between the flowing droplet and the cylinder, and there is more liquid deposited near the bottom surface of the cylinder. This may be due to the fact that wetting droplets tend to adhere closely to the surface and the thin layer of liquid could be easily pulled away by gravity. The authors¹²⁵ also found that the time for a droplet to pass through an obstacle decreased with increasing contact angle and the coalescence of daughter droplets occurred earlier than that of the non-wetting case. Kang et al.¹²⁶ found that for the case of a contact angle less than or equal to 90°, at a very small Bond number, the wet length between the droplet and the wall decreased with time until a steady shape was reached, while the Bond number exceeded a critical value, the droplets were pinched off or even entrained into the bulk. In the non-wetting case, for Bond number less than the critical value, the droplet shrunk along the wall from its static state until reaching the steady state; however, for Bond number above the critical value, the droplet completely detached from the wall. Chung et al.¹²⁷ observed that the droplet breakup also depended on the droplet size and the capillary number in the channel with cylindrical obstacles. With the increase in capillary number in their experiments, the small droplets, whose size were 1.3 times larger than the channel width, would break up, while the breakup did not occur for the longer droplets.

There are various obstacles in microchannel, and the evolution of droplets passing through the cylinder and the square obstacles is shown in Fig. 10. It can be seen that the square obstacle can cause more dead zones than the cylinder.¹²³ In addition, Protiere et al.⁹⁹ proved that as a droplet passed through the obstacle, there was a critical Ca below which the non-breakup could be achieved.

In practical applications, there are often more than one obstacle in a microchannel and the inclination of the obstacle is also random. Studying the influences of the obstacle number and the contact angle between the bubble and the obstacle on the bubble flow is of importance to perfect the understanding for structural evolution of bubble. Chen et al.¹²⁸ set two rectangular obstacles in a microchannel for studying the fluids of a low and a high viscosities, as shown in Fig. 11. They found that the bubble broke in the two narrow channels between the obstacle and the channel wall to form two gas streams. The two gas streams met at the convergence area (i.e., the front of the next obstacle), and one of the streams passed the next obstacle, resulting in retreating of the other gas stream. In their study, under the influence of liquid viscosity, the flow patterns of the bubbles could be divided into four types: total breakup, partial breakup, abnormal breakup, and annular flow patterns. The above flow patterns occurred in the different gas–liquid velocities. In the research of Lee and Son,¹²⁹ it was found that the droplet breakup depended on the inclination, the position, and the width and length of obstacles. Increasing the length of obstacles and tilting the obstacle would achieve symmetrical droplet breakup. Ma et al.¹³⁰ found that there was a critical length ratio between the droplet and the obstacle, when a droplet passed through a linear obstacle. The critical length ratio for the droplet breakup was caused by the competition between the enhancement (positive) and the hindrance (negative) factors. For the same obstacle, the droplet was apt to break at the critical length ratio.

The roughness of the wall surface has an effect on the velocity and pressure of the bubbles (droplets) as they flow and break in channels. Ansari and Nobes¹³¹ found that bubble pinning at the exit due to pore roughness created an additional drag. The pinning pressure applied at each pore roughness resulted in a positive pressure difference along the bubbles, which decelerated the bubbles.

Xiao et al.¹³² found that the difference between droplet and bubble wetting behavior and the effect of roughness varied by taking the critical Young contact angle as the turning point. Due to the coupled effect of solid–liquid affinity and rough microstructure, the roughness effect on the droplet contact angle hysteresis varies in different wetting states, while the increase in roughness decreases the apparent contact angle and contact angle hysteresis of bubbles. Marmur¹³³ theoretically investigated the influence of the size effect on the contact angle hysteresis of droplets and bubbles on heterogeneous smooth surfaces, and the results showed that the size effect had different effects on droplets and bubbles. Although the approximation of low-pitch solid surfaces or moderately heterogeneous solid surfaces allowed treating the effect of roughness and heterogeneity on contact angle in a similar way, Montes Ruiz-Cabello et al.¹³⁴ found that droplets and bubbles had different behaviors in contact with the rough surface. They¹³⁴ found that the curvature and surface characteristics were the main factors for the difference between the droplet and the bubble wetting behaviors by using sawtooth rough surface with a 3D droplet/bubble model. Zhu et al.¹³⁵ found that larger roughness would enhance the diffusion and adhesion of droplets on the hydrophilic surface, but the separation of bubbles from the hydrophilic surface would be easier. In addition, the droplet and bubble contact angle significantly decreased with increasing roughness. We might conjecture that the wall roughness would affect the wettability of bubbles (droplets) to the microchannel wall, which consequently affects their flow state and breakup mechanism.

A few researchers have studied the flow characteristics of the two-phase interface fluid, bubbles or droplets, in microchannels. The relevant achievements in this flied mainly focus on the fluid, and the analysis processes are insufficient. At the bubble scale, the numerical simulations of three-dimensional bubbles flowing at the bifurcation are lack of a detailed evolution of the two-phase interface when the structure of flowing bubble is changing. On the one hand, the studies on the bubble dynamics in microchannels lack quantitative relationships between the bubble breakup (or retention) and the relevant parameters (including bifurcation angle, branch channel sizes, and bubble parameters). On the other hand, there is still a lack of study on the dynamic interactive evolution of physical parameters, such as the bubble pressure at the bifurcation. Studying these relationships can provide an important reference for predicting and controlling the breakup or the retention behaviors of bubbles in multi-branch channels.

Compared with the symmetrical breakup, the asymmetric breakup can efficiently achieve the size reduction.¹³⁶ We can use T-junction/Y-junction channels with different widths,^30,113 lengths,¹¹³ angles,^101,137 and obstacles^125,138 to achieve asymmetric breakup of bubbles (droplets). Linnartz et al.¹³⁹ proposed a pathway for droplet formation by narrowing microfluidic porous-walled channels, and presented a method for producing uniform cylindrical droplets in porous-walled channels with internal cavity diameters down to 7 μm. Compared to the traditional method of producing droplets in bifurcated microchannels, this method could greatly improve the efficiency of droplet production. However, the existing experimental research types have not systematically included all of the beneficial aspects, and the numerical studies have not been effectively proved by experiments. If we apply the existing research results to design a practical application equipment, there is still an immense amount of work to be done. In order to design the optimal foam parameters in relevant applications, it is necessary to study dynamic structure evolution of bubble in microchannels. Based on the continuous studies on this topic, we can deeply understand the bubble application theories in oil and gas field development, soil remediation, and gas embolism treatment.^140–143 In addition, the two-phase interface flow in the microchannel is affected by the channel wall, and the liquid film is affected by the bubble shape. The flow rate around the Taylor bubble is very low, which limits the mass transfer between the bubble and the liquid phase. To overcome this defect, the mass transfer ability of liquid can be improved by changing the inner wall of microchannels.¹⁴⁴ It is also possible to design the devices with relevant structures to realize the quantitative control of the bubble generation and breakup. We can produce bubbles with specific parameters to realize the control of the fluid movement. The study on the structural dynamics of the gas–liquid interface in the branch channel can provide theoretical guidance for the new technology of tumor therapy, gas embolism therapy, and help understanding the deadly air embolism¹⁴⁵ of the vascular system, as well as the movement of bubbles in lung bronchus.¹⁴⁶ In the complex human vascular system, the arteries and capillaries of blood flow are flexible multi-branched systems. Understanding the blood flow with gas bubbles in flexible channels is essential for analyzing the associated pathologies and exploring cancer therapies. Different from ordinary rigid microchannels, the multiphase flow of bubbles in flexible complex channels involves the interaction between the bubble and the flexible wall, as well as the resulting associated structural deformation. However, most of the present studies focus on rigid channels, and the results cannot reflect the real bubble flow in flexible channels. At present, we lack a reasonable two-phase flow model of bubble flow in flexible microchannels and need to conduct an in-depth study on the bubble adhesion and transport behavior in three-dimensional flexible microchannels, which is more capable to show the real situation of blood vessels.

It should be noted that there are a few of differences in physical properties between bubbles and droplets, such as viscosity, density, and surface tension, which causes differences in the phase interfaces of bubbles and droplets in the breakup process.^5,147 The thickness of liquid film increases with decreasing surface tension of the droplet. Under the same conditions, the flow patterns and the breakup dynamics of these two-phase fluids in microchannels are also different.^27,108 The process of the droplet breakup is often accompanied by the thread breakup and the formation of satellite droplets, except for the similar early stages of the bubble breakup.¹⁴⁸ At the same time, the length of the thread and the size of the satellite droplet also increase as the droplet viscosity increases.^27,149,150 Moreover, it can be seen from the current experimental results that in microfluidic T-junctions, the process of the droplet breakup may be more complex than that of the bubble breakup.^25,103

Nowadays, the two-phase interface flow plays an important role in the development of human medicine, oil and gas field development, and soil remediation. Mastering the flow and evolution mechanism of bubbles (or droplets) in microchannels can accelerate the related technological application in the above fields.

This paper discusses the hydrodynamic characteristics of bubbles flowing in microchannels under various conditions. We have combed the experiments, the theories, and the simulations related to the liquid film thickness of bubble in circular channels, and summarized relevant formulas for predicting the film thickness. Moreover, we have sorted out the different breakup mechanisms and flow states of bubbles in the symmetric and asymmetric Y-junctions/T-junctions, and we also discussed the relevant factors affecting the bubble (or droplet) breakup in detail. It is also understood that the size of sub-bubbles (sub-droplets) can be controlled by adjusting the resistance of branch microchannels or using asymmetric T-junctions. Furthermore, we summarized the breakup mechanisms of bubbles (droplets) in microchannels with different obstacles. The breakup time, the shape, and the volume ratio of sub-bubbles (or sub-droplets) are affected by the eccentricity and the shape of obstacles. The existence of obstacles is conducive to produce the complex two-phase flow and promote the quantitative breakup.

In the existing studies, scholars have conducted extensive explorations on the basic situation of bubble (or droplet) flow. However, a gas–liquid two-phase flow within the flexible channel is still insufficiently investigated, lacking two-phase flow modeling of bubble flow in flexible microchannels, as well as the bubble flow and evolution in three-dimensional flexible microchannels, which show the real vascular situation. At the microscale, the prediction model of the gas–liquid film thickness and the breakup mechanism of the gas–liquid interface are still expected to be further studied from the three-dimensional perspective. In the relevant application, it is also necessary to continue to explore the breakup to supplement the regular model, as well as the multiple physical fields,^151,152 such as the membrane stress and the pressure fields in the multiphase fluid. In addition, we are supposed to systematically consider porous channels of different structures by combining with practical applications. With the development of micro-machining technology, we can design new micro-components based on the flow characteristics of bubbles (droplets), thereby achieving precise control on the evolution of the two-phase interface structure.

This work is supported by the Projects from National Natural Science Foundation of China (No. 51706173), the Natural Science Basic Research Program of Shaanxi Province (Nos. 2023-JC-QN-0439 and 2023-JC-QN-0029), China Postdoctoral Science Foundation (No. 2022M722556), and project (No. WJ2022G00002-001).

The authors have no conflicts to disclose.

Cheng Chen: Conceptualization (equal); Investigation (equal); Writing – original draft (equal); Writing – review & editing (equal). Zefeng Jing: Data curation (equal); Writing – original draft (equal); Writing – review & editing (equal). Chenchen Feng: Conceptualization (equal); Writing – original draft (equal); Writing – review & editing (equal). Xupeng Zou: Conceptualization (equal); Writing – original draft (equal). Mingzheng Qiao: Conceptualization (equal); Writing – original draft (equal). Donghai Xu: Data curation (equal). Shuzhong Wang: Data curation (equal).

Data sharing is not applicable to this article as no new data were created or analyzed in this study.