2 Problem statement

Let us investigate the behaviour of a solid heavy particle in the vicinity of a gas bubble rising in an incompressible viscous liquid filling a long cylindrical vertical tube of diameter $d$ which is much larger than the bubble and particle sizes, i.e.  $d\gg R_{b}\gg R_{p}$ (where $R_{b}$ and $R_{p}$ are the bubble and particle radii). The system is subjected to a gravity field and an external vibrational (ultrasound) field. The ultrasound waves propagate in the direction perpendicular to the tube axis; that is why in the modelling we consider horizontal vibrations, namely standing waves. The system is sketched in figure 1.

Figure 1.

Scheme of the set-up.

The presence of two significantly different time scales – the pulsation time $1/\unicode[STIX]{x1D714}$ and the dissipative time $2R_{b}/u$ (where $u$ is the rising velocity of the bubble) – makes it reasonable to consider the ‘fast’ (or pulsating) and ‘slow’ (or average) motions separately.

3 Pulsating motion

In this section, we consider the vibration-induced flow around the bubble and particle. We focus on the case when the particle is much smaller than the bubble, which allows neglect of the flow induced by the wave scattered on the particle. Let us consider the pulsating flow around a gas bubble suspended in a liquid filling the entire space.

We adopt the following assumptions for the problem.

1. (1) The vibration frequency is so high that the thickness of the Stokes boundary layer $\unicode[STIX]{x1D6FF}=\sqrt{2\unicode[STIX]{x1D708}/\unicode[STIX]{x1D714}}$ (where $\unicode[STIX]{x1D708}$ is the kinematic viscosity of the liquid) is much smaller than the bubble and particle radii.

2. (2) The sound wavelength at the vibration frequency is at least one order larger than the tube radius, $\unicode[STIX]{x1D706}\gg d$ .

3. (3) The amplitude of vibration $b$ is small compared with the particle radius, $b\ll R_{p}$ .

4. (4) The deviation of the bubble surface from a spherical shape is small compared with the bubble radius, which leads to the restriction $bd\ll R_{b}^{2}$ .

These assumptions permit us to describe a pulsating flow around the bubble using the inviscid fluid model (assumption 1) and to consider the imposed vibration field as spatially uniform (assumption 2). Under assumptions 3 and 4, we can linearize the Euler equation for the pulsating velocity field near the bubble. The ranges of the parameter values of practical interest are presented in table 1 (where $f$ is the vibration frequency) and are consistent with the assumptions made.

Table 1.

Reference values of real system parameters.

The pulsating velocity field includes the imposed field $\boldsymbol{U}_{vib}$ and the scattered fields of the liquid velocity cased by the monopole $\boldsymbol{U}_{m}$ and translational $\boldsymbol{U}_{tr}$ bubble motions:

(3.1) $$\begin{eqnarray}\boldsymbol{U}=\boldsymbol{U}_{vib}+\boldsymbol{U}_{m}+\boldsymbol{U}_{tr}.\end{eqnarray}$$

Here, $\boldsymbol{U}_{vib}=b\unicode[STIX]{x1D714}\boldsymbol{i}\cos \unicode[STIX]{x1D714}t$ is the velocity of the oscillatory flow at a sufficient distance from the bubble and $\boldsymbol{i}$ is the direction of vibrations.

3.1 Monopole motion

Let us first consider the liquid motion due to the monopole or radial bubble motion, which is induced by the term representing pressure forcing, $q\text{e}^{\text{i}\unicode[STIX]{x1D714}t}$ . Using the linearized Euler equation and the continuity equation for the liquid velocity in the spherical coordinate system, one can find by simple integration that the pressure and velocity amplitudes of the liquid near the bubble are as follows:

(3.2a ) $$\begin{eqnarray}\displaystyle & \displaystyle U_{m}=\frac{c_{1}}{r^{2}}, & \displaystyle\end{eqnarray}$$

(3.2b ) $$\begin{eqnarray}\displaystyle & \displaystyle p=\text{i}\unicode[STIX]{x1D714}\unicode[STIX]{x1D70C}\frac{c_{1}}{r}+q, & \displaystyle\end{eqnarray}$$

where $r$ is the radial distance from the bubble centre and the unknown constant $c_{1}$ is related to the surface deviation $\unicode[STIX]{x1D701}$ through the kinematic condition

(3.3) $$\begin{eqnarray}\text{i}\unicode[STIX]{x1D714}\unicode[STIX]{x1D701}=\frac{c_{1}}{R_{b}^{2}}.\end{eqnarray}$$

Another boundary condition is the normal stress balance condition, which in terms of gas pressure variation, $p_{g}^{\prime }$ , reads

(3.4) $$\begin{eqnarray}p-p_{g}^{\prime }=\frac{2\unicode[STIX]{x1D6FE}}{{R_{b}}^{2}}\unicode[STIX]{x1D701}.\end{eqnarray}$$

We assume that the gas behaves polytropically, so that $p_{g}V^{k}=\text{const}$ . Substitution of the pressure in a liquid ((3.2), taking into account (3.3)) in the normal stress balance condition (3.4) gives the expression for the surface deviation,

(3.5) $$\begin{eqnarray}\unicode[STIX]{x1D701}=\frac{q}{\unicode[STIX]{x1D714}^{2}\unicode[STIX]{x1D70C}R_{b}-{\displaystyle \frac{3kp_{0}}{R_{b}}}+{\displaystyle \frac{2\unicode[STIX]{x1D6FE}}{{R_{b}}^{2}}}}.\end{eqnarray}$$

Neglecting the surface tension, we come to the expression for the natural frequency,

(3.6) $$\begin{eqnarray}\unicode[STIX]{x1D714}_{0}^{2}=\frac{3p_{0}k}{\unicode[STIX]{x1D70C}{R_{b}}^{2}},\end{eqnarray}$$

which coincides with the result obtained in Plesset & Prosperetti (Reference Plesset and Prosperetti1977) and Leighton (Reference Leighton2012). Thus, the expression (3.5) can be written in a more compact form,

(3.7) $$\begin{eqnarray}\unicode[STIX]{x1D701}=\frac{q}{\unicode[STIX]{x1D70C}R_{b}(\unicode[STIX]{x1D714}^{2}-\unicode[STIX]{x1D714}_{0}^{2})}.\end{eqnarray}$$

Substituting (3.7) in (3.2), we obtain the following expressions for the velocity and pressure of monopole motion:

(3.8a ) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{U}_{m}=\left(\frac{\text{i}\unicode[STIX]{x1D714}\,qR_{b}}{(\unicode[STIX]{x1D714}^{2}-\unicode[STIX]{x1D714}_{0}^{2})\unicode[STIX]{x1D70C}}\frac{1}{r^{3}}\right)\boldsymbol{r}\text{e}^{\text{i}\unicode[STIX]{x1D714}t}, & \displaystyle\end{eqnarray}$$

(3.8b ) $$\begin{eqnarray}\displaystyle & \displaystyle p_{m}=\left(-\frac{\unicode[STIX]{x1D714}^{2}qR_{b}}{\unicode[STIX]{x1D714}^{2}-\unicode[STIX]{x1D714}_{0}^{2}}\frac{1}{r}+q\right)\text{e}^{\text{i}\unicode[STIX]{x1D714}t}, & \displaystyle\end{eqnarray}$$

where $\boldsymbol{r}$ is a radial vector along the radial coordinate $r$ .

The pressure $q$ far from the bubble can be determined by the following reasoning. Under the assumption that the sound wavelength at the vibration frequency is much larger than the tube diameter, the far-field pressure is found from the Euler equation in the form

(3.9) $$\begin{eqnarray}q\text{e}^{\text{i}\unicode[STIX]{x1D714}t}=-\text{i}\unicode[STIX]{x1D714}\unicode[STIX]{x1D70C}(\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{r}_{b}),\end{eqnarray}$$

where $\boldsymbol{r}_{b}$ is the radius vector of the bubble centre which is counted from the pressure nodes. The coordinates of the nodes are determined by the type of generation of oscillations in the liquid. We suppose that the oscillations in the liquid are generated in such a way that the pressure node is located in the tube centre. Under the assumption that the sound wavelength at the vibration frequency is much larger than the tube diameter, the acoustic wave inside the liquid is spatially uniform; thus, the action of the ultrasound can be described as the vibrations of the tube walls.

3.2 Translational motion

For the study of the translational motion of the bubble it is convenient to write the equations and boundary conditions describing the dynamics of the system in the reference frame moving with the bubble centre. In this case, we assume that the bubble motion is one-dimensional, and the direction of this motion coincides with the direction of the vibration axis $\boldsymbol{i}$ ( $U_{b}\boldsymbol{i}$ is the velocity of the centre of mass of the bubble relative to the laboratory reference frame). In this reference frame, for a fluid of density $\unicode[STIX]{x1D70C}$ the Euler and continuity equations may be written as

(3.10a ) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\boldsymbol{v}}{\unicode[STIX]{x2202}t}=-\frac{\unicode[STIX]{x1D735}p}{\unicode[STIX]{x1D70C}}-\frac{\text{d}U_{b}}{\text{d}t}\boldsymbol{i}, & \displaystyle\end{eqnarray}$$

(3.10b ) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{v}=0. & \displaystyle\end{eqnarray}$$

Far from the bubble, its effect on the fluid flow becomes insignificant; hence, the fluid velocity matches the vibration velocity, i.e.

(3.11) $$\begin{eqnarray}\boldsymbol{v}(r\rightarrow \infty )=(U_{vib}-U_{b})\boldsymbol{i}=(b\unicode[STIX]{x1D714}\cos \unicode[STIX]{x1D714}t-U_{b})\boldsymbol{i}.\end{eqnarray}$$

On the bubble surface ( $r=R_{b}$ ), the kinematic condition and the normal stress balance condition are satisfied:

(3.12a ) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D701}}{\unicode[STIX]{x2202}t}=\boldsymbol{v}\boldsymbol{\cdot }\frac{\boldsymbol{r}}{r}, & \displaystyle\end{eqnarray}$$

(3.12b ) $$\begin{eqnarray}\displaystyle & p-p_{g}=-\unicode[STIX]{x1D6FE}K, & \displaystyle\end{eqnarray}$$

where $K$ is the surface curvature, $\unicode[STIX]{x1D6FE}$ is the surface tension coefficient and $\unicode[STIX]{x1D701}$ is the surface deviation, which depends on time and on the polar angle $\unicode[STIX]{x1D717}$ measured from the direction of vibration. By solving (3.10) with boundary condition far from the bubble (3.11) we obtained the flow fields, which are convenient to write for irrotational flow in terms of a velocity potential and pressure,

(3.13a ) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D711}=\left(b\unicode[STIX]{x1D714}r+\frac{c_{2}}{r^{2}}\right)\cos \unicode[STIX]{x1D717}\cos \unicode[STIX]{x1D714}t-U_{b}r\cos \unicode[STIX]{x1D717}, & \displaystyle\end{eqnarray}$$

(3.13b ) $$\begin{eqnarray}\displaystyle & \displaystyle p=\unicode[STIX]{x1D70C}\unicode[STIX]{x1D714}\left(b\unicode[STIX]{x1D714}r+\frac{c_{2}}{r^{2}}\right)\cos \unicode[STIX]{x1D717}\sin \unicode[STIX]{x1D714}t. & \displaystyle\end{eqnarray}$$

We will now determine the surface deviation. It follows from the kinematic condition (3.12a ) that $\unicode[STIX]{x1D701}=c_{3}\cos \unicode[STIX]{x1D717}\cos \unicode[STIX]{x1D714}t$ . However, this solution describes the translational oscillations of the bubble as a whole, which is impossible in the reference frame moving with the bubble centre, so it follows that $c_{3}=0$ . Thus, in the leading order the bubble surface deviation vanishes.

The unknown constant $c_{2}$ , as well as the velocity of the bubble centre $U_{b}$ , can be obtained from boundary conditions on the bubble surface (3.12). Indeed, taking into consideration that the pressure inside the bubble is determined by the mean surface curvature $p_{g}=\unicode[STIX]{x1D6FE}K=2\unicode[STIX]{x1D6FE}/{R_{b}}^{2}$ , it can be found that $c_{2}=-b\unicode[STIX]{x1D714}{R_{b}}^{3}$ and $U_{b}=3b\unicode[STIX]{x1D714}\cos \unicode[STIX]{x1D714}t$ , which coincides with the result obtained in Lyubimov et al. (Reference Lyubimov, Lyubimova and Cherepanov2002). Substituting $c_{2}$ and $U_{b}$ in (3.13), the velocity of the liquid can be written as

(3.14) $$\begin{eqnarray}\boldsymbol{v}=b\unicode[STIX]{x1D714}\left\{\left(-2+2\frac{{R_{b}}^{3}}{r^{3}}\right)\cos \unicode[STIX]{x1D717}\boldsymbol{e}_{r}+\left(2+\frac{{R_{b}}^{3}}{r^{3}}\right)\sin \unicode[STIX]{x1D717}\boldsymbol{e}_{\unicode[STIX]{x1D717}}\right\}\cos \unicode[STIX]{x1D714}t.\end{eqnarray}$$

If we come back to the laboratory reference frame and exclude the incident wave from (3.14), the expression for the required liquid velocity caused by the translational bubble motion in the leading order reads

(3.15) $$\begin{eqnarray}\boldsymbol{U}_{tr}=b\unicode[STIX]{x1D714}\left(3\frac{{R_{b}}^{3}}{r^{5}}\left(\boldsymbol{i}\boldsymbol{\cdot }\boldsymbol{r}\right)\boldsymbol{r}-\frac{{R_{b}}^{3}}{r^{3}}\boldsymbol{i}\right)\cos \unicode[STIX]{x1D714}t.\end{eqnarray}$$

Finally, the pulsating velocity of the liquid near the gas bubble in the external pulsating field can be written as

(3.16) $$\begin{eqnarray}\boldsymbol{U}=\boldsymbol{U}_{vib}+(\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{r}_{b})\frac{\unicode[STIX]{x1D714}^{2}R_{b}}{(\unicode[STIX]{x1D714}^{2}-\unicode[STIX]{x1D714}_{0}^{2})}\frac{\boldsymbol{r}}{r^{3}}-\boldsymbol{U}_{vib}\frac{R_{b}^{3}}{r^{3}}+(\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{r})\frac{3R_{b}^{3}}{r^{4}}\frac{\boldsymbol{r}}{r}.\end{eqnarray}$$

4 Average motion

In this section, we calculate the average motion of the particle and the bubble. As the directions of vibrations and gravity are mutually perpendicular, we can separate their influences on the particle motion. The total force exerted on the particle is calculated as the sum of two contributions. The first one originates from the velocity field $\boldsymbol{V}$ caused by the rising motion of the bubble, and is decomposed as the sum of gravity force, viscous drag, Basset force and added mass. The other is the average vibration force $\boldsymbol{F}_{vib}$ originating from the oscillatory velocity field $\boldsymbol{U}$ , produced by the incident vibration field plus the monopole and dipole fields scattered by the bubble. For the range of particle sizes under consideration (see table 1), the particle radius is sufficiently small that the Faxen corrections (Faxen Reference Faxen1922) can be neglected. We also assume that the bubble rises with constant velocity in laminar flow; thus, the Saffman lift term is zero. Under the above assumptions, the momentum equation for a particle per unit volume in the coordinate system associated with the bubble centre (figure 2) reads

(4.1) $$\begin{eqnarray}\displaystyle \left(\unicode[STIX]{x1D70C}_{s}+\frac{\unicode[STIX]{x1D70C}}{2}\right)\frac{\text{d}^{2}\boldsymbol{r}}{\text{d}t^{2}} & = & \displaystyle \boldsymbol{F}_{vib}-\frac{9\unicode[STIX]{x1D702}}{2{R_{p}}^{2}}\left(\frac{\text{d}\boldsymbol{r}}{\text{d}t}-\boldsymbol{V}\right)\nonumber\\ \displaystyle & & \displaystyle -\,\frac{9\unicode[STIX]{x1D70C}}{2R_{p}}\sqrt{\frac{\unicode[STIX]{x1D708}}{\unicode[STIX]{x03C0}}}\int _{0}^{t}\frac{\text{d}}{\text{d}\unicode[STIX]{x1D701}}\left(\frac{\text{d}\boldsymbol{r}}{\text{d}\unicode[STIX]{x1D701}}-\boldsymbol{V}\right)\frac{\text{d}\unicode[STIX]{x1D701}}{\sqrt{t-\unicode[STIX]{x1D701}}}+(\unicode[STIX]{x1D70C}_{s}-\unicode[STIX]{x1D70C})\boldsymbol{g},\end{eqnarray}$$

where $\unicode[STIX]{x1D70C}_{s}$ is the particle density, and $\unicode[STIX]{x1D702}$ and $\unicode[STIX]{x1D708}$ are the dynamic and kinematic viscosities of the liquid. Assuming a creeping flow regime, $\boldsymbol{V}$ represents the velocity field caused by the rising bubble, if it were not oscillating, which according to Landau & Lufshits (Reference Landau and Lufshits1987) is given by

(4.2) $$\begin{eqnarray}\boldsymbol{V}=\frac{uR_{b}}{2r}\left[\frac{\boldsymbol{r}}{r^{2}}(\boldsymbol{k}\boldsymbol{\cdot }\boldsymbol{r})+\boldsymbol{k}\right]-u\boldsymbol{k},\end{eqnarray}$$

where $\boldsymbol{k}$ is the unit vector in the vertical direction and the velocity of the rising bubble is

(4.3) $$\begin{eqnarray}u=gR_{b}^{2}/9\unicode[STIX]{x1D708}.\end{eqnarray}$$

Figure 2.

Problem configuration.

In deriving (4.1) we neglected the rectified diffusion, which can be reasonable, as discussed in the Appendix. We consider that the bubble rises maintaining sphericity. According to Levich (Reference Levich1962), the bubble remains in spherical form if its radius is less than the critical radius, namely $R_{b}<\sqrt{\unicode[STIX]{x1D6FE}/\unicode[STIX]{x1D70C}g}\approx 3000~\unicode[STIX]{x03BC}\text{m}$ . For the range of bubble size under consideration (see table 1) this condition is satisfied.

We did not take into account the Derjaguin–Landau–Verwey–Overbeek (DLVO) forces due to the fact that our calculations became unreliable at distances less than $\times 10^{-2}R_{p}$ (the thickness of the Stokes boundary layer). As shown in Elimelech et al. (Reference Elimelech, Gregory, Jia and Williams1995), the DVLO forces are comparable to the Stokes force at distances less than $\times 10^{-5}R_{p}$ .

The average vibration force, acting on the particle from the oscillating liquid, was found in Lyubimov et al. (Reference Lyubimov, Lyubimova and Cherepanov2002) for the case of incompressible fluid under the assumption that the vibration frequency is high and the amplitude is small. These authors considered the general case of non-uniform inviscid pulsating flow and solved the equation of liquid motion in terms of potential, assuming that the particle size is small compared with the spatial scale of the pulsating field inhomogeneity. In the laboratory reference frame after averaging over the vibration period, the expression for the average force $\boldsymbol{F}_{vib}$ per unit volume was obtained in the form

(4.4) $$\begin{eqnarray}\boldsymbol{F}_{vib}=\frac{3}{4}\unicode[STIX]{x1D70C}\frac{\unicode[STIX]{x1D70C}_{s}-\unicode[STIX]{x1D70C}}{\unicode[STIX]{x1D70C}_{s}+{\textstyle \frac{1}{2}}\unicode[STIX]{x1D70C}}\unicode[STIX]{x1D735}\overline{U^{2}}.\end{eqnarray}$$

This expression can also be obtained from the general expression for the acoustic radiation force in inviscid liquid, derived in Krasil’nikov & Krylov (Reference Krasil’nikov and Krylov1984), assuming that the liquid is incompressible. In the work of Settnes & Bruus (Reference Settnes and Bruus2012), the analytical expression for the acoustic radiation force was extended to the case of viscous liquid. This expression is valid for any particle radius and boundary-layer thickness provided that both of these length scales are much smaller than the wavelength of the ultrasound wave. These authors considered separately generation by standing and travelling waves. For inviscid incompressible liquid and a standing wave, their expression for the radiation force is reduced to the ones obtained in previous works, namely (4.4). They also showed that in the case of inviscid ambient liquid, a travelling planar wave has no influence on average particle motion. For a travelling planar wave, the radiation force is proportional to the boundary-layer thickness; thus, in inviscid liquid this force is equal to zero. Despite the good agreement with other works in the limiting cases, the radiation force obtained by Settnes & Bruus (Reference Settnes and Bruus2012) does not describe the change from attraction to repulsion at small distances from the particle; this change was found experimentally (Klotsa et al. Reference Klotsa, Swift, Bowley and King2007, Reference Klotsa, Swift, Bowley and King2009) and numerically (Lyubimova, Lyubimov & Shardin Reference Lyubimova, Lyubimov and Shardin2011; Lyubimov, Baydin & Lyubimova Reference Lyubimov, Baydin and Lyubimova2013).

The vector field of this force is plotted in figure 3(a), where the vibrations are applied in the horizontal direction. As can be seen, all over the domain the vibrational force is directed towards the bubble surface; thus, it is reasonable to expect that particles will be attracted to the bubble.

Figure 3.

(a) The vector field of the average vibrational force. (b) The trajectories of particle motion for different vibration intensities: 1 – in the absence of vibrations, 2 –  $b\unicode[STIX]{x1D714}/u=1.6$ , 3 –  $b\unicode[STIX]{x1D714}/u=1.92$ , 4 –  $b\unicode[STIX]{x1D714}/u=3.2$ , 5 –  $b\unicode[STIX]{x1D714}/u=4.8$ and 6 –  $b\unicode[STIX]{x1D714}/u=6.4$ ( $\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D70C}_{s}/\unicode[STIX]{x1D70C}=8.9$ , $R_{b}=0.5~\text{mm}$ , $R_{p}=0.02~\text{mm}$ ). The vibrations are directed along the $x$ -axis and the gravitational force is directed along the $y$ -axis.

It can be shown that at large distances from the bubble, the most essential contribution to the average force is due to the interference of the external vibration field and the field scattered by the bubble, which means that consideration of only monopole oscillations of the bubble (as, for example, in Doinikov & Zavtrak Reference Doinikov and Zavtrak1996) would be inadequate.

Indeed, we need to evaluate the dependence of $\boldsymbol{F}_{vib}$ on the distance from the bubble surface. We restart from (3.1) and exhibit nearly all contributions to $\overline{U^{2}}$ one by one ( $\overline{U_{vib}^{2}}$ , $\overline{\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{U}_{m}}$ , $\overline{\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{U}_{tr}}$ , $\overline{\boldsymbol{U}_{m}\boldsymbol{\cdot }\boldsymbol{U}_{tr}}$ , $\overline{U_{m}^{2}}$ , and so on) using the full expression (3.16). The first term, $\overline{U_{vib}^{2}}$ , is spatially uniform and does not contribute to the force. The contribution of the solely monopole mode, $\overline{U_{m}^{2}}$ , decreases with distance as  $r^{-5}$ .

Now let us calculate the cross-contribution. The result of averaging the scalar product is written as

(4.5) $$\begin{eqnarray}\overline{\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{U}_{m}}=\frac{1}{2}\frac{\unicode[STIX]{x1D714}^{2}R_{b}}{(\unicode[STIX]{x1D714}^{2}-\unicode[STIX]{x1D714}_{0}^{2})}\frac{1}{r^{3}}(\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{r}_{b})(\boldsymbol{U}_{vib}\boldsymbol{\cdot }\boldsymbol{r}).\end{eqnarray}$$

Thus, for the cross-component of the velocity squared, we obtain an expression similar to that of the potential dipole field. Moreover, the expression for the contribution of the cross-component (3.7) to the force $\boldsymbol{F}_{\boldsymbol{c}\boldsymbol{r}}$ is similar to that of the electric field of the dipole, namely

(4.6) $$\begin{eqnarray}\boldsymbol{F}_{cr}\sim \frac{\boldsymbol{U}_{vib}}{r^{3}}-\frac{3\boldsymbol{r}\boldsymbol{\cdot }\boldsymbol{U}_{vib}}{r^{5}}\boldsymbol{r}.\end{eqnarray}$$

At large distances, this force decreases more slowly than the contributions of the monopole and translational modes to the force, and at distances comparable to the bubble radius all terms are of equal significance.

4.1 Trajectories of the particle motion

Equation (3.16) was solved numerically using the Euler method. For each time step, the Basset force was computed from previous time steps using an algorithm for partial differential equations, namely a Caputo-type fractional derivative (Diethelm et al. Reference Diethelm, Ford, Freed and Luchko2005).

The calculations show that in the absence of vibrations, the effect of the bubble on the particle motion becomes insignificant at a distance of approximately one bubble diameter – the particle trajectories are nearly vertical (figure 3 b, curve 1). When subjected to vibrations, the particle trajectories are gradually distorted, and at a certain value of the vibration parameter the particles hit the bubble surface (figure 3 b, curves 2–6).

The problem under consideration has some peculiarities. If the initial position of the particle lies in the plane formed by the directions of vibration and gravity, the particle does not leave this plane in its motion, i.e. the particle trajectory is two-dimensional. However, if the initial position of the particle does not belong to this plane, then the trajectory of the particle is three-dimensional (figure 4).

Figure 4.

The trajectories of the particle motion in three dimensions: $b\unicode[STIX]{x1D714}/u=3.2$ , $\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D70C}_{s}/\unicode[STIX]{x1D70C}=8.9$ , $R_{b}=0.5~\text{mm}$ , $R_{p}=0.02~\text{mm}$ , $r_{b}=0.5~\text{mm}$ . The vibrations are directed along the $x$ -axis and the gravitational force is directed along the $y$ -axis.

Figure 5.

The trajectories of the particle motion for different density ratios: (a)  $b\unicode[STIX]{x1D714}/u=3.2$ , 1 –  $\unicode[STIX]{x1D6FE}=2$ , 2 –  $\unicode[STIX]{x1D6FE}=4$ , 3 –  $\unicode[STIX]{x1D6FE}=8$ , 4 –  $\unicode[STIX]{x1D6FE}=16$ , 5 –  $\unicode[STIX]{x1D6FE}=20$ ; (b)  $b\unicode[STIX]{x1D714}/u=8$ , 1 –  $\unicode[STIX]{x1D6FE}=1.1$ , 2 –  $\unicode[STIX]{x1D6FE}=1.2$ , 3 –  $\unicode[STIX]{x1D6FE}=2$ , 4 –  $\unicode[STIX]{x1D6FE}=20$ ( $\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D70C}_{s}/\unicode[STIX]{x1D70C}=8.9$ , $R_{b}=0.5~\text{mm}$ , $R_{p}=0.02~\text{mm}$ , $r_{b}=0.5~\text{mm}$ ). The vibrations are directed along the $x$ -axis and the gravitational force is directed along the $y$ -axis.

The density ratio also influences the bubble–particle interaction (figure 5 a,b). The trajectories of particles just slightly denser than the liquid are not much different from the streamlines (figure 5 a, curves 1, 2), which means that vibrations do not have sufficient effect on such particles, so that they are carried away from the bubble surface with the flow. If the particle is much denser then the liquid, the gravity force controls the particle behaviour (figure 5 a, curve 5). However, if we slightly increase the ultrasound intensity, the vibration influence on the particle attraction grows significantly (figure 5 b); even particles 20 times heavier than the liquid will move to the bubble surface.

4.2 Collision cross-section

A quantitative estimation of the vibration effect on the bubble–particle interaction can be made by calculating the efficiency of collision,

(4.7) $$\begin{eqnarray}E_{c}=\frac{{b_{c}}^{2}}{{R_{b}}^{2}},\end{eqnarray}$$

where $b_{c}$ is the maximal radius of the flow tube, inside which all particles located inside this tube settle on the bubble surface. In the absence of vibrations the problem has axial symmetry; thus, the flow tube has the same symmetry.

According to the Flint and Howarth theory (Flint & Howarth Reference Flint and Howarth1971), one of the most important factors responsible for capturing small particles is gravitational force. The collision efficiency for inertialess particles is defined by

(4.8) $$\begin{eqnarray}E_{c}=\frac{G}{1+G},\end{eqnarray}$$

where $G$ is the terminal settling velocity of the particles in the undisturbed fluid, namely

(4.9) $$\begin{eqnarray}G=\frac{u_{g}}{u}=\frac{2}{9}\left(\frac{\unicode[STIX]{x1D70C}_{s}}{\unicode[STIX]{x1D70C}}-1\right)\frac{gR_{p}^{2}}{\unicode[STIX]{x1D708}}\div \frac{1}{9}\frac{gR_{b}^{2}}{\unicode[STIX]{x1D708}}=2\left(\frac{\unicode[STIX]{x1D70C}_{s}}{\unicode[STIX]{x1D70C}}-1\right)\frac{R_{p}^{2}}{R_{b}^{2}}.\end{eqnarray}$$

According to the assumption that the particle radius is much smaller than the bubble radius, the settling velocity $G$ is much smaller than one ( $G\ll 1$ ), which allows (4.8) to be transformed for the collision efficiency as follows:

(4.10) $$\begin{eqnarray}E_{c}=\frac{G}{1+G}=G(1-G+\cdots )\approx G.\end{eqnarray}$$

Substitution of this expression in the formula of the collision efficiency (4.7) gives the following expression for the impact distance:

(4.11) $$\begin{eqnarray}\frac{b_{c}}{R_{b}}=\sqrt{2\left(\frac{\unicode[STIX]{x1D70C}_{s}}{\unicode[STIX]{x1D70C}}-1\right)}\frac{R_{p}}{R_{b}}.\end{eqnarray}$$

Our calculations of the momentum equation for a particle (3.16), in the absence of vibrations and inertia of the particle, give the impact distance $b_{c}$ , which quantitatively agrees with the result obtained by Flint and Howarth (figure 6 a). To demonstrate the influence of vibration on impact distance, we choose the $x$ -coordinate of the collision cross-section ( $z$  is equal to 0) as an impact distance. When subjected to vibrations even of low intensity, the impact distance substantially increases (figure 6 b).

Figure 6.

(a) Impact distance versus density ratio; the dashed line is the impact distance, according to the Flint and Howarth theory. (b) Impact distance versus vibration velocity ( $\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D70C}_{s}/\unicode[STIX]{x1D70C}=2$ , $R_{b}=0.5~\text{mm}$ , $R_{p}=0.05~\text{mm}$ , $r_{b}=0.5~\text{mm}$ ).

A full understanding of the vibration effect on the efficiency of collision can be achieved by calculating the collision cross-section. For this, we consider a horizontal plane which is located at a distance from the centre of the bubble that approximately equals its diameter (figure 7 a). Analysis of the trajectories of the particles passing through this plane allows us to determine the region (collision cross-section) in which the captured particles will approach the bubble surface (figure 7 b). Normalization of the collision cross-section to that of the bubble leads to an important characteristic in the kinetic theory of flotation: the efficiency of collision. Our calculations, in the absence of vibrations and inertia, give the collision cross-section, which is axially symmetric and quantitatively agrees with the collision efficiency obtained by Flint and Howarth (figure 7 b, domain 1). When subjected to vibrations even of low intensity, the collision cross-section increases significantly (by several times, figure 7 b, domain 2), which results in a strong increase in the efficiency of the flotation process as a whole. The collision cross-section is curiously roughly circular (figure 7 b, domain 2). This form of the collision cross-section is related to the fact that for chosen parameter values, the radial motion of the bubble is dominant. If we consider a bubble located closer to the pressure node, the form of the collision cross-section (figure 8 a) changes, and for a bubble in the pressure node, for which there is translational bubble motion only, the form of the collision cross-section is not circular any more (figure 8 b).

Figure 7.

Collision cross-section: domain 1 – in the absence of vibrations; domain 2 – in the presence of vibrations, $b\unicode[STIX]{x1D714}/u=0.4$ ; domain 3 –  $b\unicode[STIX]{x1D714}/u=0.8$ ; domain 4 –  $b\unicode[STIX]{x1D714}/u=4$ . Here, $\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D70C}_{s}/\unicode[STIX]{x1D70C}=2$ , $R_{b}=0.5~\text{mm}$ , $R_{p}=0.05~\text{mm}$ , $r_{b}=0.5~\text{mm}$ . The dashed line is the projection of the cross-section of the bubble.

Figure 8.

Collision cross-section: (a)  $r_{b}=0.1~\text{mm}$ ; (b)  $r_{b}=0$ . Domain 1 – in the absence of vibrations; domain 2 – in the presence of vibrations, $b\unicode[STIX]{x1D714}/u=4$ ; domain 3 –  $b\unicode[STIX]{x1D714}/u=8$ . Here, $\unicode[STIX]{x1D6FE}=\unicode[STIX]{x1D70C}_{s}/\unicode[STIX]{x1D70C}=2$ , $R_{b}=0.5~\text{mm}$ , $R_{p}=0.05~\text{mm}$ . The dashed line is the projection of the cross-section of the bubble.