# Dynamics of two interacting acoustic bubbles at short separation distances

Gabriel Regnault, Alexander A Doinikov, Cyril Mauger, Philippe

Blanc-Benon, Claude Inserra

► To cite this version:

Gabriel Regnault, Alexander A Doinikov, Cyril Mauger, Philippe Blanc-Benon, Claude Inserra. Dynamics of two interacting acoustic bubbles at short separation distances. Physics of Fluids, 2023, 35 (3), pp.037116. 10.1063/5.0135370. hal-04829043

HAL Id: hal-04829043

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

Submitted on 10 Dec 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# Dynamics of two interacting acoustic bubbles at short separation distances

**Accepted Manuscript:** This article has been accepted for publication and undergone full peer review but has not been through the copyediting, typesetting, pagination, and proofreading process, which may lead to differences between this version and the Version of Record.

Cite as: Physics of Fluids (in press) (2022); https://doi.org/10.1063/5.0135370

Submitted: 18 November 2022 • Accepted: 19 February 2023 • Accepted Manuscript Online: 20 February 2023

ID Gabriel Regnault, ID Alexander A Doinikov, ID Cyril Mauger, et al.

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

View Online

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

Export Citation

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

CrossMark

## ARTICLES YOU MAY BE INTERESTED IN

*Flow in multi-layered vegetated compound channels with different bank slopes*

*Physics of Fluids (2023); https://doi.org/10.1063/5.0142400*

*Numerical simulations of miscible displacement in an inclined channel by lattice Boltzmann method*

*Physics of Fluids (2023); https://doi.org/10.1063/5.0135734*

*The effect of ribs on the sound radiation directivity of rectangular plates*

*Physics of Fluids 34, 122109 (2022); https://doi.org/10.1063/5.0127918*

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

**APL Machine Learning**

Machine Learning for Applied Physics
Applied Physics for Machine Learning

**First Articles**
**Now Online!**

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

AIP
Publishing
# Dynamics of two interacting acoustic bubbles at short separation distances

Gabriel Regnault,¹,² Alexander A. Doinikov,¹,² Cyril Mauger,¹ Philippe Blanc-Benon,¹ and Claude Inserra², a)

¹) Univ Lyon, CNRS, Ecole Centrale de Lyon, INSA Lyon, Univ Claude Bernard Lyon 1, LMFA, UMR5509, 69130, Ecully France

²) Univ Lyon, Université Lyon 1, Centre Léon Bérard, INSERM, LabTAU, F-69003 LYON, France

We investigate experimentally the attraction between two closely spaced, oscillating microbubbles. Above a certain value of the applied acoustic field, the bubbles jump to a new equilibrium location where they are separated by a thin fluid layer the thickness of which is much smaller than the bubbles radii. We demonstrate that this new equilibrium is caused by the sign reversal of the radiation interaction force acting between the two bubbles, attributed to the multiple rescattering of the waves emitted by the bubbles. Theoretical investigation reveals that a new stable equilibrium appears at short distances, resulting in a quasi contacting bubble pair.

**Introduction** — The dynamics of a single bubble has been widely investigated experimentally and theoretically when dealing with its spherical, translational or axisymmetric shape oscillations¹⁻⁵, the ultrasound-induced (primary) radiation force acting on it⁶,⁷, either in confined or free liquid⁸, or when investigating the bubble-induced flow motion in its vicinity⁹,¹⁰. The bridge from the single to collective bubble dynamics is still an ongoing challenge as few experimental investigations are reported on the dynamics of bubble clouds¹¹. Numerical analysis of the bubble cloud dynamics allows capturing the structure of the cavitation cloud by including an interaction force acting between bubbles¹². This force, also called secondary radiation force, takes its origin in the incoming acoustic waves from surrounding bubbles that act as secondary sound emitters⁶,¹³. When bubbles are driven below their resonance frequencies and submitted to a high-amplitude ultrasound field, they form particular cloud structures, called streamers, where intense sonochemical effects occur¹¹. In the case of bubbles driven above their resonance frequencies, they can agglomerate and form stable clusters (called bubble grapes), where the mean separation distance between bubbles is comparable to their sizes¹⁴. The formation of the bubble grapes has been theoretically attributed to the sign reversal of the radiation interaction force between bubbles, changing their motion from attraction to repulsion. The sign reversal of the interaction force has been theoretically attributed either to multiple scattering effects¹⁵⁻¹⁸ or to the existence of a particular (called transition) frequency in the two-bubble system that can equal the driving angular frequency for a given interbubble distance¹⁹. Another consequence of the sign reversal of the secondary radiation force is the theoretical existence of a small separation distance at which two bubbles can be maintained stable. The sign reversal of this force was previously observed experimentally for one freely moving bubble and the other trapped in a yield-stress fluid²⁰ or attached to a substrate²¹. Stable bubble pairs have been observed between two millimetric bubbles at an interbubble distance

of about 20 bubble radii²², between microbubbles forming lines or two-dimensional arrays in a high-frequency (1 MHz) field²³, or in confined geometries, where bubbles formed bounded bubble crystals, at an interbubble distance of about 10 bubble radii²⁴.

In this paper, we show that two attracting bubbles can be stabilized at a short distance, where they are separated by a thin fluid layer whose thickness is much smaller than the bubble radii. This effect is theoretically attributed to multiple scattering effects between the bubbles. The sign of the radiation interaction force is reversed as the bubbles move toward each other, which allows stabilizing the bubble pair at a new, nearly contacting equilibrium distance. To the best of our knowledge, this is the first experimental validation of the sign reversal of the secondary radiation force occurring at a short range between two stable bubbles, for which radial oscillations, phase shift and separation distances can be simultaneously measured.

**Experimental setup** — A detailed description of our experimental setup can be found in our previous works²⁵. Briefly, a train of gas bubbles (~100 µm in radius with a separation distance around 1 mm) is created at the tip of a thin capillary driven by a microfluidic pressure controller (Elveflow, OB1 Mk3). The bubbles are trapped at pressure nodes of a high-frequency (HF) standing-wave ultrasound field ($f_{HF} = 1$ MHz) when a HF transducer (Sofranel IDMF018), slightly immersed at the top of the water tank, is switched on. The theoretical interbubble distance at rest (without driven bubble oscillations) is therefore half the HF wavelength $λ_{HF}/2 \sim 750$ µm. Only two bubbles trapped at successive pressure nodes are kept, in order to isolate them from any external disturbances; see Fig. 1(a). To do so, the other surrounding bubbles are removed by hand using a thin wire, prior to any experimental acquisition. The effective bubble trapping location is deduced from the balance between the buoyancy force $\mathbf{B}_j$ ($j = 1, 2$) and the HF primary radiation force $\mathbf{F}_{pj}$ acting on the $j$th bubble. The buoyancy force and the primary radiation force on the $j$th bubble are given by²⁶

$$ \mathbf{B}_j = \frac{4}{3} \pi \rho g R_{0j}^3 \mathbf{e}_z, \quad (1) $$

<footnote>a) Electronic mail: claude.inserra@inserm.fr</footnote>
$$ \mathbf{F}_{pj} = \frac{\pi(P_{HF})^2 R_{0j}(\omega_j^2/\omega_{HF}^2 - 1) \sin(2k_{HF}z_j)}{\omega_{HF}\rho c [(1 - \omega_j^2/\omega_{HF}^2)^2 + (k_{HF}R_{0j})^2]} \mathbf{e}_z, \quad (2) $$

where $\omega_{HF}$ and $P_{HF}$ are the angular frequency and the pressure amplitude of the HF wave, respectively, $z_j$ is the position of the center of the jth bubble on the line joining the bubble centers, $R_{0j}$ is the equilibrium radius of the jth bubble, $\rho$ is the liquid density, $k_{HF} = \omega_{HF}/c$ is the HF wavenumber, where $c$ is the speed of sound in the liquid, $g$ is the gravity, $\mathbf{e}_z$ is the unit vector along the $z$ axis and $\omega_j = \sqrt{3\gamma P_0/\rho R_{0j}^2 + 2(3\gamma - 1)\sigma/\rho R_{0j}^3}$ is the angular resonance frequency of the jth bubble, where $P_0$ is the liquid hydrostatic pressure, $\gamma$ is the ratio of specific heats of the gas and $\sigma$ is the surface tension.

<figure><img src="image_7.png" /><figcaption>FIG. 1. Schematic representation of the bubble trap. Ideally, two bubbles can be trapped at the pressure nodes $z_j^*$ of a HF ultrasound field (dashed line), and separated by the theoretical interbubble distance $\lambda_{HF}/2$. (a) When switching the HF and LF transducers on, a new equilibrium location is reached, resulting from the balance between the buoyancy $\mathbf{B}_j$, the HF primary radiation force $\mathbf{F}_{pj}$ and the secondary radiation force $\mathbf{F}_{sj}$. (b) At a distance of $\lambda_{HF}/8$ from the theoretical pressure node, the HF primary radiation force (dashdotted line) reaches a maximum. Beyond this boundary, the radiation interaction force increases while the HF primary radiation force begins to decrease. This creates a forbidden zone in which the bubbles are supposed to encounter and coalesce.</figcaption></figure>

It is worth noting that the expression of the primary radiation force in Eq. (2) is valid for linearly oscillating bubbles that are smaller than the acoustic wavelength, $R_{0j} \ll \lambda_{HF}$. For the studied here bubbles ($\sim 150 \mu m$ in radius) trapped in the HF field with $\lambda_{HF} = 1.5 \text{ mm}$, the ratios $R_{0j}/\lambda_{HF}$ are of the order of $\sim 0.1$. In order to estimate the influence of the bubble size on the primary radiation force for arbitrary values of the ratio $R_{0j}/\lambda$, Annamalai et al.²⁷ derived a second-order mean force (the primary radiation force) exerted by a standing-

wave ultrasound field on a spherically-oscillating bubble, without restriction on the bubble size in comparison with the acoustic wavelength. Their model predicts that monopole term (corresponding to Eq. (2)) is the primary contributor to the primary radiation force up to $R_{0j}/\lambda \sim 0.02$. For higher $R_{0j}/\lambda$ ratios up to 0.1, adding higher-order terms results in the decrease of the amplitude of the radiation force by a factor of $\sim 2/3$, without modification of the spatial evolution of the force. For the sake of simplicity, the primary radiation force including only the monopole term (Eq. (2)) is therefore considered in our modeling.

Bubble oscillations are driven by a low-frequency (LF) ultrasound field ($f_{LF} = 31.9$ kHz) induced by a Langevin transducer (Sinaptec) located at the bottom of the tank. The LF wavelength is significantly larger than the interbubble distance, $\lambda_{LF} \gg \lambda_{HF}/2$, so that the bubble pair is assumed to be driven by a uniform LF pressure field. It is worth noting that assessing the pressure amplitude of the LF field is a challenging experimental task. Because of the LF frequency that requires large-sized hydrophones, measuring the pressure within the water tank would automatically induce disturbances of the acoustic field. Therefore, we can only evaluate the pressure amplitude through indirect measurements of the bubble oscillations or by solving an appropriate force balance equation, as the one that we propose in the present paper. As a result of the induced radial oscillations, the bubbles experience an interaction force that moves them to a new equilibrium location (Fig. 1(a)). This equilibrium location, where the bubbles are separated by the distance $d = z_2 - z_1$, results from the balance between the buoyancy $\mathbf{B}_j$, the HF primary radiation force $\mathbf{F}_{pj}$ and the secondary radiation force $\mathbf{F}_{sj}$. This allows quantifying the interaction force acting on the bubble pair²⁵. According to Eq. (2), $F_{pj}$ depends on $z_j$ as $\sin(2k_{HF}z_j)$, while the HF field depends on $z$ as $\sin(k_{HF}z)$ (Fig. 1(b)). As a result, for bubble 2, $F_{pj}$ reaches a maximum at a distance of $\lambda_{HF}/8$ from the HF pressure node $N_2$. Below the position corresponding to this local maximum, the amplitude of $F_{pj}$ decreases while the interaction force acting between the bubbles still increases. This results in a predominance of the mutual interaction force $F_{sj}$ over the trapping force $F_{pj}$, which forces the bubbles to collide. We call this region, where bubble coalescence is supposed to systematically occur, forbidden zone. The experimental procedure for the investigation of interacting bubble pairs is as follows. After trapping a bubble pair, we perform several successive recordings for every investigated LF pressure. Every recording is performed using a CMOS camera (Vision Research V12.1) equipped with a 12× objective lens and backlight illumination, with a frame size of 128 × 256 pixels (1 pixel $\sim 5 \mu m$) at the 130 kfps frame rate (with $\sim 4$ samples per acoustic period). From every recording, the dynamics of the bubble oscillations, $R_j(t) = R_{0j} [1 + \epsilon_j \cos(2\pi f_{LF}t + \varphi_j)]^{25}$, and the location of the bubble centres of mass are quantified through the normalized radial expansion $\epsilon_j = (R_j^{\text{max}} - R_{0j})/R_{0j}$ and
the position $z_j$ of the center of the jth bubble. Once a video recording is completed (with typically a duration of hundreds of milliseconds), the LF pressure is increased, which results in the motion (approach) and a new equilibrium position of the two bubbles. A new video is then recorded. This process is repeated until the sequence of bubble approach ends. The step-by-step approach of the bubble pairs is shown in recordings provided in Supplementary Material, corresponding to the cases shown in Figs. 2 and 5(a,b).

**The bubble jump** — Hundreds of bubble pairs have been trapped and investigated with the following general conclusions. The interaction of two bubbles whose equilibrium radii $R_{01}, R_{02}$ encircle the LF resonant radius $R_{LF}^{res} \sim 110$ µm always results in bubble repulsion. The LF resonant radius $R_{LF}^{res}$ is calculated from the formula for the bubble resonance frequency $\omega_j$ setting $\omega_j = \omega_{LF}$. Bubbles are therefore rapidly repelled away from each other, moving out of the camera field. This case will no longer be discussed. The interaction of two bubbles whose equilibrium radii are smaller than the resonance radius results in bubble attraction. Note that from the formula for $\omega_j$ it follows that a bubble radius smaller than the resonance radius means that the resonance frequency of this bubble is higher than the driving frequency. Measuring the interbubble distance as a function of the LF pressure amplitude allows quantifying the secondary radiation force that they experience$^{25}$. However, as soon as one of the bubble centres of mass reaches the boundary of the forbidden zone, coalescence systematically occurs. This phenomenon limits the achievable interbubble distance $d_{lim} = \lambda_{HF}/4$ that is roughly 3.5 times the mean equilibrium radius of the largest investigated bubbles ($\sim 100$ µm in radius). The interaction of two bubbles whose equilibrium radii are larger than the resonance radius also results in bubble attraction. However, when one of the bubbles reaches the boundary of the forbidden zone, the two bubbles can suddenly jump to a new equilibrium location. The bubbles are then so close that they are separated by a thin fluid layer whose thickness is much smaller than the bubble radius.

Figure 2 (Multimedia view) illustrates the case of two interacting bubbles of equilibrium radii $R_{01} = 170$ µm and $R_{02} = 139$ µm that jump to a new, nearly contacting equilibrium, where the fluid layer thickness is $\sim 30$ µm. A small lateral motion of the bubbles can be noticed in Fig. 2. This lateral motion comes from the secondary interaction force that acts along the line joining the centers of the two bubbles. Before applying the LF pressure the bubbles are in general not perfectly aligned: as the bubble is trapped in the nodal plane of the HF one-dimensional standing wave field, some shift from the vertical alignment can occur. However, as soon as the LF pressure is applied, the bubbles are attracted and realign. When the jump occurs, disturbances in the local pressure field can slightly shift the bubble horizontal locations. Because these horizontal motions are negligible in comparison with the vertical ones, they are disregarded

<figure><img src="image_8.png" /><figcaption>FIG. 2. Snapshots of a typical sequence of bubble approach in the case of two bubbles with radii larger than the resonant radius ($R_{01} = 170$ µm, $R_{02} = 139$ µm). When increasing the LF pressure amplitude, they attract and come closer to each other, until one of the bubble centres reaches the boundary of the forbidden zone. For this threshold (here indicated by a red vertical line), the bubbles suddenly jump to a new equilibrium location (Multimedia view).</figcaption></figure>

in the present study. In addition, the geometry of the two-bubble system is axisymmetric in respect to the bubble center-to-center axis. We make sure experimentally that the focal plane of the camera contains this symmetry axis. Therefore, the film thickness between the two bubbles for any interbubble distance can be measured without ambiguity in the focal plane of the camera.

The analysis of the phase shift $\Delta\varphi$ between the bubbles’ radial oscillations indicates that the bubble pair oscillates in phase until the jump occurs (Fig. 3(b)). After the jump, the phase shift reaches $\pi/2$. Concerning the normalized radial expansion $\epsilon_j$ (Fig. 3(a)), the parameter $\epsilon_2$ for bubble 2 increases when the interbubble distance decreases, equivalently when the applied LF pressure increases. However, for bubble 1, the parameter $\epsilon_1$ decreases from $\epsilon_1 = 0.13 \times 10^{-4}$ for $d = 670$ µm to $\epsilon_1 = 0.1 \times 10^{-4}$ for $d = 348$ µm. Therefore, the radial oscillations do not increase anymore when the interbubble distance reaches its smallest value. This effect is amplified for the cases presented in Fig. 5. This illustrates the acoustic interaction between the bubbles: two bubbles larger than the resonant size oscillate in phase with each other but out-of-phase with the external acoustic field. The resulting surrounding pressure field is therefore lower than the pressure expected in the case of uncoupled bubbles. Moreover, the radial expansion parameters $\epsilon_j$ never exceed $10^{-2}$ (Fig. 3(a)), meaning relatively small radial oscillation amplitudes.

**Sign reversal of the radiation interaction force** — The radiation interaction force between two bubbles was first derived by Bjerknes$^6$ in the form:

$$F_s = 2\pi\rho|C|^2 G_B e_z, \quad (3)$$

where $C$ is the complex amplitude of the liquid velocity potential, $e_z$ is the unit vector directed from the center of bubble 1 to the center of bubble 2 and the coefficient $G_B$ ($B$ stands for Bjerknes) is a function of the mechanical properties of the system (bubble equilibrium radii,
<figure><img src="image_9.png" /><figcaption>FIG. 3. (a) Response amplitude of the bubble radial oscillations normalized by the equilibrium radius as a function of the interbubble distance. (b) Phase shift $\Delta\varphi$ between the bubble radial oscillations as a function of the interbubble distance.</figcaption></figure>

damping coefficient), the driving angular frequency and the interbubble distance. This parameter is given by:

$$G_B = \frac{R_{01}R_{02}}{d^2} \frac{(1 - \omega_1^2/\omega^2)(1 - \omega_2^2/\omega^2) + \delta_1\delta_2}{[(1 - \omega_1^2/\omega^2)^2 + \delta_1^2][(1 - \omega_2^2/\omega^2)^2 + \delta_2^2]}, \quad (4)$$

where $\omega = \omega_{LF}$ is the acoustic angular frequency of the radial oscillations, and $\delta_j$ is the damping coefficient for bubble j (given by Eq. (A.5) in the Appendix). Equation (4) is valid when the spacing between the bubbles is much larger than their radii, or equivalently, when $d/(R_{01} + R_{02}) \gg 1$. Repulsion occurs when $\omega_{LF}$ lies between $\omega_1$ and $\omega_2$, hence when bubbles pulsate out-of-phase with each other. Moreover, Eq. (4) predicts an infinite value of the interaction force when the bubbles come sufficiently close, i.e. when $d/(R_{01} + R_{02}) \to 1$. Zabolotskaya¹⁵ derived an interaction force accounting for the influence of the scattered field of each bubble on the other one. The associated coefficient $G = G_Z$ of the interaction force is a refinement of the model of Bjerknes

up to the first order in the small parameter $R_{0i}/d$:

$$G_Z = \frac{R_{01}R_{02}}{d^2} \frac{\left(\frac{\omega_1^2}{\omega^2} - 1 + \frac{R_{01}}{d}\right) \left(\frac{\omega_2^2}{\omega^2} - 1 + \frac{R_{02}}{d}\right) + \delta_1\delta_2}{\left|\left(\frac{\omega_1^2}{\omega^2} - 1 - i\delta_1\right) \left(\frac{\omega_2^2}{\omega^2} - 1 - i\delta_2\right) - \frac{R_{01}R_{02}}{d^2}\right|^2} \quad (5)$$

For $d \to 0$, the terms inversely proportional to $d$ and $d^2$ become dominant in Eq. (5). Keeping only these terms in both numerator and denominator, one obtains $G_Z \to 1$. Therefore, the interaction force takes finite values at short separation distances and is hence more adapted to predict the interaction of two closely-spaced bubbles. Particularly, if both bubbles are larger than the resonance size, and if their resonance frequencies are close to $\omega$, the interaction force may change from attraction to repulsion. This change is related to the shift of the modal frequencies of the two-bubble system resulting from the bubble approach. The shift in the modal frequencies leads to out-of-phase bubble oscillations, hence to the sign reversal of the interaction force. However, this modeling can not capture our experimental observations as the measured phase shift between the bubble oscillations never exceeds $\pi/2$ (see Fig. 3(b)). Later on, an approach accounting for multiple scattering of sound between bubbles was developed by Doinikov and Zavtrak with no restriction on the separation distance between the bubbles¹⁶. The resulting interaction force is expressed in terms of an infinite sum of powers of the small parameter $R_{01}/d$, with the coefficient $G = G_D$ given by

$$G_D = -\mathrm{Re} \left( \sum_{n=0}^{\infty} \frac{n+1}{\eta_{1(n+1)}} a_n^{(11)*} a_{n+1}^{(11)} \right), \quad (6)$$

where Re means "the real part of", the asterisk denotes the complex conjugate, $\eta_{1(n+1)}$ is a constant coefficient and $a_n^{(11)}$ are the linear scattering coefficients of bubble 1. The numerical implementation of Eq. (6) and the expressions of the coefficients $\eta_{1(n+1)}$ and $a_n^{(11)}$ are considered in the Appendix. This model predicts that, when both bubbles are larger than the resonance size, the interaction force can change from attraction to repulsion. A stable equilibrium between bubbles is therefore expected. The complexity of the mathematical formulation of the terms in Eq. (6) does not allow for a simple interpretation of the conditions leading to the sign reversal of the interaction force (see the Appendix). The three above-mentioned models can predict the amplitude and sign of the interaction force between two bubbles through the coefficient $G$, as soon as their equilibrium radii $R_{0i}$ and their separation distance $d$ are known. Figure 4 shows the evolution of the amplitude $G$ of the secondary radiation force for the three models developed by Bjerknes, Zabolotskaya, and Doinikov and Zavtrak, for the same experimental case as in Fig. 2. The circle markers correspond to the predictions of Eqs. (4,5,6) for the measured radii $R_{0i}$ and
<figure><img src="image_10.png" /><figcaption>FIG. 4. Evolution of the radiation interaction force for $R_{01} = 170$ µm and $R_{02} = 139$ µm, through the parameter $G$, for the models of Bjerknes, Zabolotskaya, Doinikov and Zavtrak. The circle markers correspond to the predictions of Eqs. (4,5,6) for the measured radii $R_{0i}$ and measured normalized distance $L = d/(R_{01} + R_{02})$. When varying the normalized distance for $L \in [1, 4]$ and the equilibrium radii over the sequence of recordings, the predictions of Eqs. (4,5,6) are provided by dashed, dashdotted and solid lines.</figcaption></figure>

measured normalized distance $L = d/(R_{01} + R_{02})$. In order to reveal the influence of the interbubble distance on the amplitudes of the interaction force, the predictions of Eqs. (4,5,6) are shown (see dashed, dashdotted and solid lines) when the normalized interbubble distance $L$ is varied in the range $[1-4]$ for fixed values of the equilibrium radii. The equilibrium radii have been interpolated over the range of the measured $R_{0j}$ as the bubbles slightly grow (a few micrometer in radius) due to rectified diffusion. As one can see in Fig. 4, the secondary radiation force as proposed by Bjerknes rapidly goes to infinity for small interbubble distances, which we deal with in this experiment. While finite values of the force for $L \to 1$ are obtained with the model of Zabolotskaya, the interaction force always remains positive (attractive), preventing the occurrence of a stable equilibrium bound between bubbles. When accounting for multiple scattering effects (the model of Doinikov and Zavtrak), an inversion of the sign of the interaction force occurs systematically for small interbubbe distances. The values for which $G_D < 0$ always correspond to bubble equilibrium subsequent to the jump moment. Similar conclusions have been obtained for various bubble pairs, and two additional examples are provided in Fig. 5 (Multimedia view).

Bubbles location after the jump — Once a suitable the-

oretical model for the secondary radiation force has been established, one can predict the location of the bubbles after they have jumped. The bubbles are fixed in positions where the condition $\mathbf{F}_{sj} + \mathbf{F}_{pj} + \mathbf{B}_j = \mathbf{0}$ is fulfilled, using Eqs. (1,2,3) for the three applied forces and Eq. (6) for the term involved in the secondary radiation force. According to the mathematical formulation of the interaction force, the force balance equations of the bubbles are coupled. Solving these coupled equations allows determining the equilibrium locations $z_j$ and the interbubble distance $d = z_2 - z_1$. Without applying the LF acoustic field (Fig. 6(a)), the bubbles are trapped near the pressure nodes $N_j$, which are stable equilibrum locations. In order to investigate the stability of the equilibrium at the location of the pressure antinode $A$, we fix the equilibrium location of bubble 1 and evaluate the magnitude of the forces applied to a virtual bubble 2 in the vicinity of the position $A$ (the dashed bubble in Fig. 6(a)). In the case $P_{LF} = 0$, the position $A$ is unstable. With increasing LF pressure, the bubbles attract each other until the jump occurs. For two bubbles with radii $R_{01} = 170$ µm and $R_{02} = 139$ µm trapped in a HF field with a pressure amplitude of 199 kPa, the jump occurs for a LF pressure of about 1.5 kPa. At the jump moment, bubble 2 moves to the position $A$, which was initially unstable. The interbubble distance is $d = 339$ µm (Fig. 6(b)), resulting in a fluid layer 30 µm thick between the bubbles' interfaces, similarly to the experimental case shown in Fig. 2. The analysis of the total force applied to bubble 2 in the vicinity of the position $A$ reveals that this position becomes stable when the interaction force reverses. The location of bubble 2 near the position $A$ becomes stable only if the LF pressure is high enough (resulting in the sign reversal of the interaction force) and if the HF pressure is high enough too (otherwise bubble 2 collides with bubble 1 as the distance between the positions $A$ and $N_1$ is smaller than $R_{01} + R_{02}$). Considering that the governing equations depend on many parameters, the only way to establish the stability conditions for the position $A$ of bubble 2 is to make numerical calculations for varying values of the LF pressure, just as it happens in our experiment. Doing so, we obtain that for small values of the LF pressure, the shift of bubble 2 from the position $A$ produces a local force acting away from the position $A$, which means that the position $A$ is unstable. However, above a given value of the LF pressure, the shift of bubble 2 from the position $A$ produces a local force that returns bubble 2 to the position $A$, which means that the position $A$ is stable. To the best of our knowledge, this is the first experimental observation of the modification of the stability of equilibrium locations for bubbles experiencing linear oscillations in an acoustic field of moderate amplitude.
<figure><img src="image_11.png" /><figcaption>FIG. 5. Two additional cases are shown for bubble pairs with different bubble radii: (a) $R_{01}$=156 µm, $R_{02}$=142 µm (Multimedia view) and (b) $R_{01}$=160 µm, $R_{02}$=130 µm (Multimedia view). First row: Snapshots of a typical sequence of bubble approach. Second row: Response amplitude of the bubble radial oscillations normalized by the equilibrium radius as a function of the interbubble distance. Third row: Phase shift $\Delta\varphi$ between the bubble radial oscillations as a function of the interbubble distance. Fourth row: Evolution of the radiation interaction force in terms of the parameter $G$ for the models of Bjerknes, Zabolotskaya, Doinikov and Zavtrak. The circle markers correspond to the theoretical predictions for the measured radii $R_{0i}$ and measured normalized distance $L = d/(R_{01} + R_{02})$. When varying the normalized distance for $L \in [1, 4]$ and the equilibrium radii over the sequence of recordings, the predictions of Eqs. (4,5,6) are provided by dashed, dashdotted and solid lines.</figcaption></figure>
(a) $P_{LF} = 0$

(b) $P_{LF} = 1.5$ kPa

<figure><img src="image_12.png" /><figcaption>FIG. 6. Upper figures: Location of a bubble pair in a HF trapping field with an amplitude of $P_{HF} = 199$ kPa. The solid lines correspond to the sum of the primary radiation force and the buoyancy, which each bubble experiences. Lower figures: Evolution of the sum of the forces (including the secondary radiation force) that a virtual bubble (a) or the trapped bubble (b) would experience in the vicinity of the pressure antinode A. (a) For $P_{LF} = 0$ kPa, the two bubbles are located near the pressure nodes $N_i$, and the pressure antinode A is an unstable location. (b) For $P_{LF} = 1.5$ kPa, the jump occurs and the upper bubble is stabilized at the pressure antinode A, which becomes a stable location.</figcaption></figure>

bound at a very short separation distance when the interaction force between them reverses. This sign reversal is shown to be caused by multiple rescattering effects, and is responsible for a modification of stable locations that bubbles may experience in an ultrasound field. The knowledge of an appropriate modeling of bubble-bubble interactions is of great interest for the understanding of the dynamics of dense bubble clouds. In addition, the stabilization of two bubbles at short distances is a prerequisite for the investigation of more intense interactions, such as energy transfer between bubble oscillations.

## ACKNOWLEDGMENTS

This work was supported by the LabEx CeLyA of the University of Lyon (ANR-10-LABX-0060 / ANR-11-

INDEX-0007).

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Gabriel Regnault:** Data curation (lead); Formal analysis (lead); Investigation (lead); Visualization (equal); Writing - original draft preparation (equal).

**Alexander A Doinikov:** Formal analysis (equal); Methodology (lead); Software (lead); Validation (lead).

**Cyril Mauger:** Investigation (equal); Resources (equal); Supervision (equal).

**Philippe Blanc-Benon:** Funding Acquisition (lead); Resources (lead); Supervision (equal).

**Claude Inserra:** Formal analysis (equal); Investigation (equal); Project administration (lead); Supervision (lead); Writing - review & editing (lead).

### Data availability

The data that support the findings of this study are available from the corresponding author upon reasonable request.

## Appendix: Calculation of the interaction force accounting for multiple scattering

Doinikov and Zavtrak¹⁶ have derived an expression for the interaction force acting between two bubbles in the case of an arbitrary ratio between the bubble radii and the interbubble distance $d$. The interaction force is given by

$$F_s = 2\pi\rho|C|^2 G_D e_z, \quad (A.1)$$

where $\rho$ is the liquid density, $C$ is the complex amplitude of the liquid velocity potential, $e_z$ is the unit vector directed from the center of bubble 1 to the center of bubble 2, and the dimensionless quantity $G_D$ is given by

$$G_D = -\text{Re} \left( \sum_{n=0}^{\infty} \frac{n+1}{\eta_{1(n+1)}} a_n^{(11)*} a_{n+1}^{(11)} \right), \quad (A.2)$$

where Re means "the real part of", the asterisk denotes the complex conjugate, $\eta_{1(n+1)}$ is a constant coefficient and $a_n^{(11)}$ are the linear scattering coefficients of bubble 1. The coefficients $\eta_{1(n+1)}$ are calculated by

$$\eta_{1(n+1)} = \frac{1 + (n + 1)\nu_{1(n+1)}}{1 - (n + 2)\nu_{1(n+1)}}, \quad (A.3)$$

where

$$\nu_{1(n+1)} = [n(n+3) + 2\delta_{(n+1)0}] \sigma / (\rho\omega^2 R_{01}^3) + (\omega_1^2 / \omega^2 + i\delta_1) \delta_{(n+1)0}, \quad (A.4)$$
$\delta_{nm}$ is the Kronecker symbol, $\sigma$ is the surface tension, $\rho$ is the liquid density, $R_{0j}$ is the equilibrium radius of bubble $j$, $\omega = \omega_{LF}$ is the driving angular frequency, $\omega_1$ is the resonant angular frequency of bubble 1, and $\delta_j$ is the damping coefficient of bubble $j$, given by

$$ \delta_j = 3\gamma p_{j0}d_j/(\rho\omega^2 R_{0j}^2), \quad (A.5) $$

where $p_{j0} = P_0 + 2\sigma/R_{0j}$ and $P_0$ is the liquid hydrostatic pressure. The coefficient $d_j$ is calculated by using the following equations:

$$ \frac{d_j}{b_j} = 3(\gamma - 1) * \quad (A.6) $$

$$ \frac{X \sinh(X) + \sin(X) - 2(\cosh(X) - \cos(X))}{X^2(\cosh(X) - \cos(X)) + 3(\gamma - 1)X(\sinh(X) - \sin(X))} $$

$$ b_j = \left[ \left( 1 + \frac{d_j^2}{b_j^2} \right) \left( 1 + \frac{3(\gamma - 1)(\sinh(X) - \sin(X))}{X(\cosh(X) - \cos(X))} \right) \right]^{-1} \quad (A.7) $$

$$ X = R_{0j}(2\omega\rho_{gj}c_{pg}/k_g)^{1/2}, $$

where $\gamma$ is the specific heat ratio, $\rho_{gj} = \rho_{gA}p_{j0}/P_A$ is the equilibrium density of the gas inside bubble $j$, $\rho_{gA}$ is the gas density at the atmospheric pressure $P_A$, $c_{pg}$ is the specific heat of the gas at constant pressure, and $k_g$ is the thermal conductivity of the gas. The linear scattering coefficients $a_n^{(11)}$ are given by

$$ a_n^{(11)} = \sum_{m=n}^{\infty} K_{nm} \left( \frac{R_{01}}{d} \right)^m, \quad (A.9) $$

where the coefficients $K_{nm}$ are calculated by the equations:

$$ K_{nn} = \eta_{1n}\eta_{20}\left(\frac{R_{02}}{d}\right) - \eta_{10}\delta_{n0} \quad \text{for } n = m, \quad (A.10) $$

$$ K_{nm} = \eta_{1n} \sum_{i=0}^{[(m-n-1)/2]} s_{ni} K_{i(m-n-1-i)}, \quad \text{for } m > (A.11) $$

$$ s_{nm} = \sum_{k=0}^{\infty} \eta_{2j} D_{mk} D_{kn} \left(\frac{R_{02}}{d}\right)^{2k+1}, \quad (A.12) $$

$$ D_{n0} = 1, \quad (A.13) $$

$$ D_{nm} = \frac{(n+1)(n+2)(n+3)}{1 \ 2 \ 3 \dots \frac{(n+m)}{m}} \quad \text{for } m(A.04) $$

where [] denotes the omission of fractions. It follows from Eqs. (A.2) and (A.9) that the coefficient $G_D$ involves two infinite sums. However, since the terms of these sums decrease with increasing $n$ and $m$, they can be truncated in numerical computations at some finite values of $n$ and $m$ depending on the required accuracy. We have investigated the convergence rate of the magnitude of the coefficient $G_D$ for the bubble pair with $R_{01} = 170 \mu m$ and $R_{02} = 139 \mu m$ (as in Fig. 4) at the smallest interbubble distance $d = 352 \mu m$ that corresponds to the red circle at $L = 1.18$ in Fig. 4. The value of $G_D$ as a function of the number of terms $n$ in Eq. (A.2), as well as the computation time, is provided in Fig. 7.

<figure><img src="image_13.png" /><figcaption>FIG. 7. Evolution of the radiation interaction force (in terms of the parameter $G_D$) and the computation time for $R_{01} = 170 \mu m$ and $R_{02} = 139 \mu m$ as a function of the number of terms $n$ in Eq. (6). For $n = 8$, the sign of the interaction force changes from attraction to repulsion.</figcaption></figure>

As one can see, the value of the interaction force converges for $n > 15$. From $n = 1$ (corresponding to the model of Zabolotskaya) to $n = 15$, the increase of the computation time by two orders of magnitude is observed. However, the choice of the limit number of terms $n$ is important as the interaction force changes its sign from attraction to repulsion around $n = 8$. It means that reducing the number of terms in the model of Doinikov and Zavtrak for reducing the computation time is inadmissible. In the paper, we have used $n = 20$ to calculate Eq. (6). However, in order to limit the computation time when modeling a large number of interacting bubbles, a small number of terms can be kept in Eq. (A.2) at large distances between the bubbles (typically $n = 2$ for $L > 2$) and a ($n > 15$)-model can be used only when bubbles are getting closer and the condition $L < 2$ is fulfilled.

REFERENCES

¹W. Lauterborn, T. Kurz, "Physics of bubble oscillations", Reports on Progress in Physics, **73**(10), 106501 (2010).

²M.P. Brenner, D. Lohse, T.F. Dupont, "Bubble shape oscillations and the onset of sonoluminescence", Physical Review Letters, **75**(5), 954 (1995).

³M. Guédra, C. Inserra, C. Mauger, B. Gilles "Experimental evidence of nonlinear mode coupling between spherical and nonspherical oscillations of microbubbles", Physical Review E, **94**, 053115 (2016).

⁴S.J. Shaw, "Translation and oscillation of a bubble under axisymmetric deformation", Physics of Fluids, **18**, 072104 (2006).

⁵A.A. Doinikov, "Translational motion of a bubble undergoing shape oscillations", Journal of Fluid Mechanics, **501**, 1-24 (2004).

⁶V.F.K. Bjerknes, "Fields of Force". Columbia University Press (1906).
<footnote n="7">⁷A.I. Eller, "Force on a bubble in a standing acoustic wave", Journal of the Acoustic Society of America, **43**, 170-171 (1968).</footnote>

<footnote n="8">⁸O. Vincent, P. Marmottant, "On the statics and dynamics of fully confined bubbles", Journal of Fluid Mechanics, **827**, 194-224 (2017).</footnote>

<footnote n="9">⁹M.S. Longuet-Higgins, "Viscous streaming from an oscillating spherical bubble", Proceedings of the Royal Society A, **454**, 725 (1998).</footnote>

<footnote n="10">¹⁰S. Cleve, M. Guedra, C. Mauger, C. Inserra, P. Blanc-Benon, "Microstreaming induced by acoustically trapped, non-spherically oscillating microbubbles", Journal of Fluid Mechanics, **875**, 597-621 (2019).</footnote>

<footnote n="11">¹¹R. Mettin, "Bubble structures in acoustic cavitation". *Bubble and particle dynamics in acoustic fields: modern trends and applications*, 1-36 (2005).</footnote>

<footnote n="12">¹²U. Parlitz, R. Mettin, S. Luther, I.S. Akhatov, M.A. Voss, W. Lauterborn, "Spatio-temporal dynamics of acoustic cavitation bubble clouds", Philosophical Transactions of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences, **357**, 313-334 (1999).</footnote>

<footnote n="13">¹³L.A. Crum, "Bjerknes forces on bubbles in a stationary sound field", The Journal of the Acoustical Society of America, **57**, 1363 (1975).</footnote>

<footnote n="14">¹⁴P.L. Marston, E.H. Trinh, J. Depew, T.J. Asaki, "Response of bubbles to ultrasonic radiation pressure: Dynamics in low gravity and shape oscillations". In: *Bubble Dynamics and Interface Phenomena* by J.R. Blake *et al.*, Kluwer Academic, Dordrecht, 1994, p. 343-353.</footnote>

<footnote n="15">¹⁵E.A. Zabolotskaya, "Interaction of gas bubbles in a sound wave field", Sov. Phys. Acoust., **30**(5), 365-368 (1984).</footnote>

<footnote n="16">¹⁶A.A. Doinikov, S.T. Zavtrak, "On the mutual interaction of two gas bubbles in a sound field", Physics of Fluids, **7**, 1923-1930 (1995).</footnote>

<footnote n="17">¹⁷A.A. Doinikov, S.T. Zavtrak, "On the "bubble grapes" induced by a sound field", The Journal of the Acoustical Society of Amer-</footnote>

<footnote>ica, **99**, 3849-3950 (1996).</footnote>

<footnote n="18">¹⁸N.A. Pelekasis, A. Gaki, A. Doinikov, J.A. Tsamopoulos, "Secondary Bjerknes forces between two bubbles and the phenomenon of acoustic streamers", Journal of Fluid Mechanics, **500**, 313-347. (2004)</footnote>

<footnote n="19">¹⁹M. Ida, "Alternative interpretation of the sign reversal of secondary Bjerknes force acting between two pulsating gas bubbles", Physical Review E, **67**(5), 056617 (2003).</footnote>

<footnote n="20">²⁰M. Lanoy, C. Derec, A. Tourin, V. Leroy, "Manipulating bubbles with secondary Bjerknes forces", Applied Physics Letters, **107**, 214101 (2015).</footnote>

<footnote n="21">²¹K. Yoshida, T. Fujikawa, Y. Watanabe, "Experimental investigation on reversal of secondary Bjerknes force between two bubbles in ultrasonic standing wave", The Journal of the Acoustical Society of America, **130**(1), 135-144 (2011).</footnote>

<footnote n="22">²²T. Barbat, B. Ashgriz, C-S. Liu, "Dynamics of two interacting bubbles in an acoustic field", Journal of Fluid Mechanics, **389**, 137-168 (1999).</footnote>

<footnote n="23">²³D.L. Miller, "Stable arrays of resonant bubbles in a 1-MHz standing-wave acoustic field", The Journal of the Acoustical Society of America, **62**, 12 (1977).</footnote>

<footnote n="24">²⁴D. Rabaud, P. Thibault, M. Mathieu, P. Marmottant, "Acoustically bound microfluidic bubble crystals", Physical Review Letters, **106**(13), 134501 (2011).</footnote>

<footnote n="25">²⁵G. Regnault, C. Mauger, P. Blanc-Benon, C. Inserra, "Secondary radiation force between two closely spaced acoustic bubbles", Physical Review E, **102**(3), 031101 (2020).</footnote>

<footnote n="26">²⁶A.A. Doinikov, "Bjerknes forces and translational bubble dynamics", in *Bubble and Particle Dynamics in Acoustic Fields: Modern Trends and Applications* (Research Signpost, Trivandrum, Kerala, India), 95-143 (2005).</footnote>

<footnote n="27">²⁷S. Annamalai, S. Balachandar, "Mean force on a finite-sized rigid particle, droplet, or bubble in a viscous compressible medium", Physics of Fluids, **27**, 103304 (2015).</footnote>
(a) $P_{HF} \neq 0$, $P_{LF} \neq 0$

(b) Forbidden zone

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

<figure><img src="image_15.png" /></figure>
<figure><img src="image_16.png" /></figure>

Increasing LF pressure
increasing LF pressure

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

<figure><img src="image_18.png" /></figure>
<figure><img src="image_19.png" /></figure>

-O Bjerknes

-● Zabolotskaya

Doinikov
<figure><img src="image_20.png" /><figcaption>(b) $R_{01} = 160 \ \mu\text{m}, R_{02} = 130 \ \mu\text{m}$</figcaption></figure>

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

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

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

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

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

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

<figure><img src="image_27.png" /></figure>
(a) $P_{\text{LF}} = 0$

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

$F_{pi} + F_{si} + B_{i}$

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

(b) $P_{\text{LF}} = 1.5 \text{ kPa}$

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

$F_{pi} + F_{si} + B_{i}$

<figure><img src="image_31.png" /></figure>
<figure><img src="image_32.png" /></figure>