2.1. Turbulence generation and characterisation

Turbulence in a 0.37 m $^3$ water tank is generated by the convergence of eight turbulent jets created by four submerged water pumps, as sketched in figure 1(a) and described in greater detail in Ruth et al. (Reference Ruth, Vernet, Perrard and Deike2021). The flow from each pump is split into two parallel jets at a Y, with each outlet separated by 7.8 cm, with the centres of the Ys forming the vertices of a 25 cm square in the horizontal plane. Figure 1(b) presents properties of the flow as characterised in the central plane ($y=0$) of the experiment with two-dimensional, two-component particle image velocimetry (PIV). The background gives the local fluctuation velocity $u' = \sqrt {({u'_x}^2+{u'_z}^2)/2}$, where $u'_i=\sqrt {\overline {(u_i-\overline {u_i})^2}}$ and overbars denote averaging in time. Here $u'$ tends to be largest in the plane of the jets ($z \approx 0.01\,\textrm {cm}$) and in the region below their convergence zone ($x \approx y \approx 0$). PIV is performed in nine parallel planes, enabling the three-dimensional interpolation of turbulence quantities at any location within the measurement domain.

Figure 1. Turbulence generation and characterisation. (a) A sketch of the experiment (not to scale), consisting of a $0.37\,\textrm {m}^3$ tank of water in which four pumps, each split to two outlets, are arranged at the corners of a square in the horizontal plane. The turbulence is characterised with PIV performed separately in nine parallel planes, with illumination provided by a laser sheet (shown in green). (b) Properties of the turbulent flow field in the central plane of the experiment. The background shows the local value of $u'$, denoted by the colour given in the colourbar. The green dashed rectangle shows the field of view employed in the large air cavity disintegration experiments. The diameter of the black circles denotes the Hinze scale $d_{H}$ at various $x$ and $z$. The length of the cyan rectangles denotes the integral length scale $L_{int}$ at those locations.

As described in Ruth et al. (Reference Ruth, Vernet, Perrard and Deike2021), we compute the integral length scale $L_{int}$ locally at each point in the flow by integrating the spatial autocorrelation function. It changes throughout the experiment, being the shortest where the turbulence is the strongest. The cyan lines in figure 1(b) denote the value of $L_{int}$ at various locations in the central plane of the experiment: $L_{int}$ is shortest near the convergence of the jets, and grows at lower and higher depths. With $u'$ and $L_{int}$ calculated from the PIV data, we can then compute the local turbulence dissipation rate under the assumption of isotropy with $\epsilon = C_\epsilon u'^3/L_{int}$, with $C_\epsilon = 0.7$ (Sreenivasan Reference Sreenivasan1998), and the Kolmogorov microscale with $\eta = ((\mu /\rho )^3/\epsilon )^{1/4}$ (Pope Reference Pope2000). The Hinze scale $d_{H}$, calculated using (1.2), is denoted at various locations by the diameter of the black circles drawn in figure 1(b). The Hinze scale is smaller where the turbulence is more intense, meaning that more bubbles will be larger than the Hinze scale and susceptible to break-up at these locations. We refer to Ruth et al. (Reference Ruth, Vernet, Perrard and Deike2021) for more details on the structure of the turbulence field and for maps of turbulent quantities outside of the central plane.

2.2. Large cavity disintegration experiment

For the experiment on large cavity break-ups, air cavities were produced following Landel, Cossu & Caulfield (Reference Landel, Cossu and Caulfield2008) by placing a hollow hemispherical cup with $R=5\,\textrm {cm}$ underwater, sketched in figure 2(a), and bubbling a known volume of air $V_0 = {\rm \pi}d_0^3 / 6$ into it. Once bubbles in this cup have coalesced into a single air cavity, the cup is inverted by rotating it rapidly half a revolution, such that the air inside is suddenly no longer constrained by the curved cup surface. The top surface of the initial volume of air roughly conforms to the curved inner surface of the cup. The large air cavity, having been suddenly exposed to stresses from the surrounding turbulence and its buoyant rise through the water, deforms and starts a complex sequence of break-ups, leading to its disintegration. The surface of the cup rotates with a speed around 0.4–0.9 m s$^{-1}$, and we have checked that this speed does not systematically affect the early stages of the bubble size distribution. Further, similar experiments run without turbulence yield very little break-up, as evidenced in Appendix B.

Figure 2. Experiment on large cavity disintegration. (a) Schematic of the experiment. Air is bubbled into an inverted hemispherical cup located just under the convergence of the turbulent jets, and the cup is rapidly rotated to expose the air to the turbulence. The experiment is lit from behind (not shown) and filmed with a high-speed camera. (b) One representative image of a cavity breaking apart, with the cup still slightly visible at the bottom of the image. (c) The characteristic length scales $\eta$, $d_{H}$, $L_{cap}$ and $L_{int}$ taken in analysing the data, the pixel size $\varDelta x$ and the minimum bubble size considered $d_{min}$, and the diameters of the air cavities studied (circles). Distributions of the turbulence quantities in the field of view in the centre of the tank (within the green rectangle in figure 1) are also given in grey.

The turbulent flow in the region of the tank imaged in this experiment is denoted by the green rectangle in figure 1(b). The turbulence varies spatially, so to simplify the analysis, we take $u' \approx 0.2\,\textrm {m}\,\textrm {s}^{-1}$, $L_{int} \approx 1.5\,\textrm {cm}$ and $\eta \approx 37\,\mathrm {\mu }\textrm {m}$ as characteristic values, which set $d_{H} = 3.2\,\textrm {mm}$ and ${Re}_{t} = 3400$. These length scales are denoted in figure 2(c), which also gives the distribution of the length scales present in the field of view in the middle of the tank. The mean flow is downwards with $\bar {W} \approx -0.25\,\textrm {m}\,\textrm {s}^{-1}$, largely counteracting the buoyant rise speed of larger bubbles. This enables the bubble population to linger in the measurement region for a sufficient period of time to image it over multiple large-scale eddy turnover times $T_{int} = L_{int}/u' \approx 0.075\,\textrm {s}$.

The cavities range in size between $d_0/d_{H} = 2.12$ and 8.30. Data for each condition, as well as the number of runs recorded at each, are given in table 1; to build statistical data, more runs are carried out at smaller $d_0/d_{H}$, because fewer bubbles are formed in those break-ups.

Table 1. Conditions of the experiments. Characteristic values are given for each of the cavity disintegration cases. For the experiments on individual bubble break-up, the mean and standard deviation among the 162 recorded cases are given for each quantity.

One image is shown in figure 2(b). The cup is visible in the bottom of the image as it has not yet fully rotated out of the field of view. The measurement region, which spans 15.8 cm in the $x$ direction and 8.9 cm in the $z$ direction, is illuminated from the back, and the disintegration of the cavity is filmed with a high-speed camera at 500 Hz with a spatial resolution of $38\,\mathrm {\mu }$m pixel$^{-1}$. The field of view is much larger than all the bubbles considered, so it does not introduce a significant bias related to bubbles whose images extend partially outside the field of view. Bubbles are detected with an image processing method described in Appendix A.1, and their effective diameters $d$ are determined as the equivalent diameter of a circle with the same area as the projected bubble image. In analysing the data, we consider only bubbles for which $d\geq 400\,\mathrm {\mu }\textrm {m}$, for which the detection is less sensitive to the chosen image intensity threshold. Given the typical severe deformation and overlapping images of larger bubbles ($d \gtrsim 6\,\textrm {mm}$), we note that their sizes will tend to be over-estimated by this method. The air void fraction in the vicinity of the cavity is high enough that we are unable to track the dynamics of individual break-ups, so we restrict our analysis to the resulting bubble size distribution.

To account for the limited field of view in our experiments, we adjust the measured size distributions by keeping a record of bubbles which have left and entered the field of view. Those which leave are ‘locked’ into the bubble record used in computing the size distributions, whereas those that enter the field of view are excluded from the calculation of the size distribution. This process is explained in Appendix A.2 and only has a limited effect on the results reported in this paper, as we do not consider the size distribution at late times.

2.3. Individual break-up tracking experiment

In the second set of experiments on bubble break-up, we dynamically track the individual break-ups of bubbles in the turbulent region. As sketched in figure 3(a), bubbles are introduced to the bottom of the tank through a needle and rise to the turbulent region. Two cameras, which are synchronised with a function generator, film at 1000 f.p.s.. They are oriented $90^{\circ }$ from each other and their fields of view overlap in a measurement volume of approximately $200\,\textrm {cm}^3$. The cameras are calibrated by mapping their pixels to the paths of the light rays reaching the pixels, following the method presented by Machicoane et al. (Reference Machicoane, Aliseda, Volk and Bourgoin2019). Then, following a method similar to that used in Ruth et al. (Reference Ruth, Vernet, Perrard and Deike2021), we identify the three-dimensional location of the bubbles which are simultaneously captured by each camera. The spatial resolution of each camera varies with the position of the bubble, but the typical value of the two cameras are 28.9 and $57.1\,\mathrm {\mu }\textrm {m}$ pixel$^{-1}$. An approximate lower bound for the size of the smallest resolved bubble is then $d_{min} \approx 200\,\mathrm {\mu }\textrm {m}$.

Figure 3. Experiment to obtain dynamic reconstructions of individual break-up events. (a) Schematic of the experiment. Air bubbles are introduced through a needle at the bottom of the tank and rise into the turbulence created by the jets. The bubbles are filmed with two high-speed cameras, enabling the determination of the three-dimensional bubble trajectories. (b) The trajectories of parent and child bubbles involved in one break-up event, and their projections onto the horizontal ($x$–$y$) plane. The colour corresponds to the bubble's size relative to the local Hinze scale, which varies spatially with $\epsilon$ as the bubble size is fixed. The green dot denotes the first detected position of the parent bubble; the red dots denote the final detected position of the child bubbles.

The trajectories taken by the bubbles are then determined using the Python package Trackpy (Allen et al. Reference Allen, Caswell, Keim, van der Wel and Verweij2021), which implements the algorithm from Crocker & Grier (Reference Crocker and Grier1996). Such trajectories are shown in figure 3(b). Using the three-dimensional map of the turbulence statistics obtained from PIV, we compute the bubble's size relative to the local Hinze scale $d/d_{H}$ (computed with the local value of $\epsilon$) at each bubble location, which is encoded in the colour in the figure. The mean dissipation rate at the break-up locations is $\epsilon = 0.52\,\textrm {m}^2\,\textrm {s}^{-3}$, with a standard deviation of $0.21\,\textrm {m}^2\,\textrm {s}^{-3}$. The mean values and standard deviations of quantities describing the initial conditions for the break-ups studied in this experiment are given in table 1.

From the bubble trajectories, we identify each time a bubble breaks apart, which occurs when a new trajectory (or trajectories) appears in the vicinity of a previously existing bubble. As the tracking algorithm will initially link the parent bubble to only one of the child bubbles, the parent bubble trajectory is then split at this time, and both child bubbles are treated equally. These events are denoted by the grey lines connecting the ‘end’ of one bubble to the ‘beginning’ of another in figure 3(b). Given the complex deformations involved in some break-ups, the method does not always resolve the fast splitting dynamics accurately; the break-up child size distributions we report, however, are not sensitive to the order of events occurring within one break-up event.

3.1. Disintegration of cavities of increasing sizes

The break-ups of two air cavities, one with $d_0 = 0.68\,\textrm {cm}$ and one with $d_0 = 2.25\,\textrm {cm}$, are shown in figures 4 and 5, respectively. These correspond to non-dimensional sizes of $d_0/d_{H} = 2.12$ and $7.00$, $d_0/L_{int} = 0.46$ and $1.50$ and $d_0/l_{cap} = 2.51$ and $8.28$. As a reference, the constant values taken for $L_{int}$ and $d_{H}$ and the initial size of the cavity $d_0$ are denoted in the top-left corner of the first image. In both cases, the hemispherical cup which had constrained the bubble is visible at early times as it is rotated away.

Figure 4. Disintegration of an air cavity with $d_0/d_{H} = 2.12$, involving just one break-up during the interval shown.

Figure 5. Disintegration of an air cavity with $d_0/d_{H} = 7.00$.

In the disintegration of the smaller cavity, with $d_0/d_{H} = 2.12$ (shown in figure 4), the bubble emerges from the cup with a moderate deformation caused by buoyancy and the surrounding turbulence. Eventually, within approximately one integral-scale turnover time, the bubble becomes more elongated and breaks into two bubbles. One is near the parent bubble in size, and other is slightly smaller than the Hinze scale. These two bubbles persist without breaking for at least ${\sim }2$ more integral-scale turnover times, at which point the smaller of the two bubbles is advected out of the field of view by the downwards mean flow.

The deformation to the larger cavity shown in figure 5 is more severe, leading to a more complex sequence of events during its disintegration. Upon emerging from the rotating cup, the cavity is flattened due to buoyancy (as $d_0 / l_{cap} = 8.28$ for this case), and turbulent deformations to the cavity shape on the order of the cavity size itself quickly develop. By $t/T_{int} \approx 0.4$, the cavity consists of two lobes (each of which is significantly deformed), separated by a shrinking neck of air. By the time the neck has pinched apart ($t/T_{int} \approx 0.7$), the two larger child bubbles stemming from the two lobes are accompanied by much smaller child bubbles (some with $d \ll d_{H}$ and $d \ll d_0$) which were formed during the collapse of the air neck. The larger child bubbles themselves go on to further break apart in a chain of break-ups, some of which similarly involve small bubble production via the collapse of elongated air necks. Many small bubbles which are more than an order of magnitude smaller than the initial one are eventually visible. At much later times, the largest bubbles have risen out of the field of view, and the total air volume imaged is decreased significantly.

3.2. Transient evolution of the bubble size distributions

The experiment was carried out with six values of $d_0 / d_{H}$ between 2.12 and 8.30, with 10–20 runs taken at each condition, as given in table 1. Note that the largest cavities exceed the integral length scale in size, so the typical turbulent stress at their spatial scale will be saturated relative to that predicted by the Kolmogorov scaling employed in the definition of the Hinze scale. Figure 6 shows the transient evolution of $\mathcal {N}(d/d_{H}) = N(d) d_{H}$ for each condition (ensemble-averaging the 10–20 runs). Each curve is the dimensionless size distribution, with the adjustment to account for bubble advection explained in Appendix A.2, averaged over times within ${\pm }0.1 T_{int}$ of the stated time. The number of bins employed in discretising $d$ is chosen as a function of the number of bubbles present; a finer resolution is used when the distribution is based on more bubbles. Further, we show with the dashed red lines a representative noise threshold. This is defined, somewhat arbitrarily, as the distribution resulting from a total of five bubbles within each $d$ bin (at the coarsest discretisation) being imaged, each during the entirety of the averaging window, over the course of the 10–20 experimental runs at each condition. (Therefore, the noise limit can be reached for a given $d$ by the observation of five bubbles of size $d$ in one run, or by the observation of one bubble of size $d$ in five separate runs.)

Figure 6. Bubble size distributions during the disintegration of cavities with $d_0/d_{H}$ between 2.12 and 8.30 and times up to $4 T_{int}$ after the cavity is released into the turbulence. The size of the parent bubble is denoted by the dashed vertical line. Each distribution integrates to the average number of bubble observed at that condition at that time. The sub-Hinze power-law scaling exponent becomes approximately $-$3/2 once a significant number of break-ups have occurred. The dotted red line denotes a noise threshold discussed in the text.

At early times, the distributions for all $d_0/d_{H}$ exhibit a peak at $d_0/d_{H}$, denoted by the vertical dotted lines. For the two smallest cavities (given in figure 6a,b), for which no break-up was observed during many runs, only a small number of bubbles are formed over time, and the size distribution near the injection scale does not decrease appreciably with time.

Over time, as the larger cavities (given in figure 6c–f) begin to disintegrate, the size distribution for $d< d_0$ begins to be built up. Even among these larger cavities which produce a considerable number of sub-Hinze bubbles, the increase in the number of sub-Hinze bubbles is much more pronounced for the cavities that are initially larger (evidenced by comparing curves for $d_0/d_{H}=4.15$ and $d_0/d_{H} = 8.30$, for example). This suggests that there is a large separation of scales between the sub-Hinze bubbles and the parent bubbles responsible for their creation; more simply, large bubbles are needed for the production of small bubbles. For the largest cavities, the size distribution for sub-Hinze bubbles eventually follows an $\mathcal {N}(d/d_{H}) \propto (d/d_{H})^{\alpha _d}$ scaling, with $\alpha _d = -3/2$, sketched on all plots as the dashed line for reference. The final curves shown (for $t/T_{int}=4$) might constitute an under-estimation for the bubble size distribution for smaller bubbles, because some of the bubbles which may break have risen out of the field of view by this time. Although we expect the size distribution at $d \approx d_0$ to decrease with time for these large cavities, this observation is not clear due to the difficulty in sizing the largest bubbles, due to their severe deformations.

Now, we consider the size distributions averaged between $2 T_{int}$ and $4 T_{int}$. By these times, a significant number of break-ups have occurred (for larger $d_0/d_{H}$), but a significant portion of the bubbles have not yet left the field of view, and the bubble size distribution approaches a constant shape. Figure 7(a) compares the size distributions over these times for each value of $d_0/d_{H}$. For larger air cavities, the magnitude of $\mathcal {N}(d/d_{H})$ is increased, and the sub-Hinze power-law distribution steepens. The same data are shown in (b), normalised by the cavity diameter $d_0$ instead of the Hinze scale. Larger cavity sizes yield a ${\propto }d^{-3/2}$ scaling for all bubble sizes.

Figure 7. Time-averaged bubble size distributions. (a) The dimensionless bubble size distribution averaged between $t/T_{int}=2$ and $t/T_{int}=4$ for cases with varying $d_0/d_{H}$, denoted by the position of the coloured notches along the bottom axis. The dotted red line denotes a noise threshold, as in figure 6. (b) The bubble size distributions based on the diameter normalised by the initial cavity diameter $d_0$. (c) The exponent $\alpha _d$ of a power-law fit to the sub-Hinze portion of the distributions which are above the noise threshold, $\mathcal {N}(d/d_{H}) \propto (d/d_{H})^{\alpha _d}$ for $d/d_{H} < 1$, indicating that a $\mathcal {N}(d/d_{H}) \propto (d/d_{H})^{-3/2}$ scaling is realised for large values of $d_0/d_{H}$.

Figure 7(c) shows the power-law exponent fit to the sub-Hinze portion of the distributions in (a), $\mathcal {N}(d/d_{H}) \propto (d/d_{H})^{\alpha _d}$ for $d/d_{H} < 1$, for cases with size distributions above the noise threshold (which is the case for $d_0/d_{H} > 3$). An $\alpha _d \approx -3/2$ scaling, indicated by the dashed black line, is realised for all these cases. The power-law behaviour of the size distribution is affected not only by the break-up physics, but is also steepened by the rising dynamics of the bubbles: as small bubbles rise more slowly than larger bubbles, they tend to linger in the field of view for longer, increasing their concentrations (Garrett et al. Reference Garrett, Li and Farmer2000).

Integrating the transient size distributions over the bubble diameter, the temporal evolution of the number of resolved bubbles $n$ (with the minimum resolvable size $d_{min}=0.12 d_{H}$) is shown in figure 8(a). The grey shaded region denotes the times over which the bubble size distributions are averaged in figures 7 and 8(b).

Figure 8. Evolution of the number of resolved bubbles (limited to $d/d_{H} > 0.12$) with time and the initial cavity size. (a) Temporal evolution of the average number of all bubbles measured experimentally for different initial cavity sizes $d_0/d_{H}$. Circles give the values employed in § 5.3.1. (b) The total number of bubbles (black), number of sub-Hinze bubbles (orange) and number of super-Hinze bubbles (purple) averaged between $2 T_{int} \leq t < 4 T_{int}$ (the region shaded in grey in (a)) as a function of the initial cavity size.

Figure 8(b) shows the number of resolved sub-Hinze, super-Hinze and total bubbles, averaged over the time period considered in figure 7. Again, we see an increase in the number of bubbles formed with the initial size of the cavity. Further, the number of sub-Hinze bubbles increases with the parent bubble size more rapidly than the number of super-Hinze bubbles does, making sub-Hinze bubbles constitute a larger portion of the bubble size spectrum for larger $d_0/d_{H}$. This is remarkable, because as $d_0/d_{H}$ is increased, the span of bubble sizes constituting resolvable sub-Hinze bubbles ($d_{min} < d < d_{H}$) remains fixed, whereas the span of potential super-Hinze bubble sizes ($d_{H} < d < d_0$) increases.

Taken together, figures 7 and 8 are congruent with the capillary pinching mechanism for sub-Hinze bubble production proposed by Rivière et al. (Reference Rivière, Ruth, Mostert, Deike and Perrard2022). Our figures suggest that their formation relies on the break-up of cavities that are significantly larger than the Hinze scale: larger values of $d_0/d_{H}$ yield the $\mathcal {N}(d/d_{H}) \propto (d/d_{H})^{-3/2}$ power-law scaling in the sub-Hinze bubble size distribution, and the dependence on $d_0$ of the number of sub-Hinze bubbles produced (shown in figure 8b) is steeper than that of the number of super-Hinze bubbles produced. We propose in the next section an explanation of the mechanisms leading to this dependence.

3.3. Capillary splitting of ligaments prepared by the turbulence produces small bubbles

Visual observations of the large air cavities disintegrating provide clues into the mechanism responsible for the production of sub-Hinze bubbles: child bubbles much smaller than the Hinze scale are seen to originate from a Rayleigh–Plateau-like instability that occurs during the pinching apart of elongated fluid ligaments prepared by the turbulence. However, the turbulence is only able to deform bubbles that are large enough with respect to the Hinze scale that such ligaments might be created, because surface tension is effective at limiting the severity of deformations to smaller bubbles. These experimental observations parallel a recent interpretation of DNSs of bubble break-up (Rivière et al. Reference Rivière, Ruth, Mostert, Deike and Perrard2022).

Illustrative examples of bubble break-up are given in figure 9, which shows the typical break-ups of bubbles of two sizes: one is near the Hinze scale in size (figure 9a), and another is seven times larger than the Hinze scale (figure 9b). The smaller bubble, with $d_0/d_{H} = 2.12$, is initially deformed into two comparably sized lobes, and the neck separating the two splits at a single point to form two child bubbles, each of a similar scale as the parent. Here, the parent bubble is small enough that surface tension is able to prevent significant deformation during the break-up.

Figure 9. Break-up of a bubble initially close to the Hinze scale in size (a, with $d_0/d_{H} = 2.12$) and initially much bigger than the Hinze scale (b, with $d_0/d_{H} = 7.00$). The deformation to the smaller bubble produces two comparably sized lobes, which split apart to form two comparably sized child bubbles. The deformation to the larger bubble also produces two comparably sized lobes, but these are separated by a much more elongated filament of air. The unstable collapse of this filament produces the small ‘capillary’ child bubble between the two larger bubbles. The small bubbles in the lower left of the image were produced in previous break-ups.

The larger bubble, with $d_0/d_{H}=7.00$, is similarly deformed by the turbulence into two comparably sized lobes prior to pinch-off. However, the filament of air separating the two just prior to pinch-off has become significantly more elongated than the neck in the break-up of the smaller bubble. This elongation opens the door to capillary instabilities along the filament during its collapse: in the instance shown in figure 9(b), the filament pinches apart at two separate points, leaving a small child bubble (with $d \ll d_{H}$) between the two lobes.

The two examples of break-up discussed illustrate two mechanisms present in the break-up of bubbles by turbulence. The first is the deformation of the parent bubble by a turbulent structure, likely on the spatial scale of the parent bubble itself. This brings the bubble to an unstable state consisting of two lobes (which will become what we call the ‘inertial’ child bubbles) separated by a neck of air, which begins to pinch apart under capillarity. When the deformation to the bubble is severe enough, this ligament can take on an elongated, deformed shape. Its pinching can become unstable under a Rayleigh–Plateau-like mechanism, leading to the formation of small ‘capillary’ bubbles.

Figure 10 shows an additional five instances of deformed ligaments undergoing a capillary instability to produce sub-Hinze bubbles. Prior to break-up, the parent bubbles shown are highly deformed, reaching (from the case in the top row to the case in the bottom row) $10.0 d_{H}$, $8.7 d_{H}$, $8.2 d_{H}$, $7.4 d_{H}$ and $8.0 d_{H}$ in their longest axis in the images. The aspect ratios of ellipses fit to their images prior to break-up range from 2.9 to 3.8. Many cases do not solely involve one ligament separated by two well-defined lobes; turbulent deformations make the bubble shapes more irregular. However, in all instances, very small bubbles (some more than an order of magnitude smaller than $d_{H}$) are produced as an air ligament involved in the turbulent deformation collapses unstably.

Figure 10. Five cases of sub-Hinze bubble production by the unstable collapse of a deformed ligament. Each row shows four snapshots in time, spaced 10, 4, 2 and 0 ms before the time at which the sub-Hinze bubbles are first visible. The field of view in the final two columns is given by the blue square in the second column. Prior to break-up, the parent bubbles shown are highly deformed, reaching (from the case in the top row to the case in the bottom row) $10.0 d_{H}$, $8.7 d_{H}$, $8.2 d_{H}$, $7.4 d_{H}$ and $8.0 d_{H}$ in their longest axis in the images.

This description is clearly a simplified understanding of the bubble break-up process, as it does not capture the redistribution of air due to a capillary pressure difference between lobes that may be responsible for the formation of small bubbles (Andersson & Andersson Reference Andersson and Andersson2006), nor does it describe the ‘tearing off’ of very small bubbles that we observe occurring to large parent bubbles more frequently than $1/T_{turb}(d_0)$. However, the framework serves as a bridge between the inertial deformations to a bubble by the turbulence and the later-time collapse dynamics instigated by capillarity. This understanding mirrors the description of bubble pinch-off in turbulence given in Ruth et al. (Reference Ruth, Mostert, Perrard and Deike2019), in which we showed that turbulence sets an ‘initial’ deformed bubble shape before the collapse dynamics overtake the turbulent dynamics. Once the inertial collapse of the neck becomes fast enough (equivalently, once the neck becomes small enough), however, the turbulence effectively ‘freezes’ in place relative to the accelerating collapse dynamics. The end result is that the final stage of the pinching process (in this case, the production of small bubbles through the capillary instability of gas ligaments) is affected by the turbulence only insofar as the turbulence sets the ‘initial condition’ on which the remainder of the process evolves under capillary and, eventually, inertial dynamics.

4.1. Qualitative discussion of the break-up sequences

One break-up producing $m=2$ child bubbles is shown in figure 11 and one producing $m=4$ bubbles is shown in figure 12. In each, images throughout the break-up sequence are shown in panels (a)–(c), and the three-dimensional trajectories taken by the bubbles involved are shown in panel (d). At each point, the bubble's size is computed relative to the local Hinze scale, and $d/d_{H}$ is encoded in the trajectory colour. The spatial scale is given in terms of the integral length scale at the break-up location, $L_{{int},0}$, showing that the bubble trajectories are resolved over multiple integral length scales. Panel (e) shows the dimensional diameters of the bubbles involved over time.

Figure 11. One dynamically tracked bubble break-up involving the production of $m=2$ child bubbles. (a–c) Images recorded by one of the two high-speed cameras throughout the sequence. (d) The trajectories taken by the bubbles involved, with their size relative to the Hinze scale encoded in the colour. The green circle marks the first observation of the parent bubble and the red circles mark the final observation of the child bubbles. The side length of the square shown is given in terms of the integral length scale at the break-up location. (e) The ‘family tree’ for the single break-up, giving diameters of the bubbles present at each point in time. Fainter lines give the instantaneously measured diameters and straight lines give the median for each bubble, which is the quantity we consider in our analysis.

Figure 12. One dynamically tracked bubble break-up involving the production of $m=4$ child bubbles. (a–c) Images recorded by one of the two high-speed cameras throughout the sequence. (d) The trajectories taken by the bubbles involved, with their size relative to the Hinze scale encoded in the colour. The green circle marks the first observation of the parent bubble and the red circles mark the final observation of the child bubbles. The side length of the square shown is given in terms of the integral length scale at the break-up location. (e) The ‘family tree’ for the single break-up, giving diameters of the bubbles present at each point in time. Fainter lines give the instantaneously measured diameters and straight lines give the median for each bubble, which is the quantity we consider in our analysis.

In the case of binary break-up, given in figure 11, the parent bubble enters the imaged volume from the foreground, and quickly encounters a region of more intense turbulence, where $d_0/d_{H}$ increases. Eventually, having become deformed, the bubble pinches apart into two child bubbles, each of which are comparable in size to the parent. Both child bubbles persist in the field of view for at least a tenth of a second (approximately an integral-scale turnover time) without breaking.

In the more complex break-up shown in figure 12, the parent bubble similarly traverses from a region of less-intense turbulence to more-intense turbulence, increasing the value of $d_0/d_{H}$. Eventually, at $t=0.170\,\textrm {s}$ (figure 12a), the bubble becomes elongated in the vertical direction, and in a sequence of two rapid splitting events produces the three child bubbles that are visible at $t=0.192\,\textrm {s}$ (figure 12b). One is still larger than $d_{H}$, one is of the order of $d_{H}$ and the third, left between the two, is smaller than $d_{H}$. The bubble of the order of the Hinze scale is still significantly deformed, the capillary dynamics involved with the break-up not yet having relaxed. A short time later, by $t=0.201\,\textrm {s}$ (figure 12c), an additional bubble has split from it, leaving a total of four child bubbles.

4.2. Identification of break-up events

We identify bubble break-ups like those shown in figures 11 and 12 as being sequences of bubble splitting events not exceeding the eddy turnover time at the parent bubble scale, $T_{turb}(d_0) = \epsilon ^{-1/3} d_0^{2/3}$. To enforce this temporal constraint, we first construct a ‘family tree’ of all splitting events recorded in one run. Then, if any bubble is present at a time $T_{turb}(d_0)$ beyond the initial detected break-up of the first bubble (with diameter $d_0$), we truncate the family tree at that bubble and start a new family tree with the same bubble (if it later breaks apart). After doing so, we store the sizes of the parent bubble and child bubbles, as well as the turbulence characteristics spatially interpolated at the initial break-up location. To remove spurious break-ups, we discard those for which the sum of the calculated volumes of the $m$ child bubbles is less than 50 % of, or more than 200 % of, the calculated volume of the parent bubble.

In total, we captured 162 bubble break-ups with this dynamical tracking approach that fit the volume conservation criteria. Figure 13(a) shows the distributions of the break-up conditions (the Hinze scale at the break-up location and the parent bubble size) for the aggregated dataset, which we later break down by the parent bubble's size relative to the Hinze scale. The parent bubble diameter $d_0$ is typically slightly larger than the Hinze scale, as the distribution of $d_0$ (the green line) is located just to the right of that of $d_{H}$ (the dashed red line). Thus, the break-ups we capture in this experiment have $d_0/d_{H} \approx 0.4 \unicode{x2013} 3.7$. The black curve shows the distribution of the sizes of child bubbles formed during break-ups, integrating to the average number of bubbles formed per break-up event.

Figure 13. Results on individual bubble break-ups. (a) Distributions of the Hinze scale at the break-up location (red), parent bubble sizes (green) and the child bubble sizes (black). (b) The Hinze scale at the bubble's break-up position (vertical axis) as a function of the Hinze scale at the bubble's position one bubble-scale turnover time prior to break-up (horizontal axis), for the 52 % of cases in which the bubble was in the volume resolved with PIV at this time.

To gauge the effect of inhomogeneity in the generated turbulence, we consider how the local turbulence intensity experienced by the bubble (in a Lagrangian sense) varies over timescales relevant to the bubble's break-up. Ideally, a bubble would not be advected through statistically inhomogeneous turbulence during the course of its break-up. Denoting the Hinze scale at the bubble's location at time $t$ as $d_{H}(t)$, figure 13(b) shows the Hinze scale at the break-up location $d_{H}(t_0)$ as a function of the Hinze scale at the bubble's location one bubble-scale eddy turnover time prior, $d_{H}(t_0 - T_{turb}(d_0))$ for the 52 % of observed break-ups in which the bubble is inside the volume resolved by PIV (in which we are able to compute $d_{H}$) at this point in time. There is little difference between $d_{H}(t_0)$ and $d_{H}(t_0 - T_{turb}(d_0))$, suggesting that the local turbulence characteristics experienced by the bubble do not change appreciably during the break-up, and that the turbulence is homogeneous over scales relevant to the break-up.

4.3. Child size distribution

Now, we compute the dimensionless bubble child size distributions conditioned on the approximate dimensionless parent bubble size, $\mathcal {P}_d(d/d_{H};d_0/d_{H})$. The data is averaged over three ranges of $d_0/d_{H}$ (the ranges between [0.3:1.55], [1.55:1.93] and [1.93:3.70]), and results are shown in figure 14(a). As $d_0/d_{H}$ is increased, the dependence of $\mathcal {P}_d$ on $d/d_{H}$ becomes steeper. The dashed line gives the $\mathcal {P}_d(d/d_{H};d_0/d_{H}) \propto (d/d_{H})^{-3/2}$ scaling, which is approached for large $d_0/d_{H}$ due to the production of small bubbles by capillary instabilities (Rivière et al. Reference Rivière, Ruth, Mostert, Deike and Perrard2022). Qualitatively, the child size distribution for smaller parent bubbles is flatter near the Hinze scale, whereas that for larger parent bubbles increases more rapidly with decreasing bubble size as a power-law relationship. Note that the child size distribution is defined so that it integrates to the average number of child bubbles formed.

Figure 14. Dimensionless bubble break-up child size distributions for various approximate values of $d_0/d_{H}$. The value give for each curve (which is denoted by the notch on the horizontal axis) corresponds to the mean value of $d_0/d_{H}$ for that curve. (a) The distributions of child bubble diameter normalised by the Hinze scale. (b) The volumetric child size distribution, with child bubble volumes normalised by the parent bubble volume, which is approximated as the sum of the resolved child bubble volumes.

This representation of the child size distribution masks the large number of bubbles formed very close to the parent bubble size. To capture these small bubbles, we also compute the volumetric child size distribution, normalised by the volume of the parent bubble $V_0$. As the determination of the volumes of larger bubbles is difficult given their deformations, we approximate the parent bubble volume as the sum of the volumes of the child bubbles, and consider $(d/d_0)^3 \approx d^3 / \sum _{i=1}^m d_i^3$ (Vejražka et al. Reference Vejražka, Zedníková and Stanovský2018). The distribution of these dimensionless volumes is shown in figure 14(b), exhibiting a $\cup$ shape that is not strongly dependent on $d_0/d_{H}$ (though we again see increased small bubble production with larger $d_0/d_{H}$). For smaller parent bubbles, which largely undergo binary break-up, the distribution is nearly symmetric about $(d/d_0)^3 = 0.5$. The large values of these distributions near 1 suggest that in many break-up events, small bubbles are ‘torn off’ of the parent bubble, without inertial deformation producing multiple child bubbles of sizes comparable to that of the parent. We note that the resolution of our experiment (in which the smallest bubble we can detect is approximately $200\,\mathrm {\mu }\textrm {m}$ in diameter) limits the number of bubbles detected.

4.4. Small bubble production without significant inertial deformation

In many of the break-ups we observe in the large cavity disintegration and individual break-up experiments, small bubbles were seen to be ‘torn off’ from a parent bubble, without an appreciable large-scale deformation to the parent bubble. These events are reminiscent of tip-streaming (Montanero & Gañán-Calvo Reference Montanero and Gañán-Calvo2020). This phenomenon is evidenced by the right-hand side of the $\cup$-shaped child size distributions shown in figure 14(b), as a child bubble that is nearly the size of the parent is the signature of such break-ups. To understand these events, we present in figure 15 a qualitative discussion of the dynamics of individual splitting events. For each splitting event, we compare the velocity of the parent bubble at break-up $\boldsymbol {v}_{parent}$ (denoted by the grey arrow in figure 15a) with the displacement between the parent bubble's final position $\boldsymbol {x}_{parent}$ (the grey circle) and the initial positions at which the child bubbles are detected $\boldsymbol {x}_{child}$ (the black circles). The child bubble's initial detected position ahead of or behind the parent bubble, $\kappa = \boldsymbol {v}_{parent} \boldsymbol {\cdot } (\boldsymbol {x}_{child} - \boldsymbol {x}_{parent}) / (u' L_{int})$, normalised by turbulence quantities, is then computed and is plotted against the child bubble's size relative to the parent size in figure 15(b). The colour of each marker denotes the size of the splitting event's parent bubble relative to the Hinze scale. The black line shows the expected value of $\kappa$ given the normalised child bubble size. Smaller child bubbles (with $d_{child}/d_{parent} < 0.6$, below which the mean value of $\kappa$ becomes negative) tend to be left in the wake of the parent bubble ($\kappa <0$), whereas larger child bubbles tend to be produced ahead of the parent bubble $(\kappa > 0)$.

Figure 15. Statistics of the positions of bubbles after splitting events. (a) A sketch of a splitting event involving small bubble production, including the parent bubble velocity at break-up and the final and initial positions, respectively, of the parent and child bubbles. (b) The initial child bubble position relative to the parent bubble's motion, $\kappa$, for each splitting event (circles), as well as the mean value conditioned on the normalised child bubble size (black line). Here $\kappa <0$ denotes bubble production behind the parent bubble, whereas $\kappa >0$ denotes bubble production ahead of the parent bubble.

Although the conceptual picture for break-up discussed in § 3.3 describes the role of capillarity during break-ups involving large-scale deformations, it is likely that break-ups solely involving small bubble production are also regulated by capillarity: in these cases, a turbulent motion smaller than the parent bubble may succeed in producing a ligament which extends off of one side of the parent, and this ligament may pinch apart into many small bubbles in a capillary instability as it is retracted back into the bulk of the parent bubble. Specifically, figure 15 suggests that the bulk of a bubble may often be swept forward by a turbulent eddy, and the trailing ligament may become unstable as it ‘catches up’ with the rest of the parent bubble. Similar to the framework presented in § 3.3, the process is initiated by a turbulent deformation to the parent and ends with the capillary instability of a ligament involved in the deformation.

5.1. Physical ideas

The experiments presented in §§ 3 and 4, taken together with the existing literature, point to three important time scales that must be considered in developing a population balance model: the inverse of the break-up frequency, the break-up duration and the capillary pinching time.

The longest of these is the typical duration until a break-up occurs; that is, the inverse of the break-up frequency, $1/\omega (d_0)$. This time scale will control how many break-up events will occur over a given time and will be a function of $d_0/d_{H}$. The second timescale is that over which a break-up typically occurs, or the event duration (i.e. lasting from the start of the deformation until the child bubbles have all been formed), and will also be a function of $d_0/d_{H}$. The break-ups taking the longest time will be those instigated by the largest eddies capable of causing break-up, which are taken to be those at the parent bubble's scale (Luo & Svendsen Reference Luo and Svendsen1996). Thus, an upper bound and typical scale of the break-up duration is taken to be the eddy turnover time at the parent bubble's scale, $T_{turb}(d_0) = \epsilon ^{-1/3} d_0^{2/3}$, in agreement with experimental and numerical observations of the time over which bubbles are deformed prior to break-up (Risso & Fabre Reference Risso and Fabre1998; Martínez-Bazán et al. Reference Martínez-Bazán, Montañes and Lasheras1999b; Rivière et al. Reference Rivière, Mostert, Perrard and Deike2021). The final timescale we consider is that of the capillary instabilities of gas ligaments that produce a small child bubble of size $d$, which will occur over the capillary timescale of that child bubble, $T_{cap}(d) = (\rho /\gamma )^{1/2} d^{3/2} / (2 \sqrt {3})$ (Rivière et al. Reference Rivière, Ruth, Mostert, Deike and Perrard2022).

From these three relevant time scales, we define three types of events. At the shortest time, we define the individual binary splitting events. For the production of bubbles with $d \ll d_{H}$, we have $T_{cap}(d) \ll T_{turb}(d_0)$. At the eddy turnover time, we define a break-up as being a sequence composed of all the splitting events occurring in a time bounded by $\Delta T_{break\textrm {-}up} = T_{turb}(d_0)$, which permits the production of more than two bubbles in a single event (similar to the definition used for drop break-ups by Solsvik, Maaß & Jakobsen Reference Solsvik, Maaß and Jakobsen2016). Finally, following the nomenclature from Hinze (Reference Hinze1955), a disintegration is a longer-duration process involving an arbitrary number of break-ups.

These timescales are sketched in figure 16, which illustrates two break-up events that stem from a bubble of diameter $d_{A}$ encountering turbulence. The deformation to the parent bubble that instigates the break-up is assumed to happen within a time $T_{turb}(d_{A})$ before the first bubble splits from the parent. Then, within an additional time bounded by $T_{turb}(d_{A})$, subsequent splitting events occur due to capillary instabilities arising from the deformation. One such instability produces a bubble with diameter $d_{C}$, and the time over which this instability develops is set by the capillary timescale at the smaller child bubble size, $T_{cap}(d_{C})$. Later on, one of the child bubbles produced in the first break-up, with diameter $d_{B}$, itself breaks up.

Figure 16. Sketch of two bubble break-ups and the associated timescales, with $T_{turb}(d) = \epsilon ^{-1/3}d^{2/3}$ and $T_{cap}(d) = (\rho /\gamma )^{1/2} d^{3/2} / (2 \sqrt {3})$. The grey vertical lines denote the times associated with each of the two break-ups. The shaded region to the left bounds the time over which the deformation to the parent bubble is assumed to occur (the turbulent timescale at the parent bubble size) and the region to the right of the line bounds the time over which the subsequent splitting events are assumed to occur (also taken to be the same turbulent timescale). During the subsequent splitting events, the capillary timescale at the size of the smaller child bubble sets the time over which the splitting event occurs (Rivière et al. Reference Rivière, Ruth, Mostert, Deike and Perrard2022). The time between break-ups is set by the inverse of the break-up frequency $\omega$ of the bubble which is to break, which we address later in the paper.

Using these ideas, we propose a population balance model that integrates these physical elements and models the evolution of a bubble size distribution with a Boltzmann transport equation using the bubble size as an internal coordinate. The population balance model considers a break-up rate kernel $f$, constructed from child size distributions computed through a Monte Carlo approach (constrained by results from experiments and DNSs, informing the number of children and the shape of the distribution) and a parent bubble break-up frequency taken from the literature. With the kernel defined, we integrate the model in time to simulate the evolution of the size distribution during a cavity disintegration and compare to our experimental data.

5.2. Population balance modelling

In a confined region of homogeneous turbulence, the transient evolution of the absolute dimensionless volumetric bubble size distribution $\mathcal {N}_V(\tilde {V}) = N_V(V) V_{H} = \mathcal {N}(d/d_{H}) / (3 (d/d_{H})^2)$, where $N_V(V)$ is the absolute dimensional volumetric size distribution, $\tilde {V} = V/V_{H}$, and $\tilde {t} = t/T_{int}$, is described by

(5.1)\begin{equation} \frac{\partial \mathcal{N}_V(\tilde{V},\tilde{t})}{\partial \tilde{t}} ={-} \frac{\mathcal{N}_V(\tilde{V},\tilde{t})}{\langle m \rangle(\tilde{V})} \int_0^{\tilde{V}} \tilde{f}(\tilde{\delta};\tilde{V})\,\mathrm{d}{\tilde{\delta}} + \int_{\tilde{V}}^\infty \mathcal{N}_V (\tilde{\varDelta},\tilde{t}) \tilde{f}(\tilde{V};\tilde{\varDelta})\,\mathrm{d} \tilde{\varDelta}, \end{equation}

where the first term on the right-hand side gives the rate of consumption of bubbles of volume $\tilde {V}$ due to their break-ups and the second term on the right-hand side gives the rate of production of bubbles of volume $\tilde {V}$ due to the break-ups of larger bubbles (Martínez-Bazán et al. Reference Martínez-Bazán, Rodríguez-Rodríguez, Deane, Montañes and Lasheras2010). The break-up kernel $\tilde {f}(\tilde {\delta };\tilde {\varDelta }) = f(\delta ;\varDelta ) V_{H} T_{int}$ can be decomposed into a parent break-up frequency and volumetric child size distribution with $\tilde {f}(\tilde {\delta };\tilde {\varDelta }) = \tilde {\omega }(\tilde {\varDelta }) \tilde {p}(\tilde {\delta };\tilde {\varDelta })$, with the dimensionless break-up frequency $\tilde {\omega }(\tilde {\varDelta }) = \omega (d) T_{int}$ and dimensionless volumetric child size distribution $\tilde {p}(\tilde {\delta };\tilde {\varDelta }) = p(\delta ;\varDelta ) V_{H}$. Thus, we can move $\tilde {\omega }(\tilde {\varDelta })$ outside the integral in the first term on the right-hand side and invoke $\int _0^{\tilde {V}} \tilde {p}(\tilde {\delta };\tilde {V})\,\mathrm {d} \tilde {\delta } = \langle m \rangle (\tilde {V})$, with $\langle m \rangle (\tilde {V})$ the average number of bubbles formed in the break-up of a bubble of volume $\tilde {V}$, to express the bubble consumption term as simply $-\mathcal {N}_V(\tilde {V},\tilde {t}) \tilde {\omega }(\tilde {V})$. Note that we define $\tilde {p}(\tilde {\delta };\tilde {\varDelta })$ so that it integrates over $\tilde {\delta }$ to the average number of child bubbles formed by the break-up of a bubble of volume $\tilde {\varDelta }$.

5.3. Construction of the child size distributions

We develop a parameterisation of the break-up volumetric child size distribution $\tilde {p}(\tilde {\delta };\tilde {\varDelta })$ that accounts both for child bubbles produced by both the slower inertial mechanism (occurring over the eddy turnover time) and the faster capillary pinching mechanism (occurring over the capillary timescale of the small child bubbles) using a Monte Carlo approach. We consider a set of rules constrained by our experimental and numerical observations describing the outcomes of individual break-up events, then aggregate the outcomes of these events into child size distributions.

5.3.1. Statistics on the number of child bubbles formed

A key step in modelling each break-up is to constrain the distribution of the number of bubbles formed in each event. To this end, we first consider the data from our dynamical experiments given in § 4. The average number of child bubbles larger than the experimentally-resolvable minimum size $d_{min}/d_{H} \approx 0.07$, $\langle m \rangle$, is shown in figure 17(a). As the parent bubble increases in size, more child bubbles are typically produced. Given the steep dependence of $\tilde {p}(\tilde {\delta };\tilde {\varDelta })$ on $\tilde {\delta }$, we must qualify each observation of $m$ with the minimum resolved bubble size to better enable comparisons between different experiments. For compactness, however, we take all $m$ values to be the number of resolved bubbles larger than $0.07 d_{H}$ unless otherwise noted. Our experimental observations of $\langle m \rangle$, binned by $d_0/d_{H}$, are shown in the black squares and the grey region around them bounds $\pm$ one half of a standard deviation around the mean.

Figure 17. Experimental data on the number of child bubbles formed in each break-up. (a) The average number of resolved bubbles (with $d_{min}/d_{H} = 0.07$) formed in each break-up event $\langle m \rangle$ as a function of the dimensionless parent bubble size. The shaded region shows $\pm$ one half of a standard around the mean for our dynamical data. Open circles give data from the disintegration of the three largest cavities and closed circles give those data with an adjustment for the differing spatial resolution. Open stars give data from DNSs from Rivière et al. (Reference Rivière, Mostert, Perrard and Deike2021) and closed stars give those data with the spatial resolution adjustment. The open grey markers give data from experiments reported by Vejražka et al. (Reference Vejražka, Zedníková and Stanovský2018). The thick orange line is the parameterisation given in (5.2). (b) The p.d.f.s of $m'/\langle m' \rangle =(m-m_{min})/(\langle m \rangle - m_{min})$ for the experiments (squares) and DNSs (stars), along with the exponential fit employed in the Monte Carlo simulations.

Next, to consider the number of bubbles produced in the break-ups of larger bubbles, we turn to data from the disintegration of the three largest cavities presented in § 3. With the assumption that the initial splitting event happens nearly instantly after the bubble is released into the turbulence, to apply the same definition of the duration of the break-up, we define $\langle m \rangle$ for this dataset as the number of resolved bubbles present after one eddy turnover time $T_{turb}(d_0) = \epsilon ^{-1/3} d_0^{2/3}$ has elapsed after the cavity release, which are denoted by the open circles in figures 8(a) and 17(a). We invoke (C1) to apply a slight adjustment to these numbers in order to extrapolate results to the finer spatial resolution of the tracked break-up experiment, as discussed in Appendix C. The number of bubbles in the extrapolated range constitutes about $30\,\%$ of those in the observable range. These adjusted values are shown as the filled-in light blue circles in figure 17(a).

We have additionally re-analysed the DNSs of bubbles breaking in homogeneous, isotropic turbulence presented in Rivière et al. (Reference Rivière, Mostert, Perrard and Deike2021, Reference Rivière, Ruth, Mostert, Deike and Perrard2022), tracking the bubble break-up events in a similar way to what has been done on the experimental data in § 4. From these DNSs, we can compute the average number of bubbles formed per event as a function of the parent bubble size, included in figure 17(a) as the red star markers. Open stars give the original observations, for which $d_{min}/d_{H} = 0.25$, whereas the filled-in stars give the number adjusted for the spatial resolution. Note that while we consider ${We}_{c}=1$ for the experimental data, the value of $d_{H}$ for the DNS is given by ${We}_{c}=3$ (Rivière et al. Reference Rivière, Mostert, Perrard and Deike2021).

Finally, as a comparison, the open grey markers show the (un-adjusted) number of bubbles detected experimentally in break-ups by Vejražka et al. (Reference Vejražka, Zedníková and Stanovský2018), in which break-ups varied in $\epsilon$ and $d_0$ (which is denoted by the marker shape). As shown in their paper, once collapsed to $d_0/d_{H}$, the dependence on the dimensional bubble size nearly disappears.

The four datasets (our two experiments, those from Vejražka et al. (Reference Vejražka, Zedníková and Stanovský2018), and DNSs from Rivière et al. Reference Rivière, Mostert, Perrard and Deike2021) produce a coherent picture regarding the number of bubbles formed. When $d_0/d_{H}$ is small, break-ups tend to be binary, producing on average two child bubbles after $T_{turb}(d_0)$. As $d_0/d_{H}$ increases, the number of child bubbles increases. Surface tension is less effective at preventing the severe deformation of larger bubbles, leading to more complex deformed bubble shapes that yield a greater number of child bubbles. The orange curve in figure 17(a) shows a fit to the data of the form

(5.2)\begin{equation} \langle m \rangle = m_{min} + \frac{(d_0/d_{H})^{b_2}}{b_1}, \end{equation}

where $m_{min} = 2$ and the fit constants are $b_1 = 4$ and $b_2 = 2.3$.

Figure 17(b) compiles experimental and DNS data on the distribution of the number of child bubbles produced for increasing $d_0/d_{H}$. The probability density functions (p.d.f.s) of $m'/\langle m' \rangle$ are well-described by an exponential function $\exp ({-m'/\langle m' \rangle })$, with $m'=m-m_{min}$ and $\langle m' \rangle = \langle m \rangle -m_{min}$, for both the experiments (shown as the squares) and DNSs (shown as the stars). Thus, for any parent bubble size we can write the p.d.f. of $m'$ as an exponential distribution,

(5.3)\begin{equation} r(m';d_0/d_{H}) = \frac{\exp({-}m' / \langle m' \rangle)}{\langle m' \rangle},\quad m'>0, \end{equation}

with $\langle m' \rangle +m_{min}$ the mean number of children, a function of the parent bubble size.

5.3.2. A stochastic model for each break-up

The Monte Carlo approach involves running many iterations of a stochastic model and developing a statistical representation of the aggregated results. Each discrete simulation of a break-up mirrors the physical processes involved: the bubble, sketched in figure 18(a), is first deformed into two lobes, shown in (b), and then some number of capillary bubbles are created as the neck separating the lobes collapses to create the two inertial child bubbles.

Figure 18. Process of simulating one break-up with the Monte Carlo approach for a bubble of volume $\tilde {\varDelta }$, shown in (a). (b) First, the bubble is taken to be deformed into two lobes, separated by a neck of gas. (c) Next, the sizes of the $m' = m-2$ capillary bubbles are picked from a $\tilde{\delta}_{cap}^\alpha$ distribution. (d) Finally, the sizes of the two inertial bubbles $\tilde {\delta }_{{inertial},i}$ are picked from a uniform distribution over the remaining parent bubble volume (that which has not gone to the capillary bubbles).

For each iteration (i.e. one simulated breakup) at a given value of $\tilde {\varDelta }$, we first define the number of bubbles $m$ that will be produced by picking a value of $m'$ from the distribution $r(m';d_0/d_{H})$ given by (5.3), adding $m_{min} = 2$, and rounding to the nearest integer. We pick $\tilde {\delta }_{min} = 0.07^3$ in order to match the experimental dataset on which the parameterisation of $\langle m \rangle$ is based. As we show, once the p.d.f.s have been constructed for this given value of $\tilde {\delta }_{min}$, it is straightforward to extend them to lower or higher values of $\tilde {\delta }_{min}$.

For cases in which $m\geq 3$, the capillary mechanism produces $m' = m-2$ bubbles, whose sizes follow a $\propto \tilde {\delta }^{\alpha }$ distribution with $\alpha =-7/6$ (corresponding to the $\mathcal {P}_d(d/d_{H};d_0/d_{H}) \propto (d/d_{H})^{-3/2}$ scaling described by Rivière et al. (Reference Rivière, Ruth, Mostert, Deike and Perrard2022), as distributions in diameter are related to those in volume by $\mathcal {P}_d(d/d_{H};d_0/d_{H}) = 3 (d/d_{H})^2 \tilde {p}(\tilde {\delta };\tilde {\varDelta })$ (Martínez-Bazán et al. Reference Martínez-Bazán, Rodríguez-Rodríguez, Deane, Montañes and Lasheras2010; Qi et al. Reference Qi, Masuk and Ni2020)). As is sketched in figure 18(c), the volume $\tilde {\delta }_{{cap},i}$ of capillary bubble $i$ is picked from a power-law distribution with slope $\alpha$, bounded between $\tilde {\delta }_{min}$ and the maximum allowable volume for a capillary bubble given the previously produced bubbles, $\tilde {\delta }_{{cap,max},i}$. For the production of the first capillary bubble, we set $\tilde {\delta }_{{cap,max},1} = \tilde {\varDelta }$ (noting that the steep slope of $\mathcal {P}_d(d/d_{H};d_0/d_{H})$ with respect to $d/d_{H}$ makes the production of capillary bubbles this large uncommon). For the production of the remaining capillary bubbles, we set $\tilde {\delta }_{{cap,max},i} = \tilde {\varDelta } - \sum _{j=1}^{i-1} \tilde {\delta }_{{cap},j}$. At each step of the process, if $\tilde {\delta }_{{cap},i}$ is greater than $\tilde {\delta }_{{cap,max},i}/2$, we replace it with $\tilde {\delta }_{{cap,max},i}/2 - \tilde {\delta }_{{cap},i}$, such that for any splitting event, the smaller of the two produced does not further split.

Once the volumes of the $m'$ capillary bubbles are specified, we must determine the volumes of the two inertial bubbles. To that end, we first compute the portion of the parent bubble volume that has gone to the capillary bubbles, $\chi _{cap} = \sum _{i=1}^{m'} \tilde {\delta }_{{cap},i} / \tilde {\varDelta }$. The size of the first of the two inertial child bubbles $\tilde {\delta }_{{inertial},1}$ is drawn uniformly from the remaining bubble volume, $(1-\chi _{cap}) \tilde {\varDelta }$, and the second is taken as its complement, $\tilde {\delta }_{{inertial},2} = (1-\chi _{cap}) \tilde {\varDelta } - \tilde {\delta }_{\mathrm {inertial},1}$. Once this is done, the volumes of all child bubbles produced in this single break-up have been determined. The uniform distribution for $\tilde {\delta }_{{inertial},1}$ is chosen for simplicity, and future work is warranted to better describe the distribution of the sizes of the large bubbles formed by inertial deformations.

5.3.3. Aggregation of simulated break-ups into child size distributions

For a given value of $d_0/d_{H}$ (or the equivalent normalised volume $\tilde {\varDelta } = (d_0/d_{H})^3$), the process of simulating one break-up stochastically is repeated $n_{MC} = 10^5$ times.

For each $\tilde {\varDelta }$, the sizes of the bubbles produced in each of the $n_{MC}$ events are aggregated, and the distribution of all these child bubbles defines the volumetric child size distribution $\tilde {p}(\tilde {\delta };\tilde {\varDelta })$. The distribution is normalised such that $\int _0^{\tilde {\varDelta }} \tilde {p}(\tilde {\delta };\tilde {\varDelta }) = \langle m \rangle (\tilde {\varDelta })$, with $\langle m \rangle (\tilde {\varDelta })$ the average number of bubbles formed. As the size distribution is aggregated from geometrically plausible break-ups, it itself must satisfy any constraints relating to the sizes of the bubbles produced. Figure 19(a) shows the volumetric child size distributions for five values of $\tilde {\varDelta }$. When $\tilde {\varDelta }$ is small, the child size distribution is nearly uniform, as the capillary production mechanism is negligible for small bubbles; for moderate $\tilde {\varDelta }$, the child size distribution exhibits a $\tilde {p} \propto \tilde {\delta }^\alpha$ scaling for small bubbles, while remaining close to flat for bubbles near the parent bubble size. For even larger bubbles, for which the capillary production mechanism is the most effective, the entire distribution approaches a $\tilde {\delta }^\alpha$ scaling.

Figure 19. Volumetric child size distributions constructed via the Monte Carlo approach. (a) Volumetric child size distributions $\tilde {p}(\tilde {\delta };\tilde {\varDelta })$ for five values of the parent bubble size $\tilde {\varDelta }$. Distributions compiled from the Monte Carlo simulations are given by the thin lines, whereas the thick fainter lines give the fits using (5.4). The two components of the fit form of the distribution are illustrated for $\tilde {\varDelta } = 10$. (b) The average capillary fraction $\langle \chi _{cap} \rangle$ calculated from the ensemble of simulations, as a function of the parent bubble size. (c) Fit values of the exponent $\gamma (\tilde {\varDelta })$ employed in (5.4). Data for the curves in this figure and Python code to construct $\tilde {p}(\tilde {\delta };\tilde {\varDelta })$ are available online as supplementary material.

For each $\tilde {\varDelta }$, we also obtain $\langle \chi _{cap} \rangle (\tilde {\varDelta })$, shown in figure 19(b), by averaging the portion of the parent bubble volume going to the capillary child bubbles $\chi _{cap}$ over the $n_{MC}$ events. When $\tilde {\varDelta } \ll 1$, $\chi _{cap} \approx 0$, and essentially all of the parent bubble volume goes to the two inertial child bubbles. With larger $\tilde {\varDelta }$, $\chi _{cap}$ increases, reaching $\chi _{cap} = 0.1$ at $\tilde {\varDelta } = 60$. Even at $\tilde {\varDelta } = 1000$, less than half of the parent bubble volume goes to the capillary bubbles.

We then fit each volumetric child size distribution as a sum of two components, each stemming from one of the two mechanisms of child bubble production,

(5.4)\begin{equation} \tilde{p}(\tilde{\delta};\tilde{\varDelta}) = \underbrace{a(\tilde{\varDelta}) \tilde{\delta}^{\gamma (\tilde{\varDelta})}}_{inertial\ mechanism} + \underbrace{b (\tilde{\varDelta}) \tilde{\delta}^\alpha}_{capillary\ mechanism}, \end{equation}

with $\alpha =-7/6$ set by the distribution from which the capillary bubbles are picked and $\gamma (\tilde {\varDelta })$ chosen to match the aggregated Monte Carlo simulation data. The two remaining coefficients, $a(\tilde {\varDelta })$ and $b(\tilde {\varDelta })$, are constrained by the volume going to bubbles produced by each mechanism, leading to

(5.5)$$\begin{gather} a(\tilde{\varDelta}) = (1 - \langle \chi_{cap} \rangle) \left(\frac{(\gamma + 2) \tilde{\varDelta} }{\tilde{\varDelta}^{\gamma+2} - \tilde{\delta}_{min}^{\gamma+2}} \right), \end{gather}$$

(5.6)$$\begin{gather}b(\tilde{\varDelta}) = \langle \chi_{cap} \rangle \left(\frac{(\alpha+2)\tilde{\varDelta}}{\tilde{\varDelta}^{\alpha+2} - \tilde{\delta}_{min}^{\alpha+2}} \right). \end{gather}$$

The fits to each child size distribution with (5.4) are shown as the faint, thick lines in figure 19(a).

Figure 19(c) shows the evolution of the exponent $\gamma (\tilde {\varDelta })$ describing the inertial production mechanism. Values of $\langle \chi _{cap} \rangle (\tilde {\varDelta })$ and $\gamma (\tilde {\varDelta })$, which together contain all the necessary information about the child size distributions, are stored for many values of $\tilde {\varDelta }$. To implement the child size distributions in a population balance model, we interpolate $\langle \chi _{cap} \rangle (\tilde {\varDelta })$ and $\gamma (\tilde {\varDelta })$ for a given value of $\tilde {\varDelta }$. Numerical data giving the child size distributions for a range of parent bubble sizes are provided as supplementary material at https://doi.org/10.1017/jfm.2022.604, and Python code to simulate break-ups and create the child size distributions is available at https://github.com/DeikeLab/mc-child-size-distributions.

5.4. Parameterisation of the break-up frequency

The next step is to parameterise how often the break-ups will occur. It is not possible to compute a break-up frequency from our experiments which resolve break-ups dynamically, because we only record data when bubbles are observed to break. Instead, using an approach that has been successfully applied to the break-up of oil droplets in turbulent jets (Aiyer et al. Reference Aiyer, Yang, Chamecki and Meneveau2019; Aiyer & Meneveau Reference Aiyer and Meneveau2020), we integrate the effects of eddies smaller than the parent bubble size (each of dimensional diameter $d_{e}$) which contribute to break-up (Prince & Blanch Reference Prince and Blanch1990; Tsouris & Tavlarides Reference Tsouris and Tavlarides1994), yielding

(5.7)\begin{align} \tilde{\omega}(\tilde{\varDelta}) &= K \frac{d_{H}}{L_{int}} \int_0^{d_0/d_{H}} \frac{\rm \pi}{4} \left(\frac{d_0}{d_{H}} + \frac{d_{e}}{d_{H}} \right)^2 \tilde{u}_{turb}(d_{e}/d_{H}) \nonumber\\ &\quad \times\left(\frac{d_{e}}{d_{H}} \right)^{{-}4} \varOmega(d_{e}/d_{H};d_0/d_{H})\,\mathrm{d} (d_{e}/d_{H}), \end{align}

where $K$ is an order-1 constant we adjust, $\tilde {u}_{turb}(d_{e}/d_{H}) = C_2^{1/2} \epsilon ^{1/3} d_{e}^{1/3} / u' = C_\epsilon ^{1/3} (d_{e}/d_{H})^{1/3} (d_{H}/L_{int})^{1/3}$ is the dimensionless turbulent velocity scale of the eddy, $(d_{e}/L_{int})^{-4}$ is the approximate dimensionless eddy density (Solsvik et al. Reference Solsvik, Maaß and Jakobsen2016) and $\varOmega (d_{e}/d_{H};d_0/d_{H})$ is the break-up efficiency given the eddy and bubble sizes. Neglecting viscous effects (given the low viscosity of air bubbles), the break-up efficiency, which gives the probability that an eddy has sufficient energy to overcome surface tension, is taken as the inverse of the exponential of the ratio between the average change in surface energy associated with the break-up $E_\sigma (d_0)$ and the kinetic energy of the eddy $E_{eddy}(d_{e})$, $\exp (- E_{eddy}(d_{e}) / E_\sigma (d)$). The average surface energy change is given dimensionally by

(5.8)\begin{equation} E_\sigma(d_0) = \frac{\sigma {\rm \pi}}{4} \left( \int_{\delta_{min}}^{{\rm \pi} d_0^3/6} p(\delta; {\rm \pi}d_0^3/6) \delta^{2/3}\,\mathrm{d} \delta - d_0^2 \right) = \varGamma {\rm \pi}\sigma d_0^2 / 4, \end{equation}

with the proportional change in surface area due to break-up $\varGamma$ dependent on the form of the child size distribution according to

(5.9)\begin{equation} \varGamma(\tilde{\varDelta}) = \frac{\displaystyle\int_{\tilde{\delta}_{min}}^{\tilde{\varDelta}} \tilde{p}(\tilde{\delta};\tilde{\varDelta}) \tilde{\delta}^{2/3}\,\mathrm{d} \tilde{\delta}}{\tilde{\varDelta}^{2/3}} - 1. \end{equation}

The kinetic energy of the eddy is given by $E_{eddy}(d_{e}) = ({\rm \pi} /4) \rho d_{e}^3 C_2 (\epsilon d_{e})^{2/3}$. Expressed in our non-dimensional units, the break-up efficiency is then

(5.10)\begin{equation} \varOmega(d_{e}/d_{H};d_0/d_{H}) = \exp \left( - \frac{ \varGamma(\tilde{\varDelta}) (d_0/d_{H})^2}{{We}_{c} (d_{e}/d_{H})^{11/3}} \right), \end{equation}

with the critical Weber number ${We}_{c}$ necessary to link the scales of the bubble and the turbulence.

With each component specified, (5.7) is evaluated numerically and is shown in figure 20, using $K=2$ picked through a comparison with the experimental data given in § 3. The break-up rate increases as bubbles approach the Hinze scale and then plateaus due to two competing effects: while larger bubbles are susceptible to a wider range of turbulent scales that may cause break-up, they tend to break into many more bubbles than smaller bubbles, leading to a greater surface energy term in (5.7). This means that although more eddies are interacting with the parent bubble, each is less likely to have sufficient energy to cause a break-up.

Figure 20. The parent bubble break-up rate $\tilde {\omega }$ as a function of its volume $\tilde {\varDelta }$, computed using the value of $d_{H}/L_{int}$ for our dataset. The black line shows the parent bubble break-up frequency given by (5.7). The thicker grey line gives the inverse of the eddy turnover time at the parent bubble scale, which is taken to be the upper limit in the duration of each break-up event. The dotted orange line gives the inverse of the capillary timescale at the parent bubble scale.

The thicker grey line in figure 20 gives the inverse of the turbulent turnover time at the parent bubble scale, which we take to set the duration of each break-up event. The break-up frequency is thus consistent with the break-up duration, since $\tilde {\omega }(\tilde {\varDelta }) = \omega T_{int}$ being strictly less than $T_{int} / T_{turb}(d_0)$ means that the typical duration of a break-up is never longer than the typical time between such break-ups. Finally, the dotted orange line gives the inverse of the (dimensionless) capillary timescale at the parent bubble scale, $T_{int}/T_{cap}(d_0)$, showing that capillary effects happen faster than both the break-up duration and time between break-ups (up until the largest bubbles we consider). The capillary pinching events responsible for sub-Hinze bubble creation thus occur over even shorter durations, as the capillary timescales of the small child bubbles formed will be much faster than that of the parent bubble.

5.5. Summary of parameters involved in the model

To summarise, table 2 lists each parameter in the model and explains how each is determined.

Table 2. Parameters involved in the bubble break-up model, their physical origin and the experimental/numerical data constraints.

5.6. Model comparison with transient air cavity disintegration data

With $\tilde {p}(\tilde {\delta };\tilde {\varDelta })$ and $\tilde {\omega }(\tilde {\varDelta })$ now fully specifying $\tilde {f}(\tilde {\delta };\tilde {\varDelta })$, we can simulate the turbulent disintegration of cavities we studied experimentally in § 3 by picking the appropriate initial condition for each (i.e. $\mathcal {N}(d/d_{H})$ giving one bubble of size $d_0/d_{H}$) and integrating (5.1) in time. Figure 21 compares the experimental and modelled values of the dimensionless bubble size distribution $\mathcal {N}(d/d_{H})$ at $t/T_{int} = 1$ and 3 for each value of $d_0/d_{H}$, with $d_0/d_{H}=2.12$ in (a) and $d_0/d_{H}=8.30$ in figure 21( f).

Figure 21. Comparisons of the experimental and modelled values of $\mathcal {N}(d/d_{H})$ at $t/T_{int} = 1$ and $3$ for each value of $d_0/d_{H}$. The dotted vertical line gives the value of $d_0/d_{H}$ for each condition. Dotted lines give the $\mathcal {N}(d/d_{H}) \propto (d/d_{H})^{ -3/2}$ sub-Hinze scaling. Good agreement between the measured and modelled distributions is observed for the full range of $d_0/d_H$ and times.

First, the model accurately reproduces the observed magnitudes of the size distributions near the Hinze scale, both in time and in the initial cavity size. Second, an $\mathcal {N}(d/d_{H}) \propto (d/d_{H})^{-3/2}$ scaling is approached for $d/d_{H} < 1$ with larger $d_0/d_{H}$, and this scaling is adopted more rapidly with larger cavities. With $d_0/d_{H} = 2.12$ and 2.84, shown in figures 21(a) and 21(b), $\mathcal {N}(d/d_{H})$ is flat near the Hinze scale at $t/T_{int} = 1$, as the modelled child size distributions for parent bubbles of these cavity sizes are largely flat (as shown in figure 19a). With larger parent cavities, the sub-Hinze distribution steepens as the capillary mechanism contributes more significantly to the child size distributions for parent bubbles of these larger cavity sizes (as evidenced in the $\langle \chi _{cap} \rangle (\tilde {\varDelta })$ curve shown in figure 19b).