## 1. Introduction

Fluids in nature and technology are often multicomponent, with gradients in concentration when out of equilibrium (Lohse & Zhang Reference Lohse and Zhang2020). The concentration gradients give rise to both density gradients and also solutal Marangoni stresses on an interface, thus leading to the competition between Marangoni convection and gravitational convection driven by density differences. Flows with such a competition are frequently encountered in technological applications. For example, in ink-jet printing, where the ink droplets are often multicomponent (Hoath Reference Hoath2016; Lohse Reference Lohse2022); in freezing emulsions (Ghosh & Coupland Reference Ghosh and Coupland2008; Degner *et al.* Reference Degner, Chung, Schlegel, Hutkins and McClements2014; Deville Reference Deville2017; Dedovets, Monteux & Deville Reference Dedovets, Monteux and Deville2018), where next to concentration gradients, thermal gradients also come into play; and in crystal growth (Chang & Wilcox Reference Chang and Wilcox1976; Schwabe *et al.* Reference Schwabe, Scharmann, Preisser and Oeder1978; Chang, Wilcox & Lefever Reference Chang, Wilcox and Lefever1979; Chun & Wuest Reference Chun and Wuest1979; Schwabe & Scharmann Reference Schwabe and Scharmann1979; Schwabe, Scharmann & Preisser Reference Schwabe, Scharmann and Preisser1982; Preisser, Schwabe & Scharmann Reference Preisser, Schwabe and Scharmann1983; Kamotani, Ostrach & Vargas Reference Kamotani, Ostrach and Vargas1984), where the liquid column is subjected to a temperature gradient. The effect of gravity in these situations has often been ignored based on the perception that on small scales, gravitational effects are negligible owing to the small Bond number, which demonstrates the competition between gravity and capillary force (Nepomnyashchy, Legros & Simanovskii Reference Nepomnyashchy, Legros and Simanovskii2012). However, Edwards *et al.* (Reference Edwards, Atkinson, Cheung, Liang, Fairhurst and Ouali2018), Li *et al.* (Reference Li, Diddens, Lv, Wijshoff, Versluis and Lohse2019*a*) and Diddens, Li & Lohse (Reference Diddens, Li and Lohse2021) have recently discovered that gravitational effects can be important in evaporating binary droplets, in spite of the small scale reflected in small Bond numbers. It is therefore important to thoroughly understand the competing effect of gravity and Marangoni forces as some other flows where such a competition occurs may also need reconsideration.

The most basic example in which Marangoni forces compete with gravity is presumably an immiscible drop immersed in a stably stratified liquid, here made of an ethanol–water mixture. Such a simple droplet system allows for a systematic investigation of the parameter space of the concentration gradient, the drop size and the drop viscosity. Without gravitational effects (i.e. without a density stratification), the concentration gradient only leads to simple motion of the immersed drop (Lagzi *et al.* Reference Lagzi, Soh, Wesson, Browne and Grzybowski2010; Maass *et al.* Reference Maass, Krüger, Herminghaus and Bahr2016; Jin, Krüger & Maass Reference Jin, Krüger and Maass2017; Dedovets *et al.* Reference Dedovets, Monteux and Deville2018). However, for a stable stratification, owing to density differences and gravity, intriguing and counter-intuitive phenomena can emerge, such as the oscillatory Marangoni flow around a drop/bubble (Zuev & Kostarev Reference Zuev and Kostarev2006; Viviani, Kostarev & Zuev Reference Viviani, Kostarev and Zuev2008; Schwarzenberger *et al.* Reference Schwarzenberger, Aland, Domnick, Odenbach and Eckert2015; Bratsun *et al.* Reference Bratsun, Kostarev, Mizev, Aland, Mokbel, Schwarzenberger and Eckert2018) or the continuous bouncing of a drop (Blanchette & Shapiro Reference Blanchette and Shapiro2012; Li *et al.* Reference Li, Diddens, Prosperetti, Chong, Zhang and Lohse2019*b*). For low viscosity (5 cSt) drops, it has been found that the continuous bouncing is actually caused by an oscillatory instability (Li *et al.* Reference Li, Diddens, Prosperetti and Lohse2021) of the Marangoni flow, and the onset of the instability is triggered when the Marangoni advection is too strong for diffusion to restore the concentration field around the drop.

In this paper, we study the influence of the drop viscosity on the onset of the oscillatory instability in such a system by thoroughly exploring the three-dimensional parameter space spanned by the linearly stratified density gradient, the drop size and the drop viscosity. We start with very viscous silicone oil drops, for which we discover a new mechanism for the onset of the instability. We then develop a unifying scaling theory that includes both the newly discovered mechanism for high viscosity drops and the mechanism for low viscosity drops suggested by Li *et al.* (Reference Li, Diddens, Prosperetti and Lohse2021) and compare the results of this theory with the experimental results, which reveals good agreement. Finally, the instability thresholds for different viscosities are shown in one single phase diagram spanned by the Marangoni and Rayleigh numbers, thus clearly displaying the influence of the drop viscosity on the onset of the instability in the stratified flow.

## 2. Experimental procedure and methods

Linearly stratified ethanol–water mixtures were prepared in a cubic glass container (Hellma, 704.001-OG, Germany) with inner length $L=30\,\textrm {mm}$. Because the required liquid volume was very small (${<}27\,\textrm {ml}$), we used a slightly modified double-bucket method (Oster Reference Oster1965) to generate the linear stratification, see figure 1(*a*). Syringe A was filled with a lighter liquid which was rich in ethanol (with ethanol weight fraction $w_{t}$), while syringe B was filled with a heavier liquid which was rich in water. The two syringes were filled to the same level and mounted on one syringe pump (Harvard Apparatus, PHD 2000, USA) to make sure the two liquids were injected at the same rate. The injection rate was no larger than 1 ml min$^{-1}$. During injection, the two liquids were constantly mixed in syringe A by a magnetic stirrer. The rod-like magnet in syringe A had a diameter of 4 mm. The resulting mixture flowed out of syringe A through a tube which led to the bottom of the glass container. A linear gradient was formed during the process, see figure 1(*c*). Strictly speaking, only an equal-mass flux rate could lead to a perfectly linear gradient in weight fraction. The equal-volume flux rate in our case and the volume reduction during the mixing of ethanol and water both led to an imperfect linear gradient. However, these imperfections were found to be negligible, see figure 2(*e*). Especially, these imperfections were negligible when compared to the small sizes of the drops, thus we still call it linear in the following text. Notice that by this method, the linear gradient could be varied continuously. However, in the experiments, we only chose several discrete values of the concentration gradient to span the parameter space.

Before injection, the container was first filled with a liquid layer of uniform ethanol concentration $w_{t}$. After injection, a bottom layer of uniform ethanol concentration $w_{b}$ was injected with a third syringe, so that the ethanol weight fraction of the mixture $w_{e}$ increased continuously from $w_{b}$ at the bottom to $w_{t}$ at the top, see figure 1(*b*,*c*). The values of $w_{b}$ and $w_{t}$, and the depth of the linearly stratified layer were varied depending on the degree of stratification, which will be represented by $\mathrm {d}w_{e}/\mathrm {d}y$ in the following text. To avoid bubble formation during mixing, both ethanol (Boom B.V., 100 % (v/v), technical grade, the Netherlands) and Milli-Q water were degassed in a desiccator at ${\sim }$2000 Pa for 20 min before making the mixture. To reduce the disturbance to the linear gradient caused by the preferential evaporation of ethanol and subsequent Marangoni flows on the surface, the container was always covered by a lid except when releasing the drops. During injection, a lid which held the tube in place was used. Otherwise, a glass lid was used to cover the container.

The ethanol weight fraction $w_{e}(y)$ of the mixture at corresponding height $y$ was measured by laser deflection (Lin *et al.* Reference Lin, Leger, Kunkel and McCarthy2013; Li *et al.* Reference Li, Diddens, Prosperetti, Chong, Zhang and Lohse2019*b*) immediately after the stratified mixture was made. The two uniform layers $w_{b}$ and $w_{t}$ were used to increase the accuracy of the laser deflection method. The density of the mixture $\rho (y)$ was calculated from $w_{e}(y)$ using an empirical equation (Khattab *et al.* Reference Khattab, Bandarkar, Fakhree and Jouyban2012), and $y=0$ was set to the position where the density of the mixture $\rho (0)$ equalled that of the silicone oil $\rho ^\prime$, see figure 1(*d*). The ethanol weight fraction at this density-matched position $y=0$ was noted as $w_{e}^\prime$. As reported previously (Li *et al.* Reference Li, Diddens, Prosperetti and Lohse2021), the maximum extent of the flow field induced by the drop became larger for weaker density stratifications, so that the container might not be large enough for weaker density stratifications. However, because we are not interested in the effect of the finite size of the container, experiments for weak density stratifications, i.e. smaller concentration gradients $\mathrm {d}w_{e}/\mathrm {d}y$, were repeated in a larger container (Hellma, 704.003-OG, Germany) with inner length $L=50\,\textrm {mm}$. Results that were different for the two different container sizes were discarded, thus all the results presented in this paper were not influenced by the finite size of the container.

Silicone oils (Sigma-Aldrich, Germany) of different viscosities $\nu ^\prime$ were used to generate drops of different radii $R$. They were injected from a $1\,\mathrm {\mu }$l syringe (Hamilton, KH7001) through a thin needle, whose outer-diameter was 0.515 mm. The properties of the silicone oils are listed in table 1. The drops were released from the top layer of the mixture, right after the laser deflection measurement. A collimated light-emitting diode (LED; Thorlabs, MWWHL4) was used to illuminate the drops in the mixture, and a camera (Nikon D850) connected to a long working distance lens system (Thorlabs, MVL12X12Z plus 0.25X lens attachment) was used to capture the side views of the drop. The radii and the trajectories of the drops were extracted from the sideview recordings. All images were recorded at 30 frames per second.

For each concentration gradient, a large drop was first released. After it had bounced three times, the drop was carefully taken out with another thin needle with outer-diameter 0.515 mm, so as to minimize the disturbance to the stratified liquid. Then a second drop slightly smaller in size was released. In this way, we made sure that only one drop existed in the container at one time. This process was repeated until the smallest drop we could generate was released. One stratified mixture was used for no longer than 40 min. Otherwise, a new stratified liquid with the same concentration gradient was generated again to continue the process. The data were analysed after the first sweep in the drop radius $R$. If more data points were needed, the experiment was repeated at the same concentration gradient with the desired drop radii.

The flow field around the drop was revealed by adding tracer particles to the ethanol–water mixture. PSP (polyamide seeding particles, Dantec Dynamics, Denmark) particles of diameter 20$\,\mathrm {\mu }$m were added in both ethanol and water at 4 mg ml$^{-1}$ before making the linearly stratified mixture. These particles followed the flow faithfully, see Supplemental Material available at https://doi.org/10.1017/jfm.2021.983 for details. The interfacial tensions $\sigma$ between the silicone oils and the ethanol–water mixtures were measured by the pendant droplet method on a goniometer (OCA 15Pro, DataPhysics, Germany). Each measurement was repeated 6 times, and all the interfacial tension data are shown in figure 2(*c*). The solid lines are polynomial fits to the data points.

## 3. Results for 100 cSt silicone oil drops

For 100 cSt silicone oil drops, two typical states are observed after the initial sinking phase, that is, levitating or bouncing. In figure 2(*a*), we show the successive snapshots of two 100 cSt silicone oil drops in a mixture with $\mathrm {d} w_{e}/\mathrm {d}y\approx 40\,\textrm {m}^{-1}$ (see also Supplementary Movie 1). While the smaller drop ($R=54\pm 2\,\mathrm {\mu }\textrm {m}$) is levitating at $h\approx 4.8\,\textrm {mm}$, the larger drop ($R=155\pm 2\,\mathrm {\mu }\textrm {m}$) is bouncing at $h\approx 3\,\textrm {mm}$. The value $h=0$ marks the position where the density of the oil ($\rho ^\prime = 966\,\textrm {kg}\,\textrm {m}^{-3}$) equals that of the mixture (at $w_{e}^\prime =21.0\,\%$). The position of the drop $h$ as a function of time $t$ in the same stratified liquid, the ethanol weight fraction $w_{e}$ and the density of the mixture $\rho$ at the corresponding height are respectively shown in figure 2(*d*–*f*). The smallest drop ($R_1\approx 38\,\mathrm {\mu }\textrm {m}$) is levitating at $h\approx 5.9\,\textrm {mm}$. As the drop size increases, it levitates at a lower position until, above a critical radius $R_{cr}$, it starts to bounce instead of levitate. If its size is further increased, the drop bounces at a lower position (but still with $h>0$).

The smaller drop is able to levitate above the density-matched position ($h=0$) because of a stable Marangoni flow counteracting gravity, as revealed by the flow field around a levitating drop in a weaker gradient $\mathrm {d} w_{e}/\mathrm {d}y\approx 10\,\textrm {m}^{-1}$, see the particle trajectories shown in figure 2(*b*) (and also Supplementary Movie 2). The interfacial tension of the drop $\sigma$ decreases with increasing ethanol concentration of the mixture $w_{e}$, as shown in figure 2(*c*). The ethanol-concentration dependence of the interfacial tension of 20 and 50 cSt silicone oils, which will be used later in the paper, are also shown in this figure. This interfacial tension gradient at the surface of the drop pulls liquid downwards, generating a viscous force acting against gravity, which levitates the drop. The lighter liquids that are pulled down then flow sideways because of buoyancy. When the drop becomes large enough, however, the equilibrium becomes oscillatory and the drop starts to bounce between two different levels (see Supplementary Movie 3). Notice that the onset of this oscillatory instability is essentially the instability of the flow field around a drop which is not moving. Once the drop starts to move (bounce in this case), the flow field becomes drastically different because now the translational motion of the drop must also be considered. Consequently, the analysis of the bouncing behaviour is completely different from that of the instability. In this study, we are only interested in the onset of the instability, and the transition from a levitating drop to a bouncing one serves as a good indication.

Using this indication, we now explore the onset of this instability by independently varying $R$ and $\mathrm {d}w_{e}/\mathrm {d}y$. Because we cannot vary the two parameters continuously in the experiments, we use the following criterion to determine whether a drop is bouncing: if the bouncing amplitude $h_{A}$ of the drop is larger than its radius $R$, then the drop is considered to be bouncing (see Supplemental Material for more details). The results are shown in figure 3. The parameter space for $\mathrm {d}w_{e}/\mathrm {d}y\leq 40\,\textrm {m}^{-1}$ could not be explored because it was already contaminated by the finite size of the container. Different from the previously reported results of 5 cSt silicone oil drops (Li *et al.* Reference Li, Diddens, Prosperetti and Lohse2021), where there is a critical gradient $(\mathrm {d} w_{e}/\mathrm {d}y)_{cr}$ above which the Marangoni flow is always unstable, no such critical gradient is observed for 100 cSt silicone oil drops up to $\mathrm {d}w_{e}/\mathrm {d}y\approx 150\,\textrm {m}^{-1}$. On the contrary, for all the tested gradients, there is always a critical radius $R_{cr}$ below which the flow is stable. This implies a different and new mechanism to trigger the Marangoni instability of 100 cSt drops compared with that of the 5 cSt drops. In the next section, we will develop a unifying scaling theory to account for both instability mechanisms, i.e. for low and high drop viscosity $\nu ^\prime$.

## 4. A unifying scaling theory for the onset of the instability

For an immiscible drop immersed in the stably stratified ethanol–water mixture, a downwards Marangoni flow is generated on the surface of the drop, which also distorts the concentration field, see figure 4 for a sketch of the flow field and the concentration field around the drop. The system can be described in axisymmetrical spherical coordinates $(r, \theta )$ with its origin at the centre of the drop.

For the case of infinitely large solute diffusivity and zero density gradient, the analytical solution of the polar velocity $u_\theta$ outside the drop is given by (Young, Goldstein & Block Reference Young, Goldstein and Block1959)

where $\mathrm {d}\sigma /\mathrm {d}w_{e}$ is a material property (see figure 2*c*), $\mathrm {d}w_{e}/\mathrm {d}y$ the undisturbed ethanol gradient of the mixture, $\mu ^\prime$ the viscosity of the drop and $\mu$ the viscosity of the mixture. The Marangoni velocity at the equator of the drop is

In our case, the ethanol diffusivity $D$ and the density gradient both have a finite value, so the influence of advection and buoyancy cannot be neglected, see figure 2(*b*). Especially, the Marangoni advection tends to homogenize the concentration field around the drop (Li *et al.* Reference Li, Diddens, Prosperetti and Lohse2021), effectively decreasing the concentration gradient $\mathrm {d}w_{e}/\mathrm {d}y$ around the drop, thus the Marangoni flow velocity is smaller than that predicted by (4.2). Let $V_{M}$ denote the Marangoni velocity at the equator of the drop in our case, i.e. $V_{M}=u_\theta \vert _{r=R, \theta =90^{\circ }}$, it is reasonable to drop the prefactor and write (also see Appendix A)

The Marangoni flow imposes a viscous stress $\mu \nabla ^2\boldsymbol {u}\vert _{r=R}$ on the liquid close to the drop, which, in the axisymmetrical spherical coordinate, is expanded as

For the liquid close to the equator of the drop, i.e. at $\theta =90^{\circ }$, $u_\theta \vert _{r=R, \theta =90^{\circ }}=V_{M}$. Also notice that $u_r\vert _{r=R}=0$ and $\partial u_\theta /\partial \theta \vert _{\theta =90^{\circ }}=0$. The Marangoni viscous stress then becomes

Let $\delta$ denote the boundary layer thickness of the Marangoni flow outside the drop, see figure 4. Because $\partial V_{M}/\partial r=\partial u_\theta /\partial r\vert _{r=R,\theta =90^{\circ }}\sim V_{M}/\delta$ and $\partial ^2 V_{M}/\partial r^2=\partial u_\theta /\partial r\vert _{r=R,\theta =90^{\circ }}\sim -V_{M}/\delta ^2$, we obtain

Two length scales are present in (4.6), $\delta$ and $R$, and they are both in the denominator. The next question is which of them is the more relevant length scale or, in other words, which of them is smaller. This can be inferred as follows: according to (4.3), an increase in the drop viscosity $\mu ^\prime$ would lead to a decrease in the Marangoni velocity $V_{M}$. Consequently, the Reynolds number $Re=\rho V_{M}R/\mu$ becomes smaller. Because within Prandtl–Blasius–Pohlhausen boundary layer theory (Schlichting & Gersten Reference Schlichting and Gersten2016), in general, the boundary layer thickness $\delta$ is inversely proportional to the square root of the Reynolds number, then $\delta$ is larger. That is to say, an increase in $\mu ^\prime$ would lead to an increase in $\delta$. Consequently, when $\mu ^\prime$ is very small, it is possible that $\delta < R$ so that $\delta$ is the more relevant length scale. However, when $\mu ^\prime$ is very large, it is possible that $\delta >R$ so that $R$ becomes the more relevant length scale. When $\mu ^\prime$ is taking some intermediate value, it is possible that $\delta \sim R$ so that the situation is more complicated and we are entering a transitional case. The two limiting cases and the transitional case are considered subsequently in the following subsections.

### 4.1. The limiting case when $\mu ^\prime$ is very small

When $\mu ^\prime$ is very small, $\delta < R$ so that (4.6) reduces to

This is the downward viscous stress acting on the liquid close to the equator of the drop. This liquid is lighter than its surroundings as it is entrained from the top. The buoyancy force acting on it is $g{\rm \Delta} \rho =g(\rho ^*-\rho )$, where $\rho$ is its density and $\rho ^*$ is the density of the undisturbed liquid in the far field, see figure 4. This liquid is brought down by the Marangoni flow along the surface of the drop, so ${\rm \Delta} \rho \sim -R\,\mathrm {d}\rho /\mathrm {d}y$. Now let us look at the boundary layer which is thin, $\delta < R$, i.e. all the liquid in the boundary layer is close to the drop. The boundary layer experiences a downward viscous force $\mu V_{M}/\delta ^2$ and an upward buoyancy force $-gR\,\mathrm {d}\rho /\mathrm {d}y$. For a levitating drop, the two forces balance, thus

A very small $\mu ^\prime$ also leads to a very large $V_{M}$. This Marangoni advection acts to homogenize the concentration field inside the boundary layer while diffusion acts to restore it, see figure 4. When advection is too strong, the concentration field close to the drop will be smoothened, thus weakening the Marangoni flow itself, then the flow becomes unstable. The advection time scale is $\tau _{a}\sim R/V_{M}$ and the diffusion time scale is $\tau _{d}\sim \delta ^2/D$. The flow will become unstable when advection is faster than diffusion, i.e.

Because diffusion is the limiting factor that causes this instability, we call it diffusion limited. Substituting the definitions of the two time scales into (4.9), we obtain

The left-hand side of (4.10) has the form of a Péclet number, which is referred to as the Marangoni number,

where we have used (4.3) with an equal sign. The instability criterion thus is

Cancelling $\delta$ from (4.8) and (4.12), we obtain the diffusion-limited instability criterion for low viscosity $\mu ^\prime$ cases:

where

is the Rayleigh number for characteristic length $R$ and $s$ is a constant to be determined.

Equation (4.13) can also be expressed in dimensional quantities by substituting the definitions of $Ma$ and $Ra$:

Equation (4.15) actually predicts a critical concentration gradient above which the flow is always unstable.

Equations (4.13) and (4.15) were confirmed by results of 5 cSt silicone oil drops (Li *et al.* Reference Li, Diddens, Prosperetti and Lohse2021) and the instability threshold is measured to be $s=275\pm 10$.

### 4.2. The limiting case when $\mu ^\prime$ is very large

When $\mu ^\prime$ is very large, $R<\delta$ so that (4.6) reduces to

Now the boundary layer is thick ($\delta >R$) and the force balance applied to the liquid close to the drop cannot be applied to the whole boundary layer. However, the liquid closest to the drop, or in other words the most inner layer of the boundary layer, still experiences the very same two forces: the downward viscous force, as shown by (4.16), and the upward buoyancy force $-gR\,\mathrm {d}\rho /\mathrm {d}y$. However, it also experiences a third force: an upward viscous force exerted by the outer layer. This is because the Marangoni velocity decreases radially outwards. Thus, for a stable Marangoni flow, we should have $\mu V_{M}/R^2>-gR\,\mathrm {d}\rho /\mathrm {d}y$. Otherwise, the Marangoni flow will become unstable because the viscous force on the inner layer cannot overcome its buoyancy. The instability criterion thus is

Because viscosity is the limiting factor that causes this instability, we call it viscosity limited. Substituting (4.3) into (4.17) and reorganizing the left-hand side into the form of a Marangoni number, we obtain the viscosity-limited instability criterion

where $c$ is a constant to be determined.

Equation (4.18) can also be expressed in dimensional quantities by substituting the definitions of $Ma$ and $Ra$:

Equation (4.19) actually predicts a critical drop radius $R_{cr}$ above which the flow is always unstable. The fact that the diffusivity $D$ does not enter into this equation further indicates that this instability criterion is not dominated by diffusion. Note that the concentration gradient $\mathrm {d}w_{e}/\mathrm {d}y$ also does not enter into this equation, and the fluid properties $\mu$ and $\mathrm {d}\sigma /\mathrm {d}\rho$ depend on $w_{e}$ – the ethanol weight fraction at the levitation height – thus $R_{cr}$ is a function of $w_{e}$.

Before moving on, we want to stress that (4.17) is a force imbalance that only applies to the liquid closest to the drop and it does not determine the levitation height of the drop. For the levitation height of the drop, please see Appendix A. It is also found that the instability threshold $c$ should increase with $\mu ^\prime$, see Appendix B.

### 4.3. The transitional case when $\mu ^\prime$ is taking some intermediate value

In between the two limiting cases, we may also have a transitional case where $\delta \sim R$ so that the two mechanisms are intertwined. However, $R$ is an independent variable and an increase in $R$ also leads to an increase in the Reynolds number, thus decreasing $\delta$. Consequently, upon increasing the drop radius $R$ from very small to very large, the onset of the instability would change from the viscosity-limited regime to the diffusion-limited regime.

In the next section, we will compare the experimental results of drops of different viscosities ($\nu ^\prime =100$, 50 and 20 cSt) with the different cases mentioned above and measure the corresponding instability threshold $c$. Because it is found that $c$ is different for oils of different viscosities, we use $c_{100}$, $c_{50}$ and $c_{20}$ in the following text to represent the respective instability thresholds for 100, 50 and 20 cSt silicone oil drops. For general purposes, we will still use $c$, without specifying the drop viscosity.

To calculate the Marangoni and Rayleigh numbers, ethanol weight fractions at the positions where the drops levitate $w_{e}(h)$ are used to obtain the density $\rho$, viscosity $\mu$, diffusivity $D$ and the interfacial tension $\sigma$ (see Supplemental Material for the concentration dependence of $\rho$, $\mu$ and $D$). In the following, for bouncing drops, we use values corresponding to their highest positions.

## 5. Comparison with experimental results of different viscosities

The results of 100 cSt silicone oil drops, as shown in figure 3, are replotted in the $Ma$ versus $Ra$ parameter space in figure 5. As a comparison, the instability threshold of the diffusion-limited case (for 5 cSt silicone oil drops) $Ma/Ra^{1/2}\approx 275$ is shown as the red line. In the diffusion-limited case, no drop should be bouncing below the red line. However, while all the data points of 100 cSt drops are located below the red line, some of them still bounce. This indicates that diffusion is not the limiting factor to trigger the instability for the 100 cSt drops. Moreover, the blue dashed line which separates the bouncing drops from the levitating ones has a slope close to one (see (4.18)), which indicates that the instability is in the viscosity-limited regime. To further confirm this, we plot the results of 100 cSt drops in the $Ra/Ma$ versus $Ra$ parameter space in figure 6(*a*). As can be seen, there is indeed a critical value $(Ra/Ma)_{cr}$ above which the Marangoni flow is unstable, and the instability threshold is measured to be $c_{100}\approx 0.0177$ in the range $1\lesssim Ra\lesssim 5$. The critical radius $R_{cr}$ as a function of $w_{e}$ with $c_{100}=0.0177$ is calculated from (4.19) and shown as the blue curve in figure 6(*b*). The data shown in figure 6(*a*) are also replotted in figure 6(*b*) with ethanol weight fraction $w_{e}$ measured at the levitating height as the $x$-axis. As can be seen, the blue curve nicely separates the levitating drops and the bouncing ones. The dashed blue line in the range $w_{e}>77\,\textrm {wt}\%$ ($w_{e}<36\,\textrm {wt}\%$) corresponds to $Ra<1$ ($Ra>5$) and is not confirmed by experiments. We cannot confirm the region for $Ra>5$ because this would require much larger levitating drops for which container wall effects could affect the data. We also cannot confirm the range for $Ra<1$, i.e. $w_{e}>77\,\textrm {wt}\%$, because for concentration gradient $\mathrm {d}w_{e}/\mathrm {d}y$ larger than $150\,\textrm {m}^{-1}$, the drops tend to levitate close to the top layer ($w_{e}>77\,\textrm {wt}\%$), and the concentration gradient in this region is influenced by diffusion because it is in contact with the uniform top layer (see figure 1*b*). For $\mathrm {d}w_{e}/\mathrm {d}y=150\,\textrm {m}^{-1}$, the thickness of this region is ${\rm \Delta} h\approx 1.5\,\textrm {mm}$. The diffusion time scale to change the concentration in this region is $\tau \sim {\rm \Delta} h^2/D\approx 37\,\textrm {min}$. However, it normally takes ${\sim }20\,\textrm {min}$ to make sure that a small drop is indeed levitating. During this time, the concentration in this region would have changed owing to diffusion, which makes the measurement not accurate. For lower concentration gradients, the drops levitate in the middle of the linear region that is farther from the top layer. The concentration field in this region is not influenced by diffusion. Though in figure 6(*a*), only half a decade in the range of $Ra$ is confirmed by the experimental results, the confirmed range in figure 6(*b*) covers more than half of the available range: $21\,\textrm {wt}\%< w_{e}<100\,\textrm {wt}\%$ (because the drop cannot levitate below the density-matched position $w_{e}^\prime =21.0\,\textrm {wt}\%$).

The parameter spaces for 50 and 20 cSt silicone oil drops are explored experimentally in the same way. Their interfacial tensions with the ethanol–water mixture are also measured in the same way, see figure 2(*c*). Their results are plotted in figure 6(*c*,*d*) for 50 cSt oil drops and figure 6(*e*,*f*) for 20 cSt oil drops. The instability threshold for 50 cSt silicone oil is measured to be $c_{50}\approx 0.0095$, and the blue line predicted from (4.19) with $c_{50}=0.0095$ again nicely separates the levitating drops from the bouncing ones. The dashed blue line in the range $w_{e}>86\,\textrm {wt}\%$ ($w_{e}<39\,\textrm {wt}\%$) corresponds to $Ra<1$ ($Ra>2.5$) and is not confirmed by experiments owing to the same reason as above. Again, while only a small range of $Ra$ in figure 6(*c*) is confirmed by the experimental results, the confirmed range in figure 6(*d*) covers more than half of the available range: $25\,\textrm {wt}\%< w_{e}<100\,\textrm {wt}\%$.

For 20 cSt silicone oil drops, no constant instability threshold is found, see figure 6(*e*,*f*), and we are probably entering the transitional case. Still, we cannot carry out experiments for $Ra<0.35$ and $Ra>16$ because of the same reasons as above. At $Ra\approx 0.35$ ($77\,\textrm {wt}\%< w_{e}<92\,\textrm {wt}\%$), the instability threshold is measured to be $c_{20}\approx 0.003$, as shown by the black solid line in figure 6(*e*) and the blue solid line in figure 6(*f*). In the range $Ra>0.35$ ($w_{e}<77\,\textrm {wt}\%$), the instability threshold deviates from $c_{20}\approx 0.003$.

To get a better understanding of this case, results of 20 cSt silicone oil drops are re-plotted in the $Ma/Ra^{1/2}$ versus $Ra$ parameter space in figure 7. The data points confirm that this is indeed the transitional case. The instability starts from the viscosity-limited regime $Ra/Ma=0.003$ when the drops are small ($Ra<0.35$), as shown by the blue dashed line. It slowly curves towards $Ma/Ra^{1/2}=275$ as the drop radius is increased ($0.35< Ra<3$) while still maintaining the viscosity-limited feature, as shown by the blue solid line. It finally turns into the diffusion-limited regime $Ma/Ra^{1/2}=275$ for much larger drops ($Ra>3$), as indicated by the red solid line.

We also performed experiments for 10 cSt silicone oil drops, where the onset of the instability becomes diffusion limited and the instability threshold is measured to be $s=265\pm 15$. This is indistinguishable from that of the 5 cSt drops $s=275\pm 10$ with the current experimental accuracy. We speculate that the instability threshold $s$ for the diffusion-limited case does not change with $\mu ^\prime$, because the balance between advection and diffusion happens in the bulk and the boundary layer thickness $\delta$ in (4.13) has cancelled out. However, the instability threshold $c$ for the viscosity-limited case is indeed found to increase with $\mu ^\prime$ (from $c_{20}\approx 0.003$ to $c_{100}\approx 0.0177$).

## 6. Instability thresholds for different viscosities

Now with all the results presented, we can have an overview of the instability thresholds for oils of different viscosities. They are plotted in figure 8. When the oil viscosity is very low, the instability is diffusion limited, which yields $Ma/Ra^{1/2}=s$. This is the case for 5 and 10 cSt silicone oil drops, and their instability thresholds are both $s\approx 275$. When the oil viscosity is very high, the instability is viscosity limited, which yields $Ra/Ma=c$. This is the case for 100 and 50 cSt silicone oil drops, and their instability thresholds are measured to be $c_{100}\approx 0.0177$ and $c_{50}\approx 0.0095$. When the oil viscosity is taking an intermediate value, the instability changes from the viscosity-limited regime to the diffusion-limited regime as the drop radius increases. This is the case for the 20 cSt silicone oil drops, whose instability thresholds are measured to be $c_{20}\approx 0.003$ and $s\approx 275$.

## 7. Conclusions and outlook

In summary, the oscillatory instability of an immiscible oil drop immersed in a stably stratified ethanol–water mixture has been explored for oils of different viscosities. A unifying scaling theory has been developed, which predicts two different instability mechanisms depending on the two length scales in the system: the kinematic boundary layer thickness $\delta$ induced by the Marangoni flow and the drop radius $R$. (i) When $\delta < R$, the instability is triggered when advection is too strong so that diffusion cannot restore the concentration field around the drop. This is the diffusion-limited regime which yields $Ma/Ra^{1/2}=s$. (ii) When $R<\delta$, the instability is triggered when the gravitational effect is too strong so that the viscous stress cannot maintain a stable Marangoni flow. This is the viscosity-limited regime which yields $Ra/Ma=c$. An increase in the drop viscosity leads to an increase in $\delta$, which eventually moves the instability criterion from the pure diffusion-limited regime to the pure viscosity-limited regime, through a transitional case which connects both regimes. Because $\delta$ cancels out in the diffusion-limited regime, the instability threshold $s$ does not change with an increase in the drop viscosity $\mu ^\prime$, but the instability threshold $c$ in the viscosity-limited regime increases with increasing $\mu ^\prime$. The scaling theory is well supported by the experimental results.

Owing to the limitations of the scaling theory, we do not know how $\delta$ quantitatively depends on the drop viscosity $\mu ^\prime$ and drop radius $R$. Thus, we cannot predict beforehand in which regime will the instability be triggered. A quantitative understanding of $\delta$ may help to solve this problem. More detailed theoretical analysis may also be able to predict the instability thresholds beforehand. Nevertheless, our findings provide not only a deeper understanding to the simple droplet system under investigation, but also insight to other systems which allow for the direct competition between Marangoni and Rayleigh convection. For example, $Ra/Ma$ is called the dynamic Bond number in the context of Marangoni convections (Nepomnyashchy *et al.* Reference Nepomnyashchy, Legros and Simanovskii2012) $Bo_{d}=Ra/Ma=({(\mu +\mu ^\prime )}/{\mu })({\mathrm {d}\rho }/{\mathrm {d}\sigma })gR^2$. This number expresses the direct competition between Marangoni stress and gravity, as in the viscosity-limited regime. While the static Bond number $Bo_{s}={\rm \Delta} \rho gR^2/\sigma$ determines whether a drop/bubble maintains a spherical shape, the dynamic Bond number $Bo_{d}$ is more relevant for flows where the Marangoni flow is competing with a density gradient, such as in an evaporating binary sessile droplet (Edwards *et al.* Reference Edwards, Atkinson, Cheung, Liang, Fairhurst and Ouali2018; Li *et al.* Reference Li, Diddens, Lv, Wijshoff, Versluis and Lohse2019*a*; Diddens *et al.* Reference Diddens, Li and Lohse2021).

## Supplementary material and movies

Supplementary material and movies are available at https://doi.org/10.1017/jfm.2021.983.

## Acknowledgements

We thank V. Sanjay for valuable discussions. We acknowledge support from the Netherlands Center for Multiscale Catalytic Energy Conversion (MCEC), an NWO Gravitation programme funded by the Ministry of Education, Culture and Science of the government of the Netherlands, and D.L. is grateful for the ERC-Advanced Grant under project number 740479.

## Declaration of interests

The authors report no conflict of interest.

## Appendix A. The levitation heights $h$ of the drops

For a levitating drop, the viscous force imposed on it by the Marangoni flow is $F\sim \mu V_{M}R$ and the gravity of the drop is $F_{g}\sim gR^3(\rho ^*-\rho ^\prime )$. In a linear density gradient, $\rho ^*-\rho ^\prime =-h\cdot \mathrm {d}\rho /\mathrm {d}y$, where $h$ is the levitation height of the drop. Thus $F_{g}\sim -gR^3h\cdot \mathrm {d}\rho /\mathrm {d}y$. From the force balance on the drop $F\sim F_{g}$, we obtain

Substituting the definition of $V_{M}$ into (A1), we can obtain the levitation height of the drop in a linear density gradient:

For the case of infinitely large diffusivity and zero density gradient, the prefactor on the right-hand side is found to be $3/2$ (Young *et al.* Reference Young, Goldstein and Block1959). However in our case, as mentioned earlier, the Marangoni advection acts to homogenize the concentration field around the drop, which decreases the upwards viscous force $F$. Thus we expect the prefactor to be smaller than $3/2$. We write this as

with a prefactor $0<\alpha <1$, which reflects the strength of advection. A stronger advection will result in a smaller $\alpha$. Equation (A3) is confirmed by the levitation heights measured by drops of different sizes and viscosities, see figure 9. This also confirms that (4.3) is the correct form to estimate $V_{M}$. Finally, the measured $\alpha$ indeed becomes smaller for lower viscosities.

## Appendix B. The influence of the drop viscosity $\mu ^\prime$ on the viscosity-limited instability threshold $c$

The fact that $\alpha$ is of the order of one indicates that the two sides of (A1) are of the same order. To stress this, we rewrite (A1) as

At the onset of the viscosity-limited instability, see equation (4.17), we can write

where $k$ is a prefactor to be determined. Combining (B1) and (B2), one gets

Because $h_{cr}>3\,\textrm {mm}$ (see figure 9) and $R_{cr}<0.1\,\textrm {mm}$ (see figure 3), we have $k\gg 1$. Comparing equation (B2) and equations (4.17) and (4.18), one realizes that the instability threshold $c=1/k$, thus $c\ll 1$. This is consistent with the experimental results.

If the drop viscosity $\mu ^\prime$ is increased while all other quantities remain the same, the levitation height $h$ will decrease according to (A3). Then $k$ decreases according to (B3), thus $c$ increases. This explains why the instability threshold $c$ increases with increasing $\mu ^\prime$.