Interaction of cavitation bubbles with the interface of two immiscible fluids on multiple time scales

Abstract We experimentally, numerically and theoretically investigate the nonlinear interaction between a cavitation bubble and the interface of two immiscible fluids (oil and water) on multiple time scales. The underwater electric discharge method is utilized to generate a cavitation bubble near or at the interface. Both the bubble dynamics on a short time scale and the interface evolution on a much longer time scale are recorded via high-speed photography. Two mechanisms are found to contribute to the fluid mixing in our system. First, when a bubble is initiated in the oil phase or at the interface, an inertia-dominated high-speed liquid jet generated from the collapsing bubble penetrates the water–oil interface, and consequently transports fine oil droplets into the water. The critical standoff parameter for jet penetration is found to be highly dependent on the density ratio of the two fluids. Furthermore, the pinch-off of an interface jet produced long after the bubble dynamics stage is reckoned as the second mechanism, carrying water droplets into the oil bulk. The dependence of the bubble jetting behaviours and interface jet dynamics on the governing parameters is systematically studied via experiments and boundary integral simulations. Particularly, we quantitatively demonstrate the respective roles of surface tension and viscosity in interface jet dynamics. As for a bubble initiated at the interface, an extended Rayleigh–Plesset model is proposed that well predicts the asymmetric dynamics of the bubble, which accounts for a faster contraction of the bubble top and a downward liquid jet.


Introduction
Bubble dynamics is a typical multiphase flow problem, which has received much attention for many years due to its broad applications and interesting behaviours (Prosperetti 2004;Lauterborn & Kurz 2010;. In most realistic circumstances, bubble dynamic behaviours are inevitably affected by boundary conditions of the flow field or external force fields such as gravity and acoustic waves (Lauterborn & Kurz 2010;Supponen et al. 2016). Up to now many experimental, numerical and theoretical studies have been carried out of bubble dynamics near different boundaries, the most commonly seen of which are rigid walls (Blake, Taib & Doherty 1986;Zhang, Duncan & Chahine 1993;Hsiao et al. 2014;Wang 2014;Beig, Aboulhasanzadeh & Johnsen 2018), free surfaces (Chahine 1977;Wang et al. 1996;Quah et al. 2018;Kang & Cho 2019), elastic membranes (Sankin, Yuan & Zhong 2010), suspended structures (Goh et al. 2017;Wu et al. 2017;Li et al. 2019a), adjacent bubbles (Tomita, Shima & Ohno 1984;Bremond et al. 2006;Ochiai & Ishimoto 2017), etc. A significant diversity of bubble collapse patterns and jetting behaviours has been revealed. For instance, a high-speed liquid jet, an important destructive mechanism in cavitation and underwater explosions, is directed towards a rigid wall or away from a free surface (Blake & Gibson 1987;Philipp & Lauterborn 1998;Kim et al. 2014;Kiyama et al. 2021). We also benefit from bubble jets in some other applications, such as ultrasonic cleaning (Chahine et al. 2016;Reuter & Mettin 2016), sonoporation (Ohl et al. 2006b;Kooiman et al. 2011), printing (Turkoz et al. 2018), etc. In the majority of published literature, the flow field surrounding bubbles merely consists of a single type of fluid. The interaction between an oscillating bubble and the interface of two immiscible fluids is far from well understood, which has applications in ultrasonic emulsification (Canselier et al. 2002), pharmacy (Freitas et al. 2006), the food industry (Nishinari et al. 2014), sediment transport and dredging (Nielsen, Bach & Bollwerk 2015), ocean engineering (Xu et al. 2020), etc.
There have been a few experimental observations of bubble dynamics near a fluid-fluid interface. Chahine & Bovis (1980) performed experiments for spark-generated bubbles near the interface of two immiscible liquids. They discussed the dependence of the bubble jet direction on the standoff parameter and the Froude number since their motivation was cavitation damage reduction. Thereafter, this problem received little attention until very recently. Perdih, Zupanc & Dular (2019) found rich phenomena occurred near the liquid-liquid interface which was exposed to ultrasonic cavitation, including cavitation bubble oscillation near the interface, penetration of a water jet into the bulk oil phase, breakup of oil droplets, etc. Yamamoto, Matsutaka & Komarov (2021) experimentally studied the dynamics of acoustic cavitation bubbles near a gallium droplet interface. They partially demonstrated that a high-speed liquid jet from the bubble is the prime cause of liquid emulsification and a large amplitude of bubble oscillation is required to trigger a liquid jet. However, due to the limitations of spatio-temporal resolutions of such small-scale experiments (∼100 μm), the microscopic phenomena during the bubble-interface interaction are difficult to observe clearly. To better understand the fundamental ultrasonic emulsification process, Orthaber et al. (2020) studied the jetting behaviour of a laser-induced cavitation bubble (∼1 mm) near a liquid-liquid interface. They found that the direction of the bubble jet is always from the lighter liquid to the denser liquid. The dependence of the bubble dynamics on an anisotropy parameter ) was also discussed. A similar bubble jet behaviour was also reported in the work of Yin et al. (2020). In the present work, the electric discharge method (Fong et al. 2009;Cui et al. 2018) is used to generate centimetre-scale cavitation bubbles, which allows us to achieve a higher spatio-temporal resolution of bubble dynamics and interface evolution than in earlier works. Additionally, three types of oils with different densities and viscosities are used in our experiments, aiming to provide new physical insights for bubble dynamics near a fluid-fluid interface.
There are also a few numerical studies of the interaction between a cavitation bubble and a fluid-fluid interface. Klaseboer & Khoo (2004b) established a boundary integral (BI) model for bubble-interface interactions based on the potential flow theory. The effects of the density ratio α between two fluids were studied therein, including two limiting situations (α → 0 and α = ∞). The interface acts like a free surface and a rigid wall in these two situations, respectively. Curtiss et al. (2013) extended this model to study the interaction between a single ultrasound contrast agent type of bubble and a tissue layer. They found the inertial bubble provides an efficient way for removing polluted material layers by literally lifting them off an attached substrate. Rowlatt & Lind (2017) adopted the spectral element marker particle method to study the bubble collapse near a fluid-fluid interface with applications in bioengineering. Liu et al. (2019) proposed a volume of fluid model implemented in the finite element method to study the behaviour of a bubble generated at the interface of two different liquids. They revealed the role of gravity and density ratio of the two fluids in bubble migration and jet direction. The same method was also applied to study the bubble-seabed interaction in shallow water (Xu et al. 2020). Yamamoto & Komarov (2020) numerically studied the jet dynamics of acoustic cavitation bubbles near a gallium droplet using a commercial software. They found that the jet velocity of the bubble is maximized at a moderate initial bubble-interface distance.
Most of the aforementioned experimental and numerical studies were restricted to bubble dynamics on a short time scale; however, we have little knowledge of the residual flow after bubble collapse, i.e. the dynamics of the interface jet on a much longer time scale. To fill this knowledge gap, both the bubble dynamics and the interface jet evolution on multiple time scales are investigated in this study. Remarkably, besides the penetration of a high-speed bubble jet into the fluid-fluid interface, the pinch-off of an interface jet also leads to the mixing of fluids. The dependence of the associated fluid dynamics on the governing parameters is systematically studied via experiments, BI simulations and a scaling analysis. This paper is structured as follows. In § 2, we present our experimental set-up and numerical model. In § 3, the general physical phenomena are discussed for bubble initiations in different fluids or at the interface. In § 4, the experimental observations are compared with BI simulations or theoretical results from an extended Rayleigh-Plesset model. In § 5, a quantitative study is presented for bubble jet dynamics. In § 6, the dependence of the interface jet dynamics on governing parameters is discussed. Finally, this work is summarized and conclusions are drawn in § 7.

Experimental set-up
An underwater electric discharge method Cui et al. 2018) is adopted to generate cavitation bubbles in our experiments. Two copper-alloy wires with a diameter of ∼0.2 mm cross and touch at a point, which is the initial centre of the bubble. At first, a capacitor is charged to 500 V. Upon discharge, strong Joule heating at the crossing point vaporizes the surrounding water or oil, and a centimetre-scale bubble is thus generated. The maximum radius of the bubble is around 15 mm. Details of the bubble generator can be found in our previous studies (Cui et al. 2018 of local heating and boiling, it is still called a 'cavitation bubble' in the community of bubble dynamics since the basic mechanics of cavitation and boiling is similar and the main content is vapour (Kling & Hammitt 1972;Blake & Gibson 1987;Brennen 1995). Though cold boiled water and oils are used in the experiment, we still find that a small amount of non-condensable gas slowly enters the bubble due to diffusion. Experiments are performed in a water tank of 300 mm × 300 mm × 600 mm under atmospheric pressure and at room temperature (∼20 • C). The tank is firstly filled with water to 300 mm in depth and then the oil is added gently and slowly to 280 mm in depth. An oil-water interface is thus formed between the two immiscible fluids. To gain better insight into the interaction between bubbles and the oil-water interface, we use three types of oils in the experiments and their physical properties are given in table 1. Besides, the density and viscosity of water in our experiment are 0.999 g cm −3 and 0.001 Pa s, respectively. The initial centre of a cavitation bubble is arranged near or at the interface and the strong interaction is experimentally studied in a systematic manner. A continuous light source provides illumination from the back. A high-speed camera (Phantom V711) is triggered at the same time with a discharge switch. Both the transient bubble behaviours and the interface evolutions on a longer time scale are well captured by the camera working at 16 000-21 000 frames per second with an exposure time of 10 μs. The temporal resolution (47.6-62.5 μs) is within 2 % of the first period of bubble oscillation and 0.1 % of the interface jet evolution process. The uncertainty of the length measurement can be estimated as one pixel of the image (0.1 to 0.2 mm), which is around 1 % of the maximum bubble radius.
The bubble-interface interactions in our experiments are inertia-controlled on a small time scale and the viscosity effects are negligible (see further explanation in § 3.1). For micrometre-sized ultrasonic bubbles in some practical applications, the key physical process is usually accompanied by an energetic collapse and a high-speed liquid jet. For instance, a large amplitude of bubble oscillation is an essential condition for emulsification (Yamamoto et al. 2021); the jet impact plays a vital role in removing unwanted material layers on a surface (Ohl et al. 2006a;Curtiss et al. 2013). The associated Reynolds numbers are still much larger than 1. Therefore, it is convincing that our experimental data can shed light on the behaviour of micrometre-sized cavitation bubbles. For a much longer time scale, the interface evolution in our system is determined by the competing effects of inertia, surface tension, gravity and viscosity, which may also provide insights into the essential physics of problems at the micrometre scale.

Numerical model
A sketch of the physical problem is shown in figure 1. A relatively heavy fluid (fluid 1, water) and a second fluid (fluid 2, oil) are separated by a sharp fluid-fluid interface. An axisymmetric BI method is used to simulate the interaction between a cavitation bubble and the fluid-fluid interface. We define a cylindrical coordinate O-rθz with the origin O located at the fluid-fluid interface, vertically above or below an initially spherical bubble. The distance between the initial bubble centre (0, 0, z b ) and the interface is denoted by d b .
In the following introduction of the model, we suppose the bubble is initiated in fluid 1. One can easily extend the model to the situation of bubble initiation in fluid 2. According to the work of Klaseboer & Khoo (2004a), Gordillo et al. (2005) and RodríGuez-RodríGuez, Gordillo & Martínez-Bazán (2006), fluid 1 and fluid 2 are considered as inviscid and incompressible. The Laplace equation is valid in both fluids and the velocity potentials ϕ 1 and ϕ 2 satisfy the BI equation, written as where c denotes the solid angle, r and q the control and source points, respectively, and ∂/∂n the normal derivative. Here S refers to the fluid-fluid interface and the bubble surface when i = 1 (flow domain 1), while it refers to the fluid-fluid interface only when i = 2 (flow domain 2). The dynamic boundary conditions on the bubble surface and at the fluid-fluid interface can be written as where P ∞ is the hydrostatic pressure at z = 0, P L the liquid pressure on the bubble surface, P 1 and P 2 the pressures just below and above the interface, respectively, ρ the density of the fluid and g the gravitational acceleration. With the material derivative D/Dt = ∂/∂t + ∇ϕ 1 · ∇, (2.3), (2.4) and (2.5) transform into Taking the surface tension into account, the relation between the internal pressure of the bubble P b and the liquid pressure P L on the bubble surface satisfies where the subscript '0' represents initial quantities and V is the bubble volume, λ the ratio of the specific heats, σ the surface tension coefficient and κ the local curvature. For simplicity, here we use the adiabatic approximation to model the gas pressure inside the bubble, as suggested by , Lee, Klaseboer & Khoo (2007) and others. For the bubble dynamics in our experiments, the associated Péclet number (defined as R 2 m /T osi D, where R m is the maximum bubble radius, T osi the bubble period and D the thermal diffusivity) can be estimated as O(10 3 ), which justifies the adiabatic assumption.
Substituting (2.9) into (2.6), the dynamic boundary condition on the bubble surface is thus obtained as follows: (2.10) At the water-oil interface, the condition on the normal stresses is given by where μ 2 is the oil viscosity and the last term denotes the normal viscous stress. Here we ignore the water viscosity since its value is much smaller than that of the oil. The theory of viscous potential flow works well when vorticity is restricted to a thin layer near the boundary (Joseph & Wang 2004;Klaseboer et al. 2011). Some justification is given in § 4.1. From (2.7), (2.8) and (2.11), the dynamic boundary condition at the interface is obtained as follows: 12) where α = ρ 2 /ρ 1 is defined as the density ratio.
The kinematic boundary condition on all surfaces is given by If a cavitation bubble is initiated in oil (fluid 2), (2.12) can be deduced in the same manner with the material velocity u 2 = ∇ϕ 2 .
Although the above derivation process is similar to that in Klaseboer & Khoo (2004a,b), the interaction between a toroidal bubble (after the jet impact) and the interface was not considered with therein. In this study, the interaction between a toroidal bubble and the interface is numerically investigated using a vortex ring model (Wang et al. 1996). The velocity potential ϕ is decomposed into two parts, i.e. the potential due to the circulation of the vortex ring ϕ vor , which is obtained using a semi-analytical method (Zhang, Li & Cui 2015), and the remnant potential ϕ res . The induced velocity of the vortex ring is calculated using the Biot-Savart law, while the remnant velocity is calculated from the BI equation (2.2). For details of the vortex ring model for simulating the toroidal bubble motion, the reader is referred to the work of Wang et al. (1996) and Zhang et al. (2015).

Non-dimensionalization and initialization
In the present study, numerical calculations are performed in a dimensionless system. The equivalent maximum radius of the bubble R m , the hydrostatic pressure at the initial interface P ∞ and the density of the heavier fluid ρ 1 are used as three basic quantities to convert other parameters into dimensionless quantities. Four dimensionless variables are introduced as follows: where γ is defined as the standoff parameter, α the density ratio, ε the strength parameter and δ the buoyancy parameter (equivalent to the inverse Froude number). For convenience, we denote γ w and γ o as the standoff parameter for bubble initiations in water (denser fluid) and oil (lighter fluid), respectively. To match the experiment, the initial conditions should be set properly in numerical simulations. In the experiment, the discharge process lasts for about 0.3 ms (about one-tenth of the first cycle of the bubble); thus the electric energy is not transferred to the bubble immediately. Up until now, to the best knowledge of the authors, no successful attempt has been made to model the early stage of an electric discharge bubble. Following published literature (Tong et al. 1999;Klaseboer et al. 2005;Hsiao et al. 2013;Zeng et al. 2018), we use a simplified method to initialize the bubble in simulations, namely the bubble is set as an initially stationary high-pressure gas bubble. The Rayleigh-Plesset equation (Plesset 1949) in a non-dimensional form, is used to fit the free-field experiment by adjusting the initial bubble radius R 0 and the strength parameter ε. Since the dimensionless maximum bubble radius is 1, the relationship between R 0 and ε  can be derived from the energy conservation law: Since the content of the bubble is vapour and a small amount of diffused air, the ratio of the specific heats λ should be less than 1.4. However, it is difficult to measure the exact value. Fortunately, the bubble dynamics is not sensitive to the choice of λ; for example, the bubble oscillation period T osi only varies ∼3 % when λ ranges from 1.2 to 1.4 (not shown here). Hence, we set λ = 1.25 in this study according to Lee et al. ε to fit the experimental data and the initial bubble radius is calculated from (2.16).
The experimental data of the bubble radius evolution in a free field are compared with theoretical predications from (2.15) with different ε, as shown in figure 2. Four experiments are presented, including bubble initiations in sunflower oil and water. In the simulations, four different ε are chosen, namely 50, 100, 200 and 400. Discrepancies between the experimental and theoretical results are noted in the expansion phase, which are mainly due to the fact that the bubble energy is gradually increased during the discharge process in the experiments. Nevertheless, after the bubble reaches the maximum volume, the experimental data start to follow the theoretical results when ε is set between 100 and 200. The ε = 50 simulation overestimates T osi while the ε = 400 simulation underestimates T osi . The inset also illustrates the dependence of T osi on ε. A satisfactory result can be obtained if ε is chosen within the red dotted lines. Therefore, the initial non-dimensional pressure and radius of the bubble are set as ε = 100 and R 0 = 0.1485 in this study. Dozens of experiments have been performed in a free field and the difference in T osi between experimental and theoretical results is within 2.3 %, which also proves the good reproducibility of the present experimental set-up. Finally, we emphasize that the time evolutions of the scaled bubble radius nearly collapse together, indicating that the viscosity of the oil plays a minor role in bubble oscillations of our present system.

General physical phenomena
Firstly, we present and discuss several representative experimental results for bubbles initiated in water, in oil and at the interface. The dependence of the overall physical phenomena on the standoff parameter is qualitatively studied. In this section, only the sunflower oil is used.

Bubble initiation in water
In figure 3, four representative experiments are shown with a decreasing standoff parameter γ w . This allows us to anticipate an increase in the interaction between the bubble and the water-oil interface. In the first experiment (see figure 3a), the standoff parameter is γ w = 1.32 and the corresponding bubble-interface interaction is relatively weak. The bubble expands rapidly after inception (frames 1-3) and the interface elevates slightly with bubble expansion. The maximum amplitude of the interface motion is reached when the bubble attains the maximum volume (frame 3). The interface descends with the contraction of the bubble (frame 4) and almost recovers its initial shape when the bubble attains its minimum volume (frame 5). We note that the bubble keeps a spherical shape during the first oscillation cycle, indicating that the bubble is little influenced by the presence of the upper interface. In frame 6, the bubble rebounds to its maximum volume during the second cycle and we can hardly observe perturbations at the water-oil interface. The remaining energy of the bubble in the second cycle is only 16 % of that of the first cycle, which is estimated from the maximum bubble volume (Lee et al. 2007). Thereafter, the bubble oscillates for several cycles with a damping amplitude, accompanied by a downward migration of the bubble centroid (frame 7). On a much longer time scale, the interface rises very slowly (frames 8-9) and the maximum height of the interface is reached at t = 101.92 (frame 10). When γ w 1.32, this residual flow at the interface cannot be observed in our experiments. It is noted that the cavitation bubble eventually turns into some non-condensable gas bubbles due to diffusion (Moreno Soto et al. 2018). The average volume of the non-condensable gas is about 0.33 ml. Subsequently, these rising gas bubbles pass through the water-oil interface, which is beyond the scope of this study.
In the second experiment (see figure 3b), γ w is decreased to 0.91. The first cycle of the bubble is not given here (the reader is referred to Appendix B). Frame 1 shows the rebounding bubble in the early second cycle. A thin liquid jet can be seen inside the toroidal bubble, which implies that the bubble-interface interaction becomes stronger compared to the previous case. Additionally, the downward migration of the bubble is faster and the interface shows a simple smooth hump in the later stage (frames 3-5). It takes longer for the interface to reach its maximum height (t = 117.78, frame 5). In the third experiment (see figure 3c), γ w = 0.58. A downward liquid jet forms around the moment of the minimum bubble volume (not shown here). An annular neck can be seen on the toroidal bubble surface during the rebound phase (frame 1). After the bubble migrates away from the interface, the interface rises quickly (frame 2) and the interface jet grows much higher than in the previous two cases (frames 3-5). The maximum dimensionless height of the interface jet h m is around 1.32. In the fourth experiment, γ w is further reduced to 0.42. The toroidal bubble is more elongated along the axis of symmetry during the rebound phase (frame 1), indicating that the downward liquid jet of the bubble is more energetic. On a longer time scale, the interface jet is no longer a smooth hump; instead it shows an annular neck above the half-height position (frames 2-4). However, this neck disappears in the later stage (frame 5).
To understand the governing factors that influence the bubble dynamics and interface jet evolution, we first look at the associated dimensionless numbers. For the bubble dynamics on a short time scale (∼T osc ), the associated Reynolds number and Weber number can be defined as where the characteristic velocity is taken as U = √ P ∞ /ρ ≈ 10 m s −1 . Here we use the water viscosity μ 1 = 10 −3 Pa s and surface tension σ = 0.073 N m in (3.1a,b). This illustrates that the liquid viscosity and surface tension play a minor role during the  bubble-interface interaction on a short time scale. This statement will be further confirmed by our numerical simulation. However, for the long-time evolution of the interface jet, the characteristic velocity should be replaced by a new one, namely the average rising velocity of the interface jetŪ = h m /T m , where T m is the time for the interface to reach its maximum height. In the four experiments discussed above, 0.01 m s −1 <Ū < 0.14 m s −1 .
The characteristic length remains the same since the interface jet has a size comparable with that of the bubble. The oil viscosity μ 2 = 62 × 10 −3 Pa s is used. Thus we can obtain another Reynolds number and Weber number: This implies that both the viscosity and surface tension come into play during the long-time evolution of the interface jet. Additionally, we use the Bond number to estimate the ratio of buoyancy to capillarity. In the four experiments discussed above, Bo ≈ 7; thus both gravity and surface tension are important in the growth and descent of the interface jet.
In the experiments discussed above, we find that the interface evolution is closely related to the bubble motion during the first cycle of the bubble. More specifically, the interface is pushed away by the expanding bubble and attracted by the collapsing bubble, which is different from that in bubble interactions with a water-air interface (Blake & Gibson 1987;Koukouvinis et al. 2016;Kang & Cho 2019). The inertia of the air is much smaller than that of water and the pressure at the interface nearly remains at atmospheric pressure. For small standoff parameters, the water-air interface is continuously pushed upward during the whole bubble life, generating a pronounced interface spike (Chahine 1977;Wang et al. 1996;Koukouvinis et al. 2016;Kang & Cho 2019). Additionally, bubble bursting would occur if the bubble is very close to the water-air interface (i.e. γ f 0.5) (Li et al. 2019b). In the present system, however, despite the small dimensionless standoff parameter in the fourth case (γ w = 0.42), there is always a thin film between the bubble and the interface, and different interface jet dynamics are thus obtained.
Nevertheless, previous work on bubble interactions with a water-air interface stimulates us to give a similar mechanical explanation of the interface jet in the context of the conservation of linear momentum via a quantity of Kelvin impulse (Blake & Gibson 1987;Blake, Leppinen & Wang 2015;Kang & Cho 2019). As shown in figure 3, a long-lasting downward-moving bubble can be observed after the first collapse phase, which induces a downward fluid motion of a 'virtual/added mass' (Benjamin & Ellis 1966;Philipp & Lauterborn 1998). Thus, a portion of fluid is expected to move in the opposite direction to bubble migration. One can also readily understand this using Newton's third law. Specifically, the bubble propels itself downward by pushing a certain amount of fluid upwards in the form of an interface jet. The above argument is justified by a scaling analysis in § 6 and a more quantitative discussion is given.

Bubble initiation at the water-oil interface
In this section, we discuss a special case, i.e. the initiation point of the bubble is positioned at the interface (γ = 0). As shown in figure 4, both the upper half and lower half of the bubble expand hemispherically (frames 1-2), while the different densities of the two fluids lead to different velocities of the upper and lower parts of the bubble (frame 3). Finally, a downward jet forms due to the faster contraction of the upper part, carrying some oil into the water through the interface (frame 4), which is a direct consequence of the lower inertia of the oil. Yamamoto et al. (2021) supposed that such a bubble jet is an important mechanism in ultrasonic emulsification; however, no clear evidence was provided therein due to limitations of spatio-temporal resolutions. As shown in frame 5, a subsequent fast migration of the toroidal bubble in the rebound phase is observed.   Presumably, the whole bubble is fully submerged in water at this moment. Owing to the significant influence of the bubble on the interface motion, the interface jet forms earlier and rises more quickly compared with the four cases considered in figure 3. The interface evolution on a longer time scale is presented in frames 6-10. A mushroom-shaped jet is gradually formed (frames 6-7). Subsequently, its cap grows and forms a stretching water film; meanwhile, the jet becomes thinner and thinner (frame 8). Afterwards, the interface jet splits into two parts, i.e. the cap and the lower part (frames 9-10). The very distorted cap turns into an ellipsoidal one due to surface tension (frame 10). It takes much longer for the water droplets to fall down to the interface (not given here). This finding may be the second mechanism of fluid mixing in ultrasonic emulsification. Different from the first mechanism, the second one transports the heavier fluid into the lighter fluid. Additionally, the pinch-off droplet is much larger than the oil droplets carried by the bubble jet; thus these two mechanisms may play a role in different stages of emulsification.
A similar experimental observation can be found in Yin et al. (2020) (see figure 9 therein). The difference is that their laser-induced cavitation bubble is one order of magnitude smaller than ours and the associated Bo is smaller than 1. Thus the effect of surface tension dominates over the effect of gravity in their system. Later we discuss the dependence of the interface jet dynamics on Bo via numerical simulations.

Bubble initiation in oil
In figure 5, we present and discuss four representative experiments in which bubbles are generated in the sunflower oil, also starting with a large standoff parameter γ o . In the first experiment (see figure 5a), γ o = 1.27. The bubble oscillates spherically in the first cycle and no evident liquid jet can be seen (the reader is referred to Appendix B for the bubble dynamics in the first cycle). Only a tiny protrusion can be observed at the bottom of  the bubble surface during the rebound phase (frame 1). Frame 2 shows the maximum volume of the bubble in the second cycle. Subsequently, the bubble oscillates for several cycles and migrates towards the interface very slowly (frames 3-4). The bubble only causes some deformation of the interface and no penetration occurs. In frame 5, non-condensable gas bubbles are rising above the water-oil interface. Below this standoff parameter γ o ≈ 1.27, we find that the bubble can penetrate the water-oil interface.
In the second experiment shown in figure 5(b), γ o is reduced to 1.2. A sharp downward protrusion forms at the bottom of the bubble (frame 1), indicating the formation of a liquid jet. However, the jet tip cannot reach the position of the interface before disintegrating (frame 2). Though the bubble energy dissipates much after multiple oscillations, the downward-migrating bubble eventually passes through the water-oil interface (frames 3-5), resulting in a strong mass transport between the oil and water. No evident interface jet can be observed in this case.
In the third experiment shown in figure 5(c), γ o = 0.8. The downward bubble jet directly impacts and penetrates the interface (frames 1-2). More discussion of bubble jet dynamics is given in § 5. Thereafter, the whole bubble passes through the water-oil interface and further migrates downward (frame 3). Remarkably, a pronounced interface jet is generated but no necking or pinch-off occurs (frames 4-5).
In the fourth experiment shown in figure 5(d), the bubble is initiated very close to the interface, i.e. γ o = 0.4. The bubble jet penetrates deeper into the water compared to the third experiment (frame 1). When the bubble rebounds to the maximum size (frame 2), more than half of the bubble enters the water. As the bubble re-collapses and continually migrates downward, it quickly becomes fully submerged in the water bulk (frame 3). As expected, an interface jet is produced afterwards (frame 4). The jet tip assumes a mushroom shape and the cap finally separates from the main jet (frame 5). We can compare this case with that in figure 3(d). Despite a similar standoff parameter in these two experiments, the interface jet in the current experiment case grows faster and finally the pinch-off occurs; while in the former case, the interface jet grows more slowly and only a necking phenomenon is observed. We explain this as follows. When a bubble is initiated in water, the bubble migrates away from the interface and thus the subsequent bubble-interface interaction is weakened. For a bubble initiated in oil, the bubble migrates towards the interface and gradually passes through the interface, leading to a stronger bubble-interface interaction.

Comparison of experiments with simulations
In this section, we compare the observed bubble-interface interactions in three typical experiments with simulation results. For bubble initiations in water and oil ( § § 4.1 and 4.2), we use a BI method to reproduce the experimental observations; however, the present BI model cannot be applied to an extreme case in which the bubble is initiated at the water-oil interface. Instead, we propose an extended Rayleigh-Plesset equation to model the bubble dynamics ( § 4.3). Figure 6 shows a comparison of the experiment in figure 3(c) with our BI simulations. Both the inviscid BI (without normal viscous stress) and viscous BI (with normal viscous stress) simulation results are plotted. Frames 1-4 show the bubble-interface interaction during the first cycle of the bubble. In this stage, both the simulation results agree well with the experimental observations, which implies that this transient process of bubble-interface interaction is inertia-dominated and the viscosity plays a minor role. We also turn off the surface tension and gravity in our simulation and the results are slightly altered (not shown here).

Bubble initiation in water
We notice that the characteristic time scale of the interface jet evolution is much longer than that of the bubble oscillation. It is not an easy task to simulate the whole process due to the multiple time scales. Thus we adopt a simplified model to simulate the subsequent interface jet evolution, i.e. removing the bubble from the simulation when it reaches its minimum volume. The justification is as follows. We notice that the bubble-interface  interaction mainly occurs during the first bubble cycle. The bubble is repelled by the interface and thus the bubble-interface distance keeps increasing after the first collapse phase. More importantly, the energy loss of the bubble is significant (over 80 %) during the first collapse and rebound phase. The corresponding mechanism is very complex, which may be associated with acoustic radiation, heat transfer, condensation and so on (Keller & Kolodner 1956;Lee et al. 2007;Wang 2016;Zeng et al. 2018). Therefore, the bubble motion after the first cycle has a much smaller effect on the subsequent interface evolution. A similar numerical procedure can be found in previous literature (Fong et al. 2009;Borkent et al. 2009;Dadvand et al. 2012;Peters et al. 2013;Han, Zhang & Li 2014).
Frames 5-10 in figure 6 show the dynamics of the interface jet on a long time scale. At the early stage of the interface growth (frames 5-6), the difference between the two simulations is indistinguishable, indicating a minor effect of the viscosity within such a short time. As the interface jet develops, the difference between the two models gradually becomes evident (frames 7-9). Remarkably, the results from the viscous BI simulation (denoted by the black dashed lines) generally follow the experimental observation except for a slight overestimation of the height of the interface jet. However, in the inviscid BI simulation (denoted by the red solid lines), both the jet shape and height have a striking difference from the experiment. This implies that the viscosity plays an important role in interface evolutions. Figure 7 shows a quantitative comparison of the time evolution of the interface height between experimental data and numerical simulations. Apparently, there are two distinct phases of the interface dynamics. In the first phase, its motion is related to the growth and collapse of the bubble. Both the inviscid BI and viscous BI results reproduce the experiment quite well. In the second phase, the interface rises and descends slowly on a longer time scale. Both the simulation results follow the experimental data at an early stage. However, the maximum height of the interface jet is significantly overestimated by the inviscid BI while the viscous BI works much better. This again confirms the importance of the viscous effect in interface jet dynamics. More comparisons are made in § 6. 932 A8-15 We can estimate the vorticity-affected region by √ νt, where ν is the kinematic viscosity of the liquid. As for the bubble dynamics, the characteristic time is taken as the first cycle of the bubble (t = 3 ms), yielding √ ν 1 t ≈ 55 μm, which is much smaller than the bubble size. Hence, the viscous effect can be safely ignored in the stage of bubble dynamics. As for the longer-time-scale interface evolution, we take the moment when the interface reaches the maximum height as the characteristic time (t = 187 ms), yielding h m /6 < √ ν 2 t = 3.6 mm < h m /5. This justifies the viscous potential flow model in most of the growing process of the interface. Considering the simplifications, the quantitative agreement between the viscous BI simulation results and the experiment is surprising and encouraging.

Bubble initiation in oil
The experiment in figure 5(c) is compared against the numerical results, as shown in figure 8. In each frame, the experimental observation is given in the left-hand half and the simulation result is given in the right-hand half, together with the pressure field. The bubble expands spherically and the interface is pushed downward with a quite uniform pressure field surrounding the bubble (frame 1). During the collapse phase, a high-pressure region is gradually formed above the bubble (frames 2-3). A downward jet is directed towards the interface at the minimum volume moment (frame 3). Afterwards, the bubble rebounds with a continuous migration towards the interface (frame 4) and the protrusion causes a small downward bulge of the interface (frames 5-6). The high-pressure region at the jet tip moves downward with the bubble migration and finally acts on the interface. In this case, both bubble and interface profiles are accurately captured by the simulation, which again demonstrates the validity of the present model.
Although the liquid jet is not strong enough to penetrate the interface in the first cycle, the bubble migrates all the way down and finally passes through the interface in the second cycle, followed by the slow motion of the interface jet, as shown in figure 5(c). The 'bubble removal' treatment in § 4.1 cannot be applied in this case because the bubble migration towards the interface in the subsequent oscillations has a considerable effect on the interface jet evolution. In § 5, this numerical model is used to study the bubble jetting behaviours for bubbles initiated in oil.

Bubble initiation at the water-oil interface
In this section, we consider an extreme case of bubble initiation at the water-oil interface, as previously discussed in § 3.2. Unfortunately, our present BI model cannot be applied directly in this extreme situation. To better understand the bubble dynamics, we propose a simplified theoretical model to reproduce the motion of the bubble top and bottom. Since both the upper half and lower half of the bubble expand and collapse hemispherically during most of the first cycle, we extend the classic Rayleigh-Plesset equation to describe the gross motion of the bubble, given by where R 1 and R 2 denote the radii of the lower and upper hemispherical bubbles, respectively. The two hemispherical bubbles are connected; thus the gas pressures inside the two bubbles are the same, and we have (4. 2) The comparison of the experimental data and theoretical results obtained from (4.1) is given in figure 9. As can be seen, the maximum size of the upper hemispherical bubble is larger than that of the lower one. However, the oscillation period of the upper one is shorter, indicating that the bubble top oscillates faster than the bubble bottom. This is a direct consequence of the lower inertia of the oil. Remarkably, this simplified model reproduces the experimental observations very well, accounting for the faster contraction of the bubble top and a downward liquid jet.

Bubble jet dynamics
Since the bubble jet penetration into the water-oil interface is an important mechanism of emulsification or fluid mixing, especially the formation of very fine droplets, here we quantitatively study the jet dynamics for bubbles initiated in oil. Firstly, all the experimental data are divided into two regimes, namely whether the bubble jet penetrates the interface or not. The critical standoff parameter and jet impact velocity are discussed in § 5.1. Thereafter, we obtain the bubble jet dynamics as a function of governing parameters via numerical simulations, which are discussed in § 5.2. Figure 10 shows the variation of the jet impact velocity V jet versus the standoff parameter γ o , in which the magenta triangles and the blue circles denote whether the downward bubble jet can penetrate the water-oil interface or not, respectively. Figures 10(a)-10(c) represent different oils used in the experiments, namely the three types of oil given in table 1, respectively. Here the jet impact velocity is estimated from two or three adjacent frames before and after the jet penetration moment. Therefore, the experimental data of V jet are average values of the velocity of the bubble surface, which can be treated as lower bounds of the real values. Due to the limitations of the camera, unfortunately, the maximum relative error of V jet can be estimated as 50 %. The numerical results from BI simulations are also added in figure 10(a) for comparison. As expected, the numerical results are higher than the experimental data. The BI results are not given in figure 10(b,c) because the bubble-interface interaction is so weak that the liquid jet is extremely thin and hardly resolved in simulations. As shown in figure 10(a), the maximum value of V jet reaches 200 m s −1 around γ o = 0.55. For very small standoff parameters (γ o 0.3), V jet decreases to 120-130 m s −1 and is almost independent of γ o . As γ o increases from 0.55, V jet gradually decreases to ∼95 m s −1 . No evident jet formation can be found when γ o > 1.4. The bubble-interface interactions for bubbles initiated in oil are divided into two regimes. First, when the bubble is generated close to the water-oil interface, the bubble jet can penetrate the water-oil interface (denoted by the magenta triangles) and consequently the oil is transported to the water bulk by the jet (see a typical case in figure 5c). Second, when the bubble is far from the interface, the thickness of the oil layer between the lower bubble surface and the interface increases with γ o and the jet velocity is greatly reduced when travelling across the oil layer (see figure 5b); thus the jet cannot penetrate the interface (denoted by the blue circles). A sharp transition between the two regimes is present. The critical value of the standoff parameter (denoted by γ oc ) is between 0.89 and 0.91. Figure 10(b) presents the dependence of V jet on γ o for type-II oil. The density ratio increases to 0.971. Consequently, the bubble jet dynamics is greatly affected. The overall magnitude of V jet is smaller than that in figure 10(a) and the maximum value of V jet is about 170 m s −1 . Additionally, the critical value of γ oc is greatly reduced (0.54-0.61) and the jet formation disappears at γ o ≈ 0.78. Therefore, we can conclude that the bubble jet dynamics is mainly governed by the standoff parameter γ o and the density ratio α. A larger α leads to a smaller γ oc .

Penetration of the bubble jet into the interface
To better understand the influence of oil viscosity on the jet dynamics, more experiments are carried out using the type-III oil (more viscous than the type-II oil). As shown in figure 10(c), the overall magnitude of V jet is slightly lower than that in figure 10(b). This can be explained as follows. After the jet impacts the lower surface of the bubble, the surface of the toroidal bubble becomes unstable and not smooth, and thus the protrusion of the bubble receives a larger drag force with increasing viscosity. Since we use the protrusion location to estimate the jet velocity, the experimentally obtained V jet decreases. We note that the critical value of γ oc ranges from 0.58 to 0.61 and the bubble jet disappears at γ c ≈ 0.76. These two special characteristics are very similar to those in figure 10(b), which confirms the minor role of the viscosity in bubble jetting behaviours.

Bubble jet dynamics as a function of parameters
In this section, we numerically study the dependence of bubble jet dynamics on the standoff parameter γ o and density ratio α. The standoff parameter γ o is set in the range [0.3, 0.6]; thus a jet penetration into the interface is expected in most cases. Since the densities of common industrial oils are above 0.7 g cm −3 , the density ratio α varies from 0.7 to 0.95 in our simulations. Figure 11 depicts a comparison of the flow field at the bubble jet impact moment for different density ratios while other parameters are kept the same. Figure 11(a) shows the numerical results for α = 0.7. The liquid layer between the bubble surface and the interface is quite thin and an energetic downward jet is about to penetrate the lower interface. The contours in the left-hand half and right-hand half present the velocity and pressure fields, respectively. As can be seen, there exists a localized high-pressure region above the bubble surface, which is the mechanism for the emergence of a high-speed liquid jet (Blake et al. 1986;Tong et al. 1999). Only the liquid within the bubble jet has a relatively high speed at this moment, which is expected to penetrate the fluid-fluid interface. In figure 11(b), the density ratio α is increased to 0.95. The jet impact occurs around the moment of the minimum bubble volume; thus the magnitude of the pressure surrounding the bubble is quite high and the jet speed reaches about 23. However, the volume of the liquid jet is much smaller than that in figure 11(a). Additionally, the thickness of the liquid layer between the bubble and interface is larger, and thus the liquid jet will be weakened before penetrating the interface. Now we quantitatively study the volume of the liquid jet and the associated kinetic energy. Figure 12(a) shows the dependence of the jet volume on γ o and α, with a sketch of the jet volume in the inset. Surprisingly, the jet volume is not sensitive to the variation of γ o in the range 0.3 γ o 0.6. Due to the small γ o , similar bubble collapse patterns and jet dynamics are observed. Additionally, the jet impact velocity V jet is almost independent of γ o in the range 0.3 γ o 0.6 and the difference is within 15 % (not shown here). This finding is consistent with the experimental results in figure 10. The density ratio α is the main factor that influences the jet volume. Specifically, the jet volume deceases linearly as α increases from 0.7 to 0.95. However, the jet impact velocity V jet increases with α. To better understand the dependence of the bubble jet dynamics on γ o and α, we introduce the kinetic energy of the liquid jet: van der Meer 2020), as denoted by the green solid line in the inset of figure 12(a). A parameter map of the dependence of E jet on the density ratio α and the standoff parameter γ o is given in figure 12(b). Energy E jet decreases with increasing α or increasing γ o .
Obviously, E jet is more sensitive to α. This finding may provide a reference for practical operations.
6. Interface jet dynamics 6.1. Dependence of the interface height on γ As discussed in the previous sections, the dynamics of the water-oil interface is highly dependent on the standoff parameter γ . In this section, we first discuss the dependence of the maximum interface height h m on γ w . In the sunflower oil experiments, h m is very small when γ w 1.4, and thus we only discuss the situations in which γ w 1.4. Figure 13(a) shows the variation of h m versus γ w . Numerical results obtained from the inviscid BI and viscous BI simulations are also added for comparison. When γ w < 0.2, h m is about 4 due to the strong interaction between the bubble and interface. As γ w increases, the bubble-interface interaction becomes weaker and h m decreases accordingly. Remarkably, our viscous BI simulation results (the green line) show excellent agreement with the experimental data. As expected, h m is overestimated if the viscous effect is absent in the simulation (see the black line), especially for relatively small γ w .
In figure 13(b), the relationship between h m and γ w is presented on a doubly logarithmic scale. For a large standoff parameter γ w > 1, the plot reveals a power law: where θ is the power-law exponent. A numerical fit to the experimental data gives a power-law exponent θ = −3.1, with a minimum root-mean-square error of 0.036. Both the simulation results give θ = −3.2, which is very close to the experimental value.
To reveal the underlying mechanism of this power law, we perform a scaling analysis of the interface jet evolution with the help of the Kelvin impulse theory (Blake & Gibson 1987;Kang & Cho 2019). More details can be found in Appendix A. The following scaling law is obtained: Given the same α and δ in the experiments, we have h m ∼ γ −4 w . The small difference between the theoretical result and the fitting results from experiments can be explained as follows. First, the classic Kelvin impulse theory is derived from the spherical bubble equation and a linearized boundary condition. However, we use the experimental data in the range 1 < γ w 1.4 to fit θ, which may not be large enough to satisfy the linearized boundary condition. More importantly, the surface tension and viscous effects are not included in the scaling analysis. To examine this point, we turn off the surface tension and viscosity in the simulation, and the power-law exponent θ is approximately −4. More discussion of the dependence of θ on surface tension and oil viscosity is given in the following sections.
As for a bubble initiated in the oil, the relationship between h m and γ o is plotted in figure 14(a). The standoff parameter γ o ranges from 0.1 to 1.2. Generally speaking, the interface jet achieves a larger height when the bubble is initiated in oil than in water. Although the first bubble cycle has a dominant effect on the interface evolution, the subsequent multi-oscillations apparently have a considerable influence due to fact that the bubble migrates towards the interface and even penetrates into the lower dense fluid. Considering some uncertainties of the experiment, all the data shown in figure 14(a) roughly reveal a linear dependency of h m on γ o in the range 0.1 < γ o < 1.2 with the slope being −4.38.
It is interesting to discuss the dependence of h m on γ o for large γ o and compare it with the bubble-in-water cases. We present the experimental data in a doubly logarithmic plot, as shown in figure 14(b). A fitted straight line with a slope of −5.5 is obtained using 30 experimental data with γ o > 0.8. The absolute value of the slope is larger than that in bubble-in-water cases, which again demonstrates the stronger bubble-interface interaction for bubbles initiated in oil.

Effects of the oil viscosity
It is demonstrated in figure 13 that the viscosity of oil plays a role in the evolution of the interface jet. In this section, the effects of the oil viscosity are further investigated using purposely conducted experiments. As given in   Figure 15 compares the physical process between two representative experiments of bubbles initiated at the interface (γ = 0). The lighter liquid used in figure 15(a) is the type-II oil in table 1. The bubble dynamics in the first cycle is similar to that in figure 4. Specifically, the bubble keeps a quasi-spherical shape before the final stage of collapse (frames 1-3). A downward liquid jet forms near the moment of the minimum bubble volume (frame 4). Subsequently, the bubble rebounds (frame 5) and migrates downward (frame 6). The interface jet gradually develops on a longer time scale (frames 7-9), followed by the descent process (frame 10). To explore the effects of the oil viscosity, the third type of oil with a higher viscosity is used. In figure 15(b), the dimensionless times of each frame are set very close to those in figure 15(a). We notice that the overall bubble dynamics is slightly altered (frames 1-5) and only the maximum height of the interface jet becomes lower. The interface jets in both experiments assume a smooth-hump shape. Neither necking phenomenon nor pinch-off of the interface jet occurs. Therefore, one can expect only one mechanism of fluid mixing in these two series of experiments, i.e. the bubble jet penetration. This is quite different from that in the sunflower oil experiments (see figures 3-5). There are two reasons for this phenomenon. First, the density ratio α increases, and thus the asymmetry of the flow field decreases. Second, the viscosity of the oil increases, which hinders the growth of the interface. The associated Reynolds numbers of the interface jets are 15 and 3 in the above two cases, respectively. Figure 16 shows the dependence of h m on γ w and γ o in these two series of experiments. Numerical results are also added in figure 16(a) for comparison. It seems that the oil viscosity has little effect on h m for large γ w . As γ w decreases, the difference in h m between the two series of experiments becomes evident. The red line denotes the numerical results without viscous effect, which is the upper bound of h m for the specific surface tension and density ratio in the experiments. The interface height is smaller than the uncertainty of measurement when γ w > 1, and thus a power law of h m and γ w is studied via numerical simulations in the following discussion. Figure 16(b) shows the variation of h m with γ o . The distribution of the experimental data is more dispersed than in figure 16(a), indicating that the interface jet is more unstable for bubbles initiated in oil. The penetration of the bubble into the interface is mainly responsible for the dispersion.
To quantitatively investigate the effect of the viscosity on the power law exponent θ, we perform numerical simulations with five different oil viscosities, namely μ = 0, 0.05, 0.108, 0.2 and 0.35 Pa s. The surface tension is turned off in the simulations (Bo = ∞) and other parameters are set according to the experiment: α = 0.971, δ = 0.038. As shown in figure 17, h m decreases with increasing viscosity of the oil. The variation of h m with γ w can be roughly divided into two regimes, namely γ w 1 and γ w > 1. A numerical fit to the data within the green dashed rectangle (γ w 1.1) gives the slopes of the curves, i.e. the power-law exponents in (6.1a,b), as shown in the inset. In the inviscid simulations, we obtain a slope of −3.968, which is very close to the theoretical value of −4. As the oil viscosity increases, the slope increases slowly (the absolute value decreases), indicating that the oil viscosity hinders the growth of the interface jet. When μ = 0.05 Pa s, the associated Re * within the green dashed rectangle can be estimated as O(10 −1 -10 0 ). As the oil viscosity increases, Re * is expected to be smaller in the rest of the simulations and the oil viscosity has a considerable effect on the interface jet evolution.
6.3. Effects of the surface tension In this section, we numerically investigate the effects of surface tension. We perform five series of numerical simulations with different surface tension coefficients, namely σ = 0, 0.01, 0.02, 0.03 and 0.04 N m, respectively. The corresponding Bond numbers are ∞, 6.2, 3.1, 2.1 and 1.5, respectively. The viscous effect is turned off and other parameters are set as α = 0.971, δ = 0.038, γ w = 0.4 and ε = 100 in the simulations. Figure 18 shows a comparison of the interface jet morphology for different σ . At the early stage of the interface jet evolution (figure 18a), the profiles of the interface nearly coincide. The surface tension plays a minor role at this stage because the interface is slightly deformed. In figure 18(b), we can find some difference between the five cases. When σ = 0 (denoted by the blue line), the interface exhibits some local instabilities and an annular pinch-off occurs (see the inset). Therefore, the simulation stops at t ≈ 40 in this case. In the simulations with non-zero σ , the interface jet is overall a smooth hump and its height decreases with σ . In figure 18(c), the interface jet reaches its maximum height in the σ = 0.04 N m case, while in the other three cases with smaller σ the interfaces are still rising. Particularly, the interface jet in the σ = 0.01 N m case exhibits a necking phenomenon. Subsequently, this necking disappears as the interface jet further develops (see figure 18d). At t ≈ 250, the interface jet in the σ = 0.01 N m case finally reaches the maximum height, while the interface jets in the other cases are descending due to surface tension and gravity. The dependence of h m on γ w for different σ (or Bo) is presented in figure 19, showing that h m decreases with increasing σ . A numerical fit to the data when γ w 1.1 (within the green dashed rectangle) is performed with slopes plotted in the inset. As σ increases from zero, the slope increases slowly from −4 (the absolute value decreases), indicating that the surface tension also hinders the growth of the interface jet.
The associated Bo is above 1 in these simulations with σ 0.01 N m, indicating that gravity is also an important factor that decelerates the interface jet. To estimate the respective contribution of surface tension and gravity to the deceleration of the interface jet, we calculate and compare the associated surface tension energy E s and gravitational potential energy E p . Figure 20(a) shows the time evolution of E s for the non-zero σ cases in figure 18. It is noted that both the maximum value of E s and the slope of the curves increase with σ . Figure 20(b) depicts the time evolution of E p . For a smaller σ , the interface jet achieves a higher height, thus the maximum value of E p increases with decreasing σ . The magnitudes of E s and E p are comparable, so that both surface tension and gravity are important in the present system (centimetre-sized scale). As can be expected, the surface tension would dominate the interface jet dynamics for a smaller-scale problem and the gravity effect would gradually become negligible. Hence, it is interesting to study the effects of Bo on the interface jet dynamics in a large parameter range. Here we set 10 −4 Bo 10 2 with σ = 0.035 N m in the simulations, which covers most of the bubble scales encountered in practical applications (60 μm R m 0.1 m). Other parameters are γ w = 0.4 and Re = ∞. Figure 21 shows max(E s )/max(E p ) as a function of Bo for different density ratios. The blue circles and red rectangles seem to fall nicely on two parallel straight lines, respectively, which suggest a power law: A fit to the numerical results gives the power-law exponent θ * = −0.93, which is found to be nearly independent of α. We also examined other standoff parameters and found that θ * slightly varies with γ w , which is not discussed further. Here we emphasize that the surface tension greatly hinders the interface evolution when Bo < 1. In practical applications, to enhance emulsification, one need to reduce the surface tension, which is often realized by increasing the fluid temperature or adding surfactants (Canselier et al. 2002;Modarres-Gheisari et al. 2019).

Effects of the density ratio α
In this section, the effect of the density ratio α on the interface jet evolution is investigated. The simulations cover the range 0.5 γ w 2.0 and 0.7 α 0.95. Other dimensional parameters are set according to our experiments with the third type of oil: R m = 15 mm, μ 2 = 0.35 Pa s and σ = 0.035 N m. Figure 22(a) shows the time evolution of the interface height for different α when γ w = 1. During the bubble dynamics stage (t 2), the interface is pushed upward by the expanding bubble and attracted downward by the  collapsing bubble. The interface height slightly increases with decreasing α in this stage. This transient process is mainly inertia-dominated and one can easily find that the interface evolution is closely related to the bubble motion. After the bubble collapses, the interface rises again but with a much slower velocity. Apparently, the interface jet rises higher for a smaller α. The inset reveals a linear dependency of h m on (1 − α)/(1 + α) 2 in the range of 0.7 < α < 0.95 with the slope being 8.2. However, this slope depends on γ w .
Since our scaling analysis (see Appendix A) gives h m = kδ −2 (1 − α)/(1 + α) 2 /γ 4 w (k is the prefactor and δ = 0.038), we plot the relationships between h m and (1 − α)/(1 + α) 2 /γ 4 w for different α and γ w in figure 22(b). At first sight, the lines collapse together for 932 A8-29 a large standoff parameter γ w 1. However, since the surface tension and viscosity come into play, the slope of each curve (kδ −2 ) is found to increase from 8.5 to 13.1 when γ w increases from 1 to 2, as shown in the inset.

Summary and conclusions
In summary, we have studied the complicated interaction between cavitation bubbles and a water-oil interface on multiple time scales via experimental, numerical and theoretical means, including the bubble dynamics on a short time scale ∼10 0 ms and the water-oil interface evolution on a much longer time scale ∼10 2 ms. In the experiments, the electric discharge method is used to generate centimetre-sized oscillating bubbles, enabling us to obtain a higher spatio-temporal resolution of high-speed recordings than previously available. To better understand the underlying mechanisms, BI simulations are performed and the results are found to be in good agreement with the experimental observations. Some specific physical phenomena are explained with a scaling analysis and theoretical modellings. The main conclusions are drawn as follows.
During the first bubble cycle, the interface evolution is closely related to the bubble motion. Specifically, the interface is pushed away by the expanding bubble and attracted by the collapsing bubble. Regardless of the bubble initiation position, a downward migration of the bubble can be observed in the experiments, accompanied by a very thin liquid jet with the maximum velocity exceeding 200 m s −1 . It is demonstrated that this transient process is mainly inertia-dominated and can be well reproduced by our BI simulations. For a bubble initiated at the interface, the upper half and lower half of the bubble behave like two hemispherical bubbles before the final stage of collapse. An extended Rayleigh-Plesset model is proposed that well predicts the asymmetric dynamics of the bubble, which accounts for a faster contraction of the bubble top and a downward liquid jet.
For a bubble initiated in the oil phase or at the water-oil interface, the bubble jet can penetrate the interface when the standoff parameter is less than a critical value γ oc and consequently cause fluid mixing. It is experimentally demonstrated that γ oc highly depends on the density ratio α. For instance, γ oc is found around 0.9 when α = 0.919 and 0.6 when α = 0.971. Additionally, α is also the most important parameter that governs the volume and kinetic energy of the bubble jet. Surprisingly, within the parameter space in our simulations (0.7 α 0.95, 0.3 γ o 0.6), the jet velocity, volume and kinetic energy are nearly independent of γ o .
After the bubble dynamics phase, a pronounced interface jet gradually forms on a much longer time scale, the dynamics of which has been demonstrated to be affected by gravity, surface tension and viscosity simultaneously. Particularly, a pinch-off of the interface jet occurs in the sunflower oil experiments (α = 0.919) with small standoff parameters. This microscopic mechanism is reckoned as the second mechanism in ultrasonic emulsification. However, the pinch-off disappears when the density ratio increases to α = 0.971 (silicone oil). A modified BI model that incorporates the normal viscous stress can reproduce the main features of the interface jet. Further, we theoretically and numerically demonstrated the linear relationship between h m and (1 − α)/(1 + α) 2 with a varying prefactor mainly depending on γ w .
For bubbles initiated in water, a scaling analysis, without considering the surface tension and viscosity, suggests that the maximum interface jet height follows a power law h m ∝ γ −4 w when γ w > 1. However, the power-law exponent obtained from experimental data is −3.2. We use BI simulations to demonstrate the respective role of surface tension and viscosity in interface jet dynamics, which explains the deviation between theory and experiment. The dependence of the interface jet dynamics on Bo (the ratio between gravity and surface tension) is also investigated in a large parameter range (10 −4 Bo 10 2 ). The surface tension greatly hinders the interface evolution when Bo 1. Thus, the surface tension should be lowered to enhance emulsification, which is consistent with practical operations. For a bubble initiated in oil, it migrates towards the interface and subsequently passes through the interface, leading to a stronger bubble-interface interaction. As a result, the interface jet grows higher and is easier to pinch-off, compared with that in the initiation-in-water case. Therefore, considering both the mechanisms that contribute to fluid mixing/emulsification, bubbles initiated in the lighter fluid play a more important role than those in the denser one.    · ρ 1 − ρ 2 ρ 1 + ρ 2 · R 5 m P The dimensionless form of (A14) can be written as For convenience, the dimensionless maximum height of the interface jet is still denoted by h m in the main body of the paper.

Appendix B. Experimental observations of bubble dynamics in the first cycle
In figure 3, we only present the bubble dynamics in the first cycle for the first case (see figure 3a). Figure 23 shows the results of the other three cases of figure 3. As shown in figure 23(b), the bubble keeps a spherical shape during most of the first cycle. In figure 23(c), a downward liquid jet forms around the moment of the minimum bubble volume. In figure 23(d), the jet forms during the final collapse phase of the bubble. Figure 24 shows the bubble dynamics in the first cycle for the same cases as in figure 5. The bubble keeps a spherical shape during the expansion phase in all the cases and the bubble becomes elongated during the collapse phase. The jet forms during the rebound phase of the bubble (see figure 5).