Viscous dynamics of drops and bubbles in Hele-Shaw cells: drainage, drag friction, coalescence, and bursting

In this review article, we discuss recent studies on drops and bubbles in Hele-Shaw cells, focusing on how scaling laws exhibit crossovers from the three-dimensional counterparts and focusing on topics in which viscosity plays an important role. By virtue of progresses in analytical theory and high-speed imaging, dynamics of drops and bubbles have actively been studied with the aid of scaling arguments. However, compared with three dimensional problems, studies on the corresponding problems in Hele-Shaw cells are still limited. This review demonstrates that the effect of confinement in the Hele-Shaw cell introduces new physics allowing different scaling regimes to appear. For this purpose, we discuss various examples that are potentially important for industrial applications handling drops and bubbles in confined spaces by showing agreement between experiments and scaling theories. As a result, this review provides a collection of problems in hydrodynamics that may be analytically solved or that may be worth studying numerically in the near future.

When a bubble is injected in a glass of milk, the bubble rises up in the liquid and comes up at the liquid-air interface making a hemispherical air bubble encapsulated by a film of the liquid. The bubble stays at the interface with keeping the hemispherical shape but with decreasing the thickness of the liquid film due to drainage. After some waiting time, a hole is created and grows in the thin film, which leads to disappearance of the bubble. Similar processes are observed with a highly viscous liquid even with lavas during volcanic eruptions and even in the absence of surfactants.
These simple phenomena familiar to everyone involve several elementary processes, such as creation of a bubble, rising of a bubble, film drainage, nucleation of a hole, and bursting of a thin film. We will explore these phenomena and how dynamics change when the bubble is confined in between two plates of a Hele-Shaw cell. Throughout this review, we limit ourselves to the case in which there are no effects of surfactants. We also discuss the case in which the air bubble is replaced with a fluid drop and the case of coalescence of a drop with a bath of the same liquid in a Hele-Shaw cell.
For the readers' convenience, we here introduce some notations used in this review although they are re-introduced in the main text below as we proceed: • η 1 and ρ 1 : viscosity and density of a fluid drop (a bubble or a liquid drop).
• η 2 and η 2 : viscosity and density of another fluid (liquid or air) surrounding the fluid drop.
• ∆ρ = |ρ 1 − ρ 2 | > 0 the density difference of the two fluids; in discussing a bubble this is practically equal to the density of the surrounding fluid ρ 2 .
• γ 12 interfacial energy between the fluid drop and the surrounding fluid; in discussing a bubble, this is equal to the surface tension of the surrounding fluid γ 2 (in discussing a liquid drop in air in Sec.VI, this is equal to the surface tension of the liquid drop γ 1 ).
• κ −1 12 = γ 12 /(∆ρg) capillary length for the interface between the fluid drop and the surrounding liquid; in discussing a bubble, this is practically equal to the capillary length for the surrounding fluid κ −1 2 = γ 2 /(ρ 2 g)

II. LIFE TIME OF A BUBBLE
When a bubble is created in a ultra-viscous liquid (viscosity η 2 and density ρ 2 ), it raises up to make a hemispherical bubble at the liquid-air interface as shown in Fig. 1, and stays there for a considerably long time, even though there are no effects of surfactants, before it breaks up by bursting of the spherical thin film encapsulating the bubble. In such a case, the life time of the bubble is governed by the drainage of liquid from the thin film. With taking the θ coordinate as the angle from the bubble top along the bubble surface (θ = 0 at the top and θ ≃ π/2 at the liquid-air interface), the thinning dynamics of the film of thickness h is known [1] to be given by for the initial thickness h 0 at the top of the bubble, i.e., at θ = 0, with the life time τ and the function f given by with the radius of the bubble R and the gravitational acceleration g.
This scaling law for the life time can be understood as follows. In the present case, the velocity v satisfies the following lubricant equation: where v = v(θ)e θ with e θ the unit vector in the θ direction and the gravitational acceleration vector g pointing in the direction θ = π. This equation can be solved in the form v(θ) ≃ (ρ 2 gR 2 /η 2 ) sin θ, which results from the fact that the flow is a plug flow and does not change in the direction of the film thickness (i.e., the direction of the r coordinate for which r = R at the bubble surface). The conservation equation, can be solved through a separation of variables, which gives Eq. (1) with Eqs. (2) and (3), by noting that h = h(t, θ) and ∇ · (hv) =∂(sin θhv)/∂θ/(R sin θ).
When a bubble is sandwiched by two cell plates separated by the distance D in a Hele-Shaw cell, the boundary condition changes as illustrated in Fig. 2. Because of this, Eq.

III. DRAG FRICTION ACTING ON A BUBBLE
When a bubble is created in a liquid bath (viscosity η 2 and density ρ 2 ), the bubble rises up in the liquid. When the rising velocity V is small, we can assume that the shape of the rising bubble is practically a spherical shape of radius R. In such a case, the drag friction acting on the bubble is given by Note that Stokes' friction law for a solid sphere with radius R moving with a velocity V in a viscous liquid with viscosity η is given by F = 6πηV R.
When a bubble is sandwiched by the two cell walls separated by the distance D of a Hele-Shaw cell, the drag force is found [3] to be given by The experimental set up to confirm this is shown in Fig. 3(a), together with a related set up discussed later. The air drop is slightly elongated, which is characterized by the transverse and longitudinal lengths R T and R L as defined in Fig. 4(a).
The scaling law given in Eq. (9) is explained by considering the slight deformation of the bubble shape. When viscous dissipation is dominant for the velocity gradient V /D developed near the bubble between the cell plates separated by D, the balance between the viscous dissipation and gravitational energy gain (per unit time) can be expressed as This results in which is discussed in [4] and [5] as further explained below. If Eq. (10) is expressed in the force-balance form η 2 V R 2 T /D ≃ ρ 2 gR T R L D, we immediately know that the drag friction acting on the bubble from the viscous medium is given by Eq. (9).
In order to observe this scaling law, dissipation associated with the velocity gradient V /D wins over other viscous dissipations, such as those associated with velocity gradient V /R T , and V /R L , from which the condition D < R T , R L is required. An important fact to be noted here is that in the present experiment, since the liquid surrounding the bubble completely wet the cell plates, there are thin layers of the liquid between the cell plates and both surfaces of the disk-shaped bubble (see the left side view in Fig. 3). The creation of velocity gradient ≃ V /h with h the thickness of layers is energetically unfavorable due to the relation h ≪ D, and because of this the dissipation associated with the gradient V /D becomes dominant in the present case. However, this is not always the case, as shown in Sec. V.
The scaling law in Eq. (9) is confirmed in Fig. 5. In (a), the rising velocity V is given as a function of the cell thickness D for different parameters, R T , R L , and ν 2 = η 2 /ρ 2 . In (b), the same data are plotted with renormalized axes, on the basis of the dimensionless form of Eq. (9), (R T /R L )Ca = D/κ −1 2 , with the capillary number Ca = η 2 V /γ 2 and the capillary length κ −1 2 = (γ 2 /(ρ 2 g)) 1/2 , to show a clear data collapse, confirming the validity of Eq. (9). Although the values of the factor R T /R L for the data used in Fig. 5 are relatively close to one, this factor is important for the clear data collapse; it is confirmed in [3] that if this factor is set to one the quality of the collapse is clearly deteriorated.
The rising bubble in a Hele-Shaw cell discussed above was theoretically studied by Taylor and Saffman [4] (earlier than the Bretherton's paper on bubbles in tubes [6,7]) and by Tanveer [5]. Other theoretical works on fluid drops in the Hele-Shaw cell geometry includes the topological transition associated with droplet breakup [8][9][10][11].
Experimentally, a number of researchers have investigated the rising motion of a bubble in a Hele-Shaw cell [12][13][14][15], together with studies performed in different geometries, which includes studies on inertial regimes (Crossovers between viscous regimes treated in this review and inertial and other regimes are fruitful future challenges) [16][17][18]. However, they have mostly concerned with the case in which there is a forced flow in the outer fluid phase and/or the case in which the cell is strongly inclined nearly to a horizontal position.
One of the features of the present case is that thin liquid films surrounding a bubble plays a crucial role as mentioned above (see the left side view in Fig. 3): in many previous works, the existence of such thin films is not considered. In this respect, the present problem is closely related to the dynamics governed by thin film dissipation such as the imbibition of textured surfaces [19][20][21][22][23][24], as further discussed in Sec. V. In this sense, our problem is quasi two-dimensional, although the geometry of the Hele-Shaw cell is often associated with a purely two-dimensional problem.

IV. DRAG FRICTION ACTING ON A LIQUID DROP
Even if an air bubble is replaced by a liquid drop (density ρ 1 ) whose viscosity η 1 is smaller than that of the surrounding bath liquid η 2 , the drag friction and the cor-responding velocity are still given respectively by Eqs. (9) and (11): Note that, in the second equation, ρ 2 is replaced with ∆ρ = |ρ 2 − ρ 1 |, which is practically equal to ρ 2 in the previous section with a bubble of density ρ 1 .
The experimental setup for confirming this law is illustrated in Fig. 3(b). In this case, the fluid drop goes down, because the density of the drop is larger than that of the surrounding drop, and the drop shape is considerably similar to that of an air bubble, as shown in Fig. 4 However, when the drop viscosity η 1 is larger than that of the surrounding liquid η 2 as in Fig. 4(c), the drag friction is found to be given by the following scaling law [3]: This scaling law is expected because the Eq. (10) is replaced by In this case, the shape of the drop is significantly elongated compared with the previous case.
In summary, two scaling regimes of the viscous drag friction acting on a fluid drop in a Hele-Shaw cell (R L , R T > D) are identified as follows: with the corresponding velocity laws: The crucial condition for these scaling laws to appear is that the formation of thin layers of the surrounding liquid between the cell plates and both surfaces of the disk-shaped droplet (see the right side view in Fig. 3). Otherwise, the drop would directly contact with the cell plates making contact lines, and the dynamics would be affected by the contact angle hysteresis (CAH). In general, physical understanding of CAH is very limited and CAH is sensitive to even a slight contamination of the surface. Thus, if a special care is properly taken for the surface of cell walls, the study of drag friction under the existence of the contact line would be a challenging future problem.
It is worth mentioning that in the three dimensional case, in the limit of highly viscous fluid drop, η 2 ≪ η 1 , which corresponds to the viscous friction acting on a solid sphere, is given by Eq. (8), in which η 1 is absent in contrast with Eq. (12). This is again because of the existence of the thin layers of the surrounding fluid whose thickness is h. However, when η 1 is extremely larger than η 2 , depending on the thickness of the layers, dissipation associated with the viscous stress η 2 V /h becomes smaller than that with η 1 V /D. This case is discussed in Sec. V.
In Fig. 6, Eq. (13) for η 2 < η 1 is directly confirmed, and thus Eq. (12) for η 2 < η 1 are indirectly confirmed, for the data represented by the filled symbols. However, the data represented by open symbols do not follow the scaling law, which is the case mentioned in the previous paragraph and discussed in the next section. Note that (as clear from theoretical arguments above) Eqs. (12) and (13) for η 2 > η 1 are physically the same as Eqs. (9) and (11), respectively, which is shown experimentally in [3].

V. NON-LINEAR VISCOUS DRAG FRICTION
When the dissipation in the thin film is dominant, the velocity law is given by and corresponding friction law is non-linear in V : with the capillary length κ −1 12 = γ 12 /∆ρg defined for the drop-bath interface whose interfacial energy is γ 12 . For D < κ −1 12 , the velocity and force are replaced by the following expressions, respectively: η 2 V ≃ γ 12 (D/κ −1 12 ) 6 and F ≃ (η 2 V /γ 12 ) 1/3 γ 12 R T R L /D. We note here that the viscous friction forces including this nonlinear friction and that given in Eq.
Equations (14) and (15) are derived because Eq. (10) is replaced with with h the thickness of the thin layers of the surrounding liquid, and because, on the basis of the scaling arguments, [31] the thickness of the film is expected to be given by The value numerical coefficient k is estimated as 0.94 in a more rigorous calculation in Refs. [32,33]. With Eq. (14), Eq. (17) can be expressed as For D < κ −1 12 , the thickness is given by the expression [34] h ≃ D(η 2 V /γ 12 ) 2/3 ≃ D(D/κ −1 12 ) 4 , which leads to the above-mentioned velocity and drag force for D < κ −1 12 . The scaling law in Eq. (14) is directly confirmed, and thus Eq. (15) is indirectly confirmed, in Fig. 7(a), as explained in the caption. The data used in Fig. 7 The condition for the crossover between linear and non-linear viscous frictions, corresponding to that between Eq. (13) for η 2 < η 1 and Eq. (14), can be determined by comparing the magnitude of dissipations of the two regimes: This condition can be written as This condition is confirmed in Fig. 7(b), as explained in the caption, by using the data in Fig. 6. From this phase diagram, we expect, for example, the crossover between the two regimes for the series of data obtained for D ≃ 1 mm (red symbols). This crossover can be directly confirmed in Fig. 6. The left three red (open) symbols are on the horizontal line, whereas right two red (filled) symbols are on the dashed line. Note that red filled circles and red filled inverse triangles are exceptional data: these data are obtained when droplets go down nearly on the same trail more than once; as a result, mixing of the two liquids comes into play, which increases the viscosity of the thin layers sandwiched between the bubble and cell plates and makes the velocity gradient in the thin films unfavorable (see Ref. [35] for the details).
Knowing that there are different scaling regimes for viscous drag friction in the Hele-Shaw cell, a natural question to ask would be how the friction law will change if the surrounding fluids are replaced with granular particles in the cell. Such a study was performed by using the setup shown in Fig. 8. The cell filled with granular particles (diameter d and density ρ) can be moved at a constant velocity V and the drag force F acting on the obstacle of radius R can be monitored by the force gauge by virtue of a non-extensible strong fishing line connecting the obstacle and the force gauge. As a result, the drag force F is found [36,37] to be given (1) by the dynamic part scaling as V 2 governed by the momentum transfer associated with collisions between the obstacle and cluster of grains ahead of it and (2) by the static part F 0 governed by the friction force (friction coefficient µ) acting on the cluster from the bottom plate of the cell: where φ is the packing fraction and φ c is a critical density at which the drag force diverges.
This divergence, related to dynamic jamming of grains ahead of the obstacle [38,39], is clearly confirmed in Fig. 9. Interestingly, the value of φ c ≃ 0.841 obtained from the experiment exactly matches the static jamming point in two dimension [40,41]. The background of this study is briefly summarized as follows. The dynamic jamming has actively been explored [42][43][44][45]. Experiments on soft colloids [29,30] have demonstrated good agreements with a phenomenological theory [46]. However, as for granular systems, simulations and a scaling phenomenology revealed critical behaviors different from those of soft colloids [47][48][49]. Experimentally, at high velocities ( 100 mm/s) and high densities, a force component scaling with velocity squared has been reported through impact experiments [50][51][52], while different velocity dependences have been reported in particular at much slower velocities [53][54][55][56][57][58][59][60][61]. The divergence of the drag force at high velocities discussed above has been examined in recent numerical studies [62,63].

VI. DROP COALESCENCE
When a rain drop falling on the surface of a puddle, the rain drop coalesces into the puddle. This is because the area of the liquid-air interface is made smaller by the process. Figure 10(a) shows a series of snapshots of such a coalescence, but with minimizing the "falling velocity": the coalescence is here initiated by touching the two droplets that are held at the tips of pipets with a speed negligible for the coalescence dynamics [64]. From the figure, we see the formation of a neck between the tips of the droplets that grows with time t.
Theoretical prediction was made analytically by Eggers and co-workers in 1999 [65], and here we reproduce the result at the level of scaling laws. The driving force is surface tension, which is opposed by inertial and viscous forces. We estimate changes of the corresponding three energies per unit time. The change of surface energy, dE c /dt ≃ d(γ 1 r 2 )/dt ≃ γ 12 r 2 /t, with the neck size r and the surface energy (per unit area) γ 1 , which can be expressed as γ 12 (the interfacial energy between the drop and the surrounding fluid, i.e., air). This is because when a neck with radius r is formed, the total area the two droplets is decreased by an amount of order of r 2 . This energy will be lost by viscous dissipation and/or changed into kinetic energy.
In the small r limit, viscous dissipation dominates because the volume of moving region is small. In such a case, the two droplets contact with each other nearly by a point, implying that there is only a single length scale for the moving region, which is r. As a result, viscous dissipation per time is given as dE v /dt ≃ η 1 (V /r) 2 r 3 , where η 1 is viscosity of the liquid drop and V ≃ r/t is the growth rate of the neck size. This is because viscous dissipation per unit time and per unit volume scales as η 1 (V /r) 2 and the volume of dissipation scales as r 3 .
When r gets larger, kinetic energy comes into play. Since the length of neck scales as r 2 /R from a geometrical consideration, kinetic-energy increase per unit time scales as dE i /dt ≃ ρ 1 r 2 (r 2 /R)V 2 /t.

From the above arguments, at short times surface energy is lost in viscous dissipation:
dE c /dt ≃ dE v /dt, that is, r ≃ γ 12 t/η 1 . In fact, Eggers and co-workers predicted an extra correction term proportional to log t, which has not been confirmed experimentally. In the same way, at later times, surface energy is changed into kinetic energy: dE c /dt ≃ dE i /dt, that is, r ≃ (γ 12 R/ρ 1 ) 1/4 t 1/2 . The two scaling relations, one for the initial viscous regime and another for the later inertial regime can be expressed in the following dimensionless forms: with τ v = η 1 R/γ 12 , τ i = (ρ 1 R 3 /γ 12 ) 1/2 , and τ = τ 2 v /τ i . An example of experimental confirmation of Eq. (22) is given in Fig. 10(b) and (c), in which r and R are replaced with R and R 0 , respectively. The top plots in (b) and (c) show the raw data in the two regimes. The bottom plots show that the data are well collapsed onto the master curves by rescaling of the axes, which confirms Eq. (22). Fig. 11 [66], in which (r, γ 12 , η 1 , ρ 1 , t) is replaced with (r min , S, η 0 , ρ 0 , τ ). Here, this crossover is confirmed in a different system. When a liquid drop of dodecane is placed on the surface of bath water, the drop takes a thin lens shape. The coalescence of two of such floating droplets are captured in Fig. 11(a). In this geometry, the coalescence is found to be purely two-dimensional, meaning that the velocity gradient in the direction of the thickness of the lens is practically negligible, for which Eq. (22) is valid. As a result, a clear scaling crossover of the coalescence dynamics is shown in Fig. 11(b). Note that various crossovers of scaling regimes for a coalescence of a fluid drop surrounded by another fluid in three-dimension are recently discussed [67].  Fig. 12). This is because the acrylic plates are completely wetted by the oil. After the formation of two phases, a droplet of the glycerol solution is created in the upper oil phase, which slowly goes down the oil phase to touch the oil-glycerol interface, at which coalescence of the droplet and liquid bath of the same liquid is initiated. Figure 13 shows snapshots at short times (initial stage) and those at large times (final stage). These photographs are obtained by a high-speed camera. The time interval between adjacent shots in Fig. 13(b) are ten times larger than that in Fig. 13(a), indicating the dynamics slows down significantly at large times. Figure 14 shows two scaling regimes characterized by slope 1 (i.e., the neck size r scales with t at short times) and slope 1/4 (i.e., the neck size r scales with t 1/4 at large times), which helps us to construct a theory. Before this, let us remark on experimental tricks introduced in the present experiment. In this experiment, as explained in the second paragraph above in which Fig. 12 is first mentioned, the oil phase is deliberately introduced to remove the effect of contact line, which is difficult to control: since the liquid disk is sandwiched by thin layers of the oil, the disk has no direct contact with the cell walls. In addition, the viscosity of oil η 2 is set to much smaller than that of the glycerol solution η 1 : the oil phase could play the role of air phase in the three-dimensional counterpart.

The crossover between the two regimes in Eq. (22) is experimentally confirmed in
Considering these points, we see that the dynamics of the initial regime could be the same as the three-dimensional initial dynamics, which is characterized by the slope 1. This is because at short times when r is much smaller than the cell thickness D the dynamics could become independent of the existence of the cell walls: the effect of confinement by the cell walls is practically negligible.
On the contrary, the dynamics at the final stage characterized by slope 1/4 reflects the effect of confinement of wall as discussed below. In the present confined case, we have good reasons to make the following replacements from the three-dimensional viscous counterpart: (1) The gain in surface energy per unit time should be changed from dE c /dt ≃ d(γ 12 r 2 )/dt to dE c /dt ≃ d(γ 12 rD)/dt. (2) The viscous dissipation per unit time should be changed from which the volume of the neck is expressed as rDd. The first replacement is rather obvious because the three-dimensional circular shape of the area ≃ r 2 that disappears as a result of coalescence is here confined by the wall to take a quasi-rectangular shape of area ≃ rD.
The second replacement is expected as long as d D, in which case the dominant viscous dissipation per unit time and unit volume scales as η 1 (V /D) 2 is localized in the neck volume ≃ rDd. In addition, from a geometrical relation, R 2 = (R − d) 2 + r 2 , the neck height d is expressed as d ≃ r 2 /R as long as R d. By balancing the capillary energy gain with the viscous dissipation thus justified, we obtain r ≃ (γ 12 RD 2 t/η 1 ) 1/4 , which leads to the desired slope 1/4.
In summary, the dynamics exhibits a dimensional crossover from the short-time threedimensional one to the long-time quasi-two-dimensional one, expressed by the following set of equations, which are confirmed in Fig. 15.
with τ I = η 1 D/γ 12 , and τ F = η 1 R/γ 12 . Some remarks are as follows. (1) Equation (23) is another dimensionless form of the viscous regime given in Eq. (22). (2) The condition, √ RD < r < R, for the final regime to be valid implies that this dynamics is well observed when D ≪ R is satisfied, which is true for the experimental data. (3) The final regime is predicted to start from r and t are larger than √ RD and τ F , respectively, which are well confirmed in Fig. 15(b).
There are a number of future problems for the quasi-two-dimensional coalescence. One is to find an inertial regime since only viscous regimes have been found. In the threedimensional liquid drop coalescence in air, the Ohnesorg number (Oh ≃ η 1 / √ ρ 1 γ 1 R) appears when we discuss the crossover between the viscous and inertial regimes given in Eq. (22): we can express the crossover by the simple expression r/R = Oh. In the quasi-two dimensional case, the corresponding number will be useful to look for different regimes including inertial regimes. Another interesting problem would be to find partial coalescence in the quasitwo-dimensional case, which has been discussed for a number of systems in three dimension [69][70][71].

VII. FILM BURSTING IN A HELE-SHAW CELL
When the coalescence experiment is performed with a different combination of liquids, glycerol solution and olive oil, the mode of coalescence significantly changes as shown in Fig.   16, because thin film (viscosity η 2 ) formed between the drop and bath of the same liquid (viscosity η 1 ) is more stable in the present case (molecules contained in olive oil may play a role similar to surfactants). This problem can be seen as a problem of film bursting.
Bursting of liquid thin film is familiar to everybody as the rupture of soap films, and is important in many fields [72][73][74]. Although the fundamental physics of the dynamics of bursting film is rather well and simply understood [31], in all previous studies the bursting hole has circular symmetry, which is lost in the present case. Here, the problem is bursting of a thin film in a confined geometry, which is relevant in many practical situations where small amount of liquids has to be manipulated (e.g., microfluidics and biological applications).
When the circular symmetry is preserved, the bursting or dewetting frequently proceeds with a constant speed and in such cases a rim is formed as in the bursting of soap film in air [75][76][77] or in the surrounding oil phase [78] and in the bursting of liquid film on a solid substrate (i.e. dewetting) [79]; when the rim is not formed the speed grows exponentially as observed in the bursting of highly viscous films suspended in air [1,80].
On the contrary, when the circular symmetry is lost, the bursting proceeds with a constant speed but no rim is formed and the velocity V is found [81] to be a constant given by when η 1 η 2 . This scaling law can be derived by the balance of energy per unit time: This scaling regime is significantly different from the three-dimensional case because the dynamics proceeds with a constant speed although the rim is not formed. spaces such as petroleum industry [15] and microfluidic applications [82], as well as for industrial and fundamental problems associated with foams and microemulsions [25,83].
The scaling laws in the Hele-Shaw geometry discussed here have been established only by agreement between experiments and scaling phenomenologies, providing future challenges for analytical and numerical studies.    Slanted view of the setup. One layer of granular particles of average diameter d = 2 mm are packed in a two-dimensional cell made of transparent acrylic plates. The diameter 2R of the obstacle is about ten times larger than d. Reproduced from [36].   show qualitative comparison of theory with experiment. In the plots, τ I and τ F in Eqs. (23) and (24) are expressed as τ i and τ f , respectively. In addition, γ 12 and η 1 are expressed by γ and η, respectively. Adapted from [68].