Home | About | Journals | Submit | Contact Us | Français |

**|**Proc Math Phys Eng Sci**|**PMC4892283

Formats

Article sections

- Abstract
- 1. Introduction
- 2. Problem formulation and solution
- 3. Numerical simulations and experimental verification
- 4. Conclusion
- References

Authors

Related links

Proc Math Phys Eng Sci. 2016 April; 472(2188): 20160031.

PMCID: PMC4892283

CNRS and Université Grenoble-Alpes, LIPhy UMR 5588, Université Grenoble-Alpes, Grenoble, F-38401, France

Received 2016 January 13; Accepted 2016 March 11.

Copyright © 2016 The Author(s)

Published by the Royal Society. All rights reserved.

The volume oscillation of a cylindrical bubble in a microfluidic channel with planar elastic walls is studied. Analytical solutions are found for the bulk scattered wave propagating in the fluid gap and the surface waves of Lamb-type propagating at the fluid–solid interfaces. This type of surface wave has not yet been described theoretically. A dispersion equation for the Lamb-type waves is derived, which allows one to evaluate the wave speed for different values of the channel height *h*. It is shown that for *h*<λ_{t}, where λ_{t} is the wavelength of the transverse wave in the walls, the speed of the Lamb-type waves decreases with decreasing *h*, while for *h* on the order of or greater than λ_{t}, their speed tends to the Scholte wave speed. The solutions for the wave fields in the elastic walls and in the fluid are derived using the Hankel transforms. Numerical simulations are carried out to study the effect of the surface waves on the dynamics of a bubble confined between two elastic walls. It is shown that its resonance frequency can be up to 50% higher than the resonance frequency of a similar bubble confined between two rigid walls.

This study is inspired by the use of bubbles in acoustically activated microfluidic systems [1]. An example of such a system is described by Rabaud *et al*. [2]: micrometre-size cylindrical bubbles are positioned in a fluid channel, 25μm in depth, confined by elastic walls which are set into oscillation by a vibrating glass rod moulded in the upper wall. Under such conditions, the bubble dynamics exhibit a number of interesting effects, such as the generation of surface waves propagating along the solid–fluid interfaces and the self-organization of bubbles in crystal-like structures where equilibrium distances between bubbles are much smaller than distances that are predicted by the theory of secondary Bjerknes forces [3]. Another phenomenon of interest is acoustic microstreaming generated by bubbles in microfluidic devices [4]. This phenomenon is used for micromixing of fluids [5] and sorting of microparticles and cells [6,7]. A variety of experimental data on the dynamics of bubbles in narrow planar channels can be found in a series of recent papers [8–13], where both single bubbles and pairs of interacting bubbles are investigated.

Available theoretical models on the dynamics of cylindrical bubbles in planar channels consider only bulk scattered waves generated by the bubble oscillation in the fluid gap [14–16]. These models are appropriate if the walls of the channel are rigid. There are no models that allow for surface waves, which are generated at the solid–fluid interfaces if the channel walls are elastic, whereas experimental observations show that surface waves can play a decisive role in the motion of the fluid layer [2].

In acoustics, several types of surface waves are encountered [17]. The surface waves relevant to our study are Rayleigh, Scholte and Lamb waves, so it is pertinent to remind their definitions. Rayleigh waves are surface waves that travel at an interface between an elastic solid and a vacuum. The existence of these waves was theoretically predicted by Lord Rayleigh in 1885 [18]. The speed of Rayleigh waves is slightly less than that of shear waves propagating in the body of an elastic solid. Scholte waves are surface waves that propagate at an interface between an elastic solid and a fluid. They are named after Scholte, who discovered them in 1947 [19]. Scholte waves are slower than Rayleigh waves. Lamb waves propagate in an elastic plate placed in a vacuum or a fluid [20]. The mathematical analysis of these waves was first developed by Lamb in 1917 [21]. Lamb waves are divided into symmetric and antisymmetric types depending on whether the motion of the elastic medium is symmetric or antisymmetric about the middle plane of the elastic plate.

Although the type of surface waves considered in our study is intimately related to the types described above, it is different from them. In a certain sense, the case under study is opposite to the case considered by Lamb because we deal with surface waves propagating in a fluid channel embedded in an elastic solid. To the best of our knowledge, this case has not been considered as yet. Our study provides a dispersion equation for this type of surface waves, which allows one to evaluate the speed of the surface waves for different values of the channel height. The ultimate purpose of our study is to develop a theoretical model for the volume oscillation of a cylindrical bubble placed between two planar elastic walls, allowing for both bulk and surface waves generated in this system. It should be emphasized that we deal with a special case of bubble confinement. Indeed, the channel height is so small with respect to the bubble size that the bubble always remains in contact with the walls, which leads to a shape close to a cylinder, as depicted in figure 1. Therefore, there is no possibility for the bubble to take the spherical shape or be at a distance from the walls. This fact is a result of the design of microfluidic devices of our interest [2], which makes our case different from those implied by Ilinskii *et al*. [15] and Leighton *et al*. [22–24].

Let us consider a cylindrical gas bubble located in a fluid channel of height *h* between two planar elastic walls (figure 1). The bubble undergoes volume oscillation in response to an imposed acoustic pressure field. The bubble oscillation is considered to be linear and dependent on time as exp(−i*ωt*). For simplicity, the time factor will be omitted from equations. Quantities related to the upper and lower walls will be denoted by subscripts 1 and 2, respectively. The action of the bubble on the fluid–solid interfaces is modelled as a normal harmonic point load. Because the problem under consideration is axially symmetric, cylindrical coordinates are used, whose origin is located in the middle of the bubble and the *z*-axis is perpendicular to the surface of the walls, as shown in figure 1. It should be emphasized that the wave field in the fluid channel consists of two components. One component is the bulk scattered wave caused by the bubble radial oscillation. The second component is caused by the surface waves that propagate in the elastic walls and which, in the process of this propagation, induce perturbations in the fluid. This component will be called surface waves in the fluid channel.

The calculation consists of the following stages. The equations for the surface waves in the elastic walls and in the fluid channel are derived in §2a and §2b, respectively. The derivation is based on Hankel transforms. The boundary conditions on the solid–fluid interfaces are obtained in §2c. The dispersion equation for the surface wave speed is considered in §2d. The resulting expressions for the surface waves in the space domain are calculated by the inverse Hankel transforms in §2e. The equations for the bulk scattered wave in the fluid channel are given in §2f. The boundary conditions on the side bubble surface are obtained in §2g. The equation of bubble oscillation is derived in §2h.

As shown in figure 1, the walls occupy the space with $\left|z\right|>h/2.$ The displacement vector of the *j*th wall (*j*=1,2) can be written as

$${\mathit{\text{u}}}_{j}=\mathbf{\nabla}{\phi}_{j}+\mathbf{\nabla}\times {\mathit{\psi}}_{j},$$

2.1

with _{j} and *ψ*_{j} being called the scalar and vector potentials, respectively. Equation (2.1) is used as the most general mathematical representation for an arbitrary vector field (e.g. [17] or [25]). In the view of axial symmetry, the dependence on the angular coordinate *θ* is absent and hence _{j} and *ψ*_{j} can be written as _{j}=_{j}(*r*,*z*) and *ψ*_{j}=*ψ*_{j}(*r*,*z*)*e*_{θ}, respectively, where *e*_{θ} is the unit azimuth vector. The function _{j} and *ψ*_{j} obey the following equations [25]:

$$\frac{1}{r}\frac{\mathrm{\partial}}{\mathrm{\partial}r}\left(r\frac{\mathrm{\partial}{\phi}_{j}}{\mathrm{\partial}r}\right)+\frac{{\mathrm{\partial}}^{2}{\phi}_{j}}{\mathrm{\partial}{z}^{2}}+{k}_{\mathrm{l}}^{2}{\phi}_{j}=0$$

2.2

and

$$\frac{1}{r}\frac{\mathrm{\partial}}{\mathrm{\partial}r}\left(r\frac{\mathrm{\partial}{\psi}_{j}}{\mathrm{\partial}r}\right)+\frac{{\mathrm{\partial}}^{2}{\psi}_{j}}{\mathrm{\partial}{z}^{2}}+\left({k}_{\mathrm{t}}^{2}-\frac{1}{{r}^{2}}\right){\psi}_{j}=0,$$

2.3

where *k*_{l}=*ω*/*c*_{l} and *k*_{t}=*ω*/*c*_{t} are the wavenumbers of the longitudinal and transverse waves in the solid, respectively. The speeds of these waves, *c*_{l} and *c*_{t}, are given by [17]

$${c}_{\mathrm{l}}=\sqrt{\frac{\mathrm{\lambda}+2\mu}{{\rho}_{\mathrm{s}}}}\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{c}_{\mathrm{t}}=\sqrt{\frac{\mu}{{\rho}_{\mathrm{s}}}},$$

2.4

where λ and *μ* are the Lamé coefficients of the walls (*μ* is also called shear modulus) and *ρ*_{s} is the wall density. The components of the displacement vector *u*_{j} are calculated by [25]

$${u}_{jr}=\frac{\mathrm{\partial}{\phi}_{j}}{\mathrm{\partial}r}-\frac{\mathrm{\partial}{\psi}_{j}}{\mathrm{\partial}z}\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{u}_{jz}=\frac{\mathrm{\partial}{\phi}_{j}}{\mathrm{\partial}z}+\frac{1}{r}\frac{\mathrm{\partial}(r{\psi}_{j})}{\mathrm{\partial}r}.$$

2.5

Solutions to equations (2.2) and (2.3) can be found by applying Hankel transforms with respect to the radial variable *r*. The Hankel transform of order *n* of a function *f*(*r*) is defined as [26]

$${f}^{{\mathrm{H}}_{n}}(k)={\int}_{0}^{\mathrm{\infty}}f(r){J}_{n}(kr)r\hspace{0.17em}\mathrm{d}r,$$

2.6

where *k* is the parameter of the Hankel transform and *J*_{n} is the Bessel function of the first kind of order *n*. The inverse transform is given by

$$f(r)={\int}_{0}^{\mathrm{\infty}}{f}^{{\mathrm{H}}_{n}}(k){J}_{n}(kr)k\hspace{0.17em}\mathrm{d}k.$$

2.7

Applying the Hankel transforms of zero and first orders to equations (2.2) and (2.3), respectively, one obtains

$$\frac{{\mathrm{d}}^{2}{\phi}_{j}^{{\mathrm{H}}_{0}}}{\mathrm{d}{z}^{2}}-{q}_{\mathrm{l}}^{2}{\phi}_{j}^{{\mathrm{H}}_{0}}=0$$

2.8

and

$$\frac{{\mathrm{d}}^{2}{\psi}_{j}^{{\mathrm{H}}_{1}}}{\mathrm{d}{z}^{2}}-{q}_{\mathrm{t}}^{2}{\psi}_{j}^{{\mathrm{H}}_{1}}=0,$$

2.9

where

$${q}_{\mathrm{l}}^{2}={k}^{2}-{k}_{\mathrm{l}}^{2}\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{q}_{\mathrm{t}}^{2}={k}^{2}-{k}_{\mathrm{t}}^{2}.$$

2.10

Solutions of equations (2.8) and (2.9) which remain finite for large values of |*z*| are

$$\begin{array}{rl}{\phi}_{1}^{{\mathrm{H}}_{0}}& ={\mathit{\Phi}}_{1}(k){\mathrm{e}}^{-{q}_{\mathrm{l}}(z-h/2)},\end{array}$$

2.11

$$\begin{array}{rl}{\psi}_{1}^{{\mathrm{H}}_{1}}& ={\mathit{\Psi}}_{1}(k){\mathrm{e}}^{-{q}_{\mathrm{t}}(z-h/2)},\end{array}$$

2.12

$$\begin{array}{rl}{\phi}_{2}^{{\mathrm{H}}_{0}}& ={\mathit{\Phi}}_{2}(k){\mathrm{e}}^{{q}_{\mathrm{l}}(z+h/2)}\end{array}$$

2.13

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\mathrm{and}\phantom{\rule{2em}{0ex}}{\psi}_{2}^{{\mathrm{H}}_{1}}& ={\mathit{\Psi}}_{2}(k){\mathrm{e}}^{{q}_{\mathrm{t}}(z+h/2)},\end{array}$$

2.14

where *Φ*_{j}(*k*) and *Ψ*_{j}(*k*) are functions to be determined from the boundary conditions on the solid–fluid interfaces at *z* = ±*h*/2.

Application of the Hankel transform to equations (2.5) yields

$${u}_{jr}^{{\mathrm{H}}_{1}}=-k{\phi}_{j}^{{\mathrm{H}}_{0}}-\frac{\mathrm{d}{\psi}_{j}^{{\mathrm{H}}_{1}}}{\mathrm{d}z}\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{u}_{jz}^{{\mathrm{H}}_{0}}=\frac{\mathrm{d}{\phi}_{j}^{{\mathrm{H}}_{0}}}{\mathrm{d}z}+k{\psi}_{j}^{{\mathrm{H}}_{1}},$$

2.15

Substituting equations (2.11)–(2.14) into equations (2.15), one finds

$$\begin{array}{rl}{u}_{1r}^{{\mathrm{H}}_{1}}& =-k{\mathit{\Phi}}_{1}(k){\mathrm{e}}^{-{q}_{\mathrm{l}}(z-h/2)}+{q}_{\mathrm{t}}{\mathit{\Psi}}_{1}(k){\mathrm{e}}^{-{q}_{\mathrm{t}}(z-h/2)},\end{array}$$

2.16

$$\begin{array}{rl}{u}_{1z}^{{\mathrm{H}}_{0}}& =-{q}_{\mathrm{l}}{\mathit{\Phi}}_{1}(k){\mathrm{e}}^{-{q}_{\mathrm{l}}(z-h/2)}+k{\mathit{\Psi}}_{1}(k){\mathrm{e}}^{-{q}_{\mathrm{t}}(z-h/2)},\end{array}$$

2.17

$$\begin{array}{rl}{u}_{2r}^{{\mathrm{H}}_{1}}& =-k{\mathit{\Phi}}_{2}(k){\mathrm{e}}^{{q}_{\mathrm{l}}(z+h/2)}-{q}_{\mathrm{t}}{\mathit{\Psi}}_{2}(k){\mathrm{e}}^{{q}_{\mathrm{t}}(z+h/2)}\end{array}$$

2.18

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{u}_{2z}^{{\mathrm{H}}_{0}}& ={q}_{\mathrm{l}}{\mathit{\Phi}}_{2}(k){\mathrm{e}}^{{q}_{\mathrm{l}}(z+h/2)}+k{\mathit{\Psi}}_{2}(k){\mathrm{e}}^{{q}_{\mathrm{t}}(z+h/2)}.\end{array}$$

2.19

The fluid area corresponds to −*h*/2≤*z*≤*h*/2. The displacement vector in the fluid can be written as

$${\mathit{\text{u}}}_{\mathrm{f}}=\mathbf{\nabla}{\phi}_{\mathrm{f}},$$

2.20

where the potential _{f} obeys the equation [27]

$$\frac{1}{r}\frac{\mathrm{\partial}}{\mathrm{\partial}r}\left(r\frac{\mathrm{\partial}{\phi}_{\mathrm{f}}}{\mathrm{\partial}r}\right)+\frac{{\mathrm{\partial}}^{2}{\phi}_{\mathrm{f}}}{\mathrm{\partial}{z}^{2}}+{k}_{\mathrm{f}}^{2}{\phi}_{\mathrm{f}}=0,$$

2.21

in which *k*_{f} = *ω*/*c*_{f} and *c*_{f} is the speed of sound in the fluid. The fluid velocity corresponding to *u*_{f} is *v*_{f}=−i*ω**u*_{f}. The components of *u*_{f} are calculated by

$${u}_{\mathrm{fr}}=\frac{\mathrm{\partial}{\phi}_{\mathrm{f}}}{\mathrm{\partial}r}\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{u}_{fz}=\frac{\mathrm{\partial}{\phi}_{\mathrm{f}}}{\mathrm{\partial}z}.$$

2.22

Applying the Hankel transform of zero order to equation (2.21), one obtains

$$\frac{{\mathrm{d}}^{2}{\phi}_{\mathrm{f}}^{{\mathrm{H}}_{0}}}{\mathrm{d}{z}^{2}}-{q}_{\mathrm{f}}^{2}{\phi}_{\mathrm{f}}^{{\mathrm{H}}_{0}}=0,$$

2.23

where

$${q}_{\mathrm{f}}^{2}={k}^{2}-{k}_{\mathrm{f}}^{2}.$$

2.24

The general solution of equation (2.23) is given by

$${\phi}_{\mathrm{f}}^{{\mathrm{H}}_{0}}={\mathit{\Phi}}_{\mathrm{f}1}(k){\mathrm{e}}^{{q}_{\mathrm{f}}(z-h/2)}+{\mathit{\Phi}}_{\mathrm{f}2}(k){\mathrm{e}}^{-{q}_{\mathrm{f}}(z+h/2)},$$

2.25

where *Φ*_{f1}(*k*) and *Φ*_{f2}(*k*) are functions to be determined from the boundary conditions at *z* = ±*h*/2.

Application of the Hankel transform to equations (2.22) yields

$${u}_{\mathrm{fr}}^{{\mathrm{H}}_{1}}=-k{\phi}_{\mathrm{f}}^{{\mathrm{H}}_{0}}\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{u}_{fz}^{{\mathrm{H}}_{0}}=\frac{\mathrm{d}{\phi}_{\mathrm{f}}^{{\mathrm{H}}_{0}}}{\mathrm{d}z}.$$

2.26

Substituting equation (2.25) into equations (2.26), one finds

$${u}_{\mathrm{fr}}^{{\mathrm{H}}_{1}}=-k\left[{\mathit{\Phi}}_{\mathrm{f}1}(k){\mathrm{e}}^{{q}_{\mathrm{f}}(z-h/2)}+{\mathit{\Phi}}_{\mathrm{f}2}(k){\mathrm{e}}^{-{q}_{\mathrm{f}}(z+h/2)}\right]$$

2.27

and

$${u}_{fz}^{{\mathrm{H}}_{0}}={q}_{\mathrm{f}}\left[{\mathit{\Phi}}_{\mathrm{f}1}(k){\mathrm{e}}^{{q}_{\mathrm{f}}(z-h/2)}-{\mathit{\Phi}}_{\mathrm{f}2}(k){\mathrm{e}}^{-{q}_{\mathrm{f}}(z+h/2)}\right].$$

2.28

The boundary conditions on the surfaces of the elastic walls are written as follows [25,28]:

$$\begin{array}{rl}{u}_{1z}& ={u}_{fz}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{h}{2},\phantom{\rule{1em}{0ex}}{u}_{2z}={u}_{fz}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{-h}{2},\end{array}$$

2.29

$$\begin{array}{rl}{T}_{1rz}& =0\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{h}{2},\phantom{\rule{1em}{0ex}}{T}_{2rz}=0\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{-h}{2}\end{array}$$

2.30

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{T}_{1zz}& =-{p}_{\mathrm{f}}-Q\frac{\delta (r)}{2\pi r}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{h}{2},\phantom{\rule{1em}{0ex}}{T}_{2zz}=-{p}_{\mathrm{f}}-Q\frac{\delta (r)}{2\pi r}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{-h}{2},\end{array}$$

2.31

where *T*_{jrz} and *T*_{jzz} are the shear and normal stresses in the walls, *p*_{f} is the fluid pressure corresponding to *u*_{f}, *Q* is the magnitude of the point load produced by the bubble at the points with *r*=0 and *z*=±*h*/2 and *δ*(*r*) is the Dirac delta function. Equations (2.29) follow from the continuity of the normal displacement at the solid–fluid interfaces, equations (2.30) mean that the shear stress vanishes at the interfaces and equations (2.31) describe the normal stress balance at the interfaces.

The expressions for *T*_{jrz}, *T*_{jzz} and *p*_{f} are given by [25,27]

$$\begin{array}{rl}{T}_{jrz}& =\mu \left(\frac{\mathrm{\partial}{u}_{jr}}{\mathrm{\partial}z}+\frac{\mathrm{\partial}{u}_{jz}}{\mathrm{\partial}r}\right),\end{array}$$

2.32

$$\begin{array}{rl}{T}_{jzz}& =(\mathrm{\lambda}+2\mu )\frac{\mathrm{\partial}{u}_{jz}}{\mathrm{\partial}z}+\frac{\mathrm{\lambda}}{r}\frac{\mathrm{\partial}(r{u}_{jr})}{\mathrm{\partial}r}\end{array}$$

2.33

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{p}_{\mathrm{f}}& ={\rho}_{\mathrm{f}}{\omega}^{2}{\phi}_{\mathrm{f}},\end{array}$$

2.34

where *ρ*_{f} is the fluid density. Applying the Hankel transforms to equations (2.29)–(2.34), one obtains

$$\begin{array}{rl}{u}_{1z}^{{\mathrm{H}}_{0}}& ={u}_{fz}^{{\mathrm{H}}_{0}}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{h}{2},\phantom{\rule{1em}{0ex}}{u}_{2z}^{{\mathrm{H}}_{0}}={u}_{fz}^{{\mathrm{H}}_{0}}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{-h}{2},\end{array}$$

2.35

$$\begin{array}{rl}{T}_{1rz}^{{\mathrm{H}}_{1}}& =0\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{h}{2},\phantom{\rule{1em}{0ex}}{T}_{2rz}^{{\mathrm{H}}_{1}}=0\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{-h}{2},\end{array}$$

2.36

$$\begin{array}{rl}{T}_{1zz}^{{\mathrm{H}}_{0}}& =-{p}_{\mathrm{f}}^{{\mathrm{H}}_{0}}-\frac{Q}{2\pi}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{h}{2},\phantom{\rule{1em}{0ex}}{T}_{2zz}^{{\mathrm{H}}_{0}}=-{p}_{\mathrm{f}}^{{\mathrm{H}}_{0}}-\frac{Q}{2\pi}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}z=\frac{-h}{2},\end{array}$$

2.37

$$\begin{array}{rl}{T}_{jrz}^{{\mathrm{H}}_{1}}& =-\mu \left[2k\frac{\mathrm{d}{\phi}_{j}^{{\mathrm{H}}_{0}}}{\mathrm{d}z}+(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\psi}_{j}^{{\mathrm{H}}_{1}}\right],\end{array}$$

2.38

$$\begin{array}{rl}{T}_{jzz}^{{\mathrm{H}}_{0}}& =\mu \left[(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\phi}_{j}^{{\mathrm{H}}_{0}}+2k\frac{\mathrm{d}{\psi}_{j}^{{\mathrm{H}}_{1}}}{\mathrm{d}z}\right]\end{array}$$

2.39

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{p}_{\mathrm{f}}^{{\mathrm{H}}_{0}}& ={\rho}_{\mathrm{f}}{\omega}^{2}{\phi}_{\mathrm{f}}^{{\mathrm{H}}_{0}}.\end{array}$$

2.40

Substitution of equations (2.11)–(2.14), (2.17), (2.19), (2.25), (2.28) and (2.38)–(2.40) into equations (2.35)–(2.37) yields the following equations for the functions *Φ*_{j}(*k*), *Ψ*_{j}(*k*), *Φ*_{f1}(*k*) and *Φ*_{f2}(*k*):

$$\begin{array}{rl}k{\mathit{\Psi}}_{1}-{q}_{\mathrm{l}}{\mathit{\Phi}}_{1}& ={q}_{\mathrm{f}}\left({\mathit{\Phi}}_{\mathrm{f}1}-{\mathit{\Phi}}_{\mathrm{f}2}{\mathrm{e}}^{-{q}_{\mathrm{f}}h}\right),\end{array}$$

2.41

$$\begin{array}{rl}k{\mathit{\Psi}}_{2}+{q}_{\mathrm{l}}{\mathit{\Phi}}_{2}& ={q}_{\mathrm{f}}\left({\mathit{\Phi}}_{\mathrm{f}1}{\mathrm{e}}^{-{q}_{\mathrm{f}}h}-{\mathit{\Phi}}_{\mathrm{f}2}\right),\end{array}$$

2.42

$$\begin{array}{rl}2k{q}_{\mathrm{l}}{\mathit{\Phi}}_{1}-(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\mathit{\Psi}}_{1}& =0,\end{array}$$

2.43

$$\begin{array}{rl}2k{q}_{\mathrm{l}}{\mathit{\Phi}}_{2}+(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\mathit{\Psi}}_{2}& =0,\end{array}$$

2.44

$$\begin{array}{rl}(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\mathit{\Phi}}_{1}-2k{q}_{\mathrm{t}}{\mathit{\Psi}}_{1}& =-\frac{{\rho}_{\mathrm{f}}{\omega}^{2}}{\mu}\left({\mathit{\Phi}}_{\mathrm{f}1}+{\mathit{\Phi}}_{\mathrm{f}2}{\mathrm{e}}^{-{q}_{\mathrm{f}}h}\right)-\frac{Q}{2\pi \mu}\end{array}$$

2.45

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\mathit{\Phi}}_{2}+2k{q}_{\mathrm{t}}{\mathit{\Psi}}_{2}& =-\frac{{\rho}_{\mathrm{f}}{\omega}^{2}}{\mu}\left({\mathit{\Phi}}_{\mathrm{f}1}{\mathrm{e}}^{-{q}_{\mathrm{f}}h}+{\mathit{\Phi}}_{\mathrm{f}2}\right)-\frac{Q}{2\pi \mu}.\end{array}$$

2.46

The waves produced by the bubble in the elastic walls are symmetric, because the walls are assumed to be identical and the action of the bubble on them is symmetric about the point *z*=0. Therefore, the following conditions should be met [17]:

$${u}_{1r}^{{\mathrm{H}}_{1}}\left(z=\frac{h}{2}\right)={u}_{2r}^{{\mathrm{H}}_{1}}\left(z=\frac{-h}{2}\right)\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{u}_{1z}^{{\mathrm{H}}_{0}}\left(z=\frac{h}{2}\right)=-{u}_{2z}^{{\mathrm{H}}_{0}}\left(z=\frac{-h}{2}\right).$$

2.47

From equations (2.16)–(2.19), it follows that conditions (2.47) are met if

$${\mathit{\Phi}}_{2}(k)={\mathit{\Phi}}_{1}(k)\phantom{\rule{1em}{0ex}}\mathrm{and}\phantom{\rule{1em}{0ex}}{\mathit{\Psi}}_{2}(k)=-{\mathit{\Psi}}_{1}(k).$$

2.48

Substitution of equations (2.48) into equations (2.41)–(2.46) shows that

$${\mathit{\Phi}}_{\mathrm{f}2}(k)={\mathit{\Phi}}_{\mathrm{f}1}(k)={\mathit{\Phi}}_{\mathrm{f}}(k),$$

2.49

and that only three of these equations remain independent. As a result, we obtain the following set:

$$\begin{array}{rl}{q}_{\mathrm{l}}{\mathit{\Phi}}_{1}-k{\mathit{\Psi}}_{1}+{q}_{\mathrm{f}}(1-{\mathrm{e}}^{-{q}_{\mathrm{f}}h}){\mathit{\Phi}}_{\mathrm{f}}& =0,\end{array}$$

2.50

$$\begin{array}{rl}2k{q}_{\mathrm{l}}{\mathit{\Phi}}_{1}-(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\mathit{\Psi}}_{1}& =0\end{array}$$

2.51

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}(2{k}^{2}-{k}_{\mathrm{t}}^{2}){\mathit{\Phi}}_{1}-2k{q}_{\mathrm{t}}{\mathit{\Psi}}_{1}+\frac{{\rho}_{\mathrm{f}}{\omega}^{2}}{\mu}(1+{\mathrm{e}}^{-{q}_{\mathrm{f}}h}){\mathit{\Phi}}_{\mathrm{f}}& =-\frac{Q}{2\pi \mu}.\end{array}$$

2.52

The solutions of this set are

$$\begin{array}{rl}{\mathit{\Phi}}_{1}(k)& =-\frac{2{k}^{2}-{k}_{\mathrm{t}}^{2}}{D(k)}\frac{Q}{2\pi \mu},\end{array}$$

2.53

$$\begin{array}{rl}{\mathit{\Psi}}_{1}(k)& =-\frac{2k{q}_{\mathrm{l}}}{D(k)}\frac{Q}{2\pi \mu}\end{array}$$

2.54

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{\mathit{\Phi}}_{\mathrm{f}}(k)& =-\frac{{k}_{\mathrm{t}}^{2}{q}_{\mathrm{l}}}{{q}_{\mathrm{f}}(1-{\mathrm{e}}^{-{q}_{\mathrm{f}}h})D(k)}\frac{Q}{2\pi \mu},\end{array}$$

2.55

where

$$D(k)={(2{k}^{2}-{k}_{\mathrm{t}}^{2})}^{2}-4{k}^{2}{q}_{\mathrm{l}}{q}_{\mathrm{t}}+\frac{{\rho}_{\mathrm{f}}{k}_{\mathrm{t}}^{4}{q}_{\mathrm{l}}}{{\rho}_{\mathrm{s}}{q}_{\mathrm{f}}}\frac{1+{\mathrm{e}}^{-{q}_{\mathrm{f}}h}}{1-{\mathrm{e}}^{-{q}_{\mathrm{f}}h}}.$$

2.56

The theoretical analysis of Rayleigh waves produced by a point load, based on the use of Hankel integral transforms, is given in a book by Achenbach [25]. An analogous consideration of Scholte waves is presented by Zhu *et al*. [28]. These works show that at large distances from the loading point, where general solutions similar to equations (2.53) and (2.54) turn into far-field expressions giving the Rayleigh and Scholte waves, the dispersion equations for these waves are identical to the denominators of the general solutions similar to equations (2.53) and (2.54). This follows from the fact that a vanishing denominator gives rise to a dominant term in the spectrum of emitted waves. Analogously, in our case the dispersion equation for the surface waves of our interest is given by

$$D(k)={(2{k}^{2}-{k}_{\mathrm{t}}^{2})}^{2}-4{k}^{2}{q}_{\mathrm{l}}{q}_{t}+\frac{{\rho}_{\mathrm{f}}{k}_{\mathrm{t}}^{4}{q}_{\mathrm{l}}}{{\rho}_{\mathrm{s}}{q}_{\mathrm{f}}}\frac{1+{\mathrm{e}}^{-{q}_{\mathrm{f}}h}}{1-{\mathrm{e}}^{-{q}_{\mathrm{f}}h}}=0.$$

2.57

The roots of this equation are the wavenumbers of Lamb-type symmetric waves that propagate in the microfluidic system under study. If the fluid density tends to zero, *ρ*_{f}→0, equation (2.57) reduces to the equation for Rayleigh waves [17,25], while for infinite channel height, $h\to \mathrm{\infty}$, it reduces to the equation for Scholte waves [17,28]. Equation (2.57) allows one to evaluate the surface wave speed as a function of *h*. Examples of such evaluations are presented in §3a. It is noteworthy that equation (2.57) is derived for the first time so it can be considered as one of the main results of this study.

Applying the inverse Hankel transform to equations (2.16), (2.17), (2.27), (2.28) and (2.40) and using equations (2.53)–(2.55), one obtains

$$\begin{array}{rl}{u}_{1r}& =\frac{Q{k}_{\mathrm{t}}}{2\pi \mu}{I}_{sr}(r,z),\end{array}$$

2.58

$$\begin{array}{rl}{u}_{1z}& =\frac{Q{k}_{\mathrm{t}}}{2\pi \mu}{I}_{sz}(r,z),\end{array}$$

2.59

$$\begin{array}{rl}{u}_{\mathrm{fr}}& =\frac{Q{k}_{\mathrm{t}}}{\pi \mu}{I}_{\mathrm{fr}}(r,z),\end{array}$$

2.60

$$\begin{array}{rl}{u}_{fz}& =-\frac{Q{k}_{\mathrm{t}}}{\pi \mu}{I}_{fz}(r,z)\end{array}$$

2.61

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{p}_{\mathrm{f}}& =-\frac{Q{\rho}_{\mathrm{f}}{\omega}^{2}}{\pi \mu}{I}_{\mathrm{p}}(r,z),\end{array}$$

2.62

where

$$\begin{array}{rl}{I}_{sr}(r,z)& ={\int}_{0}^{\mathrm{\infty}}\frac{1}{\stackrel{~}{D}(\xi )}\left[(2{\xi}^{2}-1){\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{l}}{k}_{\mathrm{t}}(z-h/2)}-2{\stackrel{~}{q}}_{\mathrm{l}}{\stackrel{~}{q}}_{\mathrm{t}}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{t}}{k}_{\mathrm{t}}(z-h/2)}\right]{J}_{1}(\xi {k}_{\mathrm{t}}r){\xi}^{2}\mathrm{d}\xi ,\end{array}$$

2.63

$$\begin{array}{rl}{I}_{sz}(r,z)& ={\int}_{0}^{\mathrm{\infty}}\frac{{\stackrel{~}{q}}_{\mathrm{l}}}{\stackrel{~}{D}(\xi )}\left[(2{\xi}^{2}-1){\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{l}}{k}_{\mathrm{t}}(z-h/2)}-2{\xi}^{2}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{t}}{k}_{\mathrm{t}}(z-h/2)}\right]{J}_{0}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi ,\end{array}$$

2.64

$$\begin{array}{rl}{I}_{\mathrm{fr}}(r,z)& ={\int}_{0}^{\mathrm{\infty}}\frac{{\stackrel{~}{q}}_{\mathrm{l}}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h/2}\mathrm{cosh}({\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}z)}{{\stackrel{~}{q}}_{\mathrm{f}}\stackrel{~}{D}(\xi )(1-{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h})}{J}_{1}(\xi {k}_{\mathrm{t}}r){\xi}^{2}\mathrm{d}\xi \end{array}$$

2.65

$$\begin{array}{rl}{I}_{fz}(r,z)& ={\int}_{0}^{\mathrm{\infty}}\frac{{\stackrel{~}{q}}_{l}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{t}h/2}\mathrm{sinh}({\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}z)}{\stackrel{~}{D}(\xi )(1-{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h})}{J}_{0}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi \end{array}$$

2.66

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{I}_{\mathrm{p}}(r,z)& ={\int}_{0}^{\mathrm{\infty}}\frac{{\stackrel{~}{q}}_{\mathrm{l}}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h/2}\mathrm{cosh}({\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}z)}{{\stackrel{~}{q}}_{\mathrm{f}}\stackrel{~}{D}(\xi )(1-{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h})}{J}_{0}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi .\end{array}$$

2.67

Equations (2.63)–(2.67) are expressed in terms of normalized quantities that are defined as follows:

$$\begin{array}{rl}\xi & =\frac{k}{{k}_{\mathrm{t}}}=\frac{{c}_{\mathrm{t}}}{c},\phantom{\rule{1em}{0ex}}{\xi}_{\mathrm{f}}=\frac{{k}_{\mathrm{f}}}{{k}_{\mathrm{t}}}=\frac{{c}_{\mathrm{t}}}{{c}_{\mathrm{f}}},\phantom{\rule{1em}{0ex}}{\xi}_{\mathrm{l}}=\frac{{k}_{\mathrm{l}}}{{k}_{\mathrm{t}}}=\frac{{c}_{\mathrm{t}}}{{c}_{\mathrm{l}}},\end{array}$$

2.68

$$\begin{array}{rl}{\stackrel{~}{q}}_{\mathrm{f}}& =\sqrt{{\xi}^{2}-{\xi}_{\mathrm{f}}^{2}},\hspace{0.17em}{\stackrel{~}{q}}_{\mathrm{l}}=\sqrt{{\xi}^{2}-{\xi}_{\mathrm{l}}^{2}},\hspace{0.17em}{\stackrel{~}{q}}_{\mathrm{t}}=\sqrt{{\xi}^{2}-1}\end{array}$$

2.69

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}\stackrel{~}{D}(\xi )& ={(2{\xi}^{2}-1)}^{2}-4{\xi}^{2}{\stackrel{~}{q}}_{\mathrm{l}}{\stackrel{~}{q}}_{\mathrm{t}}+\frac{{\rho}_{\mathrm{f}}{\stackrel{~}{q}}_{\mathrm{l}}}{{\rho}_{\mathrm{s}}{\stackrel{~}{q}}_{\mathrm{f}}}\frac{1+{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h}}{1-{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h}},\end{array}$$

2.70

where *c* stands for the speed that corresponds to the wavenumber *k*. The normalization by *k*_{t} provides the simplest form of the equations. The method of numerical calculation of integrals (2.63)–(2.67) is described in appendix A.

According to Ilinskii *et al*. [15], the bulk scattered wave generated by the bubble oscillation in the fluid channel can be represented by

$${\mathit{\text{v}}}_{\mathrm{b}}=\mathbf{\nabla}{\phi}_{\mathrm{b}},$$

2.71

where *v*_{b} is the fluid velocity and the velocity potential _{b} is given by the equation

$${\phi}_{\mathrm{b}}=A{H}_{0}^{(1)}({k}_{\mathrm{f}}r),$$

2.72

in which ${H}_{0}^{(1)}$ is the Hankel function of the first kind of order zero [29] and *A* is a constant to be determined from the boundary conditions on the bubble surface. Equations (2.71) and (2.72) satisfy the wave equation for the fluid motion in the case of cylindrical symmetry. The fluid pressure corresponding to *v*_{b} is calculated by

$${p}_{\mathrm{b}}=\mathrm{i}\omega {\rho}_{\mathrm{f}}{\phi}_{\mathrm{b}}.$$

2.73

It should be mentioned that the wavelength of the bulk wave, λ_{f}=2*π*/*k*_{f}, is large compared with the wavelength of the surface waves as well as to the dimensions of the microfluidic system under study. This means that the spatial scale of the variations of *p*_{b} is much greater than that of the perturbations produced by the surface waves. For this reason, *p*_{b} does not contribute to the above perturbations and hence does not appear in equation (2.31).

Let us assume that the time-varying bubble radius can be represented as

$$R(t)={R}_{0}(1+a{\mathrm{e}}^{-\mathrm{i}\omega t}),$$

2.74

where *R*_{0} is the equilibrium bubble radius. Considering that the net fluid velocity is the sum of the velocities of the bulk and surface waves, the boundary condition for the normal velocity at the side bubble surface can be written as

$${v}_{\mathrm{br}}+{v}_{\mathrm{fr}}=-\mathrm{i}\omega {R}_{0}a\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}r={R}_{0},$$

2.75

where *v*_{br} is the radial component of *v*_{b},*v*_{fr} is the radial component of the surface wave velocity, given by *v*_{fr}=−i*ωu*_{fr}, and the right-hand term is the amplitude of the velocity of the bubble oscillation. Assuming that *h* is small compared with the transverse wavelength λ_{t}=2*π*/*k*_{t}, the value of *u*_{fr} can be taken at *z* = 0; see equations (2.60) and (2.65). Substitution of equations (2.60), (2.71) and (2.72) into equation (2.75) yields

$$A=\frac{\mathrm{i}\omega}{{k}_{\mathrm{f}}{H}_{1}^{(1)}({k}_{\mathrm{f}}{R}_{0})}\left[{R}_{0}a-\frac{{k}_{\mathrm{t}}Q}{\pi \mu}{I}_{\mathrm{fr}}({R}_{0},0)\right].$$

2.76

The pressure balance at the bubble surface is defined as

$${p}_{\mathrm{g}}={P}_{0}+{p}_{\mathrm{ac}}+{p}_{\mathrm{st}}+({p}_{\mathrm{b}}+{p}_{\mathrm{f}}){\mathrm{e}}^{-\mathrm{i}\omega t}\phantom{\rule{1em}{0ex}}\mathrm{at}\phantom{\rule{1em}{0ex}}r={R}_{0},$$

2.77

where *p*_{g} is the gas pressure within the bubble, *P*_{0} is the hydrostatic pressure in the fluid, *p*_{ac} is the driving acoustic pressure and *p*_{st} is the surface tension pressure. Calculating the bubble volume as

$$V=\pi {R}^{2}h\approx \pi {R}_{0}^{2}h(1+2a{\mathrm{e}}^{-\mathrm{i}\omega t})={V}_{0}(1+2a{\mathrm{e}}^{-\mathrm{i}\omega t}),$$

2.78

one obtains

$${p}_{\mathrm{g}}={P}_{\mathrm{g}0}{\left(\frac{{V}_{0}}{V}\right)}^{\gamma}\approx {P}_{\mathrm{g}0}(1-2\gamma a{\mathrm{e}}^{-\mathrm{i}\omega t}),$$

2.79

where *P*_{g0} is the equilibrium gas pressure and *γ* is the ratio of specific heats of the gas. The acoustic pressure *p*_{ac} is specified by

$${p}_{\mathrm{ac}}={P}_{\mathrm{a}}{\mathrm{e}}^{-\mathrm{i}\omega t},$$

2.80

where *P*_{a} denotes the acoustic pressure amplitude. For *p*_{st}, one has

$${p}_{\mathrm{st}}=\frac{{\sigma}_{\mathrm{f}}}{R}\approx \frac{{\sigma}_{\mathrm{f}}}{{R}_{0}}(1-a{\mathrm{e}}^{-\mathrm{i}\omega t}),$$

2.81

where *σ*_{f} is the surface tension coefficient for the fluid–gas interface. Substituting equations (2.62), (2.72), (2.73) and (2.79)–(2.81) into equation (2.77) and separating constant and time-varying terms, one obtains

$${P}_{\mathrm{g}0}={P}_{0}+\frac{{\sigma}_{\mathrm{f}}}{{R}_{0}}$$

2.82

and

$$a\left(2\gamma {P}_{\mathrm{g}0}-\frac{{\sigma}_{\mathrm{f}}}{{R}_{0}}\right)=\frac{{\rho}_{\mathrm{f}}{\omega}^{2}Q}{\pi \mu}{I}_{\mathrm{p}}({R}_{0},0)-\mathrm{i}\omega {\rho}_{\mathrm{f}}A{H}_{0}^{(1)}({k}_{\mathrm{f}}{R}_{0})-{P}_{\mathrm{a}}.$$

2.83

Substituting equations (2.76) and (2.82) into equation (2.83) and using the second of equations (2.4) to eliminate *μ*, one finds the following expression for *a*:

$$a=\frac{-({P}_{\mathrm{a}}+C)}{B},$$

2.84

where

$$B=2\gamma {P}_{0}+\frac{(2\gamma -1){\sigma}_{\mathrm{f}}}{{R}_{0}}-\frac{{\rho}_{\mathrm{f}}{\omega}^{2}{R}_{0}^{2}{H}_{0}^{(1)}({\alpha}_{\mathrm{f}})}{{\alpha}_{\mathrm{f}}{H}_{1}^{(1)}({\alpha}_{\mathrm{f}})}$$

2.85

and

$$C=-\frac{{\rho}_{\mathrm{f}}{k}_{\mathrm{t}}^{2}Q}{\pi {\rho}_{\mathrm{s}}}\left[{I}_{\mathrm{p}}({R}_{0},0)-\frac{{\alpha}_{\mathrm{t}}{H}_{0}^{(1)}({\alpha}_{\mathrm{f}})}{{\alpha}_{\mathrm{f}}{H}_{1}^{(1)}({\alpha}_{\mathrm{f}})}{I}_{\mathrm{fr}}({R}_{0},0)\right],$$

2.86

*α*_{f}=*k*_{f}*R*_{0} and *α*_{t}=*k*_{t}*R*_{0}.

To derive an equation for the amplitude of the bubble oscillation *a*, we need an expression for the load magnitude *Q* in terms of *a*. Clearly, *Q* should be proportional to the variation of the gas pressure in the bubble, given by equation (2.79). It is admissible to assume that *Q* is also proportional to the area of the lateral section of the bubble, $\pi {R}_{0}^{2}.$ Therefore, *Q* can be represented as

$$Q=-2\pi \gamma {P}_{\mathrm{g}0}{R}_{0}^{2}\epsilon a,$$

2.87

where *ε* is a proportionality coefficient that allows for the fact that the effective area that describes the action of the bubble on the wall can be different from $\pi {R}_{0}^{2}.$ In addition, *ε* can be a complex quantity because a phase shift can exist between the variation of the gas pressure and the reaction of the walls to it. This occurs because the walls may behave as a viscoelastic material with complex elastic and shear moduli, which leads to a phase lag between stress and strain [30].

Substitution of equation (2.87) into equations (2.84) and (2.86) yields

$$a=\frac{-{P}_{\mathrm{a}}}{B+S},$$

2.88

where

$$S=\frac{2\gamma {P}_{\mathrm{g}0}{\rho}_{\mathrm{f}}{\alpha}_{\mathrm{t}}^{2}\epsilon}{{\rho}_{\mathrm{s}}}\left[{I}_{\mathrm{p}}({R}_{0},0)-\frac{{\alpha}_{\mathrm{t}}{H}_{0}^{(1)}({\alpha}_{\mathrm{f}})}{{\alpha}_{\mathrm{f}}{H}_{1}^{(1)}({\alpha}_{\mathrm{f}})}{I}_{\mathrm{fr}}({R}_{0},0)\right].$$

2.89

Equation (2.88) allows one to calculate the response of the bubble as a function of *ε* for different values of the material and acoustic parameters. It should be mentioned that for *ε*=0, equation (2.88) reduces to the result obtained by Ilinskii *et al*. [15] for a cylindrical bubble between two rigid walls.

The speed of the surface waves propagating at the solid–fluid interfaces is given by solutions of equation (2.57). Numerical evaluations were performed for the following values of the physical parameters: Young’s modulus *E*=1.6MPa, Poisson’s ratio *σ*=0.499,*ρ*_{s}=970kgm^{−3},*c*_{f}=1481ms^{−1} and *ρ*_{f}=998kgm^{−3}. These values are typical of the elastic walls of a microfluidic channel made of polydimethylsiloxane (PDMS) and filled with water. The transverse wave speed for these parameters is *c*_{t}=23.456ms^{−3}. This quantity sets the upper limit for the speed of surface waves [17]. The evaluations show that equation (2.57) has one real root which gives a frequency-dependent speed *c*_{s}<*c*_{t}. Figure 2 exemplifies the frequency dependence of *c*_{s} at *h*=25μm. The dependence of *c*_{s} on *h* is demonstrated by figure 3. Note that the values of *c*_{s} are normalized by *c*_{t}, and those of *h*, by λ_{t}=*c*_{t}/*f*, where *f* is the driving frequency. With this normalization, the plot shown in figure 3 remains the same for all frequencies. Figure 3 reveals that *c*_{s} increases with increasing *h* and tends to the speed of the Scholte wave, *c*_{Sch}, when *h* tends to λ_{t}. However, at small *h*, the presence of the second boundary makes *c*_{s} much smaller than *c*_{Sch}.

The theory of surface waves predicts that the behaviour of integrals (2.63)–(2.67) should coincide with the behaviour of the Hankel functions of zero and first orders at large distances from the loading point [17]. The calculation of the integrals by the method described in appendix A conforms to this expectation. As an example, figure 4 compares the behaviour of the integral *I*_{sz}, which describes the vertical displacement of wall 1, and the Hankel function ${H}_{0}^{(1)}({k}_{\mathrm{s}}r)$, where *k*_{s}=*ω*/*c*_{s} is the wavenumber of the surface wave, given by equation (2.57). The calculations were made at *f*=30kHz and *h*=25μm, the material parameters being the same as in the preceding subsection. Figure 4*a* shows Re{*I*_{sz}(*r*,*h*/2)} (solid) and $\mathrm{Re}\{\mathrm{i}{H}_{0}^{(1)}({k}_{\mathrm{s}}r)\}$ (dashed) as functions of radial distance *r* which is normalized by λ_{s}=2*π*/*k*_{s}. The phases of these functions are compared in figure 4*b*. Note that the vertical lines in this figure do nothing but show that the sign of the phase is changed by jump. The comparison reveals that the functions give the same result for distances exceeding 0.4λ_{s}. However, for smaller distances, there is a considerable difference in their behaviour. Therefore, the approximation by the Hankel function is not suitable for use in equations (2.88) and (2.89) which describe the bubble oscillation.

Comparison with experimental measurements (presented in the next subsection) shows that agreement with the theoretical model is achieved if the parameter *ε* is considered as a complex quantity whose phase is proportional to the bubble radius *R*_{0}. Therefore, examples of resonance curves were calculated assuming that

$$\epsilon ={\epsilon}_{1}\mathrm{exp}(\mathrm{i}{\epsilon}_{2}{k}_{\mathrm{s}}{R}_{0}),$$

3.1

where *ε*_{1} and *ε*_{2} serve as variation parameters.

Resonance curves shown in figure 5 were calculated by equation (2.88) for *h*=25μm and *f*=30kHz, varying *R*_{0} in the range from 16 to 83μm. The material parameters for the solid and the fluid were the same as in the preceding subsections. The parameters related to the gas behaviour and the acoustic excitation were as follows: *σ*_{f}=0.072Nm^{−1}, *γ*=1.4, *P*_{0}=101.3kPa and *P*_{a}=10kPa. Figure 5*a* shows resonance curves for increasing values of *ε*_{1}, assuming *ε*_{2}=0. It is seen that at *ε*_{1}=0 (no surface waves), the resonance peak is at *R*_{0}*f*=1.1ms^{−1}, which corresponds to the result obtained by Ilinskii *et al*. [15]. When *ε*_{1} is increased, the resonance peak shifts to larger values of *R*_{0} and increases in magnitude. At *ε*_{1}=1, the resonance peak is found to be at *R*_{0}*f*=1.5ms^{−1} and its magnitude is maximal. Note that for *ε*_{1}=1, the effective surface area by which the bubble acts on the wall is equal to the area of the lateral section of the bubble. Further increase of *ε*_{1} decreases the magnitude of the resonance peak, moving its position to larger values of *R*_{0}*f*, which, however, do not exceed 1.7ms^{−1}.

Examples of resonance curves given by the theoretical model for different values of the parameters *ε*_{1} and *ε*_{2}. (*a*) *ε*_{1}=0,0.5,0.8,1,1.5,1.75,2; *ε*_{2}=0. (*b*) *ε*_{1}=1; *ε*_{2}=0,−0.1,−0.2,−0.4,−0.8. **...**

Figure 5*b* shows the evolution of resonance curves with increasing |*ε*_{2}| at *ε*_{1}=1. Increasing |*ε*_{2}| can be related to an increase of damping in the wall material, which leads to flattening of resonance curves and decreasing resonance bubble radius. This behaviour is observed if *ε*_{2} is negative.

The experimental data shown by circles in figures 6–8 were adopted from work by Mekki-Berrada *et al*. [31], where they were acquired by a microfluidic set-up described by Rabaud *et al*. [2]. This set-up has been already mentioned in the Introduction. The experimental data show the normalized amplitude of bubble oscillation for cylindrical bubbles of different radii, placed in a fluid channel, 25μm in depth, confined by elastic walls made of PDMS. The bubble oscillation is induced by a vibrating glass rod moulded in the upper channel wall about 150μm above the solid–fluid interface. The experimental data were obtained at three values of the driving frequency: 30, 40 and 50kHz. For each frequency, the measurements were made for a certain range of bubble radii, which was dictated by experimental conditions.

(*a*,*b*) Comparison of experimental (circles) and theoretical (solid line) results for bubble oscillation at 40kHz.

(*a*,*b*) Comparison of experimental (circles) and theoretical (solid line) results for bubble oscillation at 30kHz.

(*a*,*b*) Comparison of experimental (circles) and theoretical (solid line) results for bubble oscillation at 50kHz.

The solid lines in figures figures66–8 show the fitting of the experimental points by equation (2.88) with the fitting parameter *ε* given by equation (3.1). The material parameters were taken as pointed out in the preceding subsections. The theoretical curves were calculated at the following values of *ε*_{1} and *ε*_{2}: for 30kHz, *ε*_{1}=1.76,*ε*_{2}=−0.451; for 40kHz, *ε*_{1}=2.15,*ε*_{2}=−0.465; for 50kHz, *ε*_{1}=1.95,*ε*_{2}=−0.756. The parts (*a*) of the figures show fitting in the experimental measurement range and the parts (*b*) show the view of the full theoretical curve. The normalization of the data is performed by a maximal value in a depicted range. For this reason, the normalization is different in figure 8*a*,*b*. The experimental data indicate that there is strong damping in the PDMS walls since the resonance curves are very flat. As one can see, agreement between the theoretical and experimental data is achieved if the values of the fitting parameters are changed with frequency. This result suggests that both effective surface area of the bubble–wall contact and the damping characteristics of the PDMS walls depend on frequency by a more complicated way than that specified by equation (3.1). This problem calls for further consideration.

A theoretical model has been developed that describes the volume oscillation of an ultrasound-activated cylindrical bubble situated in a fluid channel between two planar elastic walls. The model includes both the bulk scattered wave, which propagates in the fluid gap, and the surface waves of Lamb-type, which propagate at the solid–fluid interfaces. The force exerted by the bubble oscillation on the channel walls, which is the source of the Lamb-type waves, was treated as a normal harmonic point load. A dispersion equation for the above surface waves was derived, which allows one to examine the dependence of the speed of these waves on the channel height *h*. Such an equation was calculated for the first time, which makes it one of the main results of the present study. It shows that for *h*<λ_{t}, where λ_{t} is the wavelength of the transverse wave in the elastic walls, the speed of the Lamb-type waves decreases with decreasing *h*, while for *h*→λ_{t}, their speed tends to the Scholte wave speed.

The solutions for the wave fields produced by the Lamb-type waves in the elastic walls and in the fluid channel were obtained using the Hankel transforms and expressed in terms of improper integrals. A method was proposed for numerically evaluating the integral solutions. Numerical simulations were carried out to examine the effect of the surface waves on the bubble dynamics. It was shown that the resonance frequency of a cylindrical bubble confined between two identical elastic walls could be up to 50% higher than the resonance frequency of a similar bubble confined between two rigid walls.

Experimental verification of the theoretical model has been carried out using measured values of the bubble oscillation amplitude, which were obtained for bubbles of different radii at three values of the driving frequency. The fitting of experimental data was performed by using two parameters. The first parameter describes the magnitude of the load exerted on the channel walls by the variation of the bubble gas pressure. The second parameter describes the phase shift between the variation of the bubble gas pressure and the reaction of the walls. It was found that agreement between the theoretical and experimental data is achieved on condition that the fitting parameters are considered frequency-dependent quantities.

The denominator functions of the improper integrals (2.63)–(2.67) have zeros, which makes direct numerical integration impossible. This difficulty can be overcome by using contour integration and the residue theorem as proposed by Vrettos [32]. The calculation method will be demonstrated here with integrals (2.65) and (2.67) as these integrals are used in the equation of bubble oscillation.

Expressing the Bessel function *J*_{0} as a sum of the Hankel functions of the first and second kind

$$2{J}_{0}(\xi {k}_{\mathrm{t}}r)={H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)+{H}_{0}^{(2)}(\xi {k}_{\mathrm{t}}r),$$

A 1

substituting (A1) into equation (2.67) and replacing the real variable *ξ* with the complex variable *s*=*ξ*+*iτ*, one obtains

$$2{I}_{\mathrm{p}}(r,z)={I}_{1}+{I}_{2},$$

A 2

where

$${I}_{j}={\int}_{0}^{\mathrm{\infty}}F(s,z){H}_{0}^{(j)}(s{k}_{\mathrm{t}}r)s\mathrm{d}s$$

A 3

and

$$F(s,z)=\frac{{\stackrel{~}{q}}_{\mathrm{l}}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h/2}\mathrm{cosh}({\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}z)}{{\stackrel{~}{q}}_{\mathrm{f}}\stackrel{~}{D}(s)(1-{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h})}.$$

A 4

Note that the variable *ξ* is also replaced with *s* in the radical functions ${\stackrel{~}{q}}_{\mathrm{f}},{\stackrel{~}{q}}_{\mathrm{l}}$ and ${\stackrel{~}{q}}_{\mathrm{t}}$. As a result, these functions become complex valued and branch points occur at *s*=*ξ*_{f}, *ξ*_{l},1, where *ξ*_{f}<*ξ*_{l}<1 because it is assumed that *c*_{f}>*c*_{l}>*c*_{t}; see equations (2.68). The pole of the integrand *F*(*s*,*z*) is given by the root of the equation $\stackrel{~}{D}(\xi )=0$. This root, which will be denoted by *ξ*_{s}, corresponds to a Scholte-type wave and *ξ*_{s}>1 as the speed of the Scholte wave is smaller than *c*_{t}; see the first of equations (2.68). To make the integrand single valued and take into account the presence of the pole, the integration contours *Γ*_{1} and *Γ*_{2} are taken as shown in figure 9.

Integrating *I*_{1} along the contour *Γ*_{1} and using the fact that ${H}_{0}^{(1)}(s{k}_{\mathrm{t}}r)$ vanishes along the infinite arc, one obtains

$${I}_{1}={I}_{11}+{I}_{12}+{I}_{13}+{I}_{14}+2\pi {\mathrm{iRes}}_{\mathrm{p}}({\xi}_{\mathrm{s}}),$$

A 5

where

$$\begin{array}{rl}{I}_{11}& =-{\int}_{0}^{\mathrm{\infty}}{F}^{+}(\mathrm{i}\tau ,z){H}_{0}^{(1)}(\mathrm{i}\tau {k}_{\mathrm{t}}r)\tau \mathrm{d}\tau ,\end{array}$$

A 6

$$\begin{array}{rl}{I}_{12}& =-{\int}_{0}^{{\xi}_{\mathrm{f}}}[{F}^{+}(\xi ,z)-{F}^{-}(\xi ,z)]{H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi ,\end{array}$$

A 7

$$\begin{array}{rl}{I}_{13}& =-{\int}_{{\xi}_{\mathrm{f}}}^{{\xi}_{\mathrm{l}}}[{F}^{+}(\xi ,z)-{F}^{-}(\xi ,z)]{H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi \end{array}$$

A 8

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{I}_{14}& =-{\int}_{{\xi}_{\mathrm{l}}}^{1}[{F}^{+}(\xi ,z)-{F}^{-}(\xi ,z)]{H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi ,\end{array}$$

A 9

the superscripts + and − denote that the function *F*(*s*,*z*) is calculated for the positive and negative values of the imaginary part of the radicals ${\stackrel{~}{q}}_{\mathrm{f}},{\stackrel{~}{q}}_{\mathrm{l}}$ and ${\stackrel{~}{q}}_{\mathrm{t}}$, respectively, and Res_{p}(*ξ*_{s}) is the residue at the pole *ξ*_{s} which is calculated by [33]

$${\mathrm{Res}}_{\mathrm{p}}({\xi}_{\mathrm{s}})={\left[\frac{{\stackrel{~}{q}}_{\mathrm{l}}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h/2}\mathrm{cosh}({\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}z)\xi {H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)}{\mathrm{\partial}[{\stackrel{~}{q}}_{\mathrm{f}}\stackrel{~}{D}(\xi )(1-{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h})]/\mathrm{\partial}\xi}\right]}_{\xi ={\xi}_{\mathrm{s}}}.$$

A 10

Integrating *I*_{2} along the contour *Γ*_{2} and using the fact that ${H}_{0}^{(2)}(s{k}_{\mathrm{t}}r)$ vanishes along the infinite arc, one obtains

$${I}_{2}=-{\int}_{0}^{\mathrm{\infty}}{F}^{-}(-\mathrm{i}\tau ,z){H}_{0}^{(2)}(-\mathrm{i}\tau {k}_{\mathrm{t}}r)\tau \mathrm{d}\tau .$$

A 11

Equations (A 7)–(A 9) can be simplified using the fact that *F*^{−}(*ξ*,*z*)=[*F*^{+}(*ξ*,*z*)]*, where the asterisk denotes the complex conjugate. Further, the sum of the integrals *I*_{11} and *I*_{2} can be transformed using the following identities:

$$\begin{array}{rl}{F}^{-}(-\mathrm{i}\tau ,z)& ={[{F}^{+}(\mathrm{i}\tau ,z)]}^{\ast},\end{array}$$

A 12

$$\begin{array}{rl}{H}_{0}^{(1)}(\mathrm{i}x)& =\frac{2}{\pi \mathrm{i}}{K}_{0}(x)\end{array}$$

A 13

$$\begin{array}{rl}\phantom{\rule{2em}{0ex}}\text{and}\phantom{\rule{1em}{0ex}}\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{H}_{0}^{(2)}(-\mathrm{i}x)& =-\frac{2}{\pi \mathrm{i}}{K}_{0}(x),\end{array}$$

A 14

where *K*_{n}(*x*) is the modified Bessel function of the second kind of order *n* [29]. As a result of these operations, the final expression for the integral *I*_{p}(*r*,*z*) is found to be

$$\begin{array}{rl}{I}_{\mathrm{p}}(r,z)=& -\frac{2}{\pi}{\int}_{0}^{\mathrm{\infty}}\mathrm{Im}\{{F}^{+}(\mathrm{i}\tau ,z)\}{K}_{0}(\tau {k}_{\mathrm{t}}r)\tau \mathrm{d}\tau -\mathrm{i}{\int}_{0}^{{\xi}_{\mathrm{f}}}\mathrm{Im}\{{F}^{+}(\xi ,z)\}{H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi \\ & -\mathrm{i}{\int}_{{\xi}_{\mathrm{f}}}^{{\xi}_{\mathrm{l}}}\mathrm{Im}\{{F}^{+}(\xi ,z)\}{H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi -\mathrm{i}{\int}_{{\xi}_{\mathrm{l}}}^{1}\mathrm{Im}\{{F}^{+}(\xi ,z)\}{H}_{0}^{(1)}(\xi {k}_{\mathrm{t}}r)\xi \mathrm{d}\xi +\pi {\mathrm{iRes}}_{\mathrm{p}}({\xi}_{\mathrm{s}}),\end{array}$$

A 15

where Im{*x*} means the imaginary part of *x*. These integrals are nonsingular and convergent and hence they can be calculated numerically without any difficulty.

Application of the same method to integral (2.65) yields

$$\begin{array}{rl}{I}_{\mathrm{fr}}(r,z)& =-\frac{2}{\pi}{\int}_{0}^{\mathrm{\infty}}\mathrm{Im}\{{F}^{+}(\mathrm{i}\tau ,z)\}{K}_{1}(\tau {k}_{\mathrm{t}}r){\tau}^{2}\mathrm{d}\tau -\mathrm{i}{\int}_{0}^{{\xi}_{\mathrm{f}}}\mathrm{Im}\{{F}^{+}(\xi ,z)\}{H}_{1}^{(1)}(\xi {k}_{\mathrm{t}}r){\xi}^{2}\mathrm{d}\xi \\ & \phantom{\rule{1em}{0ex}}-\mathrm{i}{\int}_{{\xi}_{\mathrm{f}}}^{{\xi}_{\mathrm{l}}}\mathrm{Im}\{{F}^{+}(\xi ,z)\}{H}_{1}^{(1)}(\xi {k}_{\mathrm{t}}r){\xi}^{2}\mathrm{d}\xi -\mathrm{i}{\int}_{{\xi}_{\mathrm{l}}}^{1}\mathrm{Im}\{{F}^{+}(\xi ,z)\}{H}_{1}^{(1)}(\xi {k}_{\mathrm{t}}r){\xi}^{2}\mathrm{d}\xi +\pi {\mathrm{iRes}}_{\mathrm{fr}}({\xi}_{\mathrm{s}}),\end{array}$$

A 16

where the residue at the pole *ξ*_{s} is calculated by

$${\mathrm{Res}}_{\mathrm{fr}}({\xi}_{\mathrm{s}})={\left[\frac{{\stackrel{~}{q}}_{\mathrm{l}}{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h/2}\mathrm{cosh}({\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}z){\xi}^{2}{H}_{1}^{(1)}(\xi {k}_{\mathrm{t}}r)}{\mathrm{\partial}[{\stackrel{~}{q}}_{\mathrm{f}}\stackrel{~}{D}(\xi )(1-{\mathrm{e}}^{-{\stackrel{~}{q}}_{\mathrm{f}}{k}_{\mathrm{t}}h})]/\mathrm{\partial}\xi}\right]}_{\xi ={\xi}_{\mathrm{s}}}.$$

A 17

A.D. carried out the development of the theory and drafted the manuscript; F.M.B. carried out the experimental measurements and helped draft the manuscript; P.T. participated in the experimental measurements and the interpretation of the results and helped draft the manuscript; P.M. participated in the development of the theory and the experimental measurements, coordinated the study and helped draft the manuscript. All authors gave final approval for publication.

We have no competing interests.

This research has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007–2013)/ERC grant agreement no. 614655 ‘Bubbleboost’.

1. Friend J, Yeo L
2011.
Microscale acoustofluidics: microfluidics driven via acoustics and ultrasonics. *Rev. Mod. Phys.*
83, 647–704. (doi:10.1103/RevModPhys.83.647)

2. Rabaud D, Thibault P, Mathieu M, Marmottant P
2011.
Acoustically bound microfluidic bubble crystals. *Phys. Rev. Lett.*
106, 134501 (doi:10.1103/PhysRevLett.106.134501) [PubMed]

3. Doinikov AA.
2005.
Bjerknes forces and translational bubble dynamics. In *Bubble and particle dynamics in acoustic fields: modern trends and applications* (ed. Doinikov AA, editor. ), pp. 95–143. Kerala, India: Research Signpost.

4. Marmottant P, Raven JP, Gardeniers H, Bomer JG, Hilgenfeldt S
2006.
Microfluidics with ultrasound-driven bubbles. *J. Fluid Mech.*
568, 109–118. (doi:10.1017/S0022112006002746)

5. Ahmed D, Mao X, Shi J, Juluri BK, Huang TJ
2009.
A millisecond micromixer via single-bubble-based acoustic streaming. *Lab Chip*
9, 2738–2741. (doi:10.1039/b903687c) [PubMed]

6. Wang C, Jalikop SV, Hilgenfeldt S
2011.
Size-sensitive sorting of microparticles through control of flow geometry. *Appl. Phys. Lett.*
99, 034101 (doi:10.1063/1.3610940)

7. Franke T, Braunmüller S, Schmid L, Wixforth A, Weitz D
2010.
Surface acoustic wave actuated cell sorting (SAWACS). *Lab Chip*
10, 789–794. (doi:10.1039/b915522h) [PubMed]

8. Zwaan E, Le Gac S, Tsuji K, Ohl C-D
2007.
Controlled cavitation in microfluidic systems. *Phys. Rev. Lett.*
98, 254501 (doi:10.1103/PhysRevLett.98.254501) [PubMed]

9. Quinto-Su PA, Lim KY, Ohl C-D
2009.
Cavitation bubble dynamics in microfluidic gaps of variable height. *Phys. Rev. E*
80, 047301 (doi:10.1103/PhysRevE.80.047301) [PubMed]

10. Quinto-Su PA, Ohl C-D
2009.
Interaction between two laser-induced cavitation bubbles in a quasi-two-dimensional geometry. *J. Fluid Mech.*
633, 425–435. (doi:10.1017/S0022112009008064)

11. Sankin GN, Yuan F, Zhong P
2010.
Pulsating tandem microbubble for localized and directional single-cell membrane poration. *Phys. Rev. Lett.*
105, 078101 (doi:10.1103/PhysRevLett.105.078101) [PMC free article] [PubMed]

12. Gonzalez-Avila SR, Klaseboer E, Khoo BC, Ohl C-D
2011.
Cavitation bubble dynamics in a liquid gap of variable height. *J. Fluid Mech.*
682, 241–260. (doi:10.1017/jfm.2011.212)

13. Yuan F, Sankin G, Zhong P
2011.
Dynamics of tandem bubble interaction in a microfluidic channel. *J. Acoust. Soc. Am.*
130, 3339–3346. (doi:10.1121/1.3626134) [PubMed]

14. Prosperetti A.
2004.
Bubbles. *Phys. Fluids*
16, 1852–1865. (doi:10.1063/1.1695308)

15. Ilinskii YA, Zabolotskaya EA, Hay TA, Hamilton MF
2012.
Models of cylindrical bubble pulsation. *J. Acoust. Soc. Am.*
132, 1346–1357. (doi:10.1121/1.4730888) [PubMed]

16. Wang C, Rallabandi B, Hilgenfeldt S
2013.
Frequency dependence and frequency control of microbubble streaming flows. *Phys. Fluids*
25, 022002 (doi:10.1063/1.4790803)

17. Viktorov IA.
1967.
*Rayleigh and Lamb waves: physical theory and applications*. New York, NY: Plenum Press.

18. Rayleigh L.
1885.
On waves propagated along the plane surface of an elastic solid. *Proc. Lond. Math. Soc.*
s1–17, 4–11. (doi:10.1112/plms/s1-17.1.4)

19. Scholte JG.
1947.
The range and existence of Rayleigh and Stoneley waves. *Geophys. J. Int.*
5, 120–126. (doi:10.1111/j.1365-246X.1947.tb00347.x)

20. Ewing WM, Jardetzky WS, Press F
1957.
*Elastic waves in layered media*. New York, NY: McGraw-Hill.

21. Lamb H.
1917.
On waves in an elastic plate. *Proc. R. Soc. Lond. A*
93, 114–128. (doi:10.1098/rspa.1917.0008)

22. Leighton TG, White PR, Schneider MF
1998.
The detection and dimension of bubble entrainment and comminution. *J. Acoust. Soc. Am.*
103, 1825–1835. (doi:10.1121/1.421374)

23. Leighton TG, White PR, Morfey CL, Clarke JWL, Heald GJ, Dumbrell HA, Holland KR
2002.
The effect of reverberation on the damping of bubbles. *J. Acoust. Soc. Am.*
112, 1366–1376. (doi:10.1121/1.1501895) [PubMed]

24. Leighton TG.
2011.
The inertial terms in equations of motion for bubbles in tubular vessels or between plates. *J. Acoust. Soc. Am.*
130, 3333–3338. (doi:10.1121/1.3638132) [PubMed]

25. Achenbach JD.
1973.
*Wave propagation in elastic solids*. Amsterdam, The Netherlands: North-Holland Publishing Company.

26. Piessens R.
2000.
The Hankel transform. In *The transforms and applications handbook* (ed. Poularikas AD, editor. ), pp. 9.1–9.30. Boca Raton, FL: CRC Press LLC.

27. Landau LD, Lifshitz EM
1987.
*Fluid mechanics*. Oxford, UK: Pergamon Press.

28. Zhu J, Popovics JS, Schubert F
2004.
Leaky Rayleigh and Scholte waves at the fluid-solid interface subjected to transient point loading. *J. Acoust. Soc. Am.*
116, 2101–2110. (doi:10.1121/1.1791718)

29. Abramowitz M, Stegun IN
1965.
*Handbook of mathematical functions*. New York, NY: Dover Publications.

30. Meyers M, Chawla K
2009.
*Mechanical behavior of materials*. Cambridge, UK: Cambridge University Press.

31. Mekki-Berrada F, Thibault P, Marmottant P
2016.
Acoustic pulsation of a microbubble confined between elastic walls. *Phys. Fluids*
28, 032004 (doi:10.1063/1.4942917)

32. Vrettos C.
2008.
Green’s functions for vertical point load on an elastic half-space with depth-degrading stiffness. *Eng. Anal. Bound. Elem.*
32, 1037–1045. (doi:10.1016/j.enganabound.2007.10.017)

33. Churchill RV.
1973.
*Operational mathematics*. New York, NY: McGraw-Hill.

Articles from Proceedings. Mathematical, Physical, and Engineering Sciences are provided here courtesy of **The Royal Society**

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |