Cavitation bubble interaction with a rigid spherical particle on a microscale

Cavitation bubble collapse close to a submerged sphere on a microscale is investigated numerically using a finite volume method in order to determine the likelihood of previously suspected mechanical effects to cause bacterial cell damage, such as impact of a high speed water jet, propagation of bubble emitted shock waves, shear loads, and thermal loads. A grid convergence study and validation of the employed axisymmetric numerical model against the Gilmore’s equation is performed for a case of a single microbubble collapse due to a sudden ambient pressure increase. Numerical simulations of bubble-sphere interaction corresponding to different values of nondimensional bubble-sphere standoff distance δ and their size ratio ε are carried out. The obtained results show vastly different bubble collapse dynamics across the considered parameter space, from the development of a fast thin annular jet towards the sphere to an almost spherical bubble collapse. Although some similarities in bubble shape progression to previous studies on larger bubbles exist, it can be noticed that bubble jetting is much less likely to occur on the considered scale due to the cushioning effects of surface tension on the intensity of the collapse. Overall, the results show that the mechanical loads on a spherical particle tend to increase with a sphere-bubble size ratio ε , and decrease with their distance δ . Additionally, the results are discussed with respect to bacteria eradication by hydrodynamic cavitation. Potentially harmful mechanical effects of bubble-sphere interaction on a micro scale are identified, namely the collapse-induced shear loads with peaks of a few megapascals and propagation of bubble emitted shock waves, which could cause spatially highly variable compressive loads with peaks of a few hundred megapascals and gradients of 100 MPa/μm.


Introduction
Although cavitation at first sparked its interest in the scientific community due to its negative effects, such as material erosion and vibrations in hydraulic machinery [12], it is nowadays being utilized for a wide variety of applications. Among others, we can find cavitation in processes of surface cleaning [7], in various medical procedures [41] and wastewater treatment [11]. Recently it has also been shown, that hydrodynamic cavitation can be efficiently used for eradication of bacteria [45] and even inactivation of viruses [24] in water. Despite the progress, the key mechanisms behind the interaction between cavitation bubbles and water contaminants are still unknown.
Considering a single cavitation bubble and a bacterial cell, we can generally distinguish between their interaction on two different spatial scales, where cavitation bubbles are either significantly larger (≫ 1 μm) or of a similar size in comparison to bacterial cells (~1 μm). Due to the large size difference in the former case, the presence of nearby microorganisms in not believed to significantly affect the behavior of bubbles. This can be supported by previous studies by Teran et al. [47], who investigated the interaction of bubbles with sediment particles. On the other hand, nearby particles or microorganisms could potentially affect the dynamics of bubbles of a similar size scale.
In the present paper, we address bubble-sphere interaction on a microscale, where the presence of a nearby suspended spherical particle could potentially have an important effect on bubble dynamics. Since we are limited with observing the considered phenomenon experimentally due to small spatial (~1 μm) and temporal scales (~10 ns), we employ a purely numerical approach. Our main interest is to determine how microbubble dynamics depends on the non-dimensional bubblesphere stand-off distance and their size ratio, and how this affects the loads exerted on a spherical particle. Here, it has to be noted that particles and microorganisms in reality adopt various shapes and sizes. For example, single bacterial cells normally measure between 0.5 and 5 μm in length, with shapes ranging from spherical to more complex, like spherocylindrical, helical, filamentous, etc. [36].
By considering an interaction of bubbles with spherical objects, the present paper falls into a branch of bubble-sphere interaction studies. Perhaps one of the first, who addressed this topic were Kalumuck et al. [22]. They conducted numerical simulations of single bubbles nearby structures of various shapes and characteristics by coupling boundary element method (BEM) and finite element method (FEM) for the fluid and structure domain, respectively. Later followed Tomita et al. [50], who showed that the velocity of a micro jet, that develops during the bubble collapse, tends to increase with higher curvatures of a nearby convex rigid structure and can surpass peak jet velocities that occur in the case of a nearby flat surface.
Li et al. [33] conducted a numerical analysis of a gas bubble close to a fixed and movable rigid sphere by using a BEM approach. The obtained temporal development of bubble radius showed no difference between the two cases for the growth phase of the bubble, while the collapse time and minimum radius were slightly lower for the case of a movable sphere. Similar topic was later addressed by Li et al. [29], who proposed the use of BEM along with an auxiliary function to decouple the mutual dependence between the force on a suspended sphere and its movement. Their results show complex bubble dynamics in relation to different values of the stand-off parameter and bubble-sphere size ratio, such as annular jet development away from the sphere, toroidal bubble splitting, and jetting towards the sphere. The same model was later successfully utilized for the studies of bubble-sphere interaction near a rigid wall [31] and beneath a free surface [30].
In the recent years this topic was also addressed experimentally. Borkent et al. [4] researched particles in diameter between 30 and 150 μm that were suspend in water and served as nucleation sites in the event of tensile wave propagation. The authors noticed an interesting phenomenon of particle-bubble interaction, where a spherical particle detaches from a collapsing bubble and is accelerated away. Additionally, a secondary bubble inception was observed after collapse of a primary bubble. Poulain et al. [43] researched the motion of a spherical particle induced by a single spark-induced cavitation bubble. The authors distinguished between three phases of particle motion: an initial movement away from the expanding bubble, a change in direction due to bubble collapse, and a postcollapse phase, where the sphere was further drawn towards the oscillating bubble. Later, an extensive experimental study of interaction between a laser-induced cavitation bubble and a rigid spherical particle was conducted by Zhang et al. [53]. They identified three general cases of bubble collapse: a mushroom-shaped, a pear-shaped, and a spherical-shaped collapse, along with their dependence on the non-dimensional bubble-sphere distance.
In certain applications, deformability of a sphere is an additional important factor to be considered besides the sphere's ability to freely move within the surrounding fluid. Both were accounted for by Goh et al. [15], who numerically and experimentally studied spark-generated bubbles near an elastic sphere. Similarly to Kalumuck et al. [22], the authors also utilized a coupled BEM-FEM solver which showed a good match with the experimental results. Silicone rubber spheres with elasticity modulus = E 312 kPa, which is similar to the elasticity of bacterial cells with E~500 kPa [44], showed minimal deformations due to bubble growth and collapse. Additionally, the sphere's elasticity seems to be of a minor importance for the dynamics of a spark-generated bubble in the presented case.
All the literature mentioned above is applicable to the bubblesphere interaction on a larger, millimeter scale, where the effects of surface tension and viscosity play a minor role on bubble dynamics. The former was taken into account by Gracewski et al. [16] and later Guo et al. [17]. Using BEM, they evaluated interaction between a cavitation microbubble ( = R 1.5 μ eq m) and a red blood cell in an ultrasonic field. Elasticity of the cell's membrane was considered as surface tension, which is linearly proportional with the areal expansion, while the cell's contents were described with the Tait equation of state [32] as water. In most cases results showed formation of an axial jet away from the cell and a maximum areal expansion of a cell in the order of 0.1%.
On the other hand, the existing experimental literature shows that a bubble can cause membrane poration of a single cell and also microporation of other materials, such as biochar [1]. Le Gac et al. [27] reported the sonoporation of suspended human promyelocytic leukemia cells using a single laser-induced cavitation bubble in a microfluidic confinement. The authors attributed membrane damage to shear stresses exerted on cells due to bubble-induced flow and observed a poration probability of more than 75% for cells closer than 0.75 times the maximum bubble radius. The same value of an effective range for membrane poration was later reported by Zhou et al. [54], who acoustically excited single microbubbles in vicinity of an oocyte. Similarly, Li et al. [34] demonstrated a membrane poration of myeloma cells, although in this case the main mechanism was the impact of a microjet. The authors utilized a microfluidic chip to trap single suspended cells and a laser induced cavitation bubble, which collapsed asymmetrically due to a nearby solid boundary of the cell trapping structure.
One can notice, that all the present numerical investigations on the topic of bubble-sphere interaction utilize potential flow theory along with BEM. We employ a different approach, a finite volume method (FVM) along with the volume of fluid (VOF) method to resolve multiphase flow. This approach has already been shown to successfully resolve various cases of aspherical bubble dynamics, such as in vicinity of a rigid wall [38,23], in a gravity field [26], bubble pair interaction [18], etc. Its advantage lies in a possibility to consider nonlinear water compressibility, which is necessary to resolve shock wave emission and propagation upon bubble collapse.
Additionally, the obtained results are discussed with respect to bacteria eradication in water. Potentially harmful mechanical effects of observed bubble-sphere interaction on a micro scale are identified, and include: impact of a high speed water jet, propagation of bubble emitted shock waves, shear loads, and thermal loads.

Theoretical background and numerical model
In this section we present the considered conservation laws and equations of state that along with adequate initial and boundary conditions enclose a system of equations that describe the phenomenon of bubble dynamics. Additionally we include Gilmore's model that describes the dynamics of a single unbounded spherical bubble in a compressible liquid. Since our approach is purely numerical, Gilmore's equation serves as a means to validate the main axisymmetric FVM-VOF model in the present study.

The governing equations of bubble dynamics
In general, bubble dynamics can be mathematically described by the principles of mass, momentum, and energy conservation. Continuity equation of mass is given below. As can be seen from the right side of Eq. (1): mass transfer is not considered in the present study, therefore accumulation of mass can be only attributed to the compressibility of fluids.
Here ρ represents density and U velocity vector field of the considered fluid.
Conservation of momentum can be written as Eq. (2), where terms on the left side represent the effect of inertial forces due to local and convective acceleration, respectively. The right-hand side of the equation includes the effects of pressure, viscous, and body forces, respectively. Here we neglect the effect of gravity due to a small scale of the considered phenomenon. The effects of surface tension are included with a body force acting on the water-bubble interface (see Section 2.4 for further details).
Here p denotes pressure, f body forces, and τ the viscous stress tensor that can be written for Newtonian fluids as: where μ is dynamic viscosity and I the unit tensor. Energy balance is described by Eq. (4), which includes the effects of thermal conduction on the right-hand side.
Here k represents thermal conductivity and e total specific energy, which can be written as where h is specific enthalpy.

The considered equations of state
Since we wish to capture shock waves upon bubble collapse, a simplified Tait equation of state for water is employed, which considers its nonlinear compressibility. It is given by Eq. (6): where the bulk modulus of water K at pressure p is calculated as where R gas and M are the universal gas constant and molar mass, respectively.

Gilmore's model for a spherical bubble
The dynamics of a single unbounded spherical bubble inside a compressible liquid can be described by the Gilmore's equation [14], given in Eq. (8).
This is a second order nonlinear ordinary differential equation, whose solution represents a time dependent bubble radius R t ( ), from which we can evaluate the bubble wall velocity Ṙand acceleration R. Looking at the equation, we can find additional terms, C and H, which represent the speed of sound and the enthalpy difference at the liquid just outside the bubble wall. Both can be evaluated following Eq.
Gilmore's model also accounts for surface tension and the principal effects of viscosity, which can play an important role on bubble dynamics at smaller scales. Both effects are modeled as a pressure jump across the bubble-liquid interface, between the internal bubble pressure P i and the liquid pressure just outside the bubble wall P, both given in Eq. (10). In the present case the bubble contains both water vapor and air, which is assumed to behave adiabatically.
Here σ denotes the surface tension coefficient, p v vapor pressure, γ adiabatic index, ∞ p ,0 the initial ambient pressure, and R eq the corresponding bubble radius at equilibrium.

Numerical model
The main numerical model [2] is based on the finite volume method to solve the aforementioned system of equations in Section 2. We employ the VOF method to formulate multiphase flow, through which we track the interface between both phases (liquid, gas) by solving the continuity equation (Eq. (11)) for the volume fraction of water α w .
Volume fraction of the gaseous phase α g can be obtained through Eq. (12): which states that in each control volume the volume fractions of both phases sum to unity. Based on the known volume fraction field, we can determine the volume-averaged material properties throughout the computational domain. Examples for calculation of density and dynamic viscosity are given in Eq. (13) below.
Following this, a single momentum (Eq. (2)) and energy (Eq. (4)) equations are solved, from which the shared velocity and temperature fields are obtained based on the material properties from Eq. (13). Additionally, the effects of surface tension are included in the procedure with a body force in the momentum equation, according to the continuum surface force model [5]. The pressure jump across the liquid-gas interface is modeled with a volume force: where κ g denotes the bubble surface curvature, calculated as: Here, n is a bubble surface normal, which is obtained as a gradient of the gas volume fraction field: = ∇α n g . For all calculations the PISO pressure-velocity coupling algorithm [21] was employed, along with a first order implicit temporal discretization. Regarding the spatial discretization, we used the pressure staggering option (PRESTO!) scheme [40] for pressure interpolation and the second order upwind scheme for density, momentum, and energy interpolation. In the cells near the water-bubble interface, a Piecewise linear interface calculation (PLIC) geometric reconstruction scheme [52] was used, which assumes that the interface between two fluids has a linear slope within each cell. Boundary conditions at the end of the computational domain were set to wave non-reflecting pressure outlet, which was placed reasonably far away from the bubble ( R100 max ) to minimize its influence on the considered phenomenon. The boundary conditions are based on characteristic wave relations derived from the Euler equations [48,49], from which the primitive flow quantities can be obtained. To solve the obtained system of equations on a boundary, the amplitude of the incoming and outgoing waves must be determined. The amplitude of the latter are computed by using extrapolated values of flow derivatives from within the computational domain, whereas the amplitude of the incoming pressure and entropy waves are computed from the Linear relaxation method of Poinsot and Lelef [42] and Selle et al. [46].
Additionally, in the case of bubble-sphere interaction simulations, we utilize a diffusion based dynamic mesh smoothing method to adapt the numerical grid to the sphere's movement after each time step. This method is based on the Laplace's equation from which vertex velocities u can be obtained according to the diffusion coefficient D and velocities at the sphere's boundary. In our case, a nonuniform diffusion coefficient was utilized, defined as a reciprocal value of the normalized cell volume.

An equilibrium micro-bubble in an infinite fluid
At first, we consider a single unbounded cavitation bubble of radius Pa. In this case, these values do not represent the operating pressures of a given hydraulic system, but rather a range of pressure values that can locally occur within a cavitating flow, e.g. the ambient pressure of a microbubble increases due to a nearby or surrounding bubble cloud collapse. Similar values to ones used here were also considered by previous authors that researched material pitting and erosion from a single cavitation bubble [19,8]. Additionally, it has been shown that the local pressures that drive bubble collapse vary across spatial scales [6] and can in some cases exceed 100 MPa [37]. The considered temporal dependence of ambient pressure is given in Eq. (17).
, for 0s, , for 0s ,0 Vapor pressure, surface tension coefficient and viscosity of water were set to = p 2339 v Pa, = σ 0.0728 N/m, and = μ 0.001 w Pa · s, respectively. Looking at the obtained results from the Gilmore's model in Fig. 1, we can observe that bubbles on a scale of~1 μm behave differently compared to the larger ones. Firstly, they are well below the critical radius of asymptotic linear growth, which is the reason, that in the case of a sudden ambient pressure drop ( = ∞ p 10 3 and 10 4 Pa), the bubble oscillates to a new equilibrium state. Interestingly, the magnitude of oscillations remains relatively small, not surpassing 35% of initial radius even in the case of = ∞ p 10 3 Pa. As expected, the case with = ∞ p 10 5 Pa does not show a notable change in bubble radius, due to bubble's initial equilibrium ambient pressure of 1 atm. Moving forward, to the bubble compression region, we can see that the oscillations at higher ambient pressures happen on roughly one order of magnitude smaller temporal scales. In the case of = ∞ p 10 6 Pa, the collapse driving pressure is not high enough to cause a bubble to collapse violently, but rather induces a weak response ( ) and is therefore not believed to have an important effect on the bubble-sphere interaction. Since we are interested in potentially harmful mechanisms of bubble-particle dynamics, only the latter case ( = ∞ p 10 7 Pa) of bubble collapse is chosen for further analyses.
Following that, we determine the effects of viscosity, surface tension and consideration of vapor pressure on bubble dynamics for the selected case. A comparison between the obtained results from the Gilmore's model is given in Fig. 2. We can see that by not accounting for surface tension effects, the bubble reaches roughly 1.6-times smaller collapse radius, which leads to the overestimation of collapse pressure by the factor of 7. Viscous effects at the bubble wall seem to play a lesser role in the considered phenomena, with 1.4% shorter collapse time and 3.9% smaller collapse radius for the inviscid case. On the other hand, the neglection of vapor pressure seems to have no notable effect on bubble size evolution. This can be explained by the fact, that vapor pressure of = p 2339 v Pa presents less than 1% of the initial, equilibrium bubble pressure, which further decreases during the bubble collapse. Because of this, we consider purely air bubbles and neglect their vapor content in further calculations.

Grid convergence study
The selected bubble collapse case, that serves as a basis for further analyses, was also numerically evaluated in the axisymmetric numerical model (see Section 2.4). In order to determine the adequate computational grid resolution and assess an error that arises from spatial discretization, we performed a grid convergence study and a Richardson extrapolation of the obtained results. Similar initial pressure conditions were used as with the previously used Gilmore  convergence criteria with the values of scaled residuals for continuity and momentum equations to − 10 6 , and energy equation to − 10 9 , which on average resulted in ten to fifteen iterations per time step.
Altogether, six simulations were carried out, with grid resolution ranging from 50 and up to 566 cells across the initial bubble radius R 0 , which corresponds to 48 thousand and 1.1 million of total cells, respectively. The obtained minimal bubble radii are presented in Table 1, along with the estimated errors according to the Richardson extrapolation.
The evolution of bubble shape until the second collapse is given in Fig. 3 for a case with a grid resolution of 283 and 33 computational cells across the initial and minimum bubble radius, respectively. We chose this spatial resolution for the forthcoming simulations of bubblesphere interaction, since it offers a reasonable trade-off between the estimated minimal bubble radius error (− 3%) and the required computational times. According to Eq. (10), such an underestimation of the collapse radius would scale to roughly 14% overestimation of maximum bubble pressure P i . Fortunately, this is not the case when one does not assume uniform and adiabatic conditions within the bubble. The peak calculated pressures in the domain for the selected resolution were = p 1.29 max GPa, which only is 3.5% above the extrapolated value of 1.25 GPa.

Comparison to the Gilmore's model
Additionally, the obtained results from the axisymmetric model were compared to the solution of the Gilmore's equation. As can be seen from Fig. 4, the match between the curves until the final stage of the first collapse is very good. With = R 120 min nm, the axisymmetric model predicts notably smaller collapse radius than 160 nm of the Gilmore's equation. Since all parameters, initial, and boundary conditions between both models were kept the same, we attribute the difference in the final stage of the collapse to the fact, that the axisymmetric model does not assume uniform and adiabatic conditions within the bubble, but solves the energy equation (Eq. (4)) instead. Since heat and work flows go in opposite directions, we can expect some heat loss from the bubble into the surrounding liquid. In the Gilmore's model, this can be indirectly accounted for by reducing the polytropic index of air to a value between 1 and 1.4, that correspond to isothermal and adiabatic process, respectively. A good match in the minimum radius can be found with a polytropic index of ≈ 1.28, which can be seen in Fig. 4. According to Koukouvinis et al. [26], this parameter can be fitted, so that the results match experimental data for each case. Looking at the bubble collapse times, a deviation of <1% can be observed between both models. Even in this case the predicted rebound by the Gilmore's model is noticeably larger. This was already observed by previous researchers [13,23], the former of whom concluded that the energy radiated into the surrounding liquid at the beginning of the rebound is not correctly predicted by the Gilmore's model, which affects the magnitude of the rebound.

Bubble-sphere interaction
In this section, we address the key topic of the present study: a bubble-sphere interaction on a micro-scale. After a brief description of the numerical model setup, follow the results which cover the effects of a nearby sphere on a bubble collapse dynamics and vice versa.

Model setup
In the present study, we treat a sphere as a rigid boundary with no slip boundary condition, that can undergo rigid body motion due to external forces. These include hydrodynamic drag and a force arising from a pressure difference across a sphere, which are computed by numerical integration of pressure and shear stress over the surface of a sphere. Since we consider an axisymmetric case, only one degree of freedom needs to be resolvedtranslation of a sphere in the axial direction. Additionally, we utilize a diffusion based dynamic mesh smoothing method to adapt the numerical grid according to the sphere movement after each time step. Density of a sphere was set to = ρ 1100 s kg/m 3 , which is the average density of bacteria, such as E. coli [3].
A scheme of the considered phenomena is shown in Fig. 5. Here, two non-dimensional parameters δ and ε are defined as an initial stand-off distance and a sphere-bubble size ratio, respectively. Both, initial and boundary conditions were kept the same as in the case of a spherically symmetric bubble collapse (see Section 3.1.1), with initial bubble radius of = R 1 μ 0 m. We used similar numerical grids to the previous cases, with an exception of cells in the direct vicinity of a sphere, which were adapted to transition between curved walls of the sphere and the orthogonal cells in the rest of the computational domain.

Results
In total, we carried out nine simulations of bubble-sphere interaction, corresponding to all the value-pairs of non-dimensional bubblesphere distance = δ 0.9, 1.1, 1.3, and their size ratio = ε 0.5, 1.0, 1.5. For reference, the two most studied cases of bubble dynamics so far, an unbounded bubble and a bubble close to an infinite rigid wall correspond to = ∞ δ and = ∞ ε , respectively. We expect the influence of bubble-sphere interaction to grow with the values of δ decreasing and ε increasing.
Looking at Fig. 6, we can see the effect of both parameters on (a) the nondimensionalized bubble collapse radius = * R R R min min 0 and (b) the maximum velocity v max within the computational domain. For both, we can notice a similar trend of asymptotic decrease towards the value corresponding to the unbounded case ( = ∞ δ ). For all axisymmetric cases the equivalent minimum bubble radius R min is calculated from the minimal bubble volume by assuming the spherical shape, as is common in this type of studies. Results also seem to adhere to the initial expectations across different values of ε, which seems to mostly affect the solutions with smaller values of δ. In addition to increasing the bubble collapse radius and therefore decreasing the overall collapse intensity ( Fig. 6(a)), the presence of a spherical particle also prolongs the collapse time. Also interesting is the fact, that the obtained results already vary significantly on a relatively small parameter space considered in the present study. Maximum velocity ranges from ≈ 500 to 2100 m/s for attached bubbles, which indicates the occurrence of vastly different bubble shape development.
The latter can be seen from Fig. 7, where bubble shapes in time of collapse are given for all nine cases. Overall, the bubble shape evolution seems to congrue within each simulation set with the same δ. Therefore we can observe a longitudinal collapse shape for attached bubbles with = δ 0.9, and almost spherical collapses with = δ 1.3. On the other hand, the cases with = δ 1.1 exhibit a shape that falsely implies a postponed development of an axial jet away from the sphere, which is not the case (see Fig. 9(d) for the detailed velocity field in the beginning of the bubble rebound). Interestingly, only the case (c) with = δ 0.9 and = ε 1.5 shows a clear development of a high-speed annular jet towards the sphere's surface. One can also observe a consistently more pronounced deviation from the spherical shape with larger values of ε. In all other cases the jetting never fully develops and is diminished upon the bubble's rebound. This is consistent with the previous work of Liu et al. [35], who researched free and encapsulated bubbles in an ultrasound field and noticed that the shape oscillations are less likely to occur for micrometer sized gas bubbles since the surface tension suppresses the development of aspherical shape modes. In our case, we do not relate this to the direct effects of surface tension, which acts in a way to keep the bubbles spherical, as bubbles in the present case are still primarily driven by pressure and inertial forces ( MPa). We attribute this to the magnifying effects of surface tension on the internal bubble pressure, which causes larger bubble collapse radii R min as opposed to the case with neglected surface tension. This was already observed in Section 3.1, where R min varied by a factor of 1.6 between both cases. Because of this the ratio between the maximal and minimal bubble radii is significantly reduced and bubbles have a smaller chance to fully develop aspherical shape features, such as axial jets, etc.
We refer the reader to Figs. 8 and 9, where a bubble shape (black contour line) development and the corresponding pressure (upper half) and velocity (lower half) fields are presented for cases (c) and (e) of  In the first stage of collapse, an attached bubble predominantly shrinks in the radial direction and establishes a sharp contact angle with a sphere, as opposed to an initially blunt angle. A typical pressure field outside the farther bubble wall develops, which drives the bubble wall towards the sphere, causing a jet development shortly before the collapse ( = t 10.0 ns). Opposite to the classic case of a single axial jet, we can observe a development of a thin and very fast annular jet, which reaches the maximal velocity of 2100 m/s. Lately, Lechner et al. [28] conducted a numerical study of a bubble growth and collapse in extreme vicinity of a rigid wall ( = δ 0.048), which shows a similar development of a thin and fast (≈ 2000 m/s), although axial, jet towards the surface. After the jet hits the sphere's surface, causing a water hammer pressure in order of a few GPa, it radially flows through the toroidal part of the bubble, away from its center. Meanwhile the bubble rebounds and expands along the sphere's surface, which is later accompanied by the development of a secondary axial jet.
All three evaluated cases with = δ 1.1 show a similar bubble behavior. At first, they deform from the initial spherical into an egg-like shape ( Fig. 9(b)). Following that, we can observe a pressure and velocity field on the farther (left) half of the bubble that is similar to the previous case, which implies a start of the jet development towards the sphere. This is consistent with the recent research of Yamamoto and Komarov [51], who numerically addressed liquid jet directionality of a collapsing bubble near a gallium and silicone droplet, and showed a consistent jetting towards the gallium (higher density than water) and away from the silicone droplet (lower density than water). Nevertheless, the difference in our case lays on the other side, closer to the sphere, which also starts progressively accelerating towards the bubble center. A somewhat similar response can be seen in the previous literature [29,53], where larger bubbles close to spherical particles were experimentally and numerically studied. There we can see bubbles developing two opposite axial jets, which causes a bubble to form a toroidal shape. This is not the case in the present study, where further jet development is prevented by the previously discussed effects of surface tension ( ≈ R 0.12 μ min m). Shortly after collapse, at = t 10.2 ns, a pressure shock wave is emitted and the bubble rebounds. The magnitude of a shock wave is not uniform along the bubble's circumference, as is exhibits its maximum of 233 MPa ( = t 10.2 ns) towards the sphere. As we saw previously, a relatively small change in the distance and bubble-sphere size ratio can manifest in a vastly different micro-bubble dynamics. Now the question remains, to what degree this also affects    the loads exerted on a spherical particle and its movement. For this reason, the peak values of maximum and average quantities at sphere's surface, such as (a) pressure, (c) shear stress, and (d) temperature in relation to δ and ε are given in Fig. 10. The reason for including the peaks of the average values is to identify whether the load on a spherical particle is of a more local or global nature. In the former case, we would expect a significant difference between both values, while they should be of a similar magnitude in the latter case. Overall, the cases with attached bubbles ( = δ 0.9) show the largest and most spatially variable loads on a nearby particle. This can be best seen on subfigures (a) and (d), which show peak pressure and temperature values. From the latter we can also notice, that thermal load is only relevant for cases of attached bubbles, where a collapsing bubble is in the direct contact with the surface of a particle. The same cannot be said for compressive and shear loads. Looking at the peak average pressures, we can see that even for the least violent case with = δ 1.3 and = ε 1.5, they surpass the ambient pressure = ∞ p 10 MPa by a factor of two and additionally exhibit a peak in average shear stress of 0.36 MPa. Also convincing are peak pressure differences across a particle in Fig. 10(b), which show a very strong relation to peak maximal pressures. This further confirms the local nature of peak pressure loads, which occur when bubbleemitted shock waves propagate past them.
Maximal (a) velocity and (b) acceleration of a sphere are also included in Fig. 11. Here we can see that the order of magnitude for velocity is similar between the evaluated cases ranging from 5 to 30 m/ s, whereas the maximum acceleration differs by roughly 15-fold between the two extreme cases. Here, it is worth mentioning that the maximal force in the axial direction on a sphere grows with its size, as one might expect, but the corresponding acceleration of the particle adheres to opposite trend. The reason for this is in mass difference between the differently sized spheres, which goes from × − 5.8 10 16 kg for = ε 0.5 to × − 1.6 10 14 kg for = ε 1.5.

Discussion
In this section, the obtained results are discussed with respect to bacteria eradication in contaminated water. We focus mainly on the potentially destructive mechanical factors that arise from the collapse of a single microbubble in vicinity of a bacterial cell, although hydrodynamic cavitation might offer additional mechanisms for cell damage, such as nucleation within the cell or its envelope, expansion of internal gas bodies, and formation of free radicals, which could inflict lethal oxidative damage on bacteria [55].
At this point, it has to be noted that bacteria in reality tend to deform under external loads and adapt to their environment. Their stiffness, shape, and size can vastly differ across species and strains. Furthermore, cell's characteristics also depend on its growth phase and environment. By considering rigid suspended spherical particles in this study, we neglect their internal structure and ability to deform under load. Nevertheless, we consider that the present results are still relevant to address the likelihood of previously suspected mechanical effects to cause cell damage, such as impact of a high speed water jet, propagation of bubble emitted shock waves, shear loads, and thermal loads [55]. Firstly, looking at the previous studies of bubble-deformable sphere interaction on both, micro [16,17] and macro scale [22,15], reported relatively small areal deformations (<0.2%) of spheres with similar stiffness as herein considered bacterial cells. Secondly, the authors report a minor influence of sphere's deformability on bubble dynamics, which is further substantiated by the present study, as the considered cases with = δ 1.3 already show a minor influence of a nearby spherical particle on bubble dynamics (see Figs. 6 and 7).
As we have seen in Section 3.2.2, the peak loads exerted on a rigid spherical particle can exceed pressures of 1 GPa, shear stresses of 100 MPa, and temperatures of a few thousand Kelvins. If we compare these values against the modulus of elasticity for bacteria E. coli, E~0.5 cell MPa, and their envelopes, − E~10 100 env MPa [9,20,44], we can identify a high potential for cell's structural disruption. As expected, the highest levels of particle loads were found for the cases of attached microbubbles. Unfortunately, the existence of bacteria-attached microbubbles has not been yet researched. Therefore we decided to select two more plausible cases for further evaluation. The first case, corresponding to = δ 1.3 and = ε 0.5, is chosen as a scenario with the minimal impact on bacterium, whereas the second case, with = δ 1.1 and = ε 1.5, is selected as the most realistic scenario with maximal loads exerted on bacterial cell.
First, we take a look at the temporal development of sphere's velocity and acceleration for both selected cases ( Fig. 12(a)). Here, we can see that initially sphere starts to accelerate towards the bubble and reaches the peak velocity of ≈ 5 to 15 m/s. Upon bubble compression, the pressure field outside the bubble eventually surpasses the ambient pressure = ∞ p 10 MPa which causes the pressure gradient at the location of a sphere to change direction towards the bubble. This causes the bubble to accelerate in the opposite direction. The peak acceleration is achieved at the time of shock wave propagation past the sphere, and ranges in the order of a few 10 10 m/s 2 , which corresponds to the maximal force in the order of 100 μN.
Not surprisingly, the peak pressures and shear stresses exerted on a sphere's surface coincide with peak accelerations, which further confirms the importance of considering the nonlinear compressibility of water and bubble emitted shock waves in studies of similar type. Additionally, maximum shear loads exhibit a lesser local maximum that follows the propagation of a shock wave and is due to sphere reaching its peak velocity away from the bubble. As we can see from Fig. 12(b) and (c), the loads during the second bubble collapse exhibit a relatively small damage potential in comparison with the first one, as their peaks differ by a roughly one order of magnitude. Additionally, thermal load also seems to be irrelevant in the case of non-attached bubbles, as we can see from Fig. 12(d), that temperature change at the location of a sphere is negligibly small, and well below the threshold of cell death, i.e. 70°C for E. coli [36]. The latter mechanism seems to be only applicable for attached microbubbles or microbubbles in extreme vicinity of bacteria.
Similar conclusions can be also drawn for an impact of a high speed water jet. As we have seen in Section 3.2.2, the vicinity of a rigid spherical particle alone does not induce bubble jetting towards its surface, which would be even less likely to occur for an actual cell due to its compliant nature. An exception to this is the case with = δ 0.9 and = ε 1.5, which also exhibits other damaging mechanisms (high temperatures, pressures, and shear loads), but is not considered as very relevant in a more realistic setting. From this we can conclude, that a high-speed jet impact and related cell damage are unlikely to occur on the investigated spatial scale. However, microbubble jetting could potentially still develop in the presence of other jet drivers, e.g. shock wave induced microbubble collapse in vicinity of bacterium, which would be an interesting topic for further research.
Looking at Fig. 12(c) reveals peak shear loads of 2 and 5 MPa for both chosen cases. As mentioned before, the peaks can be attributed to shock wave-induced particle acceleration, which results in an order of magnitude higher values compared to the bubble wall-induced flow around the sphere. By comparing peak values to shear modulus of both cell membranes of E. coli, ≈ G 10 MPa (an approximation based on the work of Hwang et al. [20], who report areal stiffness of the outer and inner E. coli membrane to be 237 mN/m and 240 mN/m in their relaxed state, respectively), we can observe a high damage potential. This is especially promising due to the fact, that the outer membrane starts to exhibit a significant stress softening after it undergoes an initial areal expansions in the order of a few percent [20]. Based on this and the findings of previous researchers [27,54] we see bubble collapse induced shear loads as one of the most likely mechanisms for bacterial cell disruption, which should be further researched in the future.
As recently reported by Pandur et al. [39], hydrostatic pressure loads of more than 120 MPa are needed for destruction of liposomes, spherical vesicles that can act as a model system for bacterial lipid bilayers. Therefore the magnitude of peak pressure alone does not seem likely to play a relevant role for bacteria destruction in the case of hydrodynamic cavitation. Perhaps one of the key attributes of peak compressive loads induced by shock wave propagation is their spatial variability, which means very high pressure gradients at the cell's   11. Overall maximum values of (a) velocity and (b) acceleration of a spherical particle for the considered δε value pairs. surface, in the order of 10 14 Pa/m. The latter is expected to induce a strong response in bacterial envelope, causing high local stresses and deformations, and could also be potentially harmful to even smaller pathogens, such as viruses [24], that normally measure from 80 to 300 nm in diameter [36].
To further support our findings, shock wave propagation is showed for both selected cases in Fig. 13. The left column corresponds to the Case 1 and right to the Case 2, respectively. For both, four pressure fields are given that correspond to a similar time in bubble-sphere interaction. In the first row, we can observe a spherical wave front for the Case 1, whereas the bubble-sphere interaction seems to magnify the peak pressures (161 MPa) in the direction towards the sphere for Case 2. Looking at the second row, we can notice a significant pressure front magnification upon its impingement on the sphere's surface. In Case 2 the peak value of 181 MPa even surpasses the one at the time of shock wave emission (161 MPa). Following that, a part of the initial pressure wave is reflected back towards the rebounding bubble, which forms a secondary wave front and a low pressure area between the bubble and a sphere. Subfigure (c2) implies formation of a tensile wave following the secondary wave front. Unfortunately, the magnitude of a tension wave could not be resolved by the current solver, since it does not allow for negative absolute pressures. In the last row of Fig. 13, we can see how the main shock front propagates past the sphere, still causing high pressure gradients at its surface even though peak pressures decrease to roughly 30 MPa.
The results of the present study show increased magnitudes of bubble-bacteria interaction with larger values of ε, which implies that cavitating flow with smaller bubbles could be more aggressive towards bacteria. This can be achieved by increasing flow velocities and system pressures, which is already known from the previous material erosion studies [10]. These conclusions are of course merely a very rude estimate, since we are considering a rigid object in an interaction with a single bubble. On the other hand, past research [45,25] point to superiority of supercavitating flow in bacteria eradication. Here we suspect a completely different type of bubble-bacteria interaction, namely the bacteria passing through a large stable bubble and hence being exposed to a very rapid pressure drop and consequent pressure recuperation. This type of effect can only be modeled once the bacterium's interior structure and envelope properties are adequately considered, which is one of the goals for the future research.

Conclusions
The aim of the present paper was to research the behavior of single cavitation microbubbles ( = R 1 μ eq m) in vicinity of a submerged sphere on a similar spatial scale, and to determine how this affects the mechanical loads exerted on a sphere. Due to a small spatial (~1 μm) and temporal (~10 ns) scales of the considered phenomena a purely numerical approach was used. We performed numerical simulations by using a FVM-VOF approach with inclusion of surface tension, viscous effects, and nonlinear compressibility of water to resolve bubble emitted shock waves upon its collapse. The final model was established based on a grid convergence study and comparison to the Gilmore's model for a simpler case of a collapsing unbounded spherical bubble.
We presented results for various value-pairs of two nondimensional parameters, δ and ε, which show vastly different bubble collapse dynamics across the considered parameter space, from the development of a fast thin annular jet towards the sphere to an almost spherical bubble collapse. Although we can find some similarities in bubble shape progression to previous studies on larger bubbles, we can generally notice that bubble jetting is much less likely to occur on the considered scale due to cushioning effects of surface tension on the intensity of the collapse. Overall, the results show that the mechanical loads on a spherical particle tend to increase with a sphere-bubble size ratio ε, and decrease with their distance δ.
The obtained results were further discussed with respect to bacteria eradication in contaminated water. We identified the bubble emitted shock wave propagation past the bacterial cell to be the most likely mechanism for inflicting structural damage, as it exhibits large pressure peaks in the order of a few hundred MPa, gradients of 100 MPa/μm, and induces the peak shear loads of a few MPa.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Fig. 13. Pressure fields showing bubble emitted shock wave propagation at the time of (a) shock wave emission after bubble's collapse, (b) shock wave impact on a sphere, (c) refraction of a shock wave from the sphere's surface and formation of a tensile pressure wave, and (d) propagation of a tensile wave through a rebounding bubble. Left column corresponds to the Case 1 and right to the Case 2, respectively. Please mind a unique pressure legend range for each subfigure.