# Immiscible Rayleigh–Taylor turbulence: Implications for bacterial degradation in oil spills

Edited by Andrea Rinaldo, Ecole Polytechnique Federale de Lausanne, Lausanne, Switzerland; received July 14, 2023; accepted January 23, 2024

## Significance

The immiscible Rayleigh–Taylor turbulence is a paradigm of buoyancy-driven flows that can serve as a model for the dynamics of deep-water oil spills. Oil from such spills is degraded by marine bacteria, which colonize oil droplets and degrade their hydrocarbon content at a rate proportional to the oil–water interface area. Here, we develop an understanding of the oil–water interface dynamics based on Rayleigh–Taylor turbulence that sheds light on the efficiency of biodegradation in oil spills.

## Abstract

An unstable density stratification between two fluids mixes spontaneously under the effect of gravity, a phenomenon known as Rayleigh–Taylor (RT) turbulence. If the two fluids are immiscible, for example, oil and water, surface tension prevents intermixing at the molecular level. However, turbulence fragments one fluid into the other, generating an emulsion in which the typical droplet size decreases over time as a result of the competition between the rising kinetic energy and the surface energy density. Even though the first phenomenological theory describing this emulsification process was derived many years ago, it has remained elusive to experimental verification, hampering our ability to predict the fate of oil in applications such as deep-water spills. Here, we provide the first experimental and numerical verification of the immiscible RT turbulence theory, unveiling a unique turbulent state that originates at the oil–water interface due to the interaction of multiple capillary waves. We show that a single, non-dimensional, and time-independent parameter controls the range of validity of the theory. Our findings have wide-ranging implications for the understanding of the mixing of immiscible fluids. This includes in particular oil spills, where our work enables the prediction of the oil–water interface dynamics that ultimately determine the rate of oil biodegradation by marine bacteria.

### Sign up for PNAS alerts.

Get alerts for new articles, or get an alert when an article is cited.

The Rayleigh–Taylor (RT) instability develops at the interface between two fluids of different densities in the presence of a relative constant acceleration, such as gravity, when the heavier fluid is placed on the top of the lighter one. After an initial transitional phase, the flow becomes turbulent and is sustained by the continuous conversion of potential into kinetic energy (1). RT turbulence can be observed in oceanographic and atmospheric flows (2), astrophysics (3, 4) and in industrial settings such as combustion chambers (5). RT turbulence also occurs in all types of nuclear fusion (6–8), where it represents the main limiting factor to achieve functional inertial confinement fusion reactors (9, 10).

RT turbulence often occurs alongside other physical effects, including rotation (11, 12), a mean shear (13, 14), a variable acceleration (15) or surface tension (16). The latter arises when the unstable density stratification consists of two immiscible fluids, such as oil and water or air and water. In this “immiscible RT” case, the interface tension inhibits molecular intermixing of the two fluids. Yet, the interface is still unstable to perturbations above a critical wavelength, which is set by the dispersion relation for capillary waves (17). As the instability evolves, the flow becomes turbulent, resulting in an emulsion-like state in which droplets of one fluid are dispersed into the other (16) (Fig. 1

*A*–*C*).Fig. 1.

As lighter oil is released below heavier water, immiscible RT turbulence is important in predicting the fate of oil spilled in deep waters because it determines the distribution of oil droplet sizes in the water column (18). RT turbulence can furthermore determine the size and distribution of oil slicks on the ocean surface, as the RT instability also occurs when waves overturn the oil–water interface (19, 20). In both cases, a reliable understanding of the oil–water interface dynamics is crucial to predict the biodegradation of oil (21), the main mechanism by which anthropogenically spilled oil is removed from the ocean (22–25), whose efficiency is strongly affected by turbulence. Indeed, turbulence increases the oil–water interface available for bacterial colonization (26) by breaking up the oil mass into droplets, and the droplets’ size distribution has a major impact on the encounter rate with bacteria and thus ultimately the degradation rate (27).

A phenomenological theory of the emulsification process driven by immiscible RT turbulence has been proposed (16) which predicts the temporal scaling of the droplet size in the emulsion as well as the occurrence of a peculiar type of turbulence—called capillary-driven turbulence—at the immiscible interface. Capillary-driven turbulence is a chaotic state that involves a large number of interactingnon-linearly coupled waves. The wave–wave interaction produces a complex statistical state characterized by random fluctuations in wave amplitude and frequency (28–31), clearly distinguishable from the classical hydro-dynamical turbulence of the Kolmogorov–Obukov fashion (32). We will show below that this distinctive feature affects the dynamics of the interface area available for bacterial biodegradation.

To date, the theory of immiscible RT turbulence has not been tested through either experiments or numerical simulations. Experimentally, it has remained challenging to reliably create an initially unstable stratification (1) and to measure the properties of the three-dimensional emulsion that is created. Meanwhile, direct numerical simulations (DNS) have been restricted to two dimensions (33–35), as the computational resources to achieve the very high spatial resolution required to compute the three-dimensional flow dynamics at the fluid–fluid interface are prohibitive (17).

In this work, we present the first experimental and numerical validation of the immiscible RT turbulence theory, thereby also identifying the range of conditions over which the theory is valid. We show that a single non-dimensional and time-independent parameter controls the existence and magnitude of capillary-driven turbulence at the immiscible interface. This enables predicting the evolution of immiscible RT and its interface area to volume ratio, which is the key quantity to determine the rate of oil biodegradation.

## Results

### The RT Emulsification Process.

We generated the immiscible RT turbulence in the laboratory through a non-stationary experiment. The instability was triggered by removing a thin plastic foil that initially maintained the unstable stratification (Fig. 1

*A*). The phase field was measured on a vertical slice through laser-induced fluorescence. In parallel, computer simulations were performed by integrating the equations of motion for an immiscible two-phase flow system (Fig. 1*B*). More details about the experimental and numerical approaches can be found in*Materials and Methods*. Laboratory experiments and simulations produced a RT emulsion that is qualitatively similar (Fig. 1*C*). Blending the computer simulations with laboratory experiments, six different interface tensions were covered. Within the laboratory experiments, a low-viscosity oil–water system was employed which mimics the surface tension representative of a broad range of applications in nature, such as air–water or oil–water interfaces, while the conducted simulations span smaller surface tension values representative of, for example, an oil–water emulsion with surfactants (see*SI Appendix*for a full description of the experiments and simulations).In both experiments and simulations, after triggering the instability, the turbulence develops, and the turbulent region grows in time by eroding the surrounding quiescent fluid. Fig. 2

*A*shows the temporal evolution of the mixing layer thickness*h*, a measure of the depth of the turbulent region, normalized by its predicted value ${h}_{0}$ at the transition to fully developed turbulence (i.e., after the time ${t}_{0}$ in Eq.**14**). Both numerical and experimental data are in reasonable agreement with the phenomenological prediction $h\sim {t}^{2}$ (Eq.**10**). During the initial transitional phase, both experiments and simulations deviate from the predicted scaling, with the experiments reaching it earlier. This can be explained by the residual effect of the instability phase, which strongly differs between simulations, where the initial perturbation is weak and of short wave-length, and experiments, where the perturbation is much stronger and of longer wave-length (36, 37). To substantiate this point, we conducted an additional numerical simulation with surface tension comparable to the exp 1-2-3 and sim 5, “sim 5 HP,” in which the long-wavelength perturbation generated hydrostatic deformation of the membrane observed in the experiments was reproduced in the initial conditions. Compared to sim 5 with the flat interface initial condition, in sim 5 HP, the predicted scaling is observed earlier and is more similar to the experiments (*SI Appendix*, Fig. S4). Further analysis of the effect of perturbation on the mixing layer evolution can be found in*SI Appendix*.Fig. 2.

As the mixing layer expands, progressively smaller isolated droplets of one phase are dispersed into the other. The 99th percentile of the droplet diameter, ${d}_{99}$—which is a measure of the largest diameter of the droplets in the flow—is shown in Fig. 2

*B*normalized by its lower cut-off ${d}_{H}$, which is the smallest possible droplet size (determined by surface tension and viscosity, see Eq.**13**). Both experiments and simulations are in reasonable agreement with the phenomenological prediction (Eq.**13**) up to the cutoff time ${t}_{H}$.From our data, it emerges that, for some combinations of interface tension, density contrast, and viscosity (corresponding to simulations 1 and 2), the theoretical prediction fails in describing the droplet diameter evolution for the whole duration of the process. This happens if the time is already larger than ${t}_{H}$ when the immiscible RT turbulence starts. In these cases, inertial fragmentation of the droplets by turbulence is inhibited because the largest droplets in the mixing layer are always smaller than the smallest turbulence length-scale

*η*. A range of scales between the droplet size*D*and the turbulence micro-scale*η*is available only if the combination of density contrast, surface tension and viscosity is such that a shaded area as depicted in Fig. 1*D*exists. Analyzing the involved length scales (see*Length-Scales Analysis*in*Materials and Methods*), we obtain that the existence and extent of this shaded area are determined by$$\begin{array}{c}\hfill R{e}_{0}=\frac{h({t}_{0})\dot{h}({t}_{0})}{\nu}={\left(\frac{{D}_{0}}{{\eta}_{0}}\right)}^{4/3},\end{array}$$

[1]

where $\dot{h}=dh/dt$ is the velocity of expansion of the mixing layer. We call this dimensionless quantity

*Reynolds number*, as it quantifies the scale separation between the droplet size ${D}_{0}$ and the smallest turbulence scale ${\eta}_{0}$ at the onset of the immiscible RT turbulence. The length scales ${D}_{0}$ and ${\eta}_{0}$ can be computed using the properties of the fluids pair only, namely the kinematic viscosity*ν*, the interface tension*γ*, the density contrast $\mathrm{\Delta}\rho $, and the mean density of the two fluids*ρ*(see*Length-Scales Analysis*in*Materials and Methods*for a full derivation of these scales). By doing so, Eq.**1**becomes $R{e}_{0}={(\gamma /\rho )}^{3/4}{(Ag)}^{-1/4}{\nu}^{-1}$, where $A=\mathrm{\Delta}\rho /(2\rho )$ is the Atwood number (a dimensionless representation of the density contrast) and*g*is the gravitational acceleration. The parameter $R{e}_{0}$ controls the range of applicability of the theory, as having $R{e}_{0}>1$ automatically implies ${D}_{0}>{\eta}_{0}$ and thus ${t}_{H}>{t}_{0}$ and ${t}_{0}>{t}_{\nu}$, where ${t}_{\nu}$ is the time at transition from laminar to turbulent flow in absence of surface tension (see Fig. 1*D*for a graphical interpretation and the*Length-Scales Analysis*in the*Materials and Methods*for the full analytical derivation). Table 1 reports the values of $R{e}_{0}$ across experiments and simulations. Consistently, the only cases for which the empirical data match the theoretical prediction are the ones with a $R{e}_{0}$ larger than unity.Table 1.

Label | Re_{0} | Re_{b} | |
---|---|---|---|

sim 1 | $2.08\xb7{10}^{-1}$ | $7.37\xb7{10}^{4}$ | |

sim 2 | $1.17\xb7{10}^{0}$ | $7.37\xb7{10}^{4}$ | |

sim 3 | $6.59\xb7{10}^{0}$ | $7.37\xb7{10}^{4}$ | |

sim 4 | $3.70\xb7{10}^{1}$ | $7.37\xb7{10}^{4}$ | |

sim 5 | $2.08\xb7{10}^{2}$ | $7.37\xb7{10}^{4}$ | |

exp 1-2-3 | $2.08\xb7{10}^{2}$ | $6.17\xb7{10}^{4}$ |

$R{e}_{0}$ is the Reynolds number defined in Eq.

**1**, namely at the onset of self-similar turbulence. $R{e}_{b}$ is the Reynolds number at half of the vertical dimension of the simulation/experiment domain ${h}_{b}$; given that $h(t)Ag{t}^{2}$, we write the time scale at which the mixing layer reach half of the box size, namely ${t}_{b}$ = ${\left({h}_{b}/Ag\right)}^{1/2}$, from which we can derive the Reynolds number $R{e}_{b}={A}^{2}{g}^{2}{t}_{b}^{3}/\nu $. The full set of dimensionless and dimensional parameter characterizing simulations and experiments can be found in*SI Appendix*.### Co-Existence of Capillary-Driven and Kolmogorov Turbulence.

The parameter $R{e}_{0}$ determines the existence and extent of the scaling range bounded above by

*D*and below by ${d}_{H}$, within which a peculiar type of turbulence manifests: capillary-driven turbulence. The latter is a chaotic state of wave turbulence where the energy is injected at the scale of the droplet size and passed on to progressively smaller waves through non-linear wave interaction before it is finally dissipated viscously at the scale ${d}_{H}$. As shown in Fig. 1*D*, in between*D*and*η*, capillary-driven turbulence is expected to co-exist with the classical turbulence of the Obukov–Kolmogorov type (38, 39) (referred to as K.O. in Fig. 1*D*), which applies in ordinary single phase turbulence. Below*η*, waves prevail and constitute the only source of fluctuations. Capillary-driven turbulence is characterized by the following distinctive velocity pair-correlation function (16, 28, 40):$$\begin{array}{c}\hfill \u27e8{\mathbf{u}}_{\mathit{int}}(\mathbf{x}+\mathbf{r})\xb7{\mathbf{u}}_{\mathit{int}}(\mathbf{x})\u27e9/{(\u03f5D)}^{2/3}\sim {(r/D)}^{-1/4},\end{array}$$

[2]

where velocity pairs are evaluated at the immiscible interface and brackets indicate an ensemble average over space for a given time. Eq.

**2**, in fact, constitutes the dynamical signature of interface tension onto the turbulence structure. By substituting the scaling of the droplet diameter (Eq.**13**) into Eq.**2**, one obtains the temporal scaling of the pair-correlation function at a given distance*r*, which, for $r\to 0$, yields the temporal scaling of the kinetic energy of the interface:$$\begin{array}{c}\hfill {k}_{\mathit{int}}=\frac{1}{2}\rho \u27e8{u}_{\mathit{int}}^{2}\u27e9\sim \rho {({h}_{0}/{t}_{0})}^{2}{(t/{t}_{0})}^{3/10}.\end{array}$$

[3]

By contrast, for the following three other cases ${k}_{\mathit{int}}$ scales with ${t}^{2}$: i) miscible RT turbulence; ii) in the bulk of the immiscible RT turbulence (i.e., averaging within the mixing layer without conditioning the statistics to the interface) at all times; and iii) after the collapse of the capillary wave range in immiscible RT turbulence (i.e., after ${t}_{H}$). In the latter case, we thus expect a transition of the slope of kinetic energy as a function of time in log–log space from $3/10$ at early times to 2 at late times.

Fig. 3

*A*shows the temporal scaling of ${k}_{\mathit{int}}$ for all simulations and experiments, rendered non-dimensional with the scales ${t}_{0}$ and ${h}_{0}$. As predicted from the theory, at early times ($t<6{t}_{0}$), the interface dynamics is dominated by capillary waves, while at later times ($t>10{t}_{0}$) the quadratic scaling of miscible RT turbulence is recovered. Only the experimental data are able to cover the time interval over which the interface dynamics are clearly dominated by capillary waves and are in the self-similar turbulent regime simultaneously. The reason is that the simulations take longer to reach the self-similar turbulent regime due to the low noise level applied in the initial perturbation. Nevertheless, simulations 4 and 5 start deviating from the quadratic regime at times lower than $t\sim 10{t}_{0}$, showing a clear trend toward a smaller exponent. Possible reasons for limited coverage of the capillary regime for the simulation with larger $R{e}_{0}$ (sim 5) might be due to the nature of the perturbation, which in the experimental data is stronger (due to the tangential shear applied when removing the foil) and of longer wavelength (due to the domain size in the horizontal direction). Moreover, numerical simulations are limited by finite resolution, preventing the resolution of the entire capillary wave range (the cutoff in sim 5 is larger than the capillary wave dissipation scale, ${d}_{H}$ in Fig. 1*D*, and of the order of the Kolmogorov scale*η*). Since the scaling exponents $3/10$ and 2 are very different, we show in the*Inset*of Fig. 3*A*a zoom on the experimental data only, which highlights the agreement between the $3/10$ scaling and the data.Fig. 3.

To provide further confirmation of the presence of capillary-driven turbulence in the experiments, we measure directly the velocity pair-correlation function of Eq.

**2**for point pairs belonging to the interface as a function of the separation*r*averaged over time and realizations after rendering the data dimensionless through the droplet diameter*D*and the energy dissipation rate*ϵ*(see the*Length Scale Analysis*in the*Materials and Methods*section). As shown in Fig. 3*B*, an incipient capillary wave range consistent with the phenomenological prediction of Eq.**2**clearly manifests, spanning roughly an order of magnitude in scales. This is perfectly consistent with the expectation provided by the value of $R{e}_{0}\sim {10}^{2}$ that yields ${D}_{0}\simeq 30{\eta}_{0}$. Such scaling is not observed in any of the numerical simulations, consistently with the observation of Fig. 3*A*, i.e., for those simulations where $R{e}_{0}$ would be high enough, the transition to developed RT turbulence takes too long.Previous experimental investigations of capillary-driven turbulence have mostly been limited to stable stratified environments in which Kolmogorov turbulence is absent or sub-leading in comparison to the gravitational force. In these flows, the upper limit of the capillary wave range is set by the length scale resulting from the balance of gravity and capillary waves. In this configuration, it has been shown that a cross-over between the gravity and the capillary wave turbulence spectrum occurs (30, 41), without any direct role of the Kolmogorov turbulence cascade. On the one hand, since RT turbulence originates from an unstable stratification, gravity waves do not exist by definition; on the other hand, it is well known that RT turbulence shows almost isotropic Kolmogorov spectra in the bulk (42–44). These factors make the immiscible RT turbulence a perfect candidate for revealing the coexistence of K.O. and capillary-driven turbulence in the absence of gravity waves. Making use of the existing phenomenological theory (16), we derived a single time-independent and non-dimensional parameter that controls the existence and extent of the capillary wave range in the immiscible RT turbulence. This allows us to reveal co-existing capillary-driven and K.O. turbulence. Rendering the equation of motion non-dimensional through the length scale ${D}_{0}={h}_{0}$ and the velocity scale $\dot{h}({t}_{0})$ permits to explicitly derive their dependence on $R{e}_{0}$ and explore the regime parametrically (see

*Materials and Methods*for the full derivation). Our experiments suggest that this type of turbulence can readily manifest in natural flows involving air and water in unstable configurations, such as plunging breaking waves (19), but also in anthropogenic flows, such as deep-water oil spills (18).### Implications for the Biodegradation of Oil-Spills.

When oil is spilled into the ocean, its mass is degraded by marine bacteria at a rate that depends on a number of different factors, among which the droplet size distribution, the type and concentration of bacteria, the seawater salinity, temperature, and the amount of nutrients (27). During the degradation process, marine bacteria cannot penetrate the oil phase as they need oxygen to survive. Nevertheless, when encountering oil droplets, they colonize the oil–water interface, metabolizing petroleum (22) at a rate that is proportional to the available interface area (45, 46). In this context, understanding how immiscible RT turbulence produces interface area over time is of primary importance since it directly relates to the rate of oil biodegradation (18, 20). To understand how RT turbulence produces interface area over time, we derive a predictor for the temporal evolution of the interface area-to-volume ratio ${A}_{\mathit{int}}/V$, using the droplet diameter scaling (Eq.

**12**) and assuming a dense emulsion. Fig. 4 shows that our phenomenological prediction ${A}_{\mathit{int}}/V\sim {t}^{2/5}$ (the full derivation can be found in*Length-Scales Analysis*in*Materials and Methods*; see Eq.**15**) forecasts with reasonable accuracy the interface area production process, for both the experimental and simulations data and up to the cutoff time ${t}_{H}$. However, the available interface area and the rate of biodegradation are intricately related, as they can simultaneously influence each other because of both physical (turbulence) and biological factors (21).Fig. 4.

To better understand how an emulsion produced by RT turbulence biodegrades in the ocean, the spill scenario can be simplified by assuming that RT fragments the oil mass over a much shorter time scale than the bio-degradation process. In fact, our experiments, as well as literature observations of RT turbulence in a range of scenarios (1), indicate that the former evolves over a time scale from seconds to minutes, while the latter takes days to months (22). This implies that RT turbulence controls the initial condition for the biodegradation process, which is plausible for most applications involving RT turbulence in the ocean, e.g., deep-water oil spills from large apertures (18) or wave-induced overturning of the oil–water interface in surface spills (19).

We observed that (

*SI Appendix*, Fig. S2*A*) RT turbulence produces a complex polydispersed emulsion, which contains a large variety of coexisting droplet diameters. To understand how this emulsion is degraded by marine bacteria, we employ the Microscale Oil DEgradation Model (MODEM) (27), which allows simulating the full process of bacterial degradation accounting for the size distribution of the oil droplets, the dynamic concentration of bacteria, and the encounter rate between bacteria and droplets. For different values of $t/{t}_{0}$ (corresponding to different levels of emulsification), we simulate the oil biodegradation process of the droplets produced by RT turbulence in seawater (see*Oil Biodegradation Simulations*in*Materials and Methods*for the technical details). Fig. 5 shows that the evolution of the total oil mass (biodegradation curves) depends on the level of emulsification. For every biodegradation curve, we estimate ${T}_{10\%}$, namely the time needed for the total oil mass to reduce to the $10\%$ of its original value. The*Inset*of Fig. 5 shows that, for seawater, ${T}_{10\%}$ is inversely proportional to the non-dimensional time $t/{t}_{0}$. Furthermore, scaling the time ${t}_{0}$ reasonably collapses the data on a single master curve, proving that, for the sea conditions of our simulations, ${t}_{0}$ is the only relevant scale controlling biodegradation.Fig. 5.

Finally, we derive a simplified phenomenological theory to explain the scaling of the biodegradation time ${T}_{10\%}$, with $t/{t}_{0}$. To this purpose, we first assume that all the oil droplets in the emulsion are colonized by a uniform distribution of bacteria. In this scenario, the rate of variation of the volume of oil $dV/dt$ can be written as:

$$\begin{array}{c}\hfill \frac{1}{V}\frac{\mathit{dV}}{\mathit{dt}}=-\frac{k}{\rho}\frac{{A}_{\mathit{int}}}{V},\end{array}$$

[4]

where ${A}_{\mathit{int}}$ is the available interface area and the constant

*k*[dimension $(M{L}^{-2}{T}^{-1})$] is an effective bacteria degradation rate (45). Eq.**4**accounts for the complex interplay between the rate of variation of the oil volume due to bio-degradation and the interface area-to-volume ratio. Solving Eq.**4**for an individual isolated droplet yields a linear-in-time shrinking of the droplet’s initial diameter. Considering an equivalent droplet monodispersion and using Eq.**12**, the $10\%$ degradation time can be expressed as:$$\begin{array}{c}\hfill {T}_{10\%}={T}_{10\%}^{0}{(t/{t}_{0})}^{-2/5}\text{.}\end{array}$$

[5]

The time ${T}_{10\%}^{0}$ is the biodegradation time at the onset of RT turbulence-driven fragmentation, as it is the value of ${T}_{10\%}$ at $t/{t}_{0}=1$. As shown in detail in

*SI Appendix*, the effective degradation*k*can then be calculated from ${D}_{0}$,*ρ*, and ${T}_{10\%}^{0}$ as $k=0.27\xb7{D}_{0}\xb7\rho /{T}_{10\%}^{0}$, where ${T}_{10\%}^{0}$ is estimated by fitting the theoretical scaling of Eq.**5**to the data (see black solid line in Fig. 5). In common practice, the RT instability is only invoked to set the maximum droplet diameter, postulating that the latter matches the first unstable wavelength of the RT system without considering any turbulence time or length scales (18). This would correspond to considering the biodegradation time equal to ${T}_{10\%}^{0}$, which however neglects the effect of turbulent fragmentation. Our phenomenology shows (Eq.**5**and*Inset*of Fig. 5) that ignoring that ${T}_{10\%}$ changes depending on the level of emulsification $t/{t}_{0}$ leads to an overestimation of the degradation time by up to*∼*400%.## Conclusions

Our study of immiscible RT turbulence reveals a unique regime of capillary turbulence at the interface between the two liquids. The existence of this regime is controlled by a single non-dimensional parameter, which balances the relative strength of inertial versus viscous forces at the onset of turbulence and depends on the properties of the fluid pair only. In the case of an oil spill in the ocean, the immiscible RT turbulence can fragment the oil phase into the water, decreasing the droplet size and increasing the available interface area for a given volume. This affects the rate of biodegradation, as marine bacteria are allowed to colonize a larger interface. Our work implies that one should take into account the full emulsification process occurring in deep-water spills or induced by the wave overturning of the oil–water interface in superficial immiscible RT turbulence rather than the initial instability phase only as has been done up to now. Biodegradation times estimated with the former versus the latter approaches can lead to vastly different biodegradation times and, hence, environmental impact assessment outcomes.

## Materials and Methods

### Laboratory Experiments.

The RT instability was created in the laboratory through a non-stationary experiment. A pair of immiscible fluids were initially separated by a thin plastic foil that keep an approximately flat interface (Fig. 1

*A*.1). To trigger the instability, the foil was rapidly removed. The time needed to completely remove the foil ($\sim 0.1s$) is negligible compared to the time scale of the instability, which is an essential factor to avoid any partial and asymmetrical overturning of the stratification. Conceptually, this method is similar to the one adopted in refs. 47, 48, 49, 50, with two differences: i) Having a very thin foil, the fluid displacement in the wake of the moving membrane is less prominent so that the potential flow produced by the moving membrane is almost negligible and ii) the thin membrane can be very easily removed from the top aperture of the tank since the plastic foil can bend up and slide near the wall (see*SI Appendix*for the drawings of the set-up). Nevertheless, the foil, by moving fast, generates two boundary layers on each side, introducing a shear perturbation at the interface between the two fluids conceptually similar to the mechanisms of a bursting bubble or a flapping flag (51). This mechanism perturbs the interface at every wavelength (Fig. 1*A*.2), but the surface tension smears out the short wavelengths. The turbulence then originates from large wave-length perturbation only (Fig. 1*A*.3), while smaller scale details grow at later times, as expected in RT turbulence.In order to measure the phase field, the laser-induced fluorescence (LIF) technique was employed. For this, a coherent laser sheet less than 0.5 mm thick was used to illuminate a slice that is orthogonal to the direction of the interface removal. To discriminate the fluids, the aqueous phase is colored with a fluorescent dye, emitting red light when illuminated by the green laser. The scene is recorded by a 1 megapixel high-speed camera running at 500 frames per second, equipped with a red filter in front of the objective so that only the fluorescent aqueous phase is visible. The camera was pointed right below the initial interface in order to observe at high enough resolution the lower half of the mixing zone. LIF requires refractive index (RI) matching of the two fluids to avoid undesirable blurring, while the physics requires matching kinematic viscosity to exclude asymmetries due to viscosity contrasts. These two conditions are achieved simultaneously by using a low viscosity and low refractive index silicon oil

^{*}($\nu =2cSt$, $RI=1.3925$) in combination with a mixture of water, salt, and glycerol. We opted for a ternary mixture since it features two degrees of freedom in the $RI-\nu $ space, through which we are able to match simultaneously both the $\mathit{RI}$ and*ν*of the silicon oil. The $RI-\nu $ matching results in a surface tension of $\gamma =32.44$ mN/m (measured with the pendant droplet method) and a density contrast of $\mathrm{\Delta}\rho =234$ kg/m${}^{3}$. The densities of the ternary mixture and the oil are 1,113 and 879 kg/m${}^{3}$ respectively, resulting in a moderate Atwood number $A=({\rho}_{\mathit{mixture}}-{\rho}_{\mathit{oil}})/({\rho}_{\mathit{mixture}}+{\rho}_{\mathit{oil}})=0.12$. The viscosity of the two fluids is further measured with a viscometer, resulting in a negligible viscosity contrast, and a mean kinematic viscosity of*∼*2cSt. We successfully matched the refractive index of different fluid pairs with low viscosity (mainly hydrocarbons such as hexane versus the water-salt-glycerol mixture), obtaining low viscosity contrasts (order unity) that can be used to explore different surface tension and density contrasts. Nevertheless, in this work, these substances have not been employed, due to their danger of ignition and toxicity.Since LIF allows to access only one slice of the volume, repetitions are needed to acquire converged statistics at each time-step. Note that some of the statistical properties of the emulsion are intrinsically three dimensional, and thus difficult to access through a two dimensional slice. However, the simulations provide both the three-dimensional phase field (Fig. 1

*B*) and two-dimensional slices, allowing to quantify the reliability of the statistics calculated from the LIF data. Fig. 1 shows a comparison from the sliced phase field of the simulation (*C*.1) and the experiment (*C*.2), proving a good qualitative agreement. In order to estimate the $99th$ percentile droplet diameter from the LIF data, we considered all the oil-in-water and water-in-oil regions that were closed and disconnected from the others, and computed the equivalent two-dimensional droplet diameter. More details on the effect of slicing can be found in section 1.C of*SI Appendix*. The interface area to volume ratio is evaluated in the two-dimensional slices as the length-to-area ratio of the interface. Finally, we employed the optical flow method to measure the interface velocity. This technique is based on the brightness constancy constrain, that relates the spatial and temporal derivative of the intensity of an image sequence to the two components of the apparent velocity. The brightness constrain alone is an under determined equation, in which the two unknowns (the two velocity components) are related by one equation only (the brightness constrain). Usually, a smoothness constrain is adopted to close the problem, which assumes that near-by patterns move coherently. Nevertheless, the smoothness constrain leads to a reduction of the final optical flow resolution. In our analysis, we opted to project the velocity field along the interface normal direction, and obtain the following equation for the interface-normal velocity component:$$\begin{array}{c}\hfill {u}_{n}=\mathbf{U}\xb7\mathbf{n}=-\frac{1}{|\mathrm{\nabla}I|}\frac{\partial I}{\partial t},\end{array}$$

[6]

where $I(\mathbf{x},t)$ is the image intensity, and $\mathrm{\nabla}$ is the spatial gradient. By this approach we are able to measure the interface-normal velocity at close-to pixel resolution, at the expense of not retrieving the tangential velocity.

### Numerical Simulations.

Simulations are performed by numerically integrating the Navier–Stokes equations in the Boussinesq approximation with the VOF method. The flow field satisfies the following continuity and momentum equations:

$$\begin{array}{c}\hfill \mathrm{\nabla}\xb7\mathbf{u}=0,\end{array}$$

[7]

$$\begin{array}{c}\hfill {\partial}_{t}\mathbf{u}+\mathbf{u}\xb7\mathrm{\nabla}\mathbf{u}=-\frac{1}{{\rho}_{0}}\mathrm{\nabla}p+\nu {\mathrm{\nabla}}^{2}\mathbf{u}+{\mathbf{f}}_{\gamma}+\mathbf{b},\end{array}$$

[8]

where ${\rho}_{0}$ is a reference density equal to the average density of the two fluids ${\rho}_{\mathit{Av}}=({\rho}_{1}+{\rho}_{2})/2$ and ${\mathbf{f}}_{\gamma}=\xi \gamma /{\rho}_{0}{\delta}_{s}\mathbf{n}$ is the surface tension force. The latter is proportional to the specific interface tension $\gamma /{\rho}_{0}$ and the local interface curvature

*ξ*.**n**is the unit vector normal to the interface, and ${\delta}_{s}$ is a Dirac function that turns on at the interface only. The term**b**is the buoyancy force. It is equal to $\beta gc{\mathbf{e}}_{3}$, where*β*is a constant equal to ${{\rho}_{\mathit{Av}}}^{-1}\partial \rho /\partial c$,*g*is the gravitational acceleration, ${\mathbf{e}}_{3}$ is the vertical unit vector, and*c*is the scalar concentration, which obeys the following advection–diffusion equation (42):$$\begin{array}{c}\hfill {\partial}_{t}c+\mathbf{u}\xb7\mathrm{\nabla}c=D{\mathrm{\nabla}}^{2}c,\end{array}$$

[9]

where

*D*is the molecular diffusivity of the scalar field. In order to locate the interface, the Multi-dimensional Tangent Hyperbola Interface Capturing (MTHINC) algorithm developed by Ii et al. (52) is adopted. In this regard, an indicator function*H*is needed to distinguish between the two phases, and obeys a pure advection equation.Initial and boundary conditions to produce RT turbulence mimic the simulations from Boffetta et al. (42) of miscible RT turbulence. Periodic boundary conditions are imposed on the lateral boundaries, while a free-slip condition is imposed on the lower and upper walls. The vertical scalar concentration profile is initiated as a step function of the form $c=\mathrm{\Delta}c/2(\text{sign}(z)+1))$, where $\mathrm{\Delta}c={\beta}^{-1}\mathrm{\Delta}\rho /{\rho}_{0}$. The interface is perturbed by adding $10\%$ noise on the scalar field in the vicinity of the interface.

The numerical integration in space is performed on a staggered fixed Eulerian grid where the fluid velocity components are located on the cell faces, while density, pressure, stress, and volume of fluid at the cell center. The time integration is performed with a fractional-step method based on the second-order Adams–Bashforth scheme. Further details can be found in refs. 53 and 54.

### Length-Scales Analysis.

According to the theoretical framework proposed by Chertkov and coauthors (16), the final stage of the RT instability of a heavier fluid (density ${\rho}_{1}$) over a lighter fluid (density ${\rho}_{2}$) consists of a self-similar turbulent mixing layer whose thickness $h(t)$ grows quadratically with time

*t*due to the continuous conversion of potential into kinetic energy (1, 55–57) (pink continuous line in Fig. 1*D*), as$$\begin{array}{c}\hfill h(t)=\alpha Ag{t}^{2},\end{array}$$

[10]

where $A=\left({\rho}_{1}-{\rho}_{2}\right)/\left({\rho}_{1}+{\rho}_{2}\right)$ is the Atwood number,

*g*the gravitational acceleration and*α*is the expansion rate coefficient, which is believed to be a universal constant and ranges between $0.03$ and $0.1$. In the self-similar turbulent stage, the integral length-scale of turbulence*L*is proportional to the mixing layer thickness*h*, i.e. $L\sim h$. Since the large-scale velocity fluctuations are of the same order of magnitude as the rate of expansion of the mixing layer, they grow linearly in time as ${u}_{L}=\dot{h}=dh/dt\sim Agt$. To derive the scaling of the energy dissipation rate, we assume that the turbulence is in a state where the energy dissipation balances production locally in time. This implies that the average turbulent dissipation rate scales as $\u03f5=\nu \left(\mathrm{\nabla}\mathbf{u}:\mathrm{\nabla}\mathbf{u}\right)\sim {u}_{L}^{3}/h\sim {A}^{2}{g}^{2}t$, where*ν*is the kinematic viscosity of the fluids and $\mathrm{\nabla}$ is the gradient operator. This hypothesis, originally named “small scale adiabaticity” (58) has been validated in several DNSs studies by analyzing the temporal scaling of the turbulent energy dissipation rate (see, e.g., ref. 42) as well as in our simulations (*SI Appendix*). Provided that the energy dissipation rate scales linearly in time, the Kolmogorov length, represented by the orange continuous line in Fig. 1*D*, is given by$$\begin{array}{c}\hfill \eta \sim {\nu}^{3/4}{\left(Ag\right)}^{-1/2}{t}^{-1/4}\text{.}\end{array}$$

[11]

In the presence of a developed inertial range, the small-scale adiabaticity assumption can be extended to the energy flux at every scale, such that the velocity fluctuations within the inertial range adapt instantaneously to the large-scale evolution of the mixing layer. As a consequence, the velocity difference between two points in the mixing layer separated by a distance

*r*at a given time*t*scales as ${\delta u}^{2}\sim {(\u03f5r)}^{2/3}\sim {(Ag)}^{4/3}{t}^{2/3}{r}^{2/3}$. Although surface tension prevents molecular mixing, the RT turbulence fragments the interface between the two fluids, generating an emulsion in which droplets of one fluid are dispersed into the other and vice-versa. Extending the adiabaticity assumption to the fragmentation process, the typical droplet diameter*D*can be estimated by balancing the interface energy density and the kinetic energy at scale*D*, namely $\gamma /D\sim \rho {\delta u}^{2}(D)$. As long as*D*lies within the inertial range of turbulence, we get the temporal evolution of the droplet size (azure dot-dashed line in Fig. 1*D*):$$\begin{array}{c}\hfill D\sim {(\gamma /\rho )}^{3/5}{(Ag)}^{-4/5}{t}^{-2/5}\text{.}\end{array}$$

[12]

Even if Eq.

**12**provides the temporal scaling of the typical droplet size, it lacks proper non-dimensionalization. To find an appropriate non-dimensionalization, we rely on the classical Kolmogorov/Hinze theory, stating that the maximum droplet size is given by $\rho \gamma D/{\mu}^{2}=C{\left({\mu}^{5}\u03f5/\rho {\gamma}^{4}\right)}^{-2/5}$, where*μ*is the dynamic viscosity,*ρ*is the fluid mean density, and*C*is a constant of order unity (59–61). By substituting the scaling for*ϵ*, we obtain the temporal scaling of*D*as a function of non-dimensional variables only:$$\begin{array}{c}\hfill D/{d}_{H}={\left(t/{t}_{H}\right)}^{-2/5},\end{array}$$

[13]

where:

$$\begin{array}{c}\hfill {d}_{H}={\left(\gamma /\rho \right)}^{-1}{\nu}^{2},{t}_{H}={\left(\gamma /\rho \right)}^{4}{\left(Ag\right)}^{-2}{\nu}^{-5}\text{.}\end{array}$$

[14]

Eq.

**13**is the non-dimensional equivalent to [**12**] but is instead made non-dimensional in the Hinze/Kolmogorov fashion to represent the phenomenon through a single universal law, i.e., droplet dynamics of immiscible RT fully developed turbulence are expected to follow this law irrespective of the value of fluid viscosity, surface tension, and density. As shown in Fig. 1*D*,*D*decreases faster than*η*, so that at some time, the droplets size reaches the turbulence micro-scale, breaking the assumption that the droplet size belongs to the inertial range of turbulence ($\eta <D<L$). This happens at the time ${t}_{H}$, at which we have that $D(t={t}_{H})=\eta (t={t}_{H})={d}_{H}$ (Fig. 1*D*).It is worth noting here that droplet inertial fragmentation is a dynamic process that, inherently, takes a finite time to occur. This fact could, in principle, limit the reliability of the adiabaticity assumption when extended to inertial-driven fragmentation. In their recent study, Vela-Martín and Avila (62) found that the breakup rate of droplets in stationary turbulence is proportional to the droplet diameter. Due to this reason, the adiabaticity hypothesis is expected to be more reliable for large droplets that can break up faster (see

*SI Appendix*for a quantitative analysis of this point).Finally, we can formulate the scaling of the area-to-volume ratio assuming a dense emulsion, namely considering the inter-droplet distance to be of the same order of magnitude as the droplet diameter (Eq.

**13**). This gives the number of droplets within the emulsion, that is $N\sim {L}_{x}{L}_{y}h/{D}^{3}$, where ${L}_{x}$ and ${L}_{y}$ are the horizontal dimension of the domain. Thus, the interface area—given by the number of droplets times the individual droplet surface area (proportional to ${D}^{2}$)—scales as ${A}_{\mathit{int}}\sim N{D}^{2}\sim {L}_{x}{L}_{y}{\left(Ag\right)}^{9/5}{\left(\gamma /\rho \right)}^{-3/5}{t}^{12/5}$. Given that the intruded volume is $V\sim {L}_{x}{L}_{y}h$, the area-to-volume ratio is:$$\begin{array}{c}\hfill {A}_{\mathit{int}}/V\sim {\left(t/{t}_{H}\right)}^{2/5}/{d}_{H}.\end{array}$$

[15]

We underline that the surface area-to-volume ratio is a distinctive signature of the interface tension, as neglecting the effect of the latter, the area to volume ratio of a scalar interface would scale in time with an exponent $11/4$, which is much faster than reality. Indeed, as turbulence is self-similar between the scales

*L*and*η*($L\sim h$ in RT turbulence), scalar interfaces are described as fractal surfaces, with co-dimension equal to $7/3$ (see, e.g., refs. 63 and 64). This means that, for a scalar interface without surface tension, the area-to-volume ratio is given by $h{(h/\eta )}^{-1/3}$, from which $A/V\sim {t}^{11/4}$.To complete the phenomenological theory, we determined under which conditions the capillary wave range is expected to manifest. In miscible RT turbulence, the timescale for the onset of turbulent self-similarity ${t}_{0}$ is the viscous timescale ${t}_{\nu}={\nu}^{1/3}{(Ag)}^{-2/3}$ obtained by balancing buoyancy and viscosity. For immiscible RT turbulence, interface tension delays the onset of turbulent self-similarity for a given type and intensity of the perturbation (37). We can thus derive ${t}_{0}$ by balancing surface tension and buoyancy, as:

$$\begin{array}{c}\hfill {t}_{0}={\left(\gamma /\rho \right)}^{1/4}{\left(Ag\right)}^{-3/4}.\end{array}$$

[16]

From ${t}_{0}$ we obtain two length scales, namely the initial droplet diameter ${D}_{0}$ (which coincides with the initial mixing layer thickness ${h}_{0}$) and the initial Kolmogorov length ${\eta}_{0}$ respectively (Fig. 1

*D*):$$\begin{array}{c}\hfill {D}_{0}={h}_{0}={\left(\gamma /\rho \right)}^{1/2}{\left(Ag\right)}^{-1/2},\\ \hfill {\eta}_{0}={\left(\gamma /\rho \right)}^{-1/16}{\left(Ag\right)}^{-5/16}{\nu}^{3/4}.\end{array}$$

[17]

These length scales allow us to define two conditions for the existence of the capillary wave range, namely ${t}_{H}/{t}_{0}>>1$ and ${t}_{0}/{t}_{\nu}>>1$. These two conditions can be merged in a single criterion, namely that:

$$\begin{array}{c}\hfill \frac{{D}_{0}}{{\eta}_{0}}={\left(\gamma /\rho \right)}^{9/16}{(Ag)}^{-3/16}{\nu}^{-3/4}>>1\text{.}\end{array}$$

[18]

Indeed, ${D}_{0}/{\eta}_{0}={({t}_{0}/{t}_{\nu})}^{9/4}={({t}_{H}/{t}_{0})}^{3/20}$. By considering the mixing layer thickness scaling of Eq.

**10**, we can calculate the mixing layer expansion velocity $\dot{h}=dh/dt\sim Agt$ at the time ${t}_{0}$, from which it follows that $R{e}_{0}={({D}_{0}/{\eta}_{0})}^{4/3}={(\gamma /\rho )}^{3/4}{(Ag)}^{-1/4}{\nu}^{-1}$. This is our main theoretical finding: the capillary wave range exists only if $R{e}_{0}>1$ and its extent is determined by the value of $R{e}_{0}$. It should be noted that the initial time is not exactly ${t}_{0}$, but it is ${C}_{0}\xb7{t}_{0}$, where ${C}_{0}$ is a constant of order unity that depends on the initial perturbation intensity. The implication of this threshold Reynolds number is that inertia-driven fragmentation of droplets will occur ($R{e}_{0}>>1$) or not ($R{e}_{0}<<1$) depending on the fluids’ density, viscosity and surface tension.Although the coefficient

*α*is believed to be a universal constant regardless of perturbation, sub-leading effects (horizontal shear, surface tension, rotation etc.), differences in the initial condition (strength and wavelength of the perturbation) or flow asymmetry (considering spikes or bubbles at finite Atwood number), may influence its value. In fact, empirical data show that*α*depends on each of these factors. As the focus of our work is not on the*α*coefficient, but rather on the relation between small-scale phenomena (e.g., capillary-driven turbulence and droplet size temporal scaling) and large-scale phenomena (mixing layer thickness, large scale velocity fluctuations), we absorb the*α*coefficient into the data scalings. This means that, for instance, by evaluating the scale ${h}_{0}$ for rendering the data non dimensional, we retain*α*, and get ${h}_{0}=\alpha {(\gamma /\rho )}^{1/2}{(Ag)}^{-1/2}$. To this purpose, we estimated two different values of the*α*coefficient, namely $0.19$ in the experiments and $0.04$ in the simulations (see*SI Appendix*for more details).### Non-Dimensional Formulation of the Governing Equations.

When considering the problem of immiscible RT turbulence, it is desirable to work directly with a non-dimensional version of the Navier–Stokes equation that considers their dependence on the parameter $R{e}_{0}$ explicitly. For this purpose, we define the following set of non-dimensional variables:

$$\begin{array}{c}\hfill {\mathbf{x}}^{\ast}=\frac{\mathbf{x}}{{h}_{0}},{\mathbf{u}}^{\ast}=\frac{\mathbf{u}}{\dot{h}({t}_{0})},{t}^{\ast}=\frac{\dot{h}({t}_{0})t}{{h}_{0}},\end{array}$$

$$\begin{array}{c}\hfill {p}^{\ast}=\frac{p}{{\rho}_{0}\dot{h}{({t}_{0})}^{2}},{c}^{\ast}=\frac{2c}{\mathrm{\Delta}c}.\end{array}$$

[19]

By substituting in the continuity and momentum equations (Eqs.

**7**and**8**respectively), power counting leads to the following non-dimensional equations:$$\begin{array}{c}\hfill \mathrm{\nabla}\xb7{\mathbf{u}}^{\ast}=0,\end{array}$$

[20]

$$\begin{array}{c}\hfill {\partial}_{{t}^{\ast}}{\mathbf{u}}^{\ast}+{\mathbf{u}}^{\ast}\xb7\mathrm{\nabla}{\mathbf{u}}^{\ast}=-\mathrm{\nabla}{p}^{\ast}+\frac{1}{R{e}_{0}}{\mathrm{\nabla}}^{2}{\mathbf{u}}^{\ast}+{\mathbf{f}}_{\gamma}^{\ast}+{\mathbf{b}}^{\ast}\text{.}\end{array}$$

[21]

The non-dimensional interface tension force is ${\mathbf{f}}_{\gamma}^{\ast}={\xi}^{\ast}{\delta}_{s}({\mathbf{x}}^{\ast})\mathbf{n}$, where ${\xi}^{\ast}{\delta}_{s}({\mathbf{x}}^{\ast})=\xi {\delta}_{s}(\mathbf{x})/{{h}_{0}}^{2}$. Recalling that $\beta \mathrm{\Delta}c=\mathrm{\Delta}\rho /{\rho}_{0}$, and that $A=\mathrm{\Delta}\rho /(2{\rho}_{0})$, the non-dimensional buoyancy force becomes ${\mathbf{b}}^{\ast}=-{\mathbf{e}}_{3}{c}^{\ast}$. Using the same scaling, the advection-diffusion equation for the scalar field (Eq.

**9**) becomes$$\begin{array}{c}\hfill {\partial}_{t}^{\ast}c+{\mathbf{u}}^{\ast}\xb7\mathrm{\nabla}{c}^{\ast}=\frac{1}{R{e}_{0}Sc}{\mathrm{\nabla}}^{2}{c}^{\ast},\end{array}$$

[22]

where $Sc=\nu /D$ is the Schmidt number and represents the ratio between momentum and molecular diffusivity. It must be noted that, by rendering the equations non-dimensional with these scales, only two control parameters emerge, that are $R{e}_{0}$ and $\mathit{Sc}$. For an immiscible pair of fluid such as an oil–water system, the solubility of the two fluid is extremely low, such as $\mathit{Sc}$ can be considered to be infinite, as no mass can diffuse across the interface. $R{e}_{0}$ remains, in fact, the only parameter controlling the overall behavior of the system, i.e., the existence, extent, and duration of capillary-driven turbulence at the interface.

### Oil Biodegradation Simulations.

The oil biodegradation process is simulated for several values of $t/{t}_{0}$ using the MODEM (27). This model represents the degradation of oil droplets by bacteria attached to the oil–water interface. Based on a Monte-Carlo approach, MOMEM allows simulations of the temporal evolution of the whole droplet diameter distribution by taking into account the oil volume fraction, the background bacteria concentration, the temporal evolution of the bacteria population, and the dependence between the rate of encounter between bacteria and oil droplets with the individual size of each droplet in the simulation. In this work, we ran MODEM on the droplet size distribution generated by sim 1 to 4 for several emulsification levels, corresponding to different values of $t/{t}_{0}$. Each run simulates the evolution of ${10}^{5}$ droplets generated using the empirical probability density function (PDF) estimated from the DNS data at the five distinct times considered throughout the paper. Only simulations 1 to 4 have been considered as they contain a large enough number of droplets to estimate the empirical PDF with reasonable accuracy. Further details about the MODEM simulations can be found in

*SI Appendix*.## Data, Materials, and Software Availability

All study data are included in the article and/or

*SI Appendix*.## Acknowledgments

S.B. and M.H. were supported by the Priority Programme SPP 1881 Turbulent Superstructures of the Deutsche Forschungsgemeinschaft under the Grant No. HO5519/1-2. M.E.R. was supported by the Okinawa Institute of Science and Technology Graduate University (OIST) with subsidy funding from the Cabinet Office, Government of Japan. M.E.R. acknowledges the computer time provided by the Scientific Computing section of the Research Support Division at OIST.

### Author contributions

S.B., R.N., M.E.R., R.S., G.B., A.M., and M.H. designed research; S.B., R.N., M.E.R., G.B., and A.M. performed research; S.B. and M.E.R. analyzed data; and S.B., R.N., M.E.R., R.S., G.B., A.M., and M.H. wrote the paper.

### Competing interests

The authors declare no competing interest.

## Supporting Information

Appendix 01 (PDF)

- Download
- 925.11 KB

## References

1

G. Boffetta, A. Mazzino, Incompressible Rayleigh–Taylor turbulence.

*Annu. Rev. Fluid Mech.***49**, 119–143 (2017).2

J. S. Turner, J. S. Turner,

*Buoyancy Effects in Fluids*(Cambridge University Press, 1979).3

A. Traxler, P. Garaud, S. Stellmach, Numerically determined transport laws for fingering (“thermohaline’’) convection in astrophysics.

*Astrophys. J. Lett.***728**, L29 (2011).4

W. H. Cabot, A. W. Cook, Reynolds number effects on Rayleigh–Taylor instability with possible implications for type Ia supernovae.

*Nat. Phys.***2**, 562–568 (2006).5

K. Nagata, S. Komori, The effects of unstable stratification and mean shear on the chemical reaction in grid turbulence.

*J. Fluid Mech.***408**, 39–52 (2000).6

R. D. Petrasso, Rayleigh’s challenge endures.

*Nature***367**, 217–218 (1994).7

R. P. Taleyarkhan et al., Evidence for nuclear emissions during acoustic cavitation.

*Science***295**, 1868–1873 (2002).8

A. Burrows, Supernova explosions in the universe.

*Nature***403**, 727–733 (2000).9

J. Kilkenny et al., A review of the ablative stabilization of the Rayleigh–Taylor instability in regimes relevant to inertial confinement fusion.

*Phys. Plasmas***1**, 1379–1389 (1994).10

A. M. Angulo et al., Design of a high-resolution Rayleigh–Taylor experiment with the crystal backlighter imager on the national ignition facility.

*J. Inst.***17**, P02025 (2022).11

G. Boffetta, A. Mazzino, S. Musacchio, Rotating Rayleigh–Taylor turbulence.

*Phys. Rev. Fluids***1**, 054405 (2016).12

K. A. Baldwin, M. M. Scase, R. J. Hill, The inhibition of the Rayleigh–Taylor instability by rotation.

*Sci. Rep.***5**, 1–12 (2015).13

S. Brizzolara, J. P. Mollicone, M. van Reeuwijk, A. Mazzino, M. Holzner, Transition from shear-dominated to Rayleigh–Taylor turbulence.

*J. Fluid Mech.***924**, A10 (2021).14

B. Akula, P. Suchandra, M. Mikhaeil, D. Ranjan, Dynamics of unstably stratified free shear flows: An experimental investigation of coupled Kelvin–Helmholtz and Rayleigh–Taylor instability.

*J. Fluid Mech.***816**, 619–660 (2017).15

S. I. Abarzhi, K. R. Sreenivasan, Self-similar Rayleigh–Taylor mixing with accelerations varying in time and space.

*Proc. Natl. Acad. Sci. U.S.A.***119**, e2118589119 (2022).16

M. Chertkov, I. Kolokolov, V. Lebedev, Effects of surface tension on immiscible Rayleigh–Taylor turbulence.

*Phys. Rev. E***71**, 055301 (2005).17

A. Celani, A. Mazzino, P. Muratore-Ginanneschi, L. Vozella, Phase-field model for the Rayleigh–Taylor instability of immiscible fluids.

*J. Fluid Mech.***622**, 115–134 (2009).18

Z. Li, M. Spaulding, D. F. McCay, D. Crowley, J. R. Payne, Development of a unified oil droplet size distribution model with application to surface breaking waves and subsea blowout releases considering dispersant effects.

*Mar. Pollut. Bull.***114**, 247–257 (2017).19

Q. Gao, G. B. Deane, L. Shen, Bubble production by air filament and cavity breakup in plunging breaking wave crests.

*J. Fluid Mech.***929**, A44 (2021).20

J. M. Shaw, A microscopic view of oil slick break-up and emulsion formation in breaking waves.

*Spill Sci. Technol. Bull.***8**, 491–501 (2003).21

D. French-McCay, D. Crowley, L. McStay, Sensitivity of modeled oil fate and exposure from a subsea blowout to oil droplet sizes, depth, dispersant use, and degradation rates.

*Mar. Pollut. Bull.***146**, 779–793 (2019).22

I. M. Head, D. M. Jones, W. F. Röling, Marine microorganisms make a meal of oil.

*Nat. Rev. Microbiol.***4**, 173–182 (2006).23

R. M. Atlas, Petroleum biodegradation and oil spill bioremediation.

*Mar. Pollut. Bull.***31**, 178–182 (1995).24

S. Harayama, H. Kishira, Y. Kasai, K. Shutsubo, Petroleum biodegradation in marine environments.

*J. Mol. Microbiol. Biotechnol.***1**, 63–70 (1999).25

R. Margesin, F. Schinner, Biological decontamination of oil spills in cold environments.

*J. Chem. Technol. Biotechnol.***74**, 381–389 (1999).26

R. R. Lessard, G. DeMarco, The significance of oil spill dispersants.

*Spill Sci. Technol. Bull.***6**, 59–68 (2000).27

V. I. Fernandez, R. Stocker, G. Juarez, A tradeoff between physical encounters and consumption determines an optimal droplet size for microbial degradation of dispersed oil.

*Sci. Rep.***12**, 1–10 (2022).28

V. E. Zakharov, V. S. L’vov, G. Falkovich,

*Kolmogorov Spectra of Turbulence I: Wave Turbulence*(Springer Science & Business Media, 2012).29

A. Pushkarev, V. Zakharov, Turbulence of capillary waves.

*Phys. Rev. Lett.***76**, 3320 (1996).30

E. Falcon, N. Mordant, Experiments in surface gravity-capillary wave turbulence.

*Annu. Rev. Fluid Mech.***54**, 1–25 (2022).31

G. Giamagas, F. Zonta, A. Roccon, A. Soldati, Propagation of capillary waves in two-layer oil–water turbulent flow.

*J. Fluid Mech.***960**, A5 (2023).32

U. Frisch, A. N. Kolmogorov,

*Turbulence: The Legacy of AN Kolmogorov*(Cambridge University Press, 1995).33

H. Liang, X. Hu, X. Huang, J. Xu, Direct numerical simulations of multi-mode immiscible Rayleigh–Taylor instability with high Reynolds numbers.

*Phys. Fluids***31**, 112104 (2019).34

H. S. Tavares, L. Biferale, M. Sbragaglia, A. A. Mailybaev, Immiscible Rayleigh–Taylor turbulence using mesoscopic lattice Boltzmann algorithms.

*Phys. Rev. Fluids***6**, 054606 (2021).35

H. S. Tavares, L. Biferale, M. Sbragaglia, A. A. Mailybaev, Validation and application of the lattice Boltzmann algorithm for a turbulent immiscible Rayleigh–Taylor system.

*Philos. Trans. R. Soc. A***379**, 20200396 (2021).36

G. Dimonte et al., A comparative study of the turbulent Rayleigh–Taylor instability using high-resolution three-dimensional numerical simulations: The alpha-group collaboration.

*Phys. Fluids***16**, 1668–1693 (2004).37

Y. N. Young, F. Ham, Surface tension in incompressible Rayleigh–Taylor mixing flow.

*J. Turbulence***7**, 1–23 (2006).38

A. N. Kolmogorov, Dissipation of energy in the locally isotropic turbulence in Dokl.

*Akad. Nauk. SSSR***32**, 19–21 (1941).39

A. Obukhov, On the distribution of energy in the spectrum of turbulent flow.

*Bull. Acad. Sci. USSR, Geogr. Geophys.***5**, 453–466 (1941).40

A. Pushkarev, V. Zakharov, Turbulence of capillary waves-theory and numerical simulation.

*Physica D***135**, 98–116 (2000).41

E. Falcon, C. Laroche, S. Fauve, Observation of gravity-capillary wave turbulence.

*Phys. Rev. Lett.***98**, 094503 (2007).42

G. Boffetta, A. Mazzino, S. Musacchio, L. Vozella, Statistics of mixing in three-dimensional Rayleigh–Taylor turbulence at low Atwood number and Prandtl number one.

*Phys. Fluids***22**, 035109 (2010).43

G. Boffetta, A. Mazzino, S. Musacchio, L. Vozella, Kolmogorov scaling and intermittency in Rayleigh–Taylor turbulence.

*Phys. Rev. E***79**, 065301 (2009).44

P. Ramaprabhu, M. Andrews, Experimental investigation of Rayleigh–Taylor mixing at small Atwood numbers.

*J. Fluid Mech.***502**, 233–271 (2004).45

S. A. Socolofsky et al., The treatment of biodegradation in models of sub-surface oil spills: A review and sensitivity study.

*Mar. Pollut. Bull.***143**, 204–219 (2019).46

J. Vilcáez, L. Li, S. S. Hubbard, A new model for the biodegradation kinetics of oil droplets: Application to the deepwater horizon oil spill in the Gulf of Mexico.

*Geochem. Trans.***14**, 1–14 (2013).47

S. Dalziel, P. Linden, D. Youngs, Self-similarity and internal structure of turbulence induced by Rayleigh–Taylor instability.

*J. Fluid Mech.***399**, 1–48 (1999).48

S. B. Dalziel, Rayleigh–Taylor instability: Experiments with image analysis.

*Dyn. Atmos. Oceans***20**, 127–153 (1993).49

M. S. D. Wykes, S. B. Dalziel, Efficient mixing in stratified flows: Experimental study of a Rayleigh–Taylor unstable interface within an otherwise stable stratification.

*J. Fluid Mech.***756**, 1027–1057 (2014).50

A. G. Lawrie, S. B. Dalziel, Rayleigh–Taylor mixing in an otherwise stable stratification.

*J. Fluid Mech.***688**, 507–527 (2011).51

H. Lhuissier, E. Villermaux, Soap films burst like flapping flags.

*Phys. Rev. Lett.***103**, 054501 (2009).52

S. Ii et al., An interface capturing method with a continuous function: The THINC method with multi-dimensional reconstruction.

*J. Comput. Phys.***231**, 2328–2358 (2012).53

M. E. Rosti, F. De Vita, L. Brandt, Numerical simulations of emulsions in shear flows.

*Acta Mech.***230**, 667–682 (2019).54

M. Crialesi-Esposito, M. E. Rosti, S. Chibbaro, L. Brandt, Modulation of homogeneous and isotropic turbulence in emulsions.

*J. Fluid Mech.***940**, A19 (2022).55

E. Fermi, J. von Neumann, Taylor instability of incompressible liquids. Part 1. Taylor instability of an incompressible liquid. Part 2. Taylor instability at the boundary of two incompressible liquids (Tech. Rep., Los Alamos National Lab. (LANL), Los Alamos, NM, 1953).

56

J. Ristorcelli, T. Clark, Rayleigh–Taylor turbulence: Self-similar analysis and direct numerical simulations.

*J. Fluid Mech.***507**, 213–253 (2004).57

A. W. Cook, W. Cabot, P. L. Miller, The mixing transition in Rayleigh–Taylor instability.

*J. Fluid Mech.***511**, 333–362 (2004).58

M. Chertkov, Phenomenology of Rayleigh–Taylor turbulence.

*Phys. Rev. Lett.***91**, 115001 (2003).59

J. O. Hinze, Fundamentals of the hydrodynamic mechanism of splitting in dispersion processes.

*AIChE J.***1**, 289–295 (1955).60

A. Kolmogoroff, The breakage of droplets in turbulent streams.

*Doki. Akad. Nauk.***66**, 825–828 (1949).61

M. E. Rosti, Z. Ge, S. S. Jain, M. S. Dodd, L. Brandt, Droplets in homogeneous shear turbulence.

*J. Fluid Mech.***876**, 962–984 (2019).62

A. Vela-Martín, M. Avila, Memoryless drop breakup in turbulence.

*Sci. Adv.***8**, eabp9561 (2022).63

K. Sreenivasan, R. Ramshankar, C. Meneveau, Mixing, entrainment and fractal dimensions of surfaces in turbulent flows.

*Proc. R. Soc. Lond. A. Math. Phys. Sci.***421**, 79–108 (1989).64

C. Meneveau, K. Sreenivasan, Interface dimension in intermittent turbulence.

*Phys. Rev. A***41**, 2246 (1990).## Information & Authors

### Information

#### Published in

#### Classifications

#### Copyright

Copyright © 2024 the Author(s). Published by PNAS. This open access article is distributed under Creative Commons Attribution-NonCommercial-NoDerivatives License 4.0 (CC BY-NC-ND).

#### Data, Materials, and Software Availability

All study data are included in the article and/or

*SI Appendix*.#### Submission history

**Received**: July 14, 2023

**Accepted**: January 23, 2024

**Published online**: March 5, 2024

**Published in issue**: March 12, 2024

#### Keywords

#### Acknowledgments

S.B. and M.H. were supported by the Priority Programme SPP 1881 Turbulent Superstructures of the Deutsche Forschungsgemeinschaft under the Grant No. HO5519/1-2. M.E.R. was supported by the Okinawa Institute of Science and Technology Graduate University (OIST) with subsidy funding from the Cabinet Office, Government of Japan. M.E.R. acknowledges the computer time provided by the Scientific Computing section of the Research Support Division at OIST.

##### Author Contributions

S.B., R.N., M.E.R., R.S., G.B., A.M., and M.H. designed research; S.B., R.N., M.E.R., G.B., and A.M. performed research; S.B. and M.E.R. analyzed data; and S.B., R.N., M.E.R., R.S., G.B., A.M., and M.H. wrote the paper.

##### Competing Interests

The authors declare no competing interest.

#### Notes

This article is a PNAS Direct Submission.

### Authors

## Metrics & Citations

### Metrics

#### Citation statements

#### Altmetrics

### Citations

If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

## View Options

### View options

#### PDF format

Download this article as a PDF file

DOWNLOAD PDF### Get Access

#### Login options

Check if you have access through your login credentials or your institution to get full access on this article.

Personal login Institutional Login#### Recommend to a librarian

Recommend PNAS to a Librarian#### Purchase options

Purchase this article to get full access to it.