# Bubbles determine the amount of alcohol in Mezcal

G. Rage, Omer Atasi, M. Wilhelmus, J. Hernández-Sánchez, B. Haut, B. Scheid, Dominique Legendre, R. Zenit

► To cite this version:

G. Rage, Omer Atasi, M. Wilhelmus, J. Hernández-Sánchez, B. Haut, et al., Bubbles determine the amount of alcohol in Mezcal. Scientific Reports, 2020, 10 (1), pp.11014. 10.1038/s41598-020-67286-x. hal-04837681

HAL Id: hal-04837681

https://hal.science/hal-04837681v1

Submitted on 13 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
SCIENTIFIC
REPORTS

nature research

<figure><img src="image_2.png" /></figure>

Check for updates

OPEN

Bubbles determine the amount of
alcohol in Mezcal

G. Rage<fnref n="1" />,⁸, O. Atasi²,³,⁴,⁸, M. M. Wilhelmus⁵, J. F. Hernández-Sánchez⁶, B. Haut², B. Scheid², D. Legendre³ & R. Zenit<fnref n="1" />,⁷✉

Mezcal is a traditional Mexican spirit, obtained from the distillation of fermented agave juices. Its preparation has been conducted for centuries in an artisanal manner. The method used to determine the correct alcohol content is of particular interest: a stream of the liquor is poured into a small vessel to induce surface bubbles. These bubbles, known as pearls by the Mezcal artisans, remain stable for tenths of seconds only if the alcohol content is close to 50%. For higher or lower alcohol content, the bubbles burst rapidly. The long bubble lifetime is the result of surfactant-induced surface tension changes. However, the precise mechanism and its relation to alcohol content remain unexplained. In this investigation, the extended lifetime of pearls was studied both experimentally and numerically. It was found that changes in surface tension, density, viscosity (resulting from mixing ethanol and water), and the presence of surfactants are all relevant to extend the bubble lifetime. The dimensionless bubble lifetime was found to reach its maximum value when the Bond number was close to unity, corresponding to 2 mm Mezcal bubbles. These findings show that the traditional empirical method does work. Beyond this, the understanding of the process provides physical insight to many other natural and industrial problems for which the stability of surface bubbles is of importance, such as bio-foams, froth floatation, and volcanic flows.

The consumption of alcoholic drinks can be traced back to the appearance of early humans. While its inception may have resulted from the need to find reliable sources of potable liquids<fnref n="1" />, these drinks have become part of the cultural identity of many societies². Most cultures around the world are known for a 'local' drink. In the case of Mexico, the most popular distilled spirit is Tequila which, belongs to a wider class of distilled agave-based products that are similar in production³ but vary regionally. It was believed that the production of agave-distilled spirits began with the arrival of Europeans by the end of the 16th century⁴. However, recent findings indicate that alcohol distillation was known in Mesoamerica long before, for at least 25 centuries⁵. Our study focuses on Mezcal, which has progressively gained worldwide recognition and has distinctive cultural value in rural villages in Oaxaca⁶. Although its production and denomination are normed⁷, its preparation remains mostly artisanal⁸.

According to popular accounts (informally documented by interviews with artisans) and a few formal reports⁹, the traditional method employed to determine the alcohol content in Mezcal consists of inspecting the lifetime of bubbles that are formed by splashing a jet of the liquor into a small container (Fig. 1(a) and S1 Video in the Supplementary Information section). As a result of the continuous splashing, the liquid surface breaks and so-called pearls form. If the amount of alcohol in the liquid is correct (about 55% volume fraction of ethanol), pearls persist for up to tenths of seconds, see Fig. 1(b). The method is surprisingly accurate.

Interestingly, a similar technique has also been used to determine the alcoholic content in other spirits. Davidson¹⁰, for example, conducted experiments on the foam stability of Bourbon diluted with different amounts of water. The technique is essentially the same as that presented here and shows the same phenomenon: for volumetric contents of alcohol of about 50%, the superficial bubbles are notably more durable than in other

<footnote n="1">¹Instituto de Investigaciones en Materiales, Universidad Nacional Autónoma de México, Apdo. Postal 70-360, Coyoacan, Ciudad de Mexico, 04510, Mexico. ²TIPs (Transfers, Interfaces and Processes), Université Libre de Bruxelles, Avenue F.D, Roosevelt, 50, CP 165/67, 1050, Brussels, Belgium. ³Institut de Mécanique des Fluides de Toulouse, Université de Toulouse, CNRS, Toulouse, France. ⁴Soft Matter, Fluidics and Interfaces, MESA+ Institute for Nanotechnology, University of Twente, P.O. Box 217, 7500 AE, Enschede, The Netherlands. ⁵Department of Mechanical Engineering, University of California, Riverside, 3401 Watkins Drive, Bourns Hall A313, Riverside, CA, 92521, USA. ⁶Instituto de Ciencias Aplicadas y Tecnología, Universidad Nacional Autónoma de México, Circuito Exterior s/n, Ciudad de México, 04510, Mexico. ⁷Center for Fluid Mechanics, School of Engineering, Brown University, 184 Hope St, Providence, RI, 02912, USA. ⁸These authors contributed equally: G. Rage and O. Atasi. ✉e-mail: zenit@brown.edu</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** (a) Traditional technique to form superficial bubbles (pearls); see also the Video S1 in the Supplementary Information section. (b) Zoom of the pearls of Mezcal in a gourd cup (jicara). Images taken by L. Diaz-Damacillo, reproduced with permission.</figcaption></figure>

mixtures. Ahmed and Dickinson¹¹ conducted experiments for whiskey and found similar results. They argue that the changes in bubble lifetime are related to the changes in solubility of surfactant molecules present in this type of beverage. In contrast, Tuinier et al.¹² found an extended foam-life for ethanol volume fractions close to 10%. In none of the previous studies, the precise mechanism responsible for the extended bubble life duration has been explained.

## Lifetime of surface bubbles

Due to its relevance to many processes beyond the production of water-alcohol mixtures, e.g. bubbles on the surface of the ocean¹³,¹⁴, fish nests¹⁵, froth flotation¹⁶ and vulcanology¹⁷, the time of residence of a superficial bubble on a free surface has been extensively studied. How long a bubble remains floating on a surface depends on the drainage of the film between the bubble and the surface, which results from the balance between two competing effects (gravitational and capillary induced drainage) and viscous forces. The film reduces its thickness in time, and when it is sufficiently thin, it spontaneously pierces and breaks. Note that the precise mechanism that determines the film rupture thickness is not known¹⁸.

Two additional factors affect bubble lifetime during Mezcal testing: the alcohol content and the presence of surfactants. Ethanol and water are the main components in Mezcal, with volume fractions of the former ranging from 36 to 55%, according to the norm⁷. The effect of alcohol on bubble lifetime is not fully understood yet as evaporation can either shorten the lifetime of films¹⁹ or extend it²⁰. Indeed, some studies suggest that alcohol can increase bubble stability¹¹, while others have shown it to have a destabilizing effect²¹. While the volatility of ethanol may induce Marangoni flows¹⁹, the properties of water-ethanol mixtures show a non-monotonic behavior with ethanol content²², making it difficult to understand the physical mechanisms behind film drainage. Furthermore, it has been long recognized that surfactants significantly alter the drainage process of bubble films²³, thereby affecting the lifetime. In addition, Mezcal has many other compounds (in small fractions) that can act as surfactants²⁴, such as methanol, acetic acid, ethyl acetate, higher alcohols, esters, ketones, furanes, acetals, aldehydes, phenols, and terpenes as well as proteins, which can have a significant effect in delaying the drainage of bubble films as in the case of beer²⁵.

## Physical mechanisms that produce surface tension gradients.
Three physical mechanisms may induce gradients of surface tension: surfactant transport, heterogenous evaporation, and local changes in temperature. All of these mechanisms affect the time the film takes to drain and, therefore, bubble lifetime. Let us first consider the effect of alcohol concentration. The addition of 10% of mass of ethanol to water, for example, decreases the surface tension of the mixture more than 30% as seen in the study by Vazquez et al.²⁶, where variations in temperature were shown to have a small effect—an increase of 30°C resulted in a surface tension reduction of less than 6%.

In ethanol-water mixtures, the alcohol evaporates at a constant rate regardless of the initial concentration²⁷. Therefore, it would be expected for the evaporation rate to be uniform along the film, avoiding gradients of surface tension. However, the evaporation rate does vary with the thickness of the film as in Mezcal pearls. Such non-uniform evaporation naturally results in a decrease of bubble lifetime. Hence, we do not discard the effect of evaporation. However, given the long lifetimes observed experimentally, we consider this effect as secondary. We can thus argue that surface tension gradients result primarily from the changes in the surfactant concentration. To fully resolve this issue, a numerical analysis would have to account for both evaporation and thermal effects. Such a scheme is beyond the scope of the present work. Instead, based on the arguments above, we account only for the surfactant-induced Marangoni effect. The details are explained in the following sections.
The aim of our investigation is to demonstrate that the traditional technique to assess the ethanol content in Mezcal actually works and to explain the relationship between the lifetime of pearls and alcohol content. Although modern methods of determining alcohol content are accessible and reliable, even for small remote rural communities, the fundamental understanding of the physical mechanisms behind bubble and foam stability is relevant due to their importance to many other fields. In this study, we conducted experiments to determine the lifetime of surface bubbles of many different types of Mezcal and water-ethanol mixtures. To gain insight into the film drainage process, we also conducted numerical simulations and analyzed all results in terms of dimensionless quantities.

**Dimensional analysis.** As proposed by Barenblatt<sup>28</sup>, dimensional analysis can be used to understand a physical phenomenon more deeply. There are several ways in which the relevant dimensionless numbers can be identified. We consider two alternatives below.

**Dimensionless groups by simple inspection.** Drainage results from a balance between viscous forces and either gravitational or surface tension forces. The ratio that determines the relative importance between gravitational and surface tension effects is the Bond number,

$$Bo = \frac{\rho g D^2}{\sigma}, \qquad (1)$$

where $\rho$ is the liquid density, $g$ is the gravitational acceleration, $D$ is the bubble diameter, and $\sigma$ is the surface tension. All of these quantities are attainable via laboratory experiments. Here, $D$ is determined by analyzing images of bubbles from a top view (Methods Section). It should be noted that this metric gives an apparent bubble size rather than an equivalent bubble diameter, which measures the diameter of a spherical bubble with equivalent volume.

Considering surface tension effects, the bubble lifetime can be scaled by incorporating viscous effects leading to:

$$T_{life}^* = \frac{T_{life}\sigma}{\mu D}, \qquad (2)$$

where $\mu$ is the liquid viscosity. One can propose the following functional relationship between $T_{life}^*$ and $Bo$:

$$T_{life}^* = \Phi(Bo). \qquad (3)$$

Now, as explained above, the ever-presence of surfactants plays an important role in the film drainage process. If one defines $\Delta\sigma$ as the change in surface tension resulting from surfactant gradients along the interface, and the relative variation as

$$\Pi = \frac{\Delta\sigma}{\sigma}, \qquad (4)$$

one can therefore argue that

$$T_{life}^* = \Phi'(Bo, \Pi). \qquad (5)$$

Here, we aim to understand the effect of each group by considering that

$$T_{life}^* = \Phi(Bo) \cdot \Psi(\Pi). \qquad (6)$$

In other words, we shall consider that the effect of surfactants can be de-correlated from the drainage. This fact is not obvious and represents an important assumption of our study. However, as discussed below, the assumption appears to be reasonable. Furthermore, the analysis indicates that the surfactant-induced (Marangoni) stresses should be considered to determine an appropriate timescale.

**Dimensionless groups from the analysis of the flow type.** The paper by Champougny et al.<sup>23</sup> discusses the unsolved problem of bubble drainage in the presence of surfactants. They propose an ad-hoc model based on an extrapolation length $\lambda$ that allows describing the transition from stress-free interfaces ($\lambda = \infty$) to no-slip interfaces ($\lambda = 0$). The flow between stress-free interfaces is a plug flow at leading-order, dominated by extensional viscous stresses. Debregeas et al.<sup>29</sup> have shown that the film thickness of such a 'bare' bubble decays exponentially in time and that the film always ruptures at the apex. Contrarily, the flow between no-slip (rigid) interfaces is a Poiseuille-like flow dominated by viscous shear stresses (as described, for instance, in<sup>30</sup>). Now, surfactant-induced Marangoni stresses have been shown to rigidify, at least partially, the interfaces, such that the flow is essentially a shear flow (e.g., Champougny et al.<sup>31</sup> for film formation or Atasi et al.<sup>32</sup> for confined bubbles in microchannels, both in the presence of surfactants). Even traces of surfactants, such as impurities, have been found to induce dominant Marangoni stresses<sup>33</sup>. Additionally, local thinning mechanisms of the film have been observed in soap bubbles<sup>18</sup>, reminiscent to marginal regeneration<sup>34</sup>, and eventually leading to the rupture of the film away from the apex. Based on these arguments, we conjecture that viscous shear stresses dominate the film drainage. This hypothesis is corroborated below via numerical simulations.
In absence of a simple model available that accounts for the presence of surfactants and with the aim to find an appropriate timescale, let us consider a drainage flow dominated by viscous shear, balanced either by gravity or capillarity. Starting with a gravity-driven drainage, the momentum balance at leading order in the lubrication approximation writes

$$ \mu \partial_{yy} u = \rho g, \qquad (7) $$

where the cross-stream coordinate y in the film can be scaled with the critical thickness for rupture $h_{rup}$. The choice of this length scale is justified by the fact that the drainage dynamics is not significantly influenced by the initial condition since the main contribution to the lifetime occurs at the later stage, i.e. when the film is the thinnest, which is also when the viscous dissipation dominates. Note, however, that $h_{rup}$ is not a quantity that can be obtained directly from the experiments. By considering the speed at which the film retracts after it ruptures, as explained in the Methods Section below, an estimate of this thickness can be inferred. For simplicity, we consider that the rupture always occurs at the same thickness value. Considering that $u \sim D/t_g$, the timescale for drainage driven by gravity can be obtained as:

$$ t_g = \frac{\mu D}{\rho g h_{rup}^2}. \qquad (8) $$

Similarly, the momentum balance for a capillary-driven drainage is given by

$$ \mu \partial_{yy} u = \partial_x P, \qquad (9) $$

where the streamwise coordinate x can be scaled with D. If we assume that the pressure P results from a Laplace effect, we can write $P \sim \sigma/D$. By estimating that in this case $u \sim D/t_c$, the timescale for drainage driven by capillarity is given by

$$ t_c = \frac{\mu D^3}{\sigma h_{rup}^2}. \qquad (10) $$

Comparing these two timescales leads to the Bond number, namely $Bo = t_c/t_g$, exactly as defined in Eq. (1).

It is interesting to note that the Bond number plays two roles in this problem, a static and a dynamic one. On the one hand, it allows evaluating the static shape of the interfaces, namely a spherical bubble underneath an almost undeformed liquid surface for $Bo \ll 1$, and a deformed bubble under a deformed liquid surface for $Bo \gg 1$. On the other hand, as explained above, it provides a way to evaluate the dominant driving force for drainage, namely the capillary force for $Bo \ll 1$ and the gravity force for $Bo \gg 1$. The striking feature of the present problem is that $Bo \sim 1$ indicates a transition which, as shown below, coincides with the maximum dimensionless bubble lifetime.

A different dimensionless lifetime can, therefore, be defined in terms of the capillary timescale, $t_c$:

$$ T_{life}^{**} = \frac{h_{rup}^2}{D^2} \frac{T_{life} \sigma}{\mu D} = \varepsilon^2 T_{life}^* \qquad (11) $$

where $T_{life}^*$ is defined in Eq. (2). We also introduce $\varepsilon = h_{rup}/D \ll 1$, a measure of the slenderness of the film. Equation (11) indicates that the viscous forces act perpendicularly to the capillary/gravity forces, as reminiscent to shear-dominated flow. Conversely, forces in the case of an extensional flow would act along in the flow direction affecting the drainage; for such a case $\varepsilon$ would be close to unity. Consequently, $T_{life}^*$ would be the appropriate group for an extensional drainage. In the present problem, $\varepsilon^2 \sim 10^{-4}$, therefore $T_{life}^*$ is more appropriate in the present context.

As will be shown later, the data shows a clear transition in trend at around $Bo \sim 1$. Using the $T_{life}^{**}$ scale, we can derive the trends in both limits of large and small Bond numbers. In the limit of large Bond number the drainage is governed by gravity, therefore, using Eq. (8), one can write

$$ T_{life}^{**} \propto Bo^{-1} \quad (Bo \ll 1). \qquad (12) $$

On the other hand, when the Bond number is small the drainage is dominated by capillarity, which, considering Eq. (10), implies that

$$ T_{life}^{**} \propto 1 \quad (Bo \ll 1). \qquad (13) $$

These two trends will be contrasted against the experimental results below.

**Review of relevant models to predict the bubble lifetime.** As stated in the Introduction, many studies have addressed the issue of bubble lifetime. Here, we summarize some relevant studies. Note that we do not aim to comprehensively summarize the state-of-the-art of the subject. To draw comparisons with our experiments and among the different models, we recast the predictions in terms of the dimensionless groups discussed above and recognize that for a given liquid, varying the Bond number amounts to varying the bubble size.
Large bubbles. Bo ≫ 1

In the limit of large Bond numbers, the bubbles are large and the film near the apex is essentially uniform and the drainage is driven by gravity.

For the case of rigid interfaces, an analytical solution exists for the dimensionless lifetime as follows³⁵

$$T_{life}^{**} = \frac{3}{2}Bo^{-1}. \qquad (14)$$

Notably, this prediction underestimates most of the experimental times shown later in the Results section, which indicates that the Marangoni stress in the experiments is large enough to retard the drainage longer than what rigid interfaces would do. This is only possible if the mean flow is upwards, i.e. towards the apex, entrained by Marangoni stresses. Such gradient-tension-induced stresses would need to built up prior to the drainage phase³⁶. The simulations shown below show this type of behavior.

In the case of stress-free interfaces, Kocarkova et al.³⁷ calculated the evolution of the film thickness considering an extensional flow in the film, assumed to be uniform and axisymmetric. They observed an exponential thinning, from which the dimensionless lifetime can be shown to be:

$$T_{life}^{*} \sim Bo^{-1}. \qquad (15)$$

Champougny et al.²³ extended the analysis of Kocarcova et al.³⁷ considering certain levels of surface rigidification of the interfaces to account for the presence of surfactants. They also found an exponentially thinning of the film with time, which leads to the same functional relation of Eq. (15) but the proportionality constant was larger in the case of surfactants. Interestingly, they also found that the puncture of the bubble film changed from the apex to the foot of the bubble as the amount of surfactants increased.

Finally, Lhuissier and Villermaux¹⁸ considered the case of bubbles in water. They argued that the scaling of a viscous drainage was not appropriate for water due to the phenomenon known as marginal regeneration³⁸. In short, the balance of capillary pressure from the film curvature, the meniscus at the foot of the bubble and surface tension gradients lead to a non-uniform film thickness which causes the appearance of a localized pinching. Moreover, Lhuissier and Villermaux¹⁸ recognized that the film breakup could also be influenced by the Bénard-Marangoni convection flows within the film. Considering the fluctuations from marginal regeneration convection cells and probabilistic arguments of the puncture breakup mechanism, they found

$$T_{life}^{*} \sim Bo^{-1/4}, \qquad (16)$$

which also indicates a decreasing lifetime with an increase of Bo, as qualitatively obtained in Eq. (12).

Small bubbles. Bo ≪ 1

In the limit of small Bond numbers, no simple model exists for the bubble lifetime because in this limit the drainage is driven by capillary forces which also induce surface deformations. The model of Howell³⁹, considers this case but assumes an extensional flow, which leads to

$$T_{life}^{*} = \kappa_{1}Bo^{1/2} \quad (Bo \ll 1), \qquad (17)$$

where k₁ is a constant that depends on the initial and rupture thicknesses. This expression indicates that the rupture time increases with Bo: the larger the bubble, the longer it will take to burst. This calculation assumes that the film drains uniformly and axisymmetrically without surfactants. The small Bo limit indicates that the bubble is nearly spherical and is mostly immersed below the liquid free surface. Note that this trend is different from that obtained considering a shear flow (Eq. 13).

Transition-sized bubbles. Bo ~ 1

From the brief review shown above, an important conclusion can be reached: the trends are opposite for small and large Bond numbers. This fact indicates that there is a transition at a certain value of Bo. From our arguments, confirmed below by the experiments, the critical Bo-value is around unity, where capillary and gravity forces are of comparable magnitude to drive the drainage. And this transition corresponds to a maximum dimensionless lifetime, as suggested by the cross-over between the models given above.

## Experimental Results

**Pearl formation.** Figure 2(a) shows an image of the formation process of pearls (see also Video S2, Supplementary Information section), considering the controlled experiment described in the Methods Section. The fluid jet fragments into droplets, which continuously splash against the surface of the liquid. The continuous splashing induces the formation of air cavities that, in turn, form bubbles that eventually rise to the surface. However, only bubbles of a certain size (of approximately 2 mm in diameter) remain on the surface for longer times.

The process of air entrainment from a plunging jet has been studied in detail⁴⁰⁻⁴², due to its relevance for many flow phenomena such as aeration of the ocean and other water reservoirs. It is believed to be well understood. Hence, we did not pursue a more complete study of the formation process. Instead, we focus on the residence time of bubbles on the surface.

**Bubble lifetime.** Figure 2(b) shows a time sequence of a typical Mezcal bubble at the moment of bursting (see also Video S3, Supplementary Information section), produced with the system described in the Methods Section. It was observed that bursting initiates with a small hole rupturing the film, which opens quickly until it
<figure><img src="image_4.png" /><figcaption>(a)</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>(b)</figcaption></figure>

**Figure 2.** (a) Reproduction of the formation process in the laboratory. (b) The bursting of a surface bubble. The images were taken from the side and top, simultaneously. The time between frames is 0.4 ms. The diameter of the bubble is 1.9 mm, corresponding to $Bo = 1.06$. The liquid is Mezcal M1, from Table 1.

<table><thead><tr><td>Fluid</td><td>type or brand agave species</td><td>density ρ, kg/m³</td><td>viscosity μ, mPa s</td><td>surface tension σ, mN/m</td><td>ethanol content % vol</td></tr></thead><tbody><tr><td>M1 (●)</td><td>Los Sánchez espadín</td><td>919.0</td><td>3.07</td><td>30.6</td><td>55</td></tr><tr><td>MA1 (○)</td><td>M1 + W or E espadín</td><td>906.2 to 965.8</td><td>1.83 to 3.03</td><td>25.0 to 41.0</td><td>22 to 82</td></tr><tr><td>MA2 (○)</td><td>M1 + W/E (55/45) espadín</td><td>897.4 to 901.2</td><td>2.82 to 2.89</td><td>28.4 to 29.1</td><td>53 to 57</td></tr><tr><td>M2 (◊)</td><td>El Amate espadín</td><td>914.4</td><td>2.98</td><td>28.7</td><td>48</td></tr><tr><td>M3a (▽)</td><td>Los peregrinos (cuerpo) espadín</td><td>912.6</td><td>2.66</td><td>29.2</td><td>41</td></tr><tr><td>M3b (□)</td><td>Los peregrinos (punta) espadín</td><td>926.9</td><td>3.07</td><td>29.2</td><td>58</td></tr><tr><td>M3c (△)</td><td>Los peregrinos (cola) espadín</td><td>933.8</td><td>2.34</td><td>32.7</td><td>31</td></tr><tr><td>M4 (★)</td><td>Los peregrinos (MC) Madre Cuishe</td><td>928.3</td><td>3.08</td><td>29.9</td><td>52</td></tr><tr><td>M5 (×)</td><td>Los peregrinos (reposado) espadín</td><td>949.1</td><td>2.58</td><td>34.0</td><td>39</td></tr><tr><td>M6 (◁)</td><td>Real Minero largo</td><td>850.0 to 950.0</td><td>2.30 to 3.06</td><td>25.9 to 31.6</td><td>30 to 70</td></tr><tr><td>M7 (▷)</td><td>Matatlan tobalá</td><td>850.0 to 950.0</td><td>2.30 to 3.06</td><td>25.8 to 32.0</td><td>30 to 70</td></tr><tr><td>WE (□)</td><td>W-E mixtures N/A</td><td>791.0 to 1000.0</td><td>1.10 to 3.01</td><td>23.0 to 72.0</td><td>0 to 100</td></tr></tbody></table>

**Table 1.** Physical properties of liquids used in the investigation. M: Mezcal, MA: modified Mezcal, W: water, E: ethanol.

retracts completely. The puncture does not typically appear at the apex of the bubble; the thinning of the film is not uniform arguably due to the partial regeneration mechanism¹⁸.

Several experiments were conducted under the same nominal conditions (same fluid and same bubble diameter) to calculate the typical value of the rupture thickness, according to the scheme described in the Methods
<figure><img src="image_6.png" /><figcaption>**Figure 3.** (a) Bubble lifetime as a function of ethanol concentration for Mezcal (M1, ●); Mezcal with the addition of water or ethanol (MA1, ○); Mezcal with the addition of 55% water-ethanol mixture (MA2, ○). (b) Bubble lifetime as a function of ethanol concentration for all tested liquids. Symbols as in Table 1.</figcaption></figure>

Section. It was found that the thickness was $h_{rup} = 23.5 \pm 4.8$ µm for the image shown in Fig. 2(b). This value is relatively large, compared to what has been measured for the case of water and seawater¹⁸. This value was used as a reference for the numerical simulations and is also discussed in the scaling of the bubble lifetime.

Experiments were initially conducted using Mezcal with the 'correct' ethanol volume fraction, identified as M1 (see Table 1). The mean lifetime in this case was 28.1 ± 12.5 s, depicted in Fig. 3(a) by the filled red circle. To vary the amount of ethanol in this liquid, either pure water or ethanol was added. The measured lifetimes for this ethanol-adjusted Mezcal are shown Fig. 3(a), denoted as liquid MA1 from Table 1. The figure shows the measured bubble lifetime as a function of ethanol volume fraction. As the amount of ethanol increases or decreases, the bubble lifetime decreases. A clear sharp maximum was observed for the Mezcal sample with 55% ethanol. It should be noted that even though the amount of other components was diluted (by adding water or ethanol), increasing or decreasing the alcohol content led to a measurable variation of the lifetime. Therefore, based on these observations, we conclude that the traditional technique does work: the bubble lifetime shows a maximum value for a certain intermediate value of ethanol in Mezcal. Although the experiments were conducted in a controlled environment (standard laboratory conditions), a large variability of the measured lifetimes was found. Large scatter is often reported in experimental reports of the lifetime of bubbles, especially with water and other surface-contaminated liquids¹³,¹⁸. To our knowledge single bubble lifetime measurements for water-alcohol mixtures have not been reported to date, but similarly large scatter can be expected due to evaporation effects¹⁹. In the present case, the large uncertainty bars in our measurements may also result from the use of unfiltered artisanal Mezcal. According to the norm⁷, particulate matter is expected to be present. The presence of particles, and maybe fibers, would also explain the relatively large bursting thickness of the film.
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Bubble lifetime as a function of ethanol concentration for water/ethanol mixtures: experiments (□, D = 2.77 ± 0.83 mm), error bars denote the standard deviation of the measurements; numerical results (○, D = 2 mm) and (−○−, D = 3 mm) symbols, for a surfactant concentration of 0.025 mol/l. The dashed line shows the trend of the change of viscosity of water-ethanol mixtures, according to [22] (right axis).</figcaption></figure>

To evaluate the effect of surfactants and other components, a second set of experiments was conducted. In this set, different amounts of a 55% water-ethanol mixture were added to the sample M1, also shown in Fig. 3(a) and denoted as liquid MA2 in Table 1. In this manner, the amount of ethanol remains approximately constant but the quantity of surfactants is diluted with respect to the original sample. The result is clear: the pearl lifetime is significantly reduced compared to the original sample, even though the volumetric concentration of ethanol is nearly the same. We can, therefore, argue that the amount of surfactants is an important factor to determine the lifetime of pearls.

To begin elucidating the underlying physical mechanism of the process, experiments were performed considering mixtures of water and ethanol, varying the volume fraction of ethanol ranging from pure water to pure ethanol, as shown in Fig. 4. It was observed that even when no surfactants are present, the maximum pearl lifetime occurs at approximately 55% of ethanol fraction. Note, however, that the lifetimes for these mixtures are one order of magnitude smaller than those observed for Mezcal. Assuming that the dimensionless bubble lifetime (Eq. 2) is only a function of the Bond number (Eq. 1), we can argue that the lifetime of pearls is directly proportional to the fluid viscosity. While density and surface tension of water-ethanol mixtures have a monotonic behavior with ethanol fraction (from the pure water to the pure ethanol values), viscosity does not: its value increases from that of water, reaching a maximum at around 55% of ethanol, to then decrease to reach the value of pure ethanol²². The dotted line in Fig. 4 shows how viscosity varies (plotted on the right scale of the figure) with the volume of ethanol. Therefore, we can argue that both fluid properties and surfactant concentration affect the bubble lifetime. Consequently, it may be inferred that the maximum lifetime observed at the ethanol content at 55% is related to the fact that the liquid viscosity has a maximum value at that concentration.

More experiments were conducted with other Mezcal types, from different regions and agave species. For two cases, the alcohol content was also modified in the same manner as described before. The results are shown in Fig. 3(b), where the bubble lifetime for all liquids is shown as a function of alcohol content. The result is essentially the same, the bubble lifetime shows a clear maximum at around 55% of ethanol volume fraction.

**Numerical Results.** From the results shown above, it is clear that the bubble lifetime is both dependent on the fluid properties and the amount of surfactants. However, to elucidate the details of the phenomena, we conducted a numerical study. Details of the numerical scheme can be found in the Methods Section.

First, the residence time of bubbles was computed for water-ethanol mixtures, with the aim to validate the numerical method and to gain some insight into the drainage process. Figure 4 shows the numerically obtained bubble lifetime, considering two bubble sizes which are close to the experimental values. The numerical lifetime was calculated by following the procedure described in the Methods Section. The simulations were obtained for a small amount of surfactant, as it is suspected to be the case for laboratory-clean aqueous mixtures. In very close agreement with the experiments, the bubble lifetime increases as the ethanol volume fraction increases, reaching a maximum for concentrations close to 55% of ethanol, to then decrease for higher concentrations. Given the variability of bubble size in the experiments, the quantitative agreement between experiments and numerical results is remarkable.

From the experiments, it is clear that the amount and type of surfactants affect the bubble lifetime. We conducted numerical simulations with increasing amounts of surfactants, as shown in Fig. 5(a). The bubble lifetime increases as the ethanol volume fraction increases, reaching a maximum for concentrations close to 55% of ethanol, to then decrease for higher concentrations, in close agreement with the experiments. Most importantly, the maximum lifetime increases with the amount of surfactants; the lifetime increases sharply with the concentration
<figure><img src="image_8.png" /><figcaption>**Figure 5.** (a) Lifetime as a function of ethanol content, for different concentrations of surfactants. (b) Maximum lifetime as a function of surface elasticity, *E*; the time is normalized with the case for *E* = 0.1. For all cases, *D* = 2 mm, *c*₀ = 0.025 mol/l, alcohol content of 55%, *Bo* = 1.24.</figcaption></figure>

but appears to not increase significantly beyond a certain value of surfactants. These results are also in good qualitative agreement with the experiments. We have conducted many more simulations and analyses to further understand this effect. These results can be found in Atasi *et al.*<sup>43</sup>.

An interesting and challenging aspect of the Mezcal case is that, in addition to the amount of surfactants, the specific type of it is not known. Furthermore, there may be changes of the surfactant type depending on the species of Agave and production region for each Mezcal. To gain some general understanding in this respect, we also conducted numerical simulations varying the surface elasticity, *E*. This quantity measures the interface resistance to deformation<sup>41</sup>. Fig. 5(b) shows the bubble lifetime as a function of the surface elasticity parameter, considering a fixed value of ethanol content. Clearly, the lifetime increases significantly with the elasticity of the surface. It is important to note that this effect has not been discussed in-depth in the literature. Some initial tests and discussions can be found in Atasi *et al.*<sup>43</sup>.

The numerical simulations can also provide significant insight into the complex drainage process of bubbles in Mezcal. In addition to the effects on the lifetime, details of the flow within the film can also be analyzed. Figure 6 shows the fluid velocity (tangent to the interface) as a function of position within the film, at a certain angle from the apex ($\theta$ = 0.35 radians) at different times, with and without surfactants. An important difference between the two cases is that speed is significantly larger and uniform (across the thickness) for the case of a bubble without surfactants, in comparison to the case with surfactants. The reason for this reduction stems, first, from the immobilization of the surfaces due to surfactants. Moreover, when surfactants are present, for the angle shown, the fluid at the inner edge of the film (s = 0) can move in the opposite direction to that of the gravitational drainage, even for early times (t = 37 ms). For later times, t = 102 ms the fluid moves at a very small downward speed all across the gap. Note that the simulations consider the effect of surfactant transport but do not contemplate evaporation
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Velocity profile across the film. Tangential fluid velocity, *V*<sub>s</sub>, as a function of the normalized distance across the film, *s*, at 0.35 rad from the bubble apex (taking the bubble centroid as reference), for different times and surfactant concentrations. *s* = 0 corresponds to the bubble surface and *s* = 1 refers to the external free surface. The black and blue lines correspond to the cases with and without surfactants, respectively. The dashed line shows the case with surfactants at a long time.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 7.** Shape of bubbles obtained numerically for different values of Bo.</figcaption></figure>

or temperature gradients. Despite this simplification, the lifetimes obtained numerically are in good agreement with those found experimentally. Therefore, we can argue that surfactant transport dominantly affects the drainage process.

Additionally, thanks to the unprecedented access to other properties of the flow from the numerical simulations, it is also possible to track the deformability of the interface. We have identified that the film becomes non-uniform along the bubble during the drainage. In particular, a neck region appears for a sufficiently large Bond number, associated to a local thinning of the film as shown in Fig. 7(b,c), which are used to validate the numerical results. Also, the position of the neck changes position from the top of the bubble toward the liquid bath, as the Bond number increases. The appearance of the neck region has been documented for the case of films³⁸, but not as extensively for the case of surface bubbles¹⁸. This thinning zone leads to the rupture of the film which, as shown in Fig. 2(b), generally does not occur at the bubble apex. The results indicate that the position where the rupture starts is also a function of the Bond number.

## Discussion

Let us now recast the bubble lifetimes in dimensionless terms. Considering the scaling from simple inspection (also arising from the extensional flow analysis), the dimensionless lifetime can be calculated according to Eq. (2). Figure 8(a) shows the bubble lifetime, $T_{life}^*$, as a function of the Bond number, *Bo*, defined in Eq. (1). Despite the large dispersion of the data, several key features can be readily identified. First, the dimensionless lifetime of the unaltered Mezcal is higher than any of the other cases, when either water or ethanol was added, or even for other Mezcal samples with other ethanol content. The dimensionless lifetime appears to increase for small values of the Bond number, to reach a maximum at around *Bo* ≈ 1, and then decrease as the value of *Bo* continues to increase.

The Bond number allows us to evaluate the driving force for drainage, namely the capillary force for $Bo \ll 1$ and the gravity force for $Bo \gg 1$. The striking feature of the results in Fig. 8(a) is that $Bo \sim 1$ indicates a transition, which coincides with the maximum dimensionless bubble lifetime. The Bond number compares the two main driving forces for drainage: gravitational drainage for large *Bo* and capillary-induced drainage for small *Bo*. Interestingly, the maximum dimensionless lifetime is reached when both effects are of the same order of magnitude.

To further support the existence of a critical *Bo*, we compared our results with predictions from the literature, shown in Fig. 8(a). For $Bo < 1$ the prediction by Howell, Eq. (17), shown by the solid line, shows an increasing
<figure><img src="image_11.png" /><figcaption>**Figure 8.** (a) Dimensionless lifetime, $T_{life}^*$ as a function of Bond number. The black lines show trends from different predictions: solid line, Eq. (17); dashed line, Eq. (15); and dotted line, Eq. (16). (b) Compensated dimensionless lifetime, $T_{life}^{**}Bo^{3/2}$ as a function of Bond number, *Bo*. Symbols as in Table 1. For clarity, error bars are omitted. The plus symbols (+) are the results of numerical simulations, for a surfactant concentration of 0.025 mol/l, for two different bubble diameters and different ethanol concentrations. The blue dashed and blue dashed-dotted lines show compensated trends predicted by Eq. (13) and Eq. (12), respectively ($T_{life}^{**}Bo^{3/2} \sim Bo^{3/2}$ and $T_{life}^{**}Bo^{3/2} \sim Bo^{1/3}$, respectively.).</figcaption></figure>

lifetime with Bond number. On the other hand, for $Bo > 1$, Eqs. (15) and (16), dashed and dotted line respectively, predict a decreasing trend of the bubble lifetime with *Bo*. These contrasting trend predictions corroborate the existence of a value of *Bo* at which $T_{life}^*$ has a maximum value.

The numerical results, in dimensionless form, are also shown in Fig. 8, along with the experiments. The numerical results show an increasing trend for small *Bo*, but at around *Bo* ≈ 1, a clear change of behavior is observed. Note that all simulations, conducted for two different bubble diameters and different ethanol concentrations, collapse into a single band when presented in dimensionless terms. The lifetime of the numerical results is shorter than that in the experiments resulting from the assumptions considered in the simulations. Nevertheless, the qualitative agreement is noteworthy. The cross-over for which gravity and capillary effects are of the same order of magnitude shows a non-monotonous behavior for which a maximum dimensionless lifetime has been observed, both numerically and experimentally. Strikingly, *Bo* ≈ 1 corresponds to a bubble size of about 2 mm, i.e. the size for which the lifetime is the longest in the experiments.

One important feature of the dimensionless bubble lifetime, $T_{life}^*$, is its large value, in the order of $10^4-10^5$ time-scale units. This is an indication that the characteristic time $\mu D/\sigma$ is orders of magnitude smaller than the bubble lifetime. Hence, its magnitude cannot be used to characterize the phenomena properly. Therefore, we use $T_{life}^{**}$ instead, as defined by Eq. (11). Furthermore, if we consider a compensated dimensionless lifetime, given by the product $T_{life}^{**}Bo^{3/2}$, the data dispersion can be reduced by eliminating the dependence on the bubble diameter.
Figure 8(b) shows the compensated dimensionless lifetime as a function of the Bond number, Bo. The time scaling considers the rupture thickness of the film. As discussed above, the factors that determine the rupture thickness are not completely understood even for pure liquids¹⁸. In the present case, the presence of particles and surfactants, pose additional difficulties⁴⁵. Hence, since we cannot reasonably assess how the rupture thickness changes for different liquids, we assume that $h_{rup} = 23.5 \pm 4.8$ µm is the same for all cases. This value corresponds to the measured film thickness for the case of unaltered Mezcal. The data presented in this form shows that the dimensional lifetime increases up to Bo of order unity, which is in good agreement with the scaling argument of Eq. (13). For bubbles for $Bo > 1$, the data clearly shows that the compensated dimensionless lifetime is relatively constant. In contrast, the prediction from Eq. (12), shows a slightly increasing trend. Note also that the magnitude of the dimensionless times is of O(10) for all experimental points, which indicates that the characteristic time $t_c$ (Eq. 10) is more representative of the flow phenomena in the draining film. Again, this argument supports the idea that the flow in the film is shear-dominated.

The numerical results are also shown in Fig. 8(b). As in the previous case, the lifetimes obtained from the simulations are shorter for all cases. However, and most importantly, the numerically obtained lifetimes collapse into a single band of data indicating that the scaling is correct. Furthermore, the trend found with the numerical results matches closely with the scaling arguments: the compensated dimensionless time increases for $Bo < 1$ and is roughly constant for $Bo > 1$. One of the advantages of the numerically obtained times is that the data scattering is greatly reduced, which allows us to observe the trend of results more clearly. From these results, it is further demonstrated that there is a critical Bo number, around unity, that characterizes a long bubble lifetime.

In summary, we conducted an experimental and numerical investigation to elucidate the phenomena behind the 'pearls of Mezcal'. By observing the extended lifetime of superficial bubbles, artisanal Mezcal makers have been able to identify the correct amount of ethanol in this distilled spirit without the need for any other measuring device. Intrigued by this empirical technique, we conducted experiments to test the technique under controlled conditions. Single bubbles were generated into a small container, and the time from reaching the surface until bursting spontaneously was measured. It was found that, indeed, bubbles in samples with ethanol volumetric content close to 55% had a characteristically long lifetime. We have found that, in dimensionless terms, the lifetime of a surface bubble increases with Bond number, reaching a maximum value at $Bo \approx 1$, and then decreases for large values of Bo. We found that both an increase of the liquid viscosity and the presence of surfactants are needed to observe a long lifetime of bubbles. These factors explain why pearls in Mezcal have a particularly long life at a certain concentration of alcohol and a certain size. Clearly, the artisanal technique is the result of observation, empiricism and tradition. The insight gained in investigating this technique may help elucidate the fundamental mechanisms of bubble formation and stability in many other contexts. For example, the lifetime of surface bubbles could be used as a diagnostic tool to infer the presence of surfactants in a liquid: if the lifetime is larger than that expected for a pure/clean liquid, then the liquid is most likely contaminated. For instance, the bubble nests of some tropical fish have long lifetimes thanks to the presence of bio-surfactants¹⁵; short-lived nest could indicate environmental contamination or health problems of the fish. Similarly, the contamination levels of water bodies could be correlated to long-lived foam and bubbles, as discussed by Zenit and Rodriguez-Rodriguez ⁴⁶. Furthermore, the lifetime of surface bubbles in applications such as froth floatation or anti-bubble formation could be adjusted by choosing the 'correct' surfactant to tune particle or gas capture rate¹⁶,⁴⁷. Finally, recent reports indicate that the rupture process of surface bubbles plays a significant role in the spread of infectious diseases⁴⁸; hence, it is important to fully understand the nature of this process. The production of aerosols resulting from bursting bubbles could be reduced if the effects of surfactants are fully understood. The results presented in this investigation may shed some light into such processes of current importance.

## Materials and methods

**Informed consent.** The person who appears in Fig. 1(a) and in the S1 Video in the Supplementary Information section has given informed consent to be shown in the paper. He understands that the images may lead to identification.

**Test fluids and physical properties.** The physical properties of all the test liquids are reported in Table 1. The origin of the Mezcal samples is also given. According to the Mexican norm, they all come from the Oaxaca state region. Seven Mezcal types were considered; three cases (M1, M6 and M7) were purposely altered to change their water or alcohol content. Also, pure water/ethanol mixtures were tested.

The surface tension was measured with a tensiometer (DynaTester, SITA). The density was measured with a 25 ml pycnometer. The alcohol content was inferred from the density measurement, considering room temperature, $T_{room} = 23$ °C. The viscosity was not measured; considering the alcohol content, its value was assumed to correspond to that of a water-ethanol mixture and was obtained from tables²².

**Visualization of bubble formation process.** The traditional technique to evaluate the alcohol content is shown in S1 Video (Supporting Information) and Fig. 1(a). The *Maestro Mezcalero* issues a stream of fluid from a long reed, with a round opening of about 2 mm in diameter. The jet impinges onto a small gourd cup, of about 10 cm in diameter. As the cup fills up, a pool of fluid is formed of several centimeters in depth. To reproduce the process in a controlled manner, the traditional reed was replaced by a glass pipette of 25 ml, with approximately the same exit diameter as the reed (2 mm). The sample of Mezcal (or other test liquids) of about 20 ml was placed in the pipette and was held fixed in a vertical position with a laboratory holding bracket. Below the exit of the pipette, a 10 × 10 cm² square transparent container was placed, filled with the same liquid to a depth of 2 cm. The process was filmed with a high-speed camera (FASTCAM-APX, Photron) using diffuse back-lighting at 1500 fps.
**Measurement of lifetime.** To accurately measure the lifetime of single bubbles, a simple experimental arrangement was used. A short cylindrical glass container of 1.6 cm in diameter and 1.2 cm in height was filled with the test liquid; the rim of the container was roughened to make it slightly hydrophobic. The amount of liquid was slightly larger than the volume of the container, such that a convex meniscus was formed. The arrangement was placed inside a closed container to minimize environmental disturbances during the experiment. At the side of the container, a needle was inserted through the wall via a plastic stopper plug. Bubbles were formed by slowly pushing air through the needle with a syringe. Since the free surface was slightly convex, when a bubble reached the surface it moved to the center of the container where it could be filmed. The bubble traveled approximately 1 cm through the liquid, lasting approximately 0.03 s from its formation until it reached the surface; the bubble moved slightly on the surface quickly reaching a stationary position. The bubble lifetime was measured from the moment when the bubble no longer moved and until it burst.

Experiments were performed with needles of different gauges ranging from 159 µm to 210 µm of internal diameter, to produce bubbles with slightly different diameters. The setup was illuminated from the bottom with a LED light. The bubble diameter was measured from the image obtained from the top; this measurement overestimates the bubble diameter of an equivalent spherical bubble by approximately 13% for the highest Bond numbers⁴⁹. The process of rupture was measured with two synchronized high-speed cameras (Phantom SpeedSense), one aligned from the top and the other from the side. Different recording rates were considered: since the rupture time could take tens of seconds, an ordinary 30 fps was used to determine the bubble lifetime. To measure the rupture speed of the film, the recording rate was as high as 5,000 fps. All experiments were performed under standard laboratory conditions. The container was thoroughly rinsed with distilled water prior to each experiment.

**Thickness of the film during bursting.** The thickness of the film can be inferred from the speed at which the leading piercing rim moves as the bubble bursts. Considering the Taylor-Culick velocity⁵⁰,⁵¹:

$$V = \sqrt{\frac{2\sigma}{\rho h_{rup}}} \qquad (18)$$

where $\sigma$ and $\rho$ are the surface tension and liquid density, respectively. $h_{rup}$ is the thickness of the film. Experiments were conducted using fluid M1 and bubbles with the same diameter ($D=1.9$ mm).

**Numerical simulations.** Direct numerical simulations were conducted by solving the Navier-Stokes equations coupled with the Level-Set method. We refer the reader to Atasi et al.³² and Abadie et al.⁵² for a detailed description of the method and its validation. Briefly, Navier-Stokes equations are solved for two Newtonian and incompressible fluids using the finite volume method (second-order accurate in time and space). Continuity is ensured through a projection method, and the capillary contribution is considered through the classical Continuum Surface Force method. The interface position is tracked using the Level-Set method where the transport of the signed distance to the interface is controlled through the re-distancing techniques. The novel aspect of the numerical scheme presented here is its ability to account for the surfactant concentration both in the liquid and on the interface as given by⁵³:

$$\frac{\partial C}{\partial t} + V_c \cdot \nabla C = \nabla \cdot (D_c \nabla C) \qquad (19)$$

$$\frac{\partial \Gamma}{\partial t} + \nabla_s \cdot (V_s \Gamma) = D_s \nabla_s^2 \Gamma + S_{\Gamma} \qquad (20)$$

where $C$ is the surfactant concentration in the liquid phase, $\Gamma$ is the surfactant concentration on the gas-liquid interface, $D_c$ and $D_s$ are the diffusion coefficients of the surfactants in the liquid phase and along the interface, respectively, $V_c$ is the velocity field of the liquid phase, $V_s$ is the projection of $V_c$ on the tangent to the interface, $\nabla_s = (\bar{I} - n \times n)$. $\nabla$ is the surface gradient operator and $S_{\Gamma}$ is the flux of surfactants from the liquid phase to the interface, due to the adsorption/desorption of the surfactants, i.e., $S_{\Gamma} = (D_c \mathbf{n} \cdot \nabla C)|_{\Gamma}$, where the subscript $I$ denotes the bubble-liquid interface. It is given by⁵⁴:

$$S_{\Gamma} = k_a C_I (\Gamma_{\infty} - \Gamma) - k_d \Gamma, \qquad (21)$$

where $k_a$ and $k_d$ are the adsorption and desorption kinetic constants, respectively, and $C_I$ is the surfactant concentration in the liquid in contact with the interface.

It is assumed that the surface tension depends on the surfactant concentration on the interface according to an equation of state derived from the Langmuir adsorption isotherm:

$$\sigma = \sigma_0 \left( 1 + Elog \left( 1 - \frac{\Gamma}{\Gamma_{\infty}} \right) \right) \qquad (22)$$

where $E = RT\Gamma_{\infty}/\sigma_0$ is the surface elasticity parameter, $R$ is the ideal gas constant, $T$ is the absolute temperature, $\sigma_0$ is the surface tension of the clean interface and $\Gamma_{\infty}$ is the maximum packing concentration of surfactants on the interface. The surface elasticity, which relates the surfactant concentration on an interface to its surface tension, has a great influence on the drainage of the liquid film atop of a bubble. In particular, for given flow
<figure><img src="image_12.png" /><figcaption>**Figure 9.** (a) Mesh convergence tests. The film thickness is shown as function of time. In all cases, the conditions are the same; only the refinement is varied.; ethanol volume fraction $B_0 = 2.36$, surfactant concentration in the liquid bulk $c_0 = 0.025$ mol/l, surfactant elasticity $E = RTT_\infty/\sigma_0 = 0.1$. (b) Illustration of the procedure to determine the bubble lifetime. In both cases $D = 2$ mm, $c_0 = 0.025$ mol/l, $E = 0.1$. The Bond number is 0.87 and 1.35 for the 17% and 87% ethanol content tests, respectively.</figcaption></figure>

conditions, the elasticity parameter dictates the strength of the Marangoni stress induced by the surfactant con- centration gradient on the interface and the consequent immobilization of this interface. Typical values of the elasticity parameter for soluble surfactants in aqueous solutions are $E = 0.05–0.25^{55}$.

The numerical solution of these equations is extensively described in Atasi et al.³². The proper implementa- tion of each term, in particular, the surfactant transport on the gas-liquid interface and the computation of the resulting Marangoni stress, has been verified by adapted validation cases and a free rising bubble situation was compared with results from the literature.

**Mesh.** Simulations were performed with a non-uniform axisymmetric orthogonal mesh characterized by 400 and 200 cells in the vertical and radial direction, respectively. The mesh size was refined in the vicinity of the interface to be able to properly capture both the film drainage and its rupture. A mesh convergence test was conducted by simulating the same case with increasing number of grid cells, as shown in Fig. 9(a). Clearly, the draining behavior is captured correctly for by the mesh size chosen (400 by 200) and higher.

We have observed that the numerical rupture of the film occurs when the film thickness is about $h_{rup}^{(num)} \approx 5\Delta$ where $\Delta$ is the grid size. From the experiments (see above), the film thickness at rupture is $h_{rup} \approx 24$ µm. Due to computing limitations, the minimum grid size we could apply in the numerical simulations was $\Delta = 10$ µm, i.e.
$h_{rup}^{(num)} \approx 50\mu m$. For the case of film drainage between two rigid interfaces, and in the limit $Bo \ll 1$, the thinning at the apex follows $h \propto t^{-1/2}$. Consequently, a decrease by a factor two for the critical film rupture should imply an increase by a factor four of the lifetime. Though this correction factor is certainly a good approximation, it was not applied since the drainage dynamics with surfactants can significantly differ from the one between rigid interfaces, especially at $Bo \approx 1$ and for $Bo \gg 1$.

To reconcile the calculated numerical bubble lifetimes with those measured experimentally, an extrapolation scheme was adopted, shown schematically in Fig. 9(b). Since the evolution of the film thickness showed a power-law dependence, the normalized film thickness was fitted to $h/h_o = at^b$. The numerical lifetime was obtained by calculating the time at which the film reached a value of $h_{rup} = 24\mu m$, using the power-law fit.

**Physical properties for surfactant modeling.** A certain surfactant is considered in the simulation to represent the effect of the interface contamination. The corresponding properties are $\Gamma_\infty/\sigma_0 = 0.1/(RT)$, where R is the ideal gas constant and T is the absolute temperature, $k_a = 3 \text{ m/(mol.s)}$, $k_d = 1 \text{ s}^{-1}$ and $C_\infty = 25 \text{ mol/m}^3$, $D = 10^{-9} \text{ m}^2/\text{s}$ and $D_s = 10^{-14} \text{ m}^2/\text{s}$. The fluid properties correspond to water-alcohol mixtures²².

**Validation.** The numerical results were extensively validated in Atasi et al. ³². One additional validation was conducted for this study, considering the shape of bubbles floating on the surface as the film drains. The shape that the bubble adopts while resting on the free surface changes with the value of the Bond number as shown in Fig. 7. The predictions are in close agreement with recent experiments⁴⁹. More importantly, the simulations do capture the retarded drainage due to surface tension gradients induced by concentration gradients of surfactants along the interfaces, as discussed in the paper.

Received: 22 April 2020; Accepted: 2 June 2020;

Published online: 03 July 2020

**References**

1. Patrick, C. H. *Alcohol, culture, and society* (AMS Press, New York, 1970).

2. Gately, I. *Drink: A Cultural History of Alcohol* (Gotham Books, New York, 2008).

3. Cedeno, M. Tequila production. *Crit. Rev. Biotech.* **15**, 1–11 (1995).

4. Colunga-GarciaMarin, P. & Zizumbo-Villarreal, D. Tequila and other agave spirits from west-central mexico: current germplasm diversity, conservation and origin. *Biodivers. Conserv.* **16**, 1653–1667 (2007).

5. Goguitchaivchili, V. *et al*. Archaeomagnetic evidence of pre-hispanic origin of mezcal. *J. Archaeol. Sci. Rep.* **21**, 504–511 (2018).

6. Gross, T. Mezcal and mexicanness: the symbolic and social connotations of drinking in Oaxaca. *Folcklore* **59**, 7–28 (2014).

7. NOM-070-SCFI-1994. Bebidas alcoholicas - Mezcal - Especificaciones. Standard, Norma Oficial Mexicana, Mexico (1994).

8. Lopez-Almendros, L., Angeles-Carreno, G., Colunga-GarciaMarin, P., Cuevas-Sagardi, M. & Franco-Martinez, I. S. Mezcal, arte tradicional. In *Revista Artes de Mexico*, vol. **98** (Fondo de Cultura Económica 2010).

9. Bowen, S. *Divided Spirits: Tequila, Mezcal, and the Politics of Production* (University of California Press, Oakland California 2015).

10. Davidson, J. A. Foam stability as an historic measure of the alcohol concentration in distilled alcoholic beverages. *J. Colloid Interface Sci.* **81**, 540–542 (1981).

11. Ahmed, M. & Dickinson, E. Effect of ethanol content on foaming of alcoholic beverages. *Food Hydrocoll.* **4**, 77–82 (1990).

12. Tuinier, R., Bisperink, C. G. J., den Berg, C. V. & Prins, A. Transient foaming behavior of aqueous alcohol solutions as related to their dilational surface properties. *J. Coll. Interface Sci.* **179**, 327–334 (1996).

13. Blanchard, D. C. & Woodcock, A. H. Bubble formation and modification in the sea and its meteorological significance. *Tellus* **9**, 145–158 (1957).

14. Garrett, W. Stabilization of air bubbles at the air-sea interface by surface-active material. *Deep-Sea Res.* **14**, 661–672 (1967).

15. Jaroensutasinee, M. & Jaroensutasinee, K. Bubble nest habitat characteristics of wild siamese fighting fish. *J. Fish Biol.* **58**, 1311–1319 (2001).

16. Shean, B. J. & Cilliers, J. J. A review of froth flotation control. *Int. J. Miner. Process.* **100**, 57–71 (2011).

17. Manga, M. & Stone, H. A. Interactions between bubbles in magmas and lavas: Effects on bubble deformation. *Volcanol. Geotherm. Res.* **63**, 267–279 (1994).

18. Lhuissier, H. & Villermaux, E. Bursting bubble aerosols. *J. Fluid Mech.* **696**, 5–44 (2012).

19. Hosoi, A. E. & Bush, J. W. M. Evaporative instabilities in climbing films. *J. Fluid Mech.* **442**, 217–239 (2001).

20. Gordeeva, V. Y. & Lyushnin, A. V. Influence of Marangoni instability on evaporation of the polar liquid film. *Eur. Phys. J. Special Topics* **219**, 45–49 (2013).

21. Edmonstone, B. D. & Matar, O. K. Simultaneous thermal and surfactant-induced marangoni effects in thin liquid films. *J. Colloid Interface Sci.* **274**, 183–199 (2004).

22. Khattab, I. S., Bandarkar, F., Fakhree, M. A. A. & Jouyban, A. Density, viscosity, and surface tension of water + ethanol mixtures from 293 to 323 k. *Kor. J. Chem. Eng.* **29**, 812–817 (2012).

23. Champougny, L., Roche, M., Drenckhan, W. & Rio, E. Life and death of not so “bare” bubbles. *Soft Matter* **12**, 5276–5284 (2016).

24. De Leon-Rodriguez, A., Gonzalez-Hernandez, L., Barba de la Rosa, A. P., Escalante-Minakata, P. & Lopez, M. G. Characterization of volatile compounds of mezcal, an ethnic alcoholic beverage obtained from agave salmiana. *J. Agric. Food Chem.* **54**, 1337–1341 (2006).

25. Blasco, L., Viñas, M. & Villa, T. G. Proteins influencing foam formation in wine and beer: the role of yeast. *Int- Microbiol.* **14**, 61–71 (2011).

26. Vazquez, G., Alvarez, E. & Navaza, J. M. Surface tension of alcohol water+ water from 20 to 50 degree c. *J. Chem. Eng. Data* **40**, 611–614 (1995).

27. Innocenzi, P. *et al*. Evaporation of ethanol and ethanol-water mixtures studied by time-resolved infrared spectroscopy. *J. Phys. Chem. A* **112**, 6512–6516 (2000).

28. Barenblatt, G. *Scaling* (Cambridge University Press, Cambridge, UK 2003).

29. Debrégeas, G., de Gennes, P.-G. & Brochard-Wyart, F. The life and death of ‘bare’ viscous bubbles. *Science* **279**, 1704–1707 (1998).

30. Bhamla, M., Giacomin, C. C. B. & Fuller, G. Influence of interfacial rheology on drainage from curved surfaces. *Soft Matter* **10**, 6917 (2014).

31. Champougny, L., Scheid, B., Restagno, F., Vermant, J. & Rio, E. Life and death of not so ‘bare’ bubbles. *Soft Matter* **11**, 2758–2770 (2015).
32. Atasi, O., Haut, B., Pedrono, A., Scheid, B. & Legendre, D. Influence of soluble surfactants and deformation on the dynamics of centered bubbles in cylindrical microchannels. *Langmuir* **34**, 10048–10062 (2018).

33. Ratulowski, J. & Chang, H.-C. Marangoni effects of trace impurities on the motion of long gas bubbles in capillaries. *J. Fluid Mech.* **210**, 303–328 (1989).

34. Mysels, K. J., Frankel, S. & Shinoda, K. *Soap films: studies of their thinning and a bibliography* (Pergamon Press 1959).

35. Scheid, B., Dorbolo, S., Arriaga, L. R. & Rio, E. Antibubble dynamics: The drainage of an air film with viscous interfaces. *Phys. Rev. Lett.* **109**, 264502 (2012).

36. Vitry, Y., Dorbolo, S., Vermant, B. & Scheid, J. Controlling the lifetime of antibubbles. *Adv. Colloid Interface Sci.* **270**, 73–86 (2019).

37. Kocarkova, H., Rouyer, F. & Pigeonneau, F. Film drainage of viscous liquid on top of bare bubble: influence of the bond number. *Phys. Fluids* **25**, 022105 (2013).

38. Nierstrasz, V. A. & Frens, G. Marginal regeneration in thin vertical liquid films. *J. Colloid Interface Sci.* **207**, 209–217 (1998).

39. Howell, P. D. The draining of a two-dimensional bubble. *J. Eng. Math.* **35**, 251–272 (1999).

40. Bin, A. K. Gas entrainment by plunging liquid jet. *Chem. Eng. Sci.* **48**, 3585–3630 (1993).

41. Ohl, C. D., Oguz, H. N. & Prosperetti, A. Mechanism of air entrainment by a disturbed liquid jet. *Phys. Fluids* **12**, 1710–1714 (2000).

42. Zhu, Y., Oguz, H. N. & Prosperetti, A. On the mechanism of air entrainment by liquid jets at a free surface. *J. Fluid Mech.* **404**, 151–177 (2000).

43. Atasi, O., Legendre, D., Haut, B., Zenit, R. & Scheid, B. Lifetime of surface bubbles in surfactant solutions. *To appear*, Langmuir https://doi.org/10.1021/acs.langmuir.9b03597 (2020).

44. Georgieva, D., Schmitt, V., Leal-Calderon, F. & Langevin, D. Lifetime of surface bubbles in surfactant solutions. *Langmuir* **25**, 5565–5573 (2009).

45. Walls, P. L. L. & Bird, J. C. Enriching particles on a bubble through drainage: Measuring and modeling the concentration of microbial particles in a bubble film at rupture. *Elem. Sci. Anth.* **5**, 34 (2017).

46. Zenit, R. & Rodriguez-Rodriguez, J. The fluid mechanics of bubbly drinks. *Phys. Today* **71**, 44 (2018).

47. Bai, L. *et al*. Formation of antibubbles and multilayer antibubbles. *Colloids Surf. A Physicochem. Eng. Asp.* **509**, 334–340 (2016).

48. Poulain, S. & Bourouiba, L. Disease transmission via drops and bubbles. *Phys. Today* **72**, 70 (2019).

49. Teixeira, M. A. C., Arscott, S., Cox, S. J. & Teixeira, P. I. C. What is the shape of an air bubble on a liquid surface? *Langmuir* **31**, 13708–13717 (2015).

50. Taylor, G. I. The dynamics of thin sheets of fluid iii. desintegration of fluid sheets. *Proc. R. Soc. Lond. A* **253**, 313–321 (1959).

51. Culick, F. E. C. Comments on a ruptured soap film. *J. Appl. Phys.* **31**, 1128 (1960).

52. Abadie, T., Aubin, J. & Legendre, D. On the combined effects of surface tension force calculation and interface advection on spurious currents within volume of fluid and level set frameworks. *J. Comp. Phys.* **297**, 611–636 (2015).

53. Stone, H. A. A simple derivation of the time dependent convective diffusion equation for surfactant transport along a deforming interface. *Phys. Fluids A* **2**, 111–112 (1990).

54. Muradoglu, M. & Tryggvason, G. A front-tracking method for computation of interfacial flows with soluble surfactants. *J. Comp. Phys.* **227**, 2238–2262 (2008).

55. Takagi, S., Ogasawara, T., Fukuta, M. & Matsumoto, Y. Surfactant effect on the bubble motions and bubbly flow structures in a vertical channel. *Fluid Dyn. Res.* **41**, 065003 (2009).

## Acknowledgements

R.Z. acknowledges the support of the Fulbright-Garcia Robles foundation during his sabbatical year at Caltech, where this investigation started. O.A. and B.S. thank the F.R.S.-FNRS for financial support.

## Author contributions

G.R., M.M.W., J.F.H.-S. and R.Z. conducted the experiments and analyzed data; O.A., B.H., B.S. and D.L. conducted numerical simulations and analyzed data; all authors discussed the results; R.Z. wrote the paper; all authors reviewed, commented and approved the manuscript.

## Competing interests

The authors declare no competing interests.

## Additional information

Supplementary information is available for this paper at https://doi.org/10.1038/s41598-020-67286-x.

Correspondence and requests for materials should be addressed to R.Z.

Reprints and permissions information is available at www.nature.com/reprints.

Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_13.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2020