2.1. Experimental set-up

As shown in figure 1, the experiment is conducted in a $300\times 300\times 300\ \textrm {mm}^3$ water tank, where the cavitation bubble is generated using an electric spark generator (ESG), whose reproducibility has been validated in our previous studies (Cui et al. Reference Cui, Zhang, Wang and Khoo2018; Han et al. Reference Han, Zhang, Tan and Li2022). The tips of two copper wires are brought into contact at a certain depth below the free surface, initiating the bubble formation. Upon discharge, a current with a voltage of up to 400 V is released by the ESG, creating an electric spark at the wire tips. This spark generates intense Joule heating, vapourising the surrounding fluid and causing it to expand into a cavitation bubble with a maximum radius $R_m$ of approximately 20 mm. A high-speed camera (Phantom V2012) operating at 30 000–40 000 frames per second is used to capture the transient interaction between the bubble and the free surface.

We establish an experimental technique to generate well-controlled surface perturbations by using the meniscus rising mechanism (Clanet & Quéré Reference Clanet and Quéré2002) and the contact line pinning effect (Schäffer & Wong Reference Schäffer and Wong1998). A thin rod of 0.2–0.5 mm in radius is vertically mounted on a micro-controlled translation table. A stable meniscoid interface where capillary forces and gravity are balanced is created by piercing the thin rod into the liquid pool. The thin rod is positioned coaxially with the bubble initiation point. When the thin rod is vertically moved at a sufficiently low speed (less than $10\ \textrm {mm}\,\textrm {s}^{-1}$ ) fixed on a high-precision displacement platform, the three-phase contact line at the rod’s surface remains pinned (Hansen & Miotto Reference Hansen and Miotto1957; Elliott & Riddiford Reference Elliott and Riddiford1967; Schäffer & Wong Reference Schäffer and Wong1998). This leads to changes in the key parameters and geometry of the meniscus. This technique allows for the generation of menisci with high controllability and reproducibility. To extend the parameter range, the rod was coated with a super-hydrophilic material (Shangmeng Technology Wuxi Co. Ltd., SM-QS3500). As a result, the maximum meniscus height increases from 0.27 mm (without coating) to 0.95 mm (with coating). Additionally, a digital single-lens reflex (DSLR) camera (Nikon Z 6II) with a macro lens is used to capture the details of the meniscus above the free surface before the bubble initiation, with a spatial resolution of $5.1\pm 0.1\ \unicode{x03BC}\textrm {m}$ (one pixel of the image), which is approximately 1.4 % of the thin rod radius. One may easily notice that the radius of the thin rod is much smaller than the capillary length $l_{\textit{c}}=\sqrt {\sigma /\rho _l g}=2.7\ \textrm {mm}$ ( $\sigma =0.073\ \textrm {N}\,\textrm {m}^{-1}$ , $\rho _l=997\ \textrm {kg}\,\textrm {m}^{-3}$ , $g=9.8\ \textrm {m}\,\textrm {s}^{-2}$ ). Under such circumstances, the influence of the thin rod itself (radius and depth) can be minimised. The demonstration for the influence of the thin rod is referred to Appendix B.

Figure 1.

(a) Schematic of experimental set-up. The bubble is initiated coaxially with the thin rod covered by super-hydrophilic coating in a $300\times 300\times 300\ \textrm {mm}^3$ water tank. An image of the meniscus is shown in the inset on the top left. (b) Comparison between representative experiments of free-surface–bubble interaction (i) without and (ii) with perturbation. The standoff parameter $\gamma$ for both experiments are 1.30. The bubble expands and collapses normally in the absence of the surface perturbation. Secondary cavitation caused by the rarefaction wave reflected by the free surface can be observed in the fluid field. In contrast, a cavity forms on the initially perturbed free surface and coalesces with the primary bubble, creating a channel that enables ventilation between the atmosphere and the bubble interior. Under such circumstances, the cavitation in the fluid field cannot be observed and the bubble exhibits a misty shape upon collapse.

Two typical experiments of free-surface–bubble interactions: (i) without and (ii) with surface perturbation, are presented in figure 1(b). As can be seen, in the case without perturbation, secondary cavitation bubbles caused by the free-surface-reflected rarefaction wave can be clearly seen. However, a cavity forms at the free surface and coalesces with the bubble in the case of a surface perturbation, while secondary cavitation is invisible. The underlying mechanisms associated with this phenomenon will be examined in § 3.3. Since the maximum bubble radius $R_m$ and the characteristic height of the meniscus $h_m$ are the two crucial parameters in the experiments, the determination procedure in an experiment is given as follows. First of all, the characteristic height of the meniscoid interface $h_m$ is directly measured from the image captured by the camera, with the axis of the macro lens aligned and parallel to the free surface. Then, the bubble oscillation is recorded with the high-speed camera. For relatively large standoff parameters, the bubble remains approximately spherical when it reaches maximum radius. Therefore, $R_m$ can be calculated with $R_m=D_m/2$ , where $D_m$ is the maximum bubble diameter measured directly from the experiments. As the bubble gets closer to the free surface, it will exhibit an apparent non-spherical shape. Under such circumstances, $R_m=D_m/2$ tends to fail for insufficient accuracy. Hence, we adopt the same correction method as Li et al. (Reference Li, Khoo, Zhang and Wang2018). To ensure the reliability and reproducibility of the experimental data, the experiments are repeated three times for each configuration. By elaborately manipulating the voltage of the ESG, and the velocity and displacement of the thin rod, the reliability and reproducibility of the experiments are encouragingly good. The relative discrepancy of the maximum bubble radius $R_m$ and meniscus height $h_m$ can be controlled below 5.5 % ( $\approx 1\ \textrm {mm}$ ) and 1.2 % ( $\approx 0.01\ \textrm {mm}$ ), respectively.

Figure 2.

(a) A schematic of the numerical model of free-surface–bubble interaction. The $z$ axis is coaxial with the thin rod, while the origin $O$ is located on the free surface. The bubble is initiated at point $O_b$ located beneath the free surface at depth $H$ with a initial radius of $R_0$ . (b) Configuration of the computational domain. The cavity length $h$ is the distance from the cavity bottom to the initial undisturbed free surface. The radius and height of the cylindrical computational domain are set as 10 $R_m$ and 12.5 $R_m$ , respectively, while the depth of the fluid is set as 10 $R_m$ .

2.2. Numerical model

To simulate the highly nonlinear interaction between the bubble and the perturbed free surface, including cavity–bubble coalescence, bubble ventilation, bursting jets and surface splashing, a segregated two-phase flow solver based on finite volume method (FVM) is employed in this study. The schematic of the physical problem is shown in figure 2(a). We establish a cylindrical coordinate system $O {-}r\theta z$ with the origin $O$ located at the free surface and aligned with both the thin rod and the bubble initiation point $O_b$ . The bubble is initialised at a depth $H$ with an initial radius $R_0$ .

Following the discharge, the bubble rapidly expands, during which the Mach number is estimated to be $\textit {Ma}\sim \mathcal{\textit {O}}(10^{-3}{-}10^{-2})$ (peaking at 0.03 when the bubble reaches its minimum volume). Here, the Mach number is defined as $\textit {Ma}=U_{\textit {collapse}}/c$ , where $U_{\textit {collapse}}$ represents the maximum velocity during the bubble collapse (reaching approximately 50 $\mathrm{m\,s}^{-1}$ in our experimental measurements) and $c$ denotes the undisturbed liquid sound speed. At such low Mach numbers, liquid compressibility effects become negligible (Wang Reference Wang2014), allowing us to safely treat the liquid as incompressible in our simulations. One can refer to a more detailed verification of the influence of the liquid compressibility in Appendix C. By employing a segregated two-phase flow solver with incompressible liquid and compressible gas (Kwak, Kiris & Kim Reference Kwak, Kiris and Kim2005; Dupuy, Toutant & Bataille Reference Dupuy, Toutant and Bataille2020), we can capture the transient cavity–bubble interactions accurately at low Mach numbers while boosting computational efficiency and minimising resource demands.

We also neglect mass and heat diffusion across the bubble surface, which is justified by the high Péclet numbers (Colombet et al. Reference Colombet, Legendre, Cockx and Guiraud2013; Li et al. Reference Li, Zhao, Zhang and Han2024): typically, $\textit {Pe}_{\textit {m}}\sim \textit {O}(10^6)$ (mass diffusion) and $\textit {Pe}_{\textit {h}}\sim \textit {O}(10^5)$ (heat diffusion) for centimetre-scale bubbles in this study. Specifically, $\textit {Pe}=R_m^2/T_{\textit{osc}}D$ , where $T_{\textit{osc}}$ is the bubble oscillation period and $D$ is either the mass diffusion coefficient or the thermal diffusivity for mass or heat transfer, respectively. The gas is treated as compressible. Viscous effects and surface tension are also retained due to the presence of the thin rod and perturbations. Based on these considerations, the compressible forms of the Navier–Stokes and continuity equations are as follows:

(2.1) \begin{align} \frac {\partial \left (\rho \boldsymbol{U}\right )}{\partial t}+\boldsymbol{\nabla }\boldsymbol{\cdot }\left (\rho \boldsymbol{U}\otimes \boldsymbol{U}\right )=-\boldsymbol{\nabla }p+\boldsymbol{\nabla }\boldsymbol{\cdot }\unicode{x1D64F}+\boldsymbol{I}_\sigma , \\[-27pt] \nonumber \end{align}

(2.2) \begin{align} \frac {\partial \rho }{\partial t}+\boldsymbol{\nabla }\boldsymbol{\cdot }\left (\rho \boldsymbol{U}\right )=0, \\[0pt] \nonumber \end{align}

where $\boldsymbol{\nabla}$ denotes the gradient, $\boldsymbol{\nabla }\boldsymbol{\cdot }$ the divergence and $\otimes$ the tensorial product. Additionally, $\rho$ is the fluid density, $\boldsymbol{U}$ the velocity field and $p$ the pressure field. Here, $\unicode{x1D64F}$ is the viscous stress tensor of the Newtonian fluid which is defined as

(2.3) \begin{align} \unicode{x1D64F}=\mu \left [\boldsymbol{\nabla }\boldsymbol{U}+\left (\boldsymbol{\nabla }\boldsymbol{U}\right )^T-\frac {2}{3}\left (\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{U}\right )\unicode{x1D644}\right ] \!, \end{align}

with $\unicode{x1D644}$ the unit tensor and $\mu$ the viscosity. The last term $\boldsymbol{I}_\sigma$ in (2.1) represents the integral of the capillary force acting on the liquid/gas interface. One can refer to studies of Tryggvason et al. (Reference Tryggvason, Bunner, Esmaeeli, Juric, Al-Rawahi, Tauber, Han, Nas and Jan2001) and Koch et al. (Reference Koch, Lechner, Reuter, Köhler, Mettin and Lauterborn2016) for more details. When it comes to the incompressible liquid, the volume of fluid elements remains constant, resulting in a velocity field without divergence and a constant density, i.e. $\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{U}=0$ . Consequently, the divergence terms in (2.1)–(2.3) can be neglected, reducing the system to the classical incompressible Navier–Stokes equations with the continuity condition.

The two phases are immiscible and subject to a no-slip condition across the interface. To capture the sharp interface, the volume-of-fluid (VOF) method (Hirt & Nichols Reference Hirt and Nichols1981) combined with the high-resolution interface capturing (HRIC) scheme (Li, Duan & Zhao Reference Li, Duan and Zhao2022) is employed in this study. The pressure and velocity fields are shared by both phases, while viscosity and density are averaged using the volume fraction $\alpha _i$ (Miller et al. Reference Miller, Jasak, Boger, Paterson and Nedungadi2013). Thus, The continuity (2.2) for each phase reads:

(2.4) \begin{align} \frac {\partial (\alpha _i \rho _i )}{\partial t}+\boldsymbol{\nabla }\boldsymbol{\cdot }\left (\alpha _i\rho _i\boldsymbol{U}\right )=0,\quad i=l,g, \end{align}

where $\alpha _i$ and $\rho _i$ denote the volume fraction and density of the liquid and gas phases. This equation applies to both compressible and incompressible fluids by retaining or neglecting the divergence term. The computational domain is discretised into cells, with the sum of volume fractions in each cell equal to unity (Koch et al. Reference Koch, Lechner, Reuter, Köhler, Mettin and Lauterborn2016), i.e. $\alpha _l+\alpha _g=1$ . Thus, the overall density and viscosity in (2.1)–(2.3) are given by $\rho =\alpha _l\rho _l+\alpha _g\rho _g$ and $\mu =\alpha _l\mu _l+\alpha _g\mu _g$ . By solving (2.4), the two phases are distinguished: $\alpha _l=1$ for liquid, $\alpha _l=0$ for gas and $0 \lt \alpha _l \lt 1$ for the interface. Details of the scheme can be found from Miller et al. (Reference Miller, Jasak, Boger, Paterson and Nedungadi2013) and Koch et al. (Reference Koch, Lechner, Reuter, Köhler, Mettin and Lauterborn2016). Assuming incompressibility and non-diffusivity, the liquid phase has a fixed density of 997 $\textrm {kg}\,\textrm {m}^{-3}$ , while the gas phase is treated as a compressible and adiabatic ideal gas with the equation of state (EoS) defined as

(2.5) \begin{align} \rho \left (P_g\right )=\rho _0\left (\frac {P_g}{P_0}\right )^{({1}/{\kappa })}, \end{align}

where $P_g$ is the gas pressure at an arbitrary moment after the bubble inception, $\rho$ the gas density, which is the function of $P_g$ , and $\kappa$ the ratio of the specific heats, which is set as 1.4 according to the study of Li et al. (Reference Li, Zhao, Zhang and Han2024). The subscript ‘0’ denotes the initial properties of the gas. The $P_0$ and $\rho _0$ of the atmosphere can be set as 101 325 Pa and 1.2 $\rm {kg\,m}^{-3}$ , while the initial condition for the gas of the bubble interior will be addressed later in § 2.3 due to the complexity of bubble inception. By using (2.5) and the constant liquid density assumption, the system defined by governing equations (2.1)–(2.4) becomes closed and solvable. The validation of our numerical model with experimental data and other analytical models is presented in Appendix A.

The simulation set-up of the free-surface–bubble system is illustrated in figure 2(b). To optimise computational efficiency, an axisymmetric computational domain is established. The radius and height of the domain is set as $10R_m$ and $12.5R_m$ . Our simulation results show that the present domain set-up is sufficiently large that the boundary barely affects the bubble dynamics. Further increase of the computational domain sizes does not significantly change the bubble characteristics, with the bubble period only increasing by 0.1 % and 0.3 % when the domain size increases to $20R_m$ and $30R_m$ , respectively. One can refer to the details of the domain verification in Appendix A. The thin rod and the boundary beneath the free surface are set as no-slip boundaries to reproduce experiments in the water tank. To accurately capture the cavity and bubble dynamics, a refined mesh with a resolution of $5\ \unicode{x03BC}\textrm {m}$ is applied to the region where the cavity and bubble evolve. This local refinement approach reduces the computational demand by maintaining the total number of mesh cells below 1 500 000 while ensuring the accuracy of the simulation. One may also notice that, since the boundary condition of the thin rod is set as no-slip, it is certain that there will exist a boundary layer around the thin rod surface, which does reflect part of the real experimental situation and may pose an influence on the numerical simulation results depending on the mesh size. The thickness of the boundary layer can be estimated using $\zeta =\sqrt {\nu t}$ , thus evaluating its influence, where $\nu =1.0\ \textrm {mm}^{2}\,\textrm {s}^{-1}$ is the kinematic viscosity of the water and characteristics time $t\approx 3\ \textrm {ms}$ is taken as the first cycle of the bubble, yielding $\zeta \approx 55\ \unicode{x03BC}\textrm {m}$ . Apparently, the mesh size we choose is much smaller than the characteristic thickness of the boundary layer, sufficiently fine to capture its detailed characteristics. To further address this issue, we perform a convergence analysis using different refinement resolutions in Appendix A, with the mesh size ranging from $2$ to $20\,\unicode{x03BC}\textrm {m}$ . Comparison between these simulation results with different mesh sizes does not exhibit significant variations in cavity behaviours (see figure 22). A reasonable explanation given by Cox et al. (Reference Cox, Pearson, Blake and Otto2004) and Li et al. (Reference Li, van der Meer, Zhang, Prosperetti and Lohse2020) can be presented as follows. The Reynolds number associated with the spark-generated bubble has the order of $\textit {O}(10^5-10^6)$ and the Mach number is much smaller than one (a detailed non-dimensional analysis is deferred to § 3.1). Moreover, the scales of bubble and cavity ( $\sim \kern-2pt 10^0{-}10^1\, \textrm {mm}$ ) are much larger than the thickness of the boundary layer ( ${\sim} 10^{-2}\, \textrm {mm}$ ). Therefore, the detailed morphology of the boundary layer and motion of the contact line do not exert a remarkable influence on the cavity dynamics. With comprehensive considerations, the mesh configuration with a total mesh number of 1 500 000 and a minimum mesh size of 5 $\unicode{x03BC}\textrm {m}$ is adopted to attain a balance between calculation efficiency and accuracy. A more detailed analysis of the boundary layer effect is performed in § 5.4.

A modified axisymmetric boundary integral (BI) method based on our previous work concerning the interaction between an oscillating bubble and a two-phase interface (Han et al. Reference Han, Zhang, Tan and Li2022; Li et al. Reference Li, Zhao, Zhang and Han2024) is also adopted, where the two fluids (gas and liquid) on both sides of the interface are considered as inviscid and incompressible. Therefore, the Laplace equation is valid in both fluids, and the velocity potential $\varphi _l$ and $\varphi _g$ satisfy the BI equation, yielding

(2.6) \begin{align} {\nabla} ^2\varphi _i=0\quad (i=l,g), \\[-27pt] \nonumber \end{align}

(2.7) \begin{align} \varPi (\boldsymbol{r})\varphi _i (\boldsymbol{r})=\iint _{S}\left [\frac {\partial \varphi _i (\boldsymbol{q})}{\partial {n}} \frac {1}{|\boldsymbol{r}-\boldsymbol{q}|}-\varphi _i (\boldsymbol{q})\frac {\partial }{\partial {n}}\left (\frac {1}{|\boldsymbol{r}-\boldsymbol{q}|}\right )\right ]\,\textrm {d}S_i (\boldsymbol{q})\quad (i=l,g), \\[0pt] \nonumber \end{align}

where $\varPi$ represents the solid angle, $\boldsymbol{r}$ and $\boldsymbol{q}$ are the control and source points, respectively, and $\partial /\partial {n}$ indicates the normal derivative. Here, $S$ refers to the gas–liquid interface and the bubble surface when $i = l$ (liquid domain below the interface), while it refers to the gas–liquid interface only when $i = g$ .

Considering the surface tension and neglecting gas viscosity, the dynamic boundary conditions on the bubble surface and at the gas–liquid interface can be written as

(2.8) \begin{align} \frac {\textrm {D}\varphi _l}{\textrm {D}t}=\frac {P_{\infty }}{\rho _l}-\frac {P_0}{\rho _l}\left (\frac {V_{0}}{V}\right )^{\kappa }+\frac {\sigma \mathcal{C}}{\rho _l}+\frac {1}{2}|\boldsymbol{\nabla }\varphi _l|^{2}-g(z+z_b)\quad \textrm {on the bubble surface,} \\[-28pt] \nonumber \end{align}

(2.9) \begin{align} \frac {\textrm {D}(\varphi _l-\beta \varphi _g)}{\textrm {D}t}=\frac {1}{2}|\boldsymbol{\nabla }\varphi _l|^2+\frac {\beta }{2}|\boldsymbol{\nabla }\varphi _g|^2-\beta \boldsymbol{\nabla }\varphi _l\boldsymbol{\cdot }\boldsymbol{\nabla }\varphi _g-(1-\beta )gz+\frac {\sigma \mathcal{C}}{\rho _l}\quad \textrm {at the interface}, \\[0pt] \nonumber \end{align}

where $\textrm {D}/\textrm {D}t=\partial /\partial t+\boldsymbol{\nabla }\varphi _i\boldsymbol{\cdot }\boldsymbol{\nabla}$ is the material derivative, $P_\infty$ the hydrostatic pressure at $z = 0$ , $P_L$ the liquid pressure on the bubble surface, $P_0$ and $V_0$ the initial pressure and volume of the bubble, $\rho$ the density of the fluid, $g$ the gravitational acceleration, $\mathcal{C}$ the local curvature, and $\beta =\rho _g/\rho _l$ the density ratio.

On all surfaces, the kinematic boundary condition is given by

(2.10) \begin{align} \frac {\textrm {D}{\boldsymbol{r}}}{\textrm {D}t}=\boldsymbol{\nabla }\varphi _l. \end{align}

Using the above-mentioned equations, one can reproduce the interaction between an oscillating bubble and a two-phase interface. The details of this method can be referred to in the studies of Han et al. (Reference Han, Zhang, Tan and Li2022), Li et al. (Reference Li, Zhang and Han2023a ) and Li et al. (Reference Li, Zhao, Zhang and Han2024). One may easily notice that this method neglects the compressibility, vorticity and other aerodynamic effects of the gas, whose influence may become prominent as the air flow becomes rapid with increasing cavity velocity, thus leading to deviations in cavity dynamics. Therefore, the purpose of implementing the BI method in this study is not to precisely reproduce the cavity dynamics, but to elucidate the factors that may have a notable influence on the cavity dynamics by comparing with the FVM simulations, such as air compressibility effects, boundary layer, etc., which cannot be reproduced by BI simulations. We will further elaborate on this issue in §§ 4.2, 5.2 and Appendix C.

Figure 3.

(a) Schematic for the analytical model of the meniscus and (b) comparison between experimental and analytical results of the meniscoid interface. $h_m$ is the distance from the initial undisturbed free surface to the meniscus apex. The analytical results calculated with (2.12) for $h_m=0.95\ \textrm {mm}$ and 0.60 mm are plotted as red and yellow dashed lines, respectively. The control experiment for $h_m=0$ is also given with the free surface plotted as a dashed blue line.

Finally, we employ an analytical model for the surface perturbation, as sketched in figure 3(a). The same cylindrical coordinate as that in figure 2(a) is adopted. The meniscus profile can be calculated by solving the Young–Laplace equation:

(2.11) \begin{align} \frac {1}{R_1}+\frac {1}{R_2}=-\frac {\Delta p}{\sigma }, \end{align}

where $R_1$ and $R_2$ denote the principal radii of curvature, $\Delta p$ the Laplace pressure on the interface, and $\sigma$ the surface tension. Rewriting (2.11) in the cylindrical coordinate (Clanet & Quéré Reference Clanet and Quéré2002) yields the following equation:

(2.12) \begin{align} \frac {\mathrm{d}\phi }{\mathrm{d}s}-\frac {\cos \phi }{r}=\frac {z}{l_{\textit{c}}^2}, \end{align}

where we adopt the expression for the principal radii of curvature $1/R_1+1/R_2=-\mathrm{d}\phi /\mathrm{d}s+\cos \phi /r$ , derived by Bouasse (Reference Bouasse1924). For an arbitrary point on the meniscus surface, $\phi$ represents the tangent angle of this point relative to the axis $O{-}z$ , while $s$ represents the arc length measured from the meniscus tip. The coordinates $(r,z)$ depict the meniscus profile in the $O{-}rz$ plane. By integrating (2.12) with the boundary condition $\phi |_{r=R_s}=\phi _0$ , we can obtain the meniscus profile analytically, where $R_s$ denotes the radius of the thin rod and $\phi _0$ the static contact angle. One may notice that the condition $z|_{r=R_s}=h_m$ is also needed to solve (2.12). Therefore, we adopt a matched asymptotic expansions method to determine $h_m$ , following the approach of James (Reference James1974) and Lo (Reference Lo1983), whose studies suggest that $h_m$ and $\phi _0$ are intrinsically related. Hence, we have

(2.13) \begin{align} \frac {h_m}{R_s}=z_1|_{\hat {r}=1}\ln \xi +z_2|_{\hat {r}=1}+z_3|_{\hat {r}=1}\xi ^2\ln ^2\xi +z_4|_{\hat {r}=1}\xi ^2\ln \xi +z_5|_{\hat {r}=1}\xi ^2, \end{align}

where $\hat {r}=r/R_s$ , $\xi =R_s/l_{\textit{c}}$ and $z_i|_{\hat {r}=1}(i=1,2, {\cdots} ,5)$ are the coefficients corresponding to different orders in the asymptotic expansion depending on the contact angle $\phi _0$ . This equation suggests that once one of $h_m$ or $\phi _0$ is determined, the other one is determined accordingly. The detailed expressions of the coefficients $z_i|_{\hat {r}=1}$ are presented in Appendix D. The full derivation of this method refers to the study of Lo (Reference Lo1983).

One can also notice that in figure 3(b), ghosting effects at the outer rim of the meniscus due to refraction at the water–air interface can be observed. Such optical artefacts may introduce systematic uncertainties in the experimentally measured $h_m$ . To mitigate this effect, we combined experiments with analysis to determine the meniscus height $h_m$ and its theoretical profile. First, by inserting the thin rod into the liquid pool and carefully adjusting its position, we establish a flat free surface and mark its position (blue dashed line in figure 3 b, $h_m=0$ ). Then, the thin rod is slowly withdrawn from the liquid pool and fixed to form a meniscus with a pinned contact line. The perpendicular distance from the contact line to the free surface is the meniscus height $h_m$ . Equations (2.12) and (2.13) show that, for a given solid–liquid pair, every value of $h_m$ corresponds to a unique meniscus profile (Bouasse Reference Bouasse1924; James Reference James1974; Lo Reference Lo1983). We therefore generated the theoretical profile corresponding to each measured $h_m$ using (2.12) and (2.13), and overlaid it on the high-resolution DSLR image. The excellent agreement (shown in figure 3 b, $h_m=0.95$ and 0.60 mm) confirms that optical distortion is negligible; any residual systematic uncertainty in $h_m$ is below 0.01 mm – roughly two pixels at our maximum magnification – and is insignificant relative to the rod radius and the capillary length.

2.3. Non-dimensionalisation and initialisation

In this section, we present the non-dimensionalisation method adopted in this study. We choose the maximum bubble radius $R_m$ , the pressure of the atmosphere on the initial free surface $P_\infty$ and the density of the liquid (water in this study) $\rho _l$ as the basic quantities. Based on these characteristic scales, the non-dimensional key variables are defined as follows:

(2.14a,b) \begin{align} \gamma =\frac {H}{R_m},\quad \varepsilon =\frac {P_{0b}}{P_\infty }, \end{align}

where $\gamma$ is the standoff parameter and $\varepsilon$ the strength parameter. Here, $H$ is the depth of the bubble inception point and $P_{0b}$ is the initial pressure of the gas of the bubble interior. Another significant parameter in our study, the cavity length $h$ and its maximum value $h_{\textit{c}}$ , is defined as the distance from the cavity bottom to the free surface (see figure 2 b). We use this definition due to the refraction effect: the cavity’s upper part is distorted and obscured by the water dome, complicating the measurement of the total cavity length from the opening to the bottom. Apparently, this definition does not apply for the cases when the cavity fails to extend below the free surface. However, our experiments indicate that these cases happen only when $\gamma$ is very small. In such scenarios, the cavity either coalesces with the bubble or the bubble bursts at the free surface (see figures 6 b and 6 c), rendering the cavity length irrelevant. Similarly, $h$ , $h_{\textit{c}}$ and the aforementioned $h_m$ are scaled by $R_m$ . All the variables in this study, unless otherwise annotated, are non-dimensional.

As mentioned in § 2.2, the initial conditions of the bubble are yet to be determined, which includes the initial pressure $P_{0b}$ , density $\rho _{0b}$ and non-dimensional initial radius $R_0$ . According to the best of the authors’ knowledge, no effective models have been developed thus far to describe the initiation stage of an electric spark bubble. Due to the intense Joule heating generated by the electric spark, the surrounding liquid is vapourised, leading to complex phase changes that play a crucial role in bubble evolution, but beyond the scope of our numerical model. Hence, we adopt a simplified method to initialise the bubble, i.e. the bubble is considered as an initially stationary gas bubble with high pressure (Klaseboer et al. Reference Klaseboer, Hung, Wang, Wang, Khoo, Boyce, Debono and Charlier2005a ; Zeng et al. Reference Zeng, Gonzalez-Avila, Dijkink, Koukouvinis, Gavaises and Ohl2018a ; Han et al. Reference Han, Zhang, Tan and Li2022) and the phase change effect is ignored in this study.

Building on these discussions, we first determine the initial pressure and density of the gas inside the bubble, denoted as $P_{0b}$ and $\rho _{0b}$ , respectively. According to previous studies of Li et al. (Reference Li, Khoo, Zhang and Wang2018) and Han et al. (Reference Han, Zhang, Tan and Li2022), the maximum radius and oscillation period of the bubble exhibit relatively reasonable agreement among theoretical, numerical and experimental results when $\varepsilon$ is set between 100 and 200. Therefore, we set $\varepsilon =125$ to achieve a satisfactory match, which yields $P_{0b}=1.27\times 10^7\ \text{Pa}$ and $\rho _{0b}=37.75\ \text{kg}\,\text{m}^{-3}$ based on (2.5), where the initial atmospheric conditions are $P_0=101\,325\ \textrm {Pa}$ and $\rho _0=1.2\ \textrm {kg}\,\textrm {m}^{-3}$ . Next, the non-dimensional initial bubble radius is set as $R_0=0.1623$ to align the numerical simulation with the experimental result according to Plesset (Reference Plesset1949) and Klaseboer et al. (Reference Klaseboer, Hung, Wang, Wang, Khoo, Boyce, Debono and Charlier2005a ), using a modified Rayleigh–Plesset equation.

3.1. Non-coalescence

Figure 4.

Two typical experiments of bubble interaction with a perturbed free surface for large standoff parameters $\gamma$ = 1.80 and 1.43. The radius of the thin rod is 0.35 mm in both experiments. (a) The interaction between the bubble and the perturbed surface is mild; thus, the free surface elevates slightly with the bubble expansion and the cavity maintains a smooth olive shape. A bursting jet occurs when the bubble reaches its minimum volume and rebounds ( $H=30.6\ \textrm {mm}$ , $R_m=17.0\ \textrm {mm}$ , $h_m=0.85\ \textrm {mm}$ ). (b) The cavity is elongated and the bottom exhibits a taper shape as $\gamma$ decreases to 1.43. A brought-forward surface-seal causes a splash when the cavity rebounds ( $H=24.6\ \textrm {mm}$ , $R_m=17.2\ \textrm {mm}$ , $h_m=0.84\ \textrm {mm}$ ). The time scale of each experiment is taken as $T_{\textit {osc}}=R_m\sqrt {\rho _l/P_\infty }$ , which are 1.69 and 1.71 ms, respectively. Two close-up views of the cavity are also presented on the right. The parameters of these two experiments are similar with those of panels (a) and (b): (i) $H=31.1\ \textrm {mm}$ , $R_m=18.0\ \textrm {mm}$ , $h_m=0.84\ \textrm {mm}$ , $\gamma =1.73$ ; (ii) $H=27.1\ \textrm {mm}$ , $R_m=19.0\ \textrm {mm}$ , $h_m=0.84\ \textrm {mm}$ , $\gamma =1.43$ .

In figure 4, two representative experiments are presented, where the free-surface–bubble interaction remains weak due to the large standoff parameter. In the first experiment (see figure 4 a), the standoff parameter is $\gamma =1.80$ . As the bubble is initiated and expands, the free surface elevates mildly (frames 1–5). The fluid motion in the vicinity of the thin rod is hindered due to the no-slip boundary condition on the solid surface. This discrepancy in velocity gives rise to a small conical depression on the fluid surface around the thin rod. As the bubble radius reaches its maximum, the internal gas pressure attains its minimum and the depression grows downward due to the pressure gradient induced by the cavitation bubble (frame 6). The free surface then descends as the bubble contracts and the depression evolves into an olive-shaped cavity with an opening on the free surface (frame 7). Finally, the cavity reaches its maximum length $h_{\textit{c}}$ (frame 8), shortly before the bubble reaches its minimum volume (frame 9). At this point, the bubble’s internal pressure exceeds that of the ambient fluid, thus reversing the pressure gradient. As a result, the cavity bottom rebounds and bursts upwards, exhibiting behaviours similar to a Worthington jet (Worthington & Cole Reference Worthington and Cole1897), due to the combined effects of the pressure wave generated by the bubble when it reaches its minimum volume and the reversed pressure gradient.

One may notice that we use the term ‘pressure wave’ here to describe the pressure response in the fluid field induced by bubble oscillation, whose propagation characteristics may not be resolved experimentally or numerically in this study. However, we primarily focus on the instantaneous response of the cavity and the fluid field to the pressure wave, while neglecting its propagation hysteresis. Such simplification is feasible, since the cavity velocity ( ${\sim} 10^{1}\, \mathrm{m\,s}^{-1}$ ) is two magnitudes smaller than the local sound speed $c$ , suggesting that the propagation of the pressure wave can be regarded as instantaneous – a similar simplification adopted in the study by Peters et al. (Reference Peters, Tagawa, Oudalov, Sun, Prosperetti, Lohse and van der Meer2013b ), Klaseboer et al. (Reference Klaseboer, Turangan, Fong, Liu, Hung and Khoo2006) and Klaseboer et al. (Reference Klaseboer, Fong, Turangan, Khoo, Szeri, Calvisi, Sankin and Zhong2007), and proven effective, where the effect of the pressure wave on the meniscoid surface is modelled as an instantaneous pressure pulse of characteristic amplitude and duration. Although the pressure wave cannot be resolved experimentally, its influence on the cavity dynamics and fluid field can still be partially reproduced in numerical simulations, which will be further discussed in § 4.3.

In the second experiment (see figure 4 b), the standoff parameter is reduced to $\gamma =1.43$ , resulting in stronger free-surface–bubble interaction. As a result, the cavity expands more rapidly once the bubble reaches its maximum radius (frames 1 and 2). In contrast to the previous case, the entire cavity is elongated and the bottom exhibits a tapered shape (frame 3), thus reaching a larger length. This profound difference indicates that the flow-focusing effect caused by the contraction of the bubble may contribute to different cavity dynamics. Additionally, an elongated cavity causes an increasing amount of air influx, resulting in a lower local pressure at the cavity opening. The resulting pressure imbalance triggers the surface-seal prematurely, even before the cavity reaches its maximum length, thereby isolating it from the atmosphere. Shortly after the cavity reaches its maximum length, the cavity rebounds rapidly due to the combined effects of the pressure wave emitted by the bubble and hydrostatic pressure. During this upward motion, the cavity bottom collides with the sealed surface, causing a sealed-splash event (frame 5).

Wrinkles on the cavity surface can also be observed in figure 4(a), frames 7 and 8. Moreover, the wrinkles grow markedly more pronounced and asymmetric as $\gamma$ decreases in figure 4(b). The largest-amplitude wrinkle structures emerge on the lower cavity wall during the late-time evolution. Since the cavity details are difficult to observe in panels (a) and (b), we present close-up views of the cavity in another two experiments with similar $\gamma$ , i.e. (i) $\gamma =1.73$ , (ii) $\gamma =1.43$ , as shown in the images on the right. Such wrinkles may derive from the late-time secondary R–T instability evolution. Earlier work (Ramaprabhu, Karkhanis & Lawrie Reference Ramaprabhu, Karkhanis and Lawrie2013) has demonstrated that complex acceleration histories – especially repeated transitions between acceleration and deceleration – promote small-scale structure formation, a behaviour consistent with the cavity acceleration and rebound process discussed later in figure 7, § 4.1. Additionally, the nonlinear late-time evolution of the primary bubble itself seeds minor bubbles and wrinkles (Ramaprabhu et al. Reference Ramaprabhu, Dimonte, Woodward, Fryer, Rockefeller, Muthuraman, Lin and Jayaraj2012). These interpretations align with our observations. Reducing $\gamma$ strengthens the nonlinear interactions between the free surface and the cavitation bubble, thereby amplifying the secondary wrinkles (figures 4i and 4ii).

To provide enlightenment for subsequent investigation of the governing parameters of cavity dynamics, we hereby perform a non-dimensional analysis. Since the cavity evolution and bubble collapse stage occur simultaneously, their time scales are similar, given by $T_{\textit {osc}}=R_m\sqrt {\rho _l/P_\infty }$ . Moreover, in figure 4(a), $h_{\textit{c}}=7.0\ \textrm {mm}$ , $U_{\textit{c}}=6.3$ $\mathrm{m\,s}^{-1}$ , $R_m=17.0$ mm, suggesting that both the maximum cavity length and velocity is comparable to those of the bubble, i.e. $h_{\textit{c}}/R_m\sim \textit {O}(10^0)$ , $U_{\textit{c}}/U_b\sim \textit {O}(10^0)$ , where $U_b=\sqrt {P_\infty /\rho _l}=10.1\ \rm {m\,s}^{-1}$ denotes the characteristic velocity of the bubble oscillation and $U_{\textit{c}}$ the maximum cavity velocity. Consequently, the Reynolds number, Weber number and Froude number of the cavity evolution can be written as

(3.1) \begin{align} \textit {We} = \frac {\rho _l U_{\textit{c}}^2 h_{\textit{c}}}{\sigma } \sim \textit {O}\big (10^3\big ),\quad \textit {Re}= \frac {\rho _l U_{\textit{c}} h_{\textit{c}}}{\mu _l} \sim \textit {O}\big (10^5\big ),\quad \textit {Fr}=\frac {U_{\textit{c}}}{\sqrt {g h_{\textit{c}}}} \sim \textit {O}\big (10^1{-}10^2\big ), \end{align}

where we set liquid viscosity as $\mu _l=10^{-3}\ \textrm {Pa}\,\textrm {s}$ , density as $\rho _l=997\ \textrm {kg}\,\textrm {m}^{-3}$ and surface tension as $\sigma =0.073\ \textrm {N}\,\textrm {m}^{-1}$ . This indicates that surface tension, viscosity and gravity effects are all negligible, while the inertia force is the primary governing factor in cavity evolution. One may also notice that, even though $h_{\textit{c}}/R_m$ decreases to the order of $\textit {O}(10^{-1})$ in some cases with relatively larger $\gamma$ , the non-dimensional parameters still have the order of $\textit {We}\sim {O}(10^2)$ , $\textit {Re}\sim {O}(10^4)$ and $\textit {Fr}\sim {O}(10^1)$ , suggesting that the aforementioned simplifications are still feasible. This leads to consistent cavity dynamics over a relatively wide range of $\gamma$ , which we will discuss in § 5.2. Moreover, we notice that the cavity evolution is highly related to the bubble collapse stage, i.e. the small depression does not grow downwards at the beginning. Instead, it rises with the water hump as a whole due to the rapid bubble expansion stage, then the cavity evolution is initiated after the bubble reaches its maximum radius and starts to contract. This indicates that the pressure gradient induced by the bubble oscillation is the dominant factor that drives the cavity evolution. This dependence will be examined in greater detail in §§ 4 and 5.

3.2. Critical condition

This section examines the critical condition of cavity–bubble interactions. As $\gamma$ decreases into a transitional region, the system shifts from the non-coalescence state to one exhibiting increasingly complex and ambiguous interaction patterns. Despite these complex dynamics, we identify (i) a critical threshold of $\gamma$ , and (ii) a classification framework to differentiate the critical conditions from the coalescence and non-coalescence cases.

Figure 5.

Four representative experiments of critical condition for $\gamma$ = 1.36, 1.38, 1.36 and 1.34. The radius of the thin rod is 0.35 mm in all four experiments. (a) The cavity is highly elongated and forms a tapered shape when $\gamma =1.36$ . The cavity bottom surpasses the bubble’s upper surface and becomes invisible. After pinch-off, the cavity splits into two segments: the upper part forms an isolated bubble, while the lower part merges with the primary bubble, creating a roughened bubble surface ( $H=25.2\ \textrm {mm}$ , $R_m=18.5\ \textrm {mm}$ , $h_m=0.91\ \textrm {mm}$ ). Additionally, feathery bright spots (possibly caused by continuous discharge of the wire) can be observed in frames 1 and 2. Panels (b)–(d) show distinct cavity–bubble interaction patterns as $\gamma$ varies slightly around 1.36, with emphasis on the dynamics just before cavity pinch-off. (b) Though the cavity is elongated, it does not reach deep enough to coalesce with the bubble. The bubble retains a smooth surface and strong luminescence occurs at collapses (see the inset in figure 5 b) ( $H=25.0\ \textrm {mm}$ , $R_m=18.1\ \textrm {mm}$ , $h_m=0.93\ \textrm {mm}$ ). (c) The cavity bottom forms a bulge inside the bubble. The primary bubble merges with the lower part of the cavity after pinch-off, forming a roughened surface (see the inset in figure 5 c) ( $H=25.5\ \textrm {mm}$ , $R_m=18.7\ \textrm {mm}$ , $h_m=0.90\ \textrm {mm}$ ). (d) The cavity coalesces with the bubble, creating splashes and a microbubble cloud around the lower bubble surface (see the inset in figure 5 d) ( $H=24.6\ \textrm {mm}$ , $R_m=18.3\ \textrm {mm}$ , $h_m=0.88\ \textrm {mm}$ ). The time scales are 1.84, 1.80, 1.86 and 1.82 ms, respectively.

Four representative experiments are presented in figure 5. Since $h_{\textit{c}}$ varies slightly with $h_m$ , we focus on a subset of experiments where $h_m\approx 0.90\ \textrm {mm}$ and perform over 100 experiments, with $\gamma$ varying from 0.49 to 3.58. It turns out that there exists a certain critical area for the cavity dynamics and the critical condition occurs when $\gamma =1.36\pm 0.02$ . In figure 5(a), a typical critical case where $\gamma =1.36$ is given. The cavity is further elongated, exhibiting the same taper shape as that in figure 4(b) (frame 1). Then, the cavity bottom surpasses the upper surface of the bubble (frames 2 and 3), rendering it invisible and preventing accurate cavity length measurements. After pinch-off, the cavity splits into two segments: the upper part detaches to form an isolated bubble, while the lower part merges with the primary bubble, creating a roughened bottom surface (frames 4 and 5). Unfortunately, direct visualisation of the cavity through the bubble surface is not feasible with the current experimental set-up, preventing clear confirmation of the cavity–bubble coalescence. Nonetheless, we can still infer the interaction pattern from post-collapse bubble morphology. The classification criteria are elaborated in figure 5(b–d), focusing on the cavity–bubble interactions immediately before cavity pinch-off.

In the second experiment (see figure 5 b), $\gamma =1.38$ , the cavity bottom surpasses the horizontal plane of the bubble’s upper surface. However, there is no observable evidence that the cavity comes into contact with or merges into the bubble. Moreover, the bubble surface remains smooth until it collapses and the bubble luminescence is still visible (see the inset in figure 5 b). Therefore, the cavity–bubble interaction pattern is almost the same as that in figure 4(b).

It is worth emphasising that the light emission of a spark-generated bubble is, in fact, a convoluted phenomenon. Earlier work has ascribed it to a number of concurrent mechanisms, including extreme temperature and pressure (Huang et al. Reference Huang, Zhang, Chen, Zhu, Liu and Yan2015; Zhang et al. Reference Zhang, Zhu, Yan, Huang, Liu and Yan2017), plasma formation (Vokurka & Plocek Reference Vokurka and Plocek2013; Stelmashuk et al. Reference Stelmashuk, Prukner, Kolacek, Tuholukov, Hoffer, Straus, Frolov and Jirasek2022), among others. Similar luminescence has been reported for non-spherical bubbles collapsing near either a free surface or a rigid boundary (Ohl et al. Reference Ohl, Wu, Klaseboer and Khoo2015; Ma et al. Reference Ma, Huang, Zhao, Wang, Chang, Qiu, Fu and Wang2018). Experimental observations also reveal feathery bright spots (as shown in figure 4 a, frame 7 and figure 5 a, frames 1–3) that emanate directly from the copper-wire tips. These filaments probably originate from a sustained micro-spark between the electrodes and could still be present when the bubble reaches minimum volume, thereby contributing to the recorded luminescence. At present, however, we do not know exactly how these structures form or how much of the total luminescence they comprise. Consequently, we cannot unambiguously attribute the observed luminescence to bubble collapse alone, to the electrode discharge or to a combination of both. What we can state with confidence is that a pronounced luminescence appears in every non-coalescence case and is completely absent whenever coalescence occurs (figures 4–6). The most plausible explanation is that coalescence-driven droplets roughen the bubble wall; the resulting multiple scattering scatters most of the internal light, regardless of its physical origin. Hence, the presence or absence of luminescence provides a robust, experimentally accessible marker for distinguishing the two cavity–bubble interaction regimes.

In the third experiment (see figure 5 c), $\gamma =1.36$ , a bulge caused by the extending cavity inside the bubble can be observed. Though the existence of the bulge implies that the cavity is close to the bubble, there remains no definitive evidence of coalescence before the cavity neck pinch-off. After the pinch-off, the lower segment detaches from the cavity and merges with the bubble. The mixture of air and droplets collides on the bubble surface, forming some minor splashes (see the inset in figure 5 c).

In the fourth experiment (see figure 5 d), $\gamma =1.34$ , coalescence undoubtedly occurs before the cavity pinch-off, as evidenced by the clearly visible splashes resulting from droplet collisions on the bubble surface.

Building on the previous discussions, we propose a criterion for classifying these three scenarios: (1) if the cavity does not coalesce with the bubble throughout its entire evolution and luminescence during bubble collapse is clearly visible, the case is classified as non-coalescence (figure 5 b); (2) if the cavity does not coalesce with the bubble before cavity neck pinch-off, but coalescence occurs afterward, the case is classified as the critical condition (figure 5 c); (3) if the cavity has coalesced with the bubble before cavity neck pinch-off and luminescence is hardly visible, the case is classified as coalescence (figure 5 d).

As demonstrated in the experiments discussed previously, the cavity–bubble interaction patterns differ significantly when $\gamma$ varies by 0.02 around the critical value, which may indicate that $\gamma$ is the dominant factor governing cavity dynamics. It is worth emphasising that this drastic change between different cavity–bubble interaction patterns only occurs around the critical value, while in the rest of the $\gamma$ region, the cavity characteristics exhibit rather good consistency. We will further investigate the critical condition through numerical simulations in § 5.

3.3. Coalescence

In this section, we present and discuss three typical experiments for cavity–bubble coalescence in figure 6. A coalescence case occurring at the end of the first oscillation cycle has already been discussed in figure 5(d). Despite the onset of coalescence, the bubble undergoes a relatively regular collapse, with the absence of a roughened surface and splashes. Thus, the coalescence has a limited effect on the overall bubble collapse dynamics in this case. As $\gamma$ decreases further below 1.36, the coalescence is brought forward, and the ventilation between the bubble interior and the atmosphere is initiated.

Figure 6.

Three typical experiments of cavity–bubble coalescence scenario for $\gamma$ = 1.30, 0.71 and 0.49. The radius of the thin rod is 0.35 mm in all three experiments. (a) Cavity expands rapidly and coalesces with the bubble before the surface-seal, creating a channel that connects the bubble and the atmosphere, and triggering the ventilation, allowing the air to flow into the bubble ( $H=22.8\ \textrm {mm},\ R_m=17.5\ \textrm {mm},\ h_m=0.83\ \textrm {mm}$ ). (b) The decreased $\gamma$ gives rise to an earlier coalescence and causes the water dome to elevate rapidly. As a result, the cavity neck pinches off quickly and detaches from the bubble, which ceases the ventilation. Despite the advanced pinch-off, an adequate amount of air influx balances the pressure difference between the gas pressure of the bubble interior and the ambient hydrostatic pressure ( $H=11.9\ \textrm {mm},\ R_m=16.8\ \textrm {mm},\ h_m=0.90\ \textrm {mm}$ ). (c) The bubble is so close to the free surface that the surface perturbation has little influence on the coalescence and ventilation. The bubble bursts into the atmosphere and creates an open cavity ( $H=7.9\ \textrm {mm},\ R_m=16.0\ \textrm {mm},\ h_m=0.92\ \textrm {mm}$ ). The time scales are 1.74, 1.67 and 1.59 ms, respectively.

In the first experiment shown in figure 6(a), $\gamma =1.30$ , the interaction between the bubble and the free surface is intensified, resulting in a faster downward growth of the cavity after the bubble reaches its maximum radius (frames 4–5). Such rapid growth, combined with the decreased standoff parameter, allows the cavity to coalesce with the bubble earlier than that in figure 5. Meanwhile, the cavity opening remains unsealed, allowing ambient air to flow into the bubble at high velocity through the channel. This establishes continuous ventilation between the bubble and the atmosphere above the surface until the surface-seal finally occurs (frames 6 and 7). Although the ventilation lasts for only 0.24 ms, the intense air influx is sufficient to deliver a significant amount of air into the bubble. This helps balance the pressure difference between the gas pressure of the bubble interior and the ambient hydrostatic pressure, despite the occurrence of surface-seal, until the cavity eventually pinches off (frame 8). Consequently, the bubble collapses and rebounds with a rough and rugged appearance (frames 9 and 10), in contrast to the regular-shaped bubbles observed in figure 4. A daughter bubble is also observed between the collapsing bubble and the free surface, which comes from the segment separated from the cavity following surface-seal and pinch-off. Notably, no evidence of secondary cavitation derived from the refraction wave reflected by the free surface is observed in the fluid field, suggesting that the amplitude of the collapse-induced pressure wave may have been attenuated due to the ventilation. Given that small cavitation bubbles might not be resolved by the high-speed camera due to large exposure time, we review the complete data set and reveal the same pattern in every non-coalescence case. The cavitation cloud persists for at least two frames, i.e. longer than the interframe time 27.02 $\unicode{x03BC}$ s. This duration is an order of magnitude above the exposure time (1 $\unicode{x03BC}$ s), indicating that the camera is sufficient to resolve secondary cavitation bubbles. More detailed dynamics of such microbubbles/nanobubbles is available from Rosselló & Ohl (Reference Rosselló and Ohl2021). Though the presence or absence of cavitation alone cannot serve as definitive evidence for pressure attenuation, this observation offers valuable insight into how the cavity influences the bubble dynamics. It also motivates further investigation into the underlying cavity behaviour. The attenuation of pressure amplitude will be examined in detail in the upcoming section.

In the second experiment shown in figure 6(b), $\gamma =0.71$ , the water dome rises higher and faster due to the enhanced interaction between the bubble and the free surface (frame 1). Meanwhile, the decreased standoff parameter and enhanced cavity–bubble interaction lead to coalescence and ventilation significantly earlier than in any of the previous experiments (frame 2). However, after the onset of coalescence and ventilation, the combined effect of hydrostatic pressure, aerodynamic pressure and the continued ascent of the water dome pinches the cavity neck, detaching the cavity from the bubble, which terminates the ventilation almost immediately. Driven by inertia, the water dome continues to rise rapidly and forms a primary jet that envelops the cavity. The bubble then starts to contract and diminishes the volume shortly after the coalescence and ventilation (frame 3). Apparently, the enormous influx of air substantially alters the internal environment of the bubble, driving the pressure across the bubble wall towards equilibrium. The weakened pressure difference slows down the contraction of the bubble, resulting in a weaker collapse and attenuated pressure wave (frame 4). The bubble and cavity both exhibit a vague and misty appearance upon and after the collapse, which is completely different from the previous experiments, where phenomena such as bursting jet (figure 4 a), collapse luminescence (figure 5 b) and splashes caused by droplets (figures 5 c and 5 d) can still be recognised. As in the cases of figure 5, the majority of the cavity and the top of the bubble are invisible due to the obstruction by the water dome, so it is difficult to discuss the cavity dynamics with experimental results. Therefore, we will numerically elaborate on this case in § 4.

In the last experiment (see figure 6 c), $\gamma$ is reduced to 0.49, where the influence of perturbation on the interaction between the bubble and the free surface becomes negligible. In this case, the standoff parameter is sufficiently small ( $\gamma \lesssim 0.5$ ) that the bubble can burst directly into the atmosphere, leaving behind an open cavity. This may give rise to a series of rich dynamic processes, such as air entry and sealed splash (Wang et al. Reference Wang, Wang and Wang2024). Under such circumstances, the influence of the surface perturbation becomes insignificant, without whose existence the bubble will burst and ventilate nonetheless. Previous studies such as Wang, Duan & Wang (Reference Wang, Duan and Wang2015) and Li et al. (Reference Li, Zhang, Wang, Li and Liu2019) have thoroughly discussed this phenomenon; thus, we will not delve into this scenario in the present study.

4.1. Non-coalescence

Figure 7 presents a comparison between the experiment shown in figure 4(a) and the corresponding numerical simulation results. The pressure field (coloured contour faces) and interface (magenta lines) obtained from the FVM simulation and the BI simulation results (blue lines) are given. The FVM simulation results exhibit good agreement with the experimental observations, which further validates our numerical model. Frame 1 illustrates the bubble-free-surface interaction when the bubble reaches its maximum radius. A small conical depression can be observed as the water dome starts to rise. Once the bubble starts to contract, the cavity grows downward, driven by the pressure gradient (frame 2) and reaches its maximum length $h_{\textit{c}}$ (frame 3). Due to the large $\gamma$ , the cavity can barely exert a notable influence on the bubble. As a result, the bubble collapses with a relatively undisturbed shape, generating a pressure peak of approximately 100 in the fluid field (frame 4). Apparently, the bubble has not yet collapsed to the minimum volume when the cavity reaches its maximum length and there exists a hysteresis between these two crucial moments. In contrast, the BI simulation results exhibit a slight discrepancy with both experimental observations and numerical results in the early stage of the cavity evolution (frames 1 and 2). Nonetheless, both numerical simulations successfully reproduce the overall characteristics of the cavity.

Figure 7.

Comparison of cavity and bubble evolution between numerical simulation and experiment for the same case as in figure 4(a). The result obtained from the FVM simulation is presented in pressure contours and magenta solid lines in the right half, while the BI simulation results is plotted with blue solid lines in the left half, along with the experimental observations. The non-dimensional parameters in this simulation are set as $\gamma =1.80$ , $R_0=0.1623$ , $\varepsilon =125$ , and $\kappa =1.4$ . All the non-dimensional parameters for the numerical simulations in § 4 are the same except for $\gamma$ . The time scale is 1.69 ms.

Figure 8.

Comparison between simulation results and experimental data for the evolution of (a) cavity length $h$ and (b) cavity velocity $U$ when $\gamma =1.77$ , $\gamma =1.45$ and $\gamma =1.39$ . Solid lines represent simulation results and circular markers denote experimental data. An inflection point is observed in both the $h$ and $U$ curves for each case, indicating a deceleration and rebound of the cavity during the final stage of bubble collapse. In panel (a), the evolution of the pressure difference between the bubble and atmosphere is also plotted (dotted lines), scaled by $P_\infty$ . The sharp drop in pressure difference coinciding with the cavity rebound suggests that the pressure gradient is a key factor in the cavity dynamics. Cavity velocity $U$ is scaled by $\sqrt {P_\infty /\rho _l}$ . The experimental uncertainties of $h$ and $U$ are approximately 0.004 and 0.36 (0.076 mm and $3.6\ \mathrm{m\,s}^{-1}$ in dimensional form), respectively.

To further discuss this phenomenon, we compare three simulation results with corresponding experimental data for the evolution of non-dimensional cavity length $h$ and velocity $U$ before the bubble reaches its minimum volume, as shown in figure 8. The standoff parameter of the three cases are $\gamma =1.77$ , $\gamma =1.45$ and $\gamma =1.39$ , respectively. The cavity velocity $U$ is determined using a five-point local temporal fitting procedure: the cavity bottom position is extracted from each frame and $U$ is obtained from the local slope of a five-frame linear fit centred at each time instant. Due to the limited spatial resolution of the high-speed camera, the cavity bottom position exhibits a spatial uncertainty of one pixel (0.076 mm). Although $R_m$ varies from different experiments (ranges from 17 to 20 mm), the non-dimensional uncertainty can be estimated to ${\sim} 0.004$ . Consequently, the uncertainty in $U$ is estimated to be $3.6\,\mathrm{m\,s}^{-1}$ ( $0.36$ non-dimensional). To facilitate our discussion, the cavity velocity is taken as positive when directed towards the bubble. We use the ratio between the atmosphere–bubble pressure difference and $\gamma$ , denoted by $(1-P_b)/(\gamma -h)$ , as a qualitative proxy for the real pressure gradient, since the direction of the pressure gradient changes synchronously with the sign of $1-P_b$ , where $P_b$ is scaled by $P_\infty$ . This simplified approach helps us analyse the effect of the pressure gradient on cavity dynamics qualitatively. A quantitative analysis of this dependence will be presented in § 5. The simulation result of the cavity length (see figure 8 a) and velocity (see figure 8 b) agree well with the experimental data in the final stage of the bubble collapse, and are slightly later than the experimental data in the early stage as we expected. This discrepancy in $h$ is attributed to the continuous electric discharge discussed in § 2.3, which causes the bubble to grow more slowly than in the simulations. As a result, the conical depression around the thin rod forms with a slightly lower altitude and velocity, triggering an earlier cavity downward growth during bubble contraction.

Upon cavity initiation, $h$ and $U$ accumulate gradually, while the pressure difference $(1-P_b)/(\gamma -h)$ remains relatively low (approximately 0.5–0.6). However, as the bubble keeps contracting, the internal pressure of the bubble increases rapidly from a minimum value to exceeding $P_\infty$ within a transient period (approximately 0.5 ms according to numerical simulations). The increasing internal pressure lowers and finally reverses the pressure gradient between the bubble and the atmosphere. Consequently, this reversed gradient causes the cavity to decelerate, eventually bringing its velocity to zero and triggering a rebound. The deceleration and rebound process occurs almost simultaneously with the plunge of the pressure difference and takes place before the bubble reaches its minimum volume, which explains the aforementioned hysteresis in the cavity–bubble interaction. This hysteresis, along with the close synchronisation between the cavity evolution and the pressure gradient changes, suggests that the pressure gradient between the atmosphere and the bubble interior governs the cavity evolution.

We also notice that both consistency and inconsistency exist in the cavity evolution shown in figure 8. On the one hand, the rebound velocity of the cavity increases as $\gamma$ decreases, which aligns with the experimental observations in figure 4: when $\gamma$ decreases, the cavity–bubble interaction becomes stronger, resulting in a larger $h_{\textit{c}}$ and a slender cavity shape with a tapered bottom (see figure 4 b, frame 4). The tapered bottom of the cavity creates a high local curvature, which leads to an earlier, faster and thinner bursting jet than that in figure 4(a) (Terasaki et al. Reference Terasaki, Kiyama, Kang, Tomita and Sato2024). On the other hand, though the rebound of the cavity occurs almost synchronously with the plunge of the pressure difference in the case of $\gamma$ = 1.77 and 1.45, a slight temporal lag of cavity rebound can still be observed in the case of $\gamma$ = 1.39. This deviation may be attributed to the close distance between the cavity and the bubble at a smaller $\gamma$ . In this case, the cavity evolution is no longer governed solely by the pressure gradient; flow-focusing effects become significant during the final collapse stage, which opposes the reversed pressure gradient, stretching the cavity bottom and intensifying the nonlinear interaction with the bubble. As a result, the rebound process is postponed.

4.2. Critical condition

Figure 9 shows the comparison among the results obtained with FVM, BI simulation and the experimental observations of the critical condition case in figure 5(a). The simulation enables a more detailed examination of the cavity–bubble interaction characteristics, which we fail to have a clear sight of in the experiments. As can be seen in figure 9 (frame 2), the upper part of the bubble is depressed inward deeper due to the existence of the cavity, compared with the BI simulation results, and forms an internal bulge. As the bubble keeps contracting, the cavity elongates downwards and the bulge keeps developing (frame 3). Different from the cases in figure 7 where the cavity is far from the bubble, the cavity stretches downwards through the high-pressure zone above the bubble in the case of the critical condition. Affected by the high local pressure, the cavity starts to shrink in the middle. This explains why the cavity is squeezed and forms a thin neck near the bubble instead of elsewhere. After pinch-off, the localised high pressure near the tips initiates a rebound motion of the segments in opposite directions, producing a bursting jet reminiscent of a Worthington jet (Worthington & Cole Reference Worthington and Cole1897). The upper segment evolves into an isolated bubble that rises and bursts upwards, while the lower one merges with the primary bubble and collapses, producing a collapse pressure of $P\approx$ 100 (frame 4). The neck pinch-off also creates rapid jets pointing oppositely inside the cavity. According to numerical simulations, the maximum jet tip velocity can exceed 300 $\rm {m\,s}^{-1}$ . Owing to the roughened cavity surface caused by droplet impact and the mist-like appearance of the jet tip (frame 4), a direct, accurate measurement of the jet velocity after pinch-off is not feasible. As a practical alternative, we track the cavity bottom point between the frames immediately at (frame 3) and after (frame 4) pinch-off, yielding an experimental cavity bottom velocity of 184.5 $\mathrm{m\,s}^{-1}$ . The same procedure applied to the numerical data gives 186.0 $\mathrm{m\,s}^{-1}$ , confirming consistency. Clearly, the actual jet tip moves faster than either value and the procedure deliberately underestimates the true jet speed.

Figure 9.

Comparison of cavity and bubble evolution between numerical simulation and experiment for the same case as in figure 5(a). Simulation results obtained from FVM (magenta lines) and the BI method (blue lines) are both plotted. FVM simulation results exhibit better agreement with the experimental observations than the BI method. The standoff parameter is set as $\gamma =1.36$ . The time scale is 1.84 ms.

As discussed in § 3.2, although coalescence between the lower segment of the cavity and the bubble occurs in the critical condition, cavitation induced by the rarefaction wave can still be observed in the fluid field (see figure 5 a, frame 5). Consequently, the peak pressure produced by bubble collapse remains comparable to that in figure 7, suggesting that the coalescence between a small segment of the cavity and the bubble alone does not significantly affect the collapse intensity. Therefore, stable ventilation that sustains for a sufficiently long duration of time is essential to achieve a significant attenuation in collapse pressure, which will be further examined in the next section.

Additionally, numerical results obtained with the BI simulation are also brought into comparison (plotted with blue solid lines on the left half). Within our expectations, the BI simulation results exhibit vast deviation from both experimental observations and FVM simulation results with respect to cavity characteristics in this scenario, while the bubble behaviour has no major difference except for the internal bulge caused by the cavity growth. This discrepancy derives from several aspects. First, the BI model is incapable of simulating the boundary layer around the thin rod, whose influence may not be significant, but still can affect the cavity dynamics to a certain extent (the influence of the boundary layer is deferred to § 5.4). Second, which is the most significant, the air phase in a typical two-phase BI simulation is treated as an incompressible and inviscid ideal fluid. Such assumptions are feasible in the general free-surface–bubble interactions, since the air velocity above the free surface is not that high. However, when it comes to the critical condition case, the cavity grows rapidly and induces rapid air flow into the cavity, subject to the pressure difference. According to FVM simulation results, the velocity of the air flow inside the cavity can reach approximately 200–280 $\mathrm{m\,s}^{-1}$ ( $\textit {Ma}\approx$ 0.6–0.8), and maximises to approximately 300 $\mathrm{m\,s}^{-1}$ at the cavity opening. Due to the high air velocity, the gas density and pressure inside the cavity drop by approximately 25 %–30 % compared with those of the atmosphere. This pressure plunge creates another gradient between the cavity interior and the atmosphere, thereby pumping more air into the cavity and elongating it. Such rapid air inflow may also lead to choked flow at the cavity opening, resulting in rich aerodynamics (Lee, Longoria & Wilson Reference Lee, Longoria and Wilson1997; Gekle et al. Reference Gekle, Peters, Gordillo, van der Meer and Lohse2010). However, our numerical simulations show maximum gas velocities below the sonic threshold ( ${\sim} 280\ \mathrm{m\,s^{-1}}$ , corresponding to $\textit {Ma}\approx 0.8$ at ambient conditions), as current mesh resolution is insufficient to fully capture the aerodynamics at the moment of surface seal. All the above-mentioned evidence suggests that the compressibility effects of the air become prominent and cannot be neglected under such circumstances, which play a significant role in the cavity evolution (Peters et al. Reference Peters, Gekle, Lohse and van der Meer2013a ; Kuang et al. Reference Kuang, Han, Zhang and Li2025). However, this exceeds the capability of the BI model, thereby introducing deviations in cavity dynamics.

4.3. Coalescence

When $\gamma \lesssim 1.36$ , the coalescence and ventilation between the cavity and the bubble become inevitable, giving rise to a series of complex cavity–bubble interaction dynamics. Figure 10 presents a comparison between the simulation results and the experimental observations in figure 6(a) at several representative moments, to illustrate the entire coalescence-ventilation process.

The numerical simulation enables us to have a clearer vision of the cavity–bubble coalescence process. The rapidly growing cavity pushes aside the fluid between the cavity and the bubble. Then, the cavity bottom coalesces with the bubble, establishing a channel between the bubble and the atmosphere. At this critical moment, the pressure inside the bubble is significantly lower than that of the atmosphere. The pressure gradient pumps air into the bubble (frames 3–4). This causes prominent air compressibility effects inside the cavity, including pressure and density plunge, as we discussed in the last section. Then, the pressure difference between the cavity interior and the atmosphere triggers the surface-seal and terminates the ventilation almost immediately, shortly after the coalescence starts (frame 5). As a result, the entire ventilation process is transient and only lasts for only about 1.05 ms (frames 3–5) according to FVM simulation results. However, the transient ventilation process still leads to pressure increase inside the bubble, which can be observed from frames 2–5. Therefore, it is reasonable to suggest that the increased bubble pressure attenuates the collapse intensity, thereby eliminating secondary cavitation across the fluid field after the bubble reaches its minimum volume. After surface-sealing terminates the ventilation process, cavity–bubble coalescence persists (frame 6) until the cavity pinches off in the middle (frame 7), forming two cavity segments. Although the simulation results deviate slightly from the experiment observations in terms of cavity and bubble morphology after the cavity pinch-off, the overall dynamics of cavity–bubble interactions are well captured.

Figure 10.

Comparison of cavity and bubble evolution between numerical simulation and experiment for the same case as in figure 6(a). Frames 3, 5 and 7 show the three critical moments: (I) coalescence initiation, (II) surface-seal and (III) cavity pinch-off. The standoff parameter is set as $\gamma =1.30$ . The time scale is 1.74 ms.

As $\gamma$ keeps decreasing to 0.72, the coalescence and ventilation occur even earlier and more intensely. As shown in figure 11, the simulation results agree well with the experimental observation on the morphology of the bubble and cavity. As we expected, the coalescence occurs even earlier than the previous case due to small $\gamma$ (frame 1), before the bubble starts to contract. As a result, a more stable channel is established between the bubble and the atmosphere (frame 2), lasting much longer than the case in figure 10. The elongated ventilation allows a large amount of air and droplets influx into the bubble, resulting in a misty bubble shape during the collapse stage (frame 4), instead of a transient collapse with strong luminescence.

Figure 11.

Comparison of cavity and bubble evolution between numerical simulation and experiment for the same case as in figure 6(b). The standoff parameter is set as $\gamma =0.72$ . The time scale is 1.67 ms.

Figure 12.

Comparison of the fluid field pressure between the simulations with (red circles) and without (blue diamonds) surface perturbation. A pressure monitoring point is located horizontally at a distance $R_m$ from the bubble inception location. (a) Critical condition in figure 9. No significant attenuation can be observed. (b) An evident attenuation in the fluid pressure field is observed in the case of figure 10, suggesting that coalescence and ventilation can attenuate the collapse intensity. (c) The collapse intensity is further attenuated by the coalescence and ventilation in figure 11.

One may also notice that the maximal $P$ when the bubble reaches its minimum volume in figures 10 and 11 are attenuated to approximately 50 and 10, respectively, compared with $P\approx$ 100 in figures 7 and 9. This attenuation originates from the ventilation process based on the previous discussions. To justify this argument, we compare the numerical results with and without surface perturbations by examining the fluid field pressure of a monitoring point positioned horizontally to the bubble inception point at a distance of 2 $R_m$ . Three groups of comparisons are presented in figure 12, corresponding to the critical condition in § 4.2 and two coalescence scenarios in this section. As shown in figure 12, when only the coalescence between cavity and bubble occurs, the fluid field pressure in the critical condition is barely attenuated, differing by 2.1 % compared with the case without surface perturbation. In contrast, the peak pressure is attenuated to approximately 70 % in the two cases where ventilation between the bubble and the atmosphere takes place. These comparisons vividly demonstrate that it is the bubble-atmosphere ventilation that causes the collapse attenuation.

5.1. Overall cavity dynamics

Figure 13 depicts the dependence of different cavity behaviours on $\gamma$ and $h_m$ . The simulation results obtained with FVM is plotted in blue contours. It can be clearly seen that $h_{\textit{c}}$ is positively correlated with $h_m$ , while negatively with $\gamma$ . Different cavity behaviours in the experimental observations are plotted with circles, diamonds and triangles, representing cavity–bubble interaction scenarios of non-coalescence, critical condition and coalescence. The experimental circles share the same colour gradient as the contours. It is worth emphasising that the colour distribution of the experimental results is almost the same as the simulation results, and so is the distribution of different cavity behaviours in the parameter map, suggesting that the experimental results of the coalescence cases agree well with the numerical simulation, which further demonstrates the effectiveness of our numerical model.

Figure 13.

Dependence of cavity length on governing parameters of $\gamma$ and $h_m$ . The cavity lengths obtained with numerical simulation are presented with the blue contours. The $\gamma$ value of critical conditions are plotted with the magenta line. The experimental observations of different cavity behaviours are plotted with circles (non-coalescence), diamonds (critical condition) and triangles (coalescence). The experimental uncertainty of $h_{\textit{c}}$ is approximately 0.004 (0.076 mm in dimensional form). The corresponding experimental images for the three different cavity–bubble interaction scenarios are presented on the right.

Moreover, the critical conditions obtained with FVM simulations are presented with a magenta line. As we expected, the critical $\gamma$ is also correlated positively with $h_m$ . One may easily notice that there exists a gap between the non-coalescence contours and the critical line. This gap corresponds to the case in figure 5(b), where the cavity stretches and surpasses the bubble top yet does not coalesce with the bubble. Under such circumstances, the cavity length cannot be measured.

Figure 14.

(a) Maximum cavity length $h_{\textit{c}}$ as a function of the standoff parameter $\gamma$ . Magenta circles correspond to experimental data with $h_m$ ranging from 0.2 to 0.8 mm, while blue solid lines represent simulation results. $h_{\textit{c}}$ decreases monotonically with increasing $\gamma$ . (b) The same data in panel (a) are presented in a doubly logarithmic plot. The plot within the range of $1.5\lesssim \gamma \lesssim 3$ reveals a power law of $h_{\textit{c}}\propto \gamma ^\alpha$ , with fitted exponents $\alpha =-2.7$ (experiments) and $-2.6$ (simulations). The experimental uncertainty of $h_{\textit{c}}$ is approximately 0.004 (0.076 mm in dimensional form).

Figure 14 presents the variation of the maximum cavity length $h_{\textit{c}}$ with respect to the standoff parameter $\gamma$ . A comparison is performed between the simulation results (blue solid lines) and the experimental data (magenta circles). Overall, $h_{\textit{c}}$ decreases monotonically as $\gamma$ increases, consistent with the experimental observations discussed in § 3. Interestingly, the experimental data exhibit remarkable consistency with varying $h_m$ (from 0.8 to 0.2 mm), suggesting that the size of the surface perturbation plays a secondary role in cavity dynamics. To further verify this trend, numerical results obtained from FVM simulations with $h_m$ ranging from 0.8 to 0.2 mm are included for comparison. The simulation results show excellent agreement with the experiments, which demonstrates the reliability of the numerical approach and motivates further investigation. As discussed in § 3, nonlinear cavity–bubble interactions intensify and the flow-focusing effect becomes prominent when $\gamma \lesssim 1.5$ . Consequently, the cavity length $h_{\textit{c}}$ increases sharply and the slope of the $\gamma {-}h_{\textit{c}}$ curve steepens in this regime. A doubly logarithmic plot of the same datasets shown in figure 14(a) is presented in figure 14(b). In the range $1.5\lesssim \gamma \lesssim 3$ , the data follow a power-law relationship of the form $h_{\textit{c}}\propto \gamma ^\alpha$ , with a fitted exponent $\alpha =-2.7$ for the experimental data (magenta circles) and $\alpha =-2.6$ for the simulation results.

The intriguing consistency of $h_{\textit{c}}$ versus $\gamma$ , despite the variation in $h_m$ , motivates a further investigation into the underlying mechanism governing the cavity dynamics. As discussed in § 3, the entire life cycle of the cavity is closely coupled with the bubble’s collapse stage. Specifically, though the small depression is initiated and develops with the bubble expansion, the cavity does not start to grow downward rapidly until the bubble reaches its maximum radius. Conventional analytical approaches based on Kelvin impulse theory (Benjamin & Ellis Reference Benjamin and Ellis1966; Blake & Cerone Reference Blake and Cerone1982), which has been successfully adopted in several previous studies of a liquid jet (Kang & Cho Reference Kang and Cho2019; Han et al. Reference Han, Zhang, Tan and Li2022), are not applicable in this context, as the momentum of the cavity is difficult to estimate. Moreover, the majority of the bubble momentum cannot be transferred to the cavity due to its transient life cycle. Therefore, alternative approaches are needed to analyse the cavity dynamics.

After reviewing the relation between $h_{\textit{c}}$ and $\gamma$ , we notice that the power law of the doubly logarithmic plot in figure 14 exhibits slight deviations from the power law relationship at both ends of the parameter range of $\gamma$ . That is to say, the curve within the range of $1.5\lesssim \gamma \lesssim 3$ satisfies the fitted power law with exponent of $\alpha =-2.6$ , while the curve becomes slightly steeper when $\gamma \lesssim 1.5$ and $\gamma \gtrsim 3$ . The steeper slope when $\gamma \lesssim 1.5$ can be explained by the aforementioned flow-focusing effect in § 4.1, which becomes even more significant as $\gamma$ decreases to the vicinity of the critical value. This results in the strong nonlinear cavity–bubble interaction, manifested as a tapered cavity bottom. However, when $\gamma \gtrsim 3$ , the free-surface–bubble interactions tend to weaken. The cavity velocity and length are reduced to the order of $\textit {O}(10^{-2}{-}10^{-1})$ and $\textit {O}(10^{-1})$ , respectively. As a result, We decreases to the order of $\textit {O}(10^0{-}10^1)$ . Under such circumstances, the effects of surface tension are comparable to those of inertial forces. Meanwhile, the surface tension and viscosity of the fluid tend to resist the acceleration induced by the bubble oscillation and stabilise the free surface, thereby leading to smaller cavities. The above-mentioned aspects jointly result in a deviated $h_{\textit{c}}\rm {-}\gamma$ curve from the scaling law.

5.2. Theoretical modelling of cavity dynamics

To advance fundamental insights into cavity dynamics, the simulation results focusing on the cavity evolution are shown in figure 15. Velocity and pressure magnitude are presented in contour faces, while the details of the fluid field are given by velocity arrows and pressure contour lines. After the bubble inception, the free surface is slightly elevated and forms a water hump due to the bubble expansion (figure 15 a–d). The no-slip boundary condition of the thin rod creates a boundary layer that hinders the fluid motion, causing a bowl-shaped small depression to form atop the water hump. Though the fluid around the thin rod is hindered, the velocity is not significantly smaller than in the outer region. Thus, the depression still moves upwards with the water hump. As a result, the depression remains relatively small during the expansion stage.

Figure 15.

Flow fields in the vicinity of the cavity. Contours of velocity magnitude (left) and pressure (right) are plotted. The arrows represent the velocity direction, while the black curves represents the contour line of pressure. The free surface is plotted with magenta lines. This is the same case as in figure 4(a).

After the bubble attains its maximum radius and contracts (figure 15 e), the water hump starts to fall and triggers the cavity downward growth (figure 15 f). Slight distortion of the pressure field beneath the cavity can be observed, resulting in a higher local pressure gradient, manifested as the twisted and compressed pressure contour lines, which drag the cavity downward and trigger the accelerating process. Consequently, the growing cavity further distorts the pressure field in return. Such a positive feedback mechanism between the pressure gradient and the cavity growth drives a rapid cavity expansion, which is terminated only when the pressure gradient becomes reversed at the final stage of bubble collapse (figure 15 h).

A schematic for the cavity evolution is shown in figure 16. As we discussed in § 5.1 and previous sections, the cavity evolution can be divided into three stages: (I) a small depression forms concurrently with the rise of the water hump due to the bubble expansion (figure 16 a); (II) the cavity is initiated to expand after the bubble reaches its maximum radius (figure 16 b); (III) the cavity stretches rapidly downwards and attains its maximum length slightly before the bubble reaches its minimum volume (figure 16 c).

Figure 16.

Schematic of cavity evolution. (a) A small bowl-shaped depression forms due to the hindrance of the no-slip boundary condition around the thin rod during the rise of the water hump (blue dotted area). (b) The water hump reaches its maximum height $h_{\kern-1pt f}$ when the bubble reaches its maximum radius. The cavity grows downwards due to the relatively large local pressure gradient (red dotted area) and attains its initial velocity of $U_0$ . (c) The cavity reaches its maximum length $h_{\textit{c}}$ .

According to the definition in § 2.2, the cavity length $h$ starts to be recorded only upon the cavity’s arrival at the original horizontal surface $z=0$ . However, the growth of the cavity starts before this critical moment, which means that the cavity has an initial velocity of $U_0$ when $h=0$ (see figure 16 b). Then, the cavity starts to grow downward based on this initial velocity. Therefore, it will be convenient for the following analysis if we can obtain the analytical expression of $U_0$ . However, since the downward cavity growth and the rise of the water hump occur simultaneously during the bubble expansion stage, it is difficult to describe the exact cavity evolution during this period. Therefore, we adopt a simplified model to analyse $U_0$ theoretically.

To obtain the analytical expression of $U_0$ , several assumptions are made. First, we assume that the concurrent motion of the rising water hump and the expanding depression can be decoupled into two separated stages: the water hump first rises to the height of $h_{\kern-1pt f}$ , then the cavity is initiated from the crest of the water hump to expand downwards and arrive at the horizontal surface. Here, $h_{\kern-1pt f}$ represents the maximum height of the water hump in the first bubble cycle. This assumption is reasonable since the cavity does not start to grow rapidly until the water hump and the bubble reach their maximum height/radius. In fact, during most of the rising stage of the water hump, the depression retains a relatively small size. Second, we assume that the water hump is undisturbed, with its crest being the starting position for the cavity evolution. According to Yan et al. (Reference Yan, Liu, Kominiarczuk and Yue2009), for relatively fast cavity development ( $\textit {Fr}\gtrsim 1$ ), the wave and splash effects can be neglected. Finally, the acceleration of the cavity $a_{\textit{c}}$ has the same order as the pressure gradient $\boldsymbol{\nabla }p$ , yielding $a_{\textit{c}}\sim -\boldsymbol{\nabla }p$ . This assumption is feasible in the early stage of cavity growth since the pressure field has not been remarkably disturbed yet, as shown in figure 15(a–e).

When $\gamma \gtrsim 1.5$ , the bubble can be treated as spherical during most of the first cycle due to the weak bubble-free-surface interaction (Chahine Reference Chahine1977; Blake & Gibson Reference Blake and Gibson1981) and the influence of the cavity on the bubble in the non-coalescence scenarios can be neglected based on the discussion in § 4.2. Hence, a simplified image bubble model is adopted (Klaseboer et al. Reference Klaseboer, Khoo and Hung2005b ; Zhang et al. Reference Zhang, Li, Cui, Li and Liu2023). The bubble near the free surface can be approximated by a bubble pair symmetric about the free surface. The bubble pair can be described with the non-dimensional RP equation, which reads

(5.1) \begin{align} R_i\ddot {R_i}+\frac {3}{2}\dot {R_i}^2=\varepsilon \left (\frac {R_0}{R_i}\right )^{3\kappa }-1-\frac {2R_j\dot {R_j}^2+R_j^2\ddot {R_j}^2}{2\gamma },\quad i,j=1,2,\ i\neq j. \end{align}

The suffixes $i=1,2$ represent the original and image bubble. In the case of free surface, it is required that $R_2=R_1$ , $\dot {R}_2=-\dot {R}_1$ and $\ddot {R}_2=-\ddot {R}_1$ . The pressure of an arbitrary point $(r,z)$ in the fluid field induced by the bubble oscillation can be written as

(5.2) \begin{align} & p (r,z) \nonumber \\ & = \frac {2R_1\dot {R}_1^2+R_1^2\ddot {R_1}}{\sqrt {r^2+\left (\gamma + z\right )^2}}+\frac {2R_2\dot {R}_2^2+R_1^2\ddot {R_2}}{\sqrt {r^2+\left (\gamma -z\right )^2}}-\frac {R_1^4\dot {R}_1^2}{2\left [r^2+\left (\gamma + z\right )^2\right ]^2}-\frac {R_2^4\dot {R}_2^2}{2\left [r^2+\left (\gamma -z\right )^2\right ]^2}+1. \end{align}

The relative error introduced by the image bubble model can be controlled below 5 % if $\gamma \gt 1.5$ and 1 % if $\gamma \gt 2$ according to the previous study of Li et al. (Reference Li, van der Meer, Zhang, Prosperetti and Lohse2020).

Since the depression maintains a relatively small size in the bubble expansion stage and the initiation of the depression is transient, the influence that the depression imposes on the free surface can be neglected (Yan et al. Reference Yan, Liu, Kominiarczuk and Yue2009), and we can safely assume that the initiation process occurs when the water hump reaches its maximum height $h_{\kern-1pt f}$ when the bubble attains its maximum radius $R_m$ .

Building on the previous discussions, theoretical solutions of $h_{\kern-1pt f}$ and $U_0$ can be obtained henceforth. Velocity of the water hump crest $U_{\kern-1.5pt f}$ can be estimated by the fluid velocity at $ (z,r )=(0,0)$ induced by the real and image bubbles, which reads

(5.3) \begin{align} U_{\kern-1.5pt f}=\frac {R_1^2\dot {R}_1-R_2^2\dot {R}_2}{\gamma ^2}. \end{align}

Integrating (5.3) from 0 to $t_e$ with respect to time, where $t_e$ is the time that the bubble expands to its maximum radius, $h_{\kern-1pt f}$ can be obtained as

(5.4) \begin{align} h_{\kern-1pt f}=\frac {2\left (1-R_0^3\right )}{3\gamma ^2}. \end{align}

The acceleration of the cavity can be estimated by the $z$ component of the pressure gradient derived from (5.2) on the free surface:

(5.5) \begin{align} a_{\textit{c}}\sim \left . -\frac {\partial p}{\partial z}\right |_{(0,0)}=\frac {1}{\gamma ^2}\sum _{i=1,2}\big (2R_i\dot {R}_i^2+R_i^2\ddot {R_i}\big ). \end{align}

The third and fourth high-order terms in (5.2) are neglected in case of relatively large $\gamma$ . It is worth emphasising that this simplification will fail if the bubble is close to the free surface (Li et al. Reference Li, van der Meer, Zhang, Prosperetti and Lohse2020). However, the non-coalescence cases discussed in this section only take place when $\gamma \gtrsim 1.36$ , and for the specific regime where the fitted power law exists, $\gamma$ is large enough to adopt this simplification safely. By substituting (5.1) into (5.5), $a_{\textit{c}}$ can be written as

(5.6) \begin{align} a_{\textit{c}}\sim \frac {2}{\gamma \left (2\gamma +1\right )} \sum _{i=1,2}\left [\varepsilon \left (\frac {R_0}{R_i}\right )^{3\kappa }+\frac {1}{2}\dot {R}_i^2-1\right ]\!. \end{align}

As we mentioned already, the initiation process occurs around the bubble reaching its maximum radius. Under such circumstances, $R_i\approx R_m\gg R_0$ , $\dot {R_i}\rightarrow {0}$ . Thus, the first two terms in the square brackets can be neglected. By assuming the concavity grows with a uniform acceleration $a_{\textit{c}}$ , we can obtain the theoretical expression of $U_0$ , yielding

(5.7) \begin{align} U_0=\sqrt {2a_{\textit{c}}h_{\kern-1pt f}}=2.31\sqrt {\frac {1}{\gamma ^3\left (2\gamma +1\right )}}. \end{align}

It is obvious that (5.7) approaches the power law of $U_0\sim \gamma ^{-2}$ asymptotically as $\gamma$ increases. To justify the mentioned theoretical model, we perform a comparison of $U_0$ among experimental data, FVM simulation results, BI simulation results and theoretical results obtained from (5.7) in figure 17(a). We choose the configuration of $h_m=0.8\ \textrm {mm}$ in figure 14. The results obtained with FVM simulations, BI simulations and (5.7) reveal a similar fitted power law as that of the experimental data, with the fitted exponent being $-2.2$ , $-2$ and $-1.9$ , which are relatively close to that of the experimental data. However, the FVM and BI simulation tend to overestimate $U_0$ , while, in contrast, (5.7) underestimates. This discrepancy derives from the uniform acceleration assumption in (5.7), which is actually non-uniform and continuously increasing in the early stage of the cavity evolution. Nevertheless, the similar fitted power law relationships and the relatively reasonable agreement with the experimental data justify the effectiveness of the simplified model. Additionally, results of $h_{\kern-1pt f}$ obtained with all three methods are plotted in the inset, where the good agreement between the three lines justifies the accuracy of (5.4).

Once the cavity surpasses the horizontal plane of $z=0$ and starts to accelerate downwards, the distortion of the pressure field becomes prominent, where the pressure gradient increases, and the uniform acceleration assumption is no longer applicable. However, we can still estimate $a_{\textit{c}}$ using the bubble-induced pressure gradient at the cavity bottom with a coordinate of $(0,-h)$ , yielding

(5.8) \begin{align} a_{\textit{c}}\sim \left . -\frac {\partial p}{\partial z}\right |_{(0,-h)}\approx \frac {2R_1\dot {R}_1^2+R_1^2\ddot {R_1}}{\left (\gamma -h\right )^2}+\frac {2R_2\dot {R}_2^2+R_2^2\ddot {R_2}}{\left (\gamma +h\right )^2}. \end{align}

Figure 17.

(a) Doubly logarithmic plot for variation of $U_0$ with $\gamma$ . Experimental data, numerical results obtained from the FVM and BI method, and theoretical results obtained with (5.7) are brought into comparison. The FVM simulation results are the same as those when $h_m=0.8\ \textrm {mm}$ in figure 14. $h_{\kern-1pt f}$ obtained from (5.4), FVM simulation and BI simulation are plotted in the inset. The experimental uncertainties of $U_0$ and $h_{\textit{c}}$ are approximately 0.36 and 0.004, respectively ( $3.6\ \mathrm{m\,s}^{-1}$ and 0.076 mm in dimensional form). (b) Comparison among experimental data, FVM and BI simulation results, and theoretical results obtained with (5.7) and (5.8). The BI simulation and theoretical results converge with both experimental data and FVM simulation results as $\gamma$ increases, and reveal a power law with a fitted exponent of $\alpha =-2.2$ and $-2.1$ .

By integrating (5.8) twice over the interval of $t_e$ to $t_{\textit{c}}$ with initial conditions $h=0$ and $U=U_0$ , we obtain the theoretical results of $h_{\textit{c}}$ , where $t_{\textit{c}}$ is characteristic time of the bubble collapse stage. The results are plotted in figure 17(b) alongside experimental data, FVM simulation results and BI simulation results. The theoretical and BI simulation results exhibit a similar fitted power-law relationship, with their fitted exponent being $\alpha =-2.1$ and $-2.2$ . Interestingly, both curves converge towards the FVM simulation results and experimental data when $\gamma \to 3$ , while being distracted from them when $\gamma \to 1.5$ . The causes for this discrepancy are obvious. When $\gamma$ approaches 1.5, the nonlinear effects become significant, which are not considered in both models. Specifically, the BI method does not account for the compressibility effects of the air inflow into the cavity, viscosity and vorticity effects of the fluid, and the boundary layer on the thin rod surface. In addition to these, the theoretical model further neglects the nonlinear boundary condition of the free surface, pressure field distortion caused by cavity growth, and the non-spherical behaviours of the bubble. Despite these simplifications, both the results of the BI simulation and the theoretical model agree reasonably well with experimental data and FVM simulation results. Especially when $\gamma \to 3$ , where inertia dominates while nonlinear effects are not prominent, the consistency between the four results is particularly good, which further justifies that the pressure gradient induced by bubble oscillation plays a major role in the cavity evolution.

However, one may notice that the theoretical model we present in this section does not comprise another key parameter $h_m$ , although whose influence on the cavity evolution is not as significant as $\gamma$ , but still cannot be neglected. We will address this parameter in the next section.

Figure 18.

Comparisons between (a) numerical results obtained from FVM (blue solid lines) and (b) analytical results obtained from (5.11) (orange solid lines). The dependence between $h_m$ and $h_{\textit{c}}$ are plotted in panels (a) and (b). $\gamma$ increases from 1.4 to 3 at a 0.2 interval. Both results are normalised to (0,1) for clearer comparison. The analytical calculations exhibit a similar tendency to the numerical results. The non-dimensional parameters in the FVM simulations are set to the same as § 4. The length scale is $R_m=17.2\ \textrm {mm}$ . The relationship between $h_m$ and $\phi _0$ is presented in panel (c).

5.3. Nonlinear Rayleigh–Taylor instability analysis

Based on the previous discussions, we have already learned that cavity evolution is highly related to the collapse stage of the cavitation bubble. During this stage, the denser fluid (water) is accelerated by the lighter fluid (air). We also learned that the dominant factor in bubble dynamics is the standoff parameter $\gamma$ and the pressure gradient/acceleration field in the fluid induced by the bubble oscillation is the major cause for cavity growth. Meanwhile, as we mentioned in § 5.1, the surface perturbation can also exert a secondary, yet not negligible effect on the cavity dynamics, i.e. the cavity length increases monotonically with the increasing $h_m$ . This acceleration-driven cavity growth from the lighter fluid to the heavier one, affected by the amplitude of the initial interface perturbations, is very similar to the characteristics of the Rayleigh–Taylor instability (RTI), thus inspiring us to perform an RTI analysis in this section.

We first perform several groups of FVM simulations to have a preliminary insight into the $h_{\textit{c}} {-}h_m$ dependence. The numerical results are presented in figure 18(a) with blue solid lines. As we expected, the $h_{\textit{c}} {-}h_m$ curve exhibits a monotonically increasing relationship. If the cavity evolution is subject to the linear RTI, $h_{\textit{c}}$ should be proportional to $h_m$ following the equation of $h=h_m \cosh (\sqrt {a_{\textit{c}} kA}t)$ , where $k=2\pi /\lambda$ is the wave number, $\lambda$ the wavelength of the initial perturbation, $A=(\rho _l-\rho _g)/(\rho _l+\rho _g)$ the Atwood number and $t$ the characteristic evolution time (assuming as constant). However, the $h_{\textit{c}}{-}h_m$ curves exhibit decreasing slope, suggesting that RTI evolution of the cavity is in fact nonlinear, which is reasonable since the cavity grows into a scale much larger than that of the initial perturbation ( $h_{\textit{c}}\gg h_m$ ) within a very short period.

Therefore, a nonlinear RTI equation proposed by Mikaelian (Reference Mikaelian2010) based on the widely used axisymmetric nonlinear RTI model developed by Layzer (Reference Layzer1955) is adopted, which reads

(5.9) \begin{align} \frac {{d}^2 \hat {h}}{\mathrm{d} s^2}-\hat {k} \hat {A} \hat {h}+\frac {1}{2a_{\textit{c}}^2}\frac {\mathrm{d}a_{\textit{c}}}{\mathrm{d}t}\frac {\mathrm{d}\hat {h}}{\mathrm{d}t}=0, \end{align}

where $\hat {h}=e^{(h-h_m)k}$ , $\hat {k}=d(1+d)(1+A)k/2(1+d+dA-A)$ and $\hat {A}=2A/(1+d+dA-A)$ are the nonlinear counterparts for $h$ , $k$ and $A$ , while $d$ is the generalisation parameter for models of different dimensions, i.e. $d=1$ for 3-D, $d=2$ for 2-D. Here, $s=\int _{t_0}^{t}\sqrt {a_{\textit{c}}}\,\mathrm{d}t$ is the substituting variable, where $t_0$ and $t$ is the initial and terminal time of instability evolution, and $a_{\textit{c}}$ the driven acceleration in (5.8). We notice that the negative acceleration in the cavity evolution will lead to imaginary numbers in the $\sqrt {a_{\textit{c}}}$ terms. By introducing a modification $\mathrm{sgn}(a_{\textit{c}})\sqrt {|a_{\textit{c}}|}$ , we can take the negative acceleration into account (Ramaprabhu et al. Reference Ramaprabhu, Karkhanis and Lawrie2013). Such modification is reasonable, since the cavity acceleration is positive (pointing downwards) during the majority of the cavity evolution, while the deceleration and rebound process only comprises a minor part, as we discussed in § 4.1. Therefore, it is also safe to further neglect the $\mathrm{d}a_{\textit{c}}/\mathrm{d}t$ terms in (5.9) without significantly changing the cavity dynamics. It is also suggested by Yan et al. (Reference Yan, Liu, Kominiarczuk and Yue2009) that the flow around the cavity can be treated as two-dimensional and axisymmetric when the Fr is relatively large, yielding $d=2$ . Given the vast density difference between the two fluids (air and water) in this study, $A$ can be set as 1. Hence, (5.9) can be simplified as

(5.10) \begin{align} \frac {\mathrm{d}^2 \hat {h}}{\mathrm{d} s^2}-\hat {k}\hat {A}\hat {h}=0. \end{align}

Solving this equation with the first-order Wentzel–Kramers–Brillouin (WKB) approximation (Bender & Orszag Reference Bender and Orszag2013) yields

(5.11) \begin{align} h=h_m+\frac {1}{\hat {k}}\ln \left [\cosh \left (s\sqrt {\hat {A} \hat {k}}\right )\right ]. \end{align}

The determination of $\hat {k}$ is relatively challenging in the present study, since the meniscoid surface perturbation does not have a representative wavelength/wavenumber. Unlike the common periodic initial perturbation, where the wavenumber and wavelength are readily obtainable both experimentally and numerically (Emmons, Chang & Watson Reference Emmons, Chang and Watson1960; Herbert Reference Herbert1983; Jacobs & Catton Reference Jacobs and Catton1988). Many precedent works have adopted several methods to analyse the influence of single-/multi-mode initial perturbations, including spectrum analysis (Cohen et al. Reference Cohen, Dannevik, Dimits, Eliason, Mirin, Zhou, Porter and Woodward2002; Olson & Jacobs Reference Olson and Jacobs2009; Ramaprabhu et al. Reference Ramaprabhu, Karkhanis and Lawrie2013) and Fourier expansions (Luo et al. Reference Luo, Liang, Si and Zhai2019; Liang & Luo Reference Liang and Luo2021; Li et al. Reference Li, Gao, Guo, Zhai and Luo2025). In the present study, the initially perturbed interface contains both high-frequency (the capillary length scale meniscus) and low-frequency (the remaining flat free surface) compositions. Although the high-frequency composition only comprises a small part of the interface, it can exert a remarkable influence on the free-surface–bubble system according to the aforementioned discussions. Meanwhile, the influence of the remaining low-frequency compositions may not be prominent, but it still constitutes the majority of the interface. Moreover, the vast scale difference between the meniscus and the bubble suggests that the Fourier series needs to be expanded from low order to very high order, otherwise the interface may not be accurately depicted, not to mention the pinning point at the contact line may introduce a strong singularity. Therefore, efforts to approximate the interface analytically are challenging and costly in both time and calculation resources, which deviates from the original intention of adopting the analytical model.

Nonetheless, this nonlinear RTI model can still explain some issues by choosing a characteristic wavelength/wavenumber. Mikaelian (Reference Mikaelian2010) suggested $h_m\sim 1$ or $h_mk\gtrsim 1$ , since many nonlinear RTI processes start with a weakly nonlinear amplitude and grow from there, and rapidly enter the nonlinear regime $hk\gg 1$ , which corresponds to the transient cavity evolution in the present study. Meanwhile, Huang et al. (Reference Huang, Zou, Liu, Tan and Zhang2010) and Zhao et al. (Reference Zhao, Li, Wang, Xue, Wu, Xiao, Ye, Ding, Zhang and He2023) pointed out that the compositions with larger wavelength/smaller wavenumber in the initial perturbations tend to grow faster within the nonlinear region and stimulate the development of the large-scale RTI structures. This argument is often related to the bubble competition and merging in the multi-mode RTI evolution, where multiple instability bubbles evolve simultaneously. However, in the present study, only one primary cavity can be observed in experiments and simulations, which is actually the superposition of multiple cavities that evolve from various modes. Therefore, the dominant mode demands a smaller $k$ according to this criterion to achieve a faster growth velocity. The paradox between these two criteria indicates that the perturbation composition with the largest wavelength/smallest wavenumber that satisfies the weakly nonlinear assumption tend to dominate the cavity evolution, thus leading us to a relation of $h_mk=1$ and $\hat {k}=1.5/h_m$ .

Once $\hat {A}$ and $\hat {k}$ are determined, (5.11) can be solved by iterating $a_{\textit{c}}$ and $h$ throughout the cavity life cycle. The calculation results are presented in figure 18(b) with orange solid lines. To gain a clearer vision of the analytical results and compare the $h_m {-}h_{\textit{c}}$ tendency in both results, $h_{\textit{c}}$ is normalised to (0,1). Encouragingly, the analytical results exhibit a similar dependence as the numerical results, i.e. the $h_{\textit{c}}$ increases monotonically as $h_m$ with a decreasing slope. However, the analytical results still have deviations from those of the numerical simulations, as expected. Specifically, the ratio between the numerical and analytical results varies from the order of $\textit {O}(10^1)$ to $\textit {O}(10^0)$ , as $\gamma$ increases from 1.4 to 3. The best agreement between the two results is obtained when $\gamma \to$ 3, while they deviate from each other when $\gamma \to 1.5$ . This deviation can be explained as follows. First, the same as the analytical model that we proposed in § 5.2, the nonlinear RTI model does not account for the influence of the boundary layer, whose influence would be increasingly prominent as $\gamma$ decreases. Moreover, non-spherical bubble and flow-focusing effects both can cause strong nonlinear cavity–bubble interactions and increases $h_{\textit{c}}$ . Additionally, other factors such as viscosity, vorticity, air compressibility and surface tension are also neglected, whose influence may not be significant as the previously mentioned ones, but still can exert slight changes on the cavity dynamics. Second, since we only choose one characteristic wavelength/wavenumber for simplification and neglect the remaining compositions of the initial perturbation, it is foreseeable that this deviation occurs. Moreover, the meniscus contact angle $\phi _0$ , which characterises the meniscus geometry, is also considered. Equation (2.13) shows that the meniscus height $h_m$ decreases monotonically with increasing $\phi _0$ (see figure 18 c), indicating that $h_m$ and $\phi _0$ are intrinsically related. Consequently, once the $h_{\textit{c}}{-}h_m$ relationship is established, the corresponding $h_{\textit{c}} {-}\phi _0$ dependence is implicitly determined and $h_{\textit{c}}$ is expected to decrease with increasing $\phi _0$ . Since $h_m$ can be measured far more accurately than $\phi _0$ in the current experimental set-up, we therefore focus on $h_m$ as the primary parameter in the analysis.

One can also notice that the significant increase of $h_{\textit{c}}$ (from $0.6\times 10^{-2}$ to $5.2\times 10^{-2}$ ) does not result in equivalent changes in the magnitude of $h_m$ . This observation is consistent with the discussions in § 5.1 associated with figure 14, where the scaling law $h_{\textit{c}}\propto \gamma ^\alpha$ remains relatively identical as $h_m$ varies, thus demonstrating the secondary role of $h_m$ in the cavity dynamics.

Nonetheless, the dependence of $h_m {-}h_{\textit{c}}$ revealed by the nonlinear RTI model is similar to the numerical results, suggesting that the cavity evolution on the initially perturbed interface driven by an oscillating cavitation bubble is actually a non-uniform acceleration driven nonlinear Rayleigh–Taylor instability phenomenon, where the acceleration field induced by the bubble oscillation governs the cavity dynamics, while the amplitude of the initial perturbations plays a secondary role by varying the cavity maximum length without significantly changing the fitted power law.

5.4. Influence of boundary layer

Apparently, the analytical model in § 5.3 will fail if $h_m=0$ . However, during the experiments, we notice that if we use the contact line pinning effect by carefully moving the thin rod to eradicate the meniscoid and obtain a flat surface, i.e. $h_m\approx 0$ , the cavity still occurs, but with a significantly smaller $h_{\textit{c}}$ . If the $\gamma$ is not very large at the same time, the cavity may be pierced by the downward jet and break up into pieces. This suggests that the boundary layer around the thin rod may have a minor influence on the cavity dynamics.

A representative experiment of the free-surface–bubble interaction without the meniscus perturbation is shown in figure 19. It clearly shows that the cavity is significantly smaller than that in the regular meniscus-perturbed experiments presented in the red box. Additionally, due to a much earlier surface-seal, the cavity becomes an isolated small bubble and moves downward under the influence of the pressure gradient. At the same time, the sealed surface creates a downward jet and pierces the bubble before the cavity attains its maximum length.

We perform a numerical simulation of the surface without meniscus perturbation to study the difference between the two experiments. By setting the free surface as flat, we can obtain a non-meniscus-perturbed interface $h_m=0$ . The comparison between the numerical simulation results and experimental observation is presented in figure 19. The numerical simulation successfully reproduces the main characteristics of the isolated cavity and the piercing jet. The simulation results intuitively explain the significant difference between the cavity dynamics with and without meniscus perturbation. In the absence of the meniscus perturbation, only the boundary layer in the close proximity of the thin rod is hindered by the no-slip boundary condition. As a result, the depression on the top of the water hump is much smaller than that in the regular perturbed case when the bubble attains its maximum radius (frame 2). The smaller depression leads to a much smaller radial and longitudinal size of the cavity when the bubble starts to contract (frame 3). Moreover, the influx of air into the cavity causes the small cavity opening to seal much earlier, creating an isolated bubble and a downward jet piercing the cavity bottom (frame 4). All the aforementioned factors conjointly cease the cavity evolution, resulting in a much smaller cavity.

Figure 19.

Comparison of cavity and bubble evolution between numerical simulation and experimental observation in a representative case without meniscus perturbation ( $H=27.5\ \textrm {mm}$ , $R_m=18.5\ \textrm {mm}$ , time scale is 1.84 ms). Under such circumstances, the cavity is significantly smaller than that of the regular cases (red frame on the right, $H=28.2\ \textrm {mm}$ , $R_m=18.6\ \textrm {mm}$ , $h_m=0.94\ \textrm {mm}$ , time scale is 1.85 ms), and the cavity becomes an isolated bubble due to the much earlier surface-seal. The sealing surface creates a downwards jet piercing the bubble rapidly. Although slight deviation exists in the size and position of the isolated cavity between the two results, the numerical simulation successfully reproduces the main features.

However, a slight deviation still exists in the size and position of the cavity between simulation results and experimental observations. This deviation derives from the disturbance on the ‘flat’ free surface. Since the meniscus and the free surface is highly sensitive to external disturbances (5.3), minor distractions such as slight vibrations and breezes could still introduce deformation on the initially flat free surface, even if the macro camera did not capture any evidence of shape changes (see the meniscus image with blue marker in figure 3 b). Therefore, the ‘flat’ free surface in the experiments may not be strictly flat and the small meniscus may result in a slightly increased cavity length $h_{\textit{c}}$ .

Figure 20.

Comparison of cavity evolution across four distinct configurations: (i) free surface without perturbation (green solid lines), (ii) flat surface with no-slip BC (yellow solid lines), (iii) meniscus perturbed surface with no-slip BC (red solid lines) and (iv) meniscus perturbed surface with slip BC (blue dashed lines). Profiles of the interface at four different times are presented: (a) 0; (b) 0.62; (c) 1.17; (d) 1.86. The time scale is 1.69 ms and the length scale is $R_m=17.0\ \textrm {mm}$ . The non-dimensional parameter of the bubble is the same as that in figure 4(a).

Given our model’s limitations in quantitatively analysing the boundary layer, we qualitatively assess its effects by toggling the rod surface’s no-slip boundary condition (BC). Figure 20 compares four configurations: (i) unperturbed free surface (green solid lines); (ii) flat surface with no-slip BC (yellow solid lines); (iii) meniscus perturbed interface with no-slip BC (red solid lines); and (iv) meniscus perturbed interface with slip BC (blue dashed lines). The unperturbed free surface serves as a benchmark. Comparisons between each one of these results suggest the following arguments. (1) The slight discrepancy between configurations (iii) and (iv) indicates that turning off the no-slip BC does not significantly vary the cavity characteristics, suggesting that the boundary layer does not exert a remarkable influence on the cavity dynamics when $\gamma$ is relatively large. (2) The vast difference between configurations (ii) and (iii) indicates that the boundary layer alone is not enough for the sufficient development of the cavity. If only the meniscoid perturbation exists, the main characteristics of the cavity occur. (3) The comparison among configurations (ii), (iii) and (iv) suggests that the meniscoid perturbation plays a major role in the cavity dynamics, while the existence of the boundary layer facilitates the cavity evolution. Changing the boundary condition varies the partial cavity details, e.g. the tapered or semi-spherical cavity bottom morphology in configurations (iii) and (iv). In conclusion, the overall radial and longitudinal scale of the cavity is still determined by the geometrical characteristics of the meniscoid perturbation. The detailed discussions between the effects of air compressibility and boundary layer are presented in Appendix C.