Standard and inverse transducer-plane streaming patterns in resonant acoustofluidic devices: Experiments and simulations

Standard and inverse transducer-plane streaming patterns in resonant acoustofluidic devices: Experiments and simulations

Applied Mathematical Modelling 77 (2020) 456–468 Contents lists available at ScienceDirect Applied Mathematical Modelling journal homepage: www.else...

2MB Sizes 0 Downloads 1 Views

Applied Mathematical Modelling 77 (2020) 456–468

Contents lists available at ScienceDirect

Applied Mathematical Modelling journal homepage: www.elsevier.com/locate/apm

Standard and inverse transducer-plane streaming patterns in resonant acoustofluidic devices: Experiments and simulations Junjun Lei∗, Feng Cheng, Zhongning Guo School of Electromechanical Engineering, Guangdong University of Technology, Guangzhou 510006, China

a r t i c l e

i n f o

Article history: Received 3 January 2019 Revised 17 July 2019 Accepted 23 July 2019 Available online 30 July 2019 Keywords: Acoustic streaming Transducer-plane streaming patterns Acoustic particle manipulation Standing wave Travelling wave Acoustic field

a b s t r a c t In this paper, we report a new transducer-plane streaming pattern, which we call here inverse four-quadrant transducer-plane streaming, in an ultrasound-excited thin-layer glass capillary device, illustrate the importance of the travelling wave component and emphasize its phase difference with the standing wave component of three-dimensional acoustic fields excited in thin-layer resonant acoustofluidic particle manipulation devices on the formation of transducer-plane streaming patterns through numerical simulations. A mathematical model was created for solving the three-dimensional acoustic streaming fields from limiting velocities derived from predefined acoustic fields, capturing both the standing wave and travelling wave components and their phase relations. This model was then successfully applied to interpret the unusually reversed in-plane streaming patterns seen in a silicon-based fluid channel. It was found that, by tuning the phase difference between the standing and travelling wave components, ϕ , the transducer-plane streaming vortices could be shifted along the fluid channel and reversed streaming patterns can be excited when ϕ differ by π , which could provide insights for the control of transducerplane streaming in thin-layer acoustofluidic devices and may have promising potential for acoustofluidic applications such as nanoscale particle manipulation. © 2019 Elsevier Inc. All rights reserved.

1. Introduction Acoustic streaming is a steady fluid motion, which is driven by the absorption of acoustic momentum flux in a fluid. A premise behind the control of acoustic streaming for enhancing or suppressing the acoustic streaming effects for applications such as particle manipulation [1,2] is to understand their driving mechanisms within acoustofluidic devices. In most bulk acoustic wave (BAW) microfluidic particle manipulation systems of interest, where the fluid channels usually have dimensions comparable to the acoustic wavelength, the acoustic streaming fields formed in the microfluidic channels are generally dominated by boundary-driven streaming [3], which arises from the acoustic absorption in the acoustic boundary layer region due to the non-slip condition on the fluid channel walls. Another class of acoustic streaming, the Eckart type streaming [4], which is usually seen in surface acoustic wave (SAW) systems [5], generally requires acoustic absorption over longer distances than those typically found in BAW resonators [6]. Classical boundary-driven streaming patterns are usually referred to as Schlichting streaming [7] and Rayleigh streaming [8], which describe streaming vortices inside and outside of the acoustic boundary layer, respectively. The inner and outer streaming fields can be modelled by considering the Reynolds stresses [9], which are the volume forces representing the ∗

Corresponding author. E-mail address: [email protected] (J. Lei).

https://doi.org/10.1016/j.apm.2019.07.044 0307-904X/© 2019 Elsevier Inc. All rights reserved.

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

457

Nomenclature a F h Ii Ji kis kit l p pi p0s p0t p1s p1t Q u ui u2 uM 2 uL u1 u∗1 vL v1

v∗1

Vpp w w1 w∗1

side length of cubic mesh elements [m] Reynolds stress force [N/m3 ] model height [m] active intensity in i direction [W/m2 ] reactive intensity in i direction [W/m2 ] standing wave number in i direction [rad/m] travelling wave number in i direction [rad/m] model length [m] fluid pressure [Pa] ith -order fluid pressure [Pa] standing wave pressure amplitude [Pa] travelling wave pressure amplitude [Pa] standing wave pressure [Pa] travelling wave pressure [Pa] quality factor [1] fluid velocity [m/s] ith -order fluid velocity [m/s] acoustic streaming velocity [m/s] mass transport velocity of acoustic streaming [m/s] x-component limiting velocity [m/s] x-component acoustic velocity [m/s] conjugate of x-component acoustic velocity [m/s] y-component limiting velocity [m/s] y-component acoustic velocity [m/s] conjugate of y-component acoustic velocity [m/s] peak-to-peak voltage [V] model width [m] z-component acoustic velocity [m/s] conjugate of z-component acoustic velocity [m/s]

Greek letters ρ fluid density [kg/m3 ] ρi ith -order fluid density [kg/m3 ] μ dynamic viscosity [Pa s] μb bulk viscosity [Pa s] ϕ phase difference between the standing and travelling wave components [rad] ω angular frequency [rad/s] Abbreviations 3D three-dimensional ASV average streaming velocity BAW bulk acoustic wave LHS left-hand-side LV limiting velocity LVM limiting velocity method MSV maximum streaming velocity PIV particle image velocimetry RHS right-hand-side RSF Reynolds stress force SAW surface acoustic wave SW standing wave TPS transducer-plane streaming TW travelling wave

non-zero time-averaged acoustic momentum flux resulting from the acoustic dissipation near the no-slip walls [6]. Alternatively, the outer streaming fields can be effectively modelled using the limiting velocity method (LVM) [10,11], which ignores the acoustic boundary layer and replaces its non-slip boundary condition with a slip velocity boundary condition on the walls where limiting velocities (LVs) derived from the acoustic field are applied to.

458

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

While Rayleigh streaming patterns, whose orientations are generally perpendicular to the driving boundaries, have been extensively studied within the field of ultrasonic particle manipulation in last decades [12–21], there are acoustic streaming patterns observed experimentally in the bulk of the fluid that cannot be explained by the classical boundary-driven streaming theory [22–29]. Particularly, acoustic streaming patterns whose orientations are parallel to the driving boundaries and the transducer radiating surface have been observed in thin-layer resonant acoustofluidic particle manipulation devices, which are adverse to the conditions found in Rayleigh streaming and thus are referred to as transducer-plane streaming (TPS) [25]. Recently, we have modelled and explored the mechanisms behind a four-quadrant TPS pattern [23,25] and an eight-octant TPS pattern [26], of which each streaming vortex respectively occupies approximately a quadrant and an octant of the visualisation plane and is parallel to the driving boundaries (which is parallel to the transducer radiating surface). It was found that the LV fields, which drive the TPS patterns, are closely related to the active sound intensity field (the real part of the complex sound intensity) on the driving boundaries (which tends to circulate about pressure nodal points [30]). The analyses of resonant BAW microfluidic particle manipulation devices of interest mainly focus on the standing wave (SW) modes, where the acoustic radiation force drives particles exposed to the fluid to the acoustic pressure nodes or antinodes [31], while less attention has been paid to the travelling wave (TW) component established in the fluid channel of such systems. We have recently found that a TW component is required, although typically at a smaller magnitude compared to the SW component, for the formation of TPS patterns in thin-layer resonant acoustofluidic manipulation devices [26]. However, it is still unclear how in reality the TW component in resonant acoustofluidic channel affects the resulting acoustic streaming patterns as only a very limited variety of TPS patterns that have been observed in acoustofluidic manipulation devices. In this work, motived by a new TPS pattern observed for the first time in our thin-layer glass capillary device, which is typically an inverse pattern to that we have discussed before [23,25], we investigate further the importance of the TW component of three-dimensional (3D) acoustic wave fields established in the fluid channel of thin-layer resonant acoustofluidic devices on the formation of TPS patterns. Importantly, the influence of the phase difference between the SW and TW components established in the fluid channel on the resulting TPS patterns are elucidated in detail. We first introduce the new TPS pattern observed in a thin-layer glass capillary device, which has not been discussed or shown experimentally before, then investigate the underlying physics of four-quadrant TPS patterns with a numerical model, which is further applied to analyse unusual in-plane streaming patterns observed in other acoustofluidic manipulation devices, followed by a conclusion and discussion of potential applications of this work. 2. Governing equations for acoustic streaming In the following, bold and normal-emphasis fonts are used to represent vector and scalar quantities, respectively. The derivation of basic acoustic streaming equations has been extensively presented in the literature. Here, we assume a homogeneous isotropic fluid, where the continuity and momentum equations for the fluid motion are [32]

∂ρ + ∇ · ( ρ u ) = 0, ∂t     1 ∂u ρ + u · ∇ u = −∇ p + μ∇ 2 u + μb + μ ∇ ∇ · u, ∂t 3

(1)

(2)

where ρ is the fluid density, t is time, u is the fluid velocity, p is the pressure, and μ and μb are respectively the dynamic and bulk viscosity coefficients of the fluid. The left-hand-side (LHS) of Eq. (2) represents the inertia force per unit volume, where the two terms in the bracket are respectively the unsteady acceleration and convective acceleration of a fluid particle and the right-hand-side (RHS) indicates the divergence of stresses. The two main numerical methods for modelling boundary-driven streaming in acoustofluidic particle manipulation devices have recently been compared [33]. Both methods are based on the perturbation theory [34–36] and assume that the steady acoustic streaming velocity field is superposed on the harmonic acoustic velocity field. Following this theory, the fluid density ρ , pressure p, and velocity u are expressed as:

ρ = ρ0 + ρ1 + ρ2 + · · · ,

(3)

p = p0 + p1 + p2 + · · · ,

(4)

u = u1 + u2 + · · · ,

(5)

where the subscripts 0, 1 and 2 represent the static (absence of sound), first-order and second-order quantities, respectively. Substituting Eqs. (3)–(5) into Eqs. (1) and (2) and taking the first-order into account, Eqs. (1) and (2) take the form,

∂ ρ1 + ρ0 ∇ · u 1 = 0 , ∂t

(6)

ρ0

(7)

  ∂ u1 1 = −∇ p1 + μ∇ 2 u1 + μb + μ ∇∇ · u1 . ∂t 3

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

459

Fig. 1. (Colour online) Characterization of a new transducer-plane streaming pattern in a thin-layer acoustofluidic manipulation device: (a) sketch of the experimental device; (b) yz cross-section, materials and dimensions of the device; (c) visualization plane; and (d) PIV results of the acoustic streaming velocity field measured at plane z = 0 at the resonant frequency of 0.8018 MHz and voltage of 28.6 Vpp , where the reference arrow at the right-bottom corner shows a velocity of 30 μm/s.

Repeating the above procedure, and taking the first- and second-order into account and the time average of Eqs. (1) and (2), the continuity and momentum equations are turned into

∇ · ρ 1 u 1 + ρ0 ∇ · u 2 = 0 , −F = −∇ p2 + μ∇ 2 u2 +

(8)



 1 μb + μ ∇∇ · u2 , 3

(9)

where the upper bar ¯· indicates the time average of the quantity below and F = −ρ0 u1 ∇ · u1 + u1 · ∇ u1 is the so-called Reynolds stress force (RSF) [6]. The divergence free velocity uM = u2 + ρ1 u1 /ρ0 , obtained from the LHS of Eq. (8), is the 2 average mass transport velocity, which is generally closer to the velocity of tracer particles experimentally measured in an acoustic streaming flow than u2 . [37] 3. Experimental setup and results The experiments were conducted in a simple transducer-capillary device, sketched in Fig. 1(a), which is of interest as it has been widely used for acousto-microfluidic applications, such as blood and bacterial capture [1,2,38], sonoporation [39,40], cell deformation [41], and cartilage tissue engineering [42]. Fig. 1(b) depicts the cross-section of our experimental device, in which properties and dimensions of device layers are presented. The glass capillary (Vitrocom) has inner dimensions of 1 mm × 15 mm, wall thickness of 1 mm, and length of 100 mm. The PZT transducer (Ferroperm) has dimensions of 5 mm × 5 mm × 1 mm (x × y × z), which was glued to approximately the centre of the bottom surface of glass capillary. Fluidic connections to the capillary were made via PTFE tubing (inner diameter of 1 mm) attached via heat-shrink sleeving. A function generator (Tektronix, AFG1062) drives a RF amplifier (JUNTEK, DPA2698) that powers the transducer with signal monitored by an oscilloscope (Tektronix, TBS1102B). An Olympus fluorescent microscope with a dual-frame CCD camera was used to image the measurements. The half-wave resonant frequency (corresponding to a half wavelength SW in the z direction of the fluid channel) was identified from impedance measurements, which was found at 0.8016 MHz. It was observed that 10 μm polystyrene particles (Fluoresbrite microspheres, Polysciences Inc.) were firstly levitated at the half-height of the fluid channel (z = 0) and then

460

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

Fig. 2. The simplified reduced-fluid model considered: (a) the 3D model regime (unit: mm), where the dash-dot lines show the symmetric planes; and (b) the mesh distribution used in simulations, where a is the side length of mesh elements.

slowly lined up. This relates to the predominant acoustic resonance in the z direction and relative weak acoustic resonance in the y-direction of the fluid channel (see modelling below). Meanwhile, steady, in-plane circulations due to the acoustic streaming effects were observed if 1 μm polystyrene particles (Fluoresbrite microspheres, Polysciences Inc.) were introduced. The acoustic streaming measurements were performed within xy horizontal planes and the investigation areas were right above the transducer radiating surface (see Fig. 1(c)). The observed in-plane streaming pattern at the resonant frequency was shown in Fig. 1(d), where particle-image-velocimetry (PIV) results of the streaming field is presented. It can be seen that a four-quadrant acoustic streaming pattern, of which each vortex occupies approximately one quadrant of the viewed xy horizontal plane, is formed. The plane of these four-quadrant streaming vortices is parallel to the transducer radiating surface (i.e. perpendicular to the axis of the main SW in the z direction), as is the case for the four-quadrant TPS we have demonstrated previously in another resonant acoustofluidic manipulation device [25] but with a different orientation on each streaming vortex. We thus call it here inverse four-quadrant TPS.

4. Numerical modelling To understand the genesis of this new TPS pattern and to explore what makes it different to the four-quadrant streaming pattern discussed previously [25], in this section, we present numerical simulations of acoustic streaming fields in the thinlayer resonant acoustofluidic device shown in Fig. 1. The finite element software COMSOL Multiphysics 5.2 [43] was used to solve the governing equations and all the numerical simulations were performed on a HP EliteBook 820G4 running Windows 10 (64-bit), equipped with 8 GB RAM and Intel (R) Core (TM) i7-7500 U processor of clock frequency 2.7 GHz.

4.1. Numerical model Fig. 4(a) shows the 3D model, where only the third-quadrant of the fluid channel was considered as the fluid channel and the acoustic streaming pattern are symmetric to planes x = 0 and y = 0 (see Fig. 1). Our previous work have shown that a reduced-fluid model can be used to effectively predict 3D acoustic streaming patterns in layered acoustofluidic devices [23]. For the computational mesh, as shown in Fig. 2(b), a uniform distribution of cubic mesh elements was used in simulations. The choice of mesh size (side length of a, shown in Fig. 2(b)) was based on the following factors: (1) the outer acoustic streaming fields are modelled from the LVM, and thus no boundary-layer mesh is required; (2) in previous modelling work [25], where the acoustic pressure field was solved from the linearized acoustic wave equations and boundary vibrations, a mesh dependency study was applied to optimize the computation mesh, i.e. to ensure that a further refining of mesh size does not significantly change the numerical solution, and a pretty fine mesh with five elements in the z-direction of the fluid channel was used, which resulted in a mesh-induced numerical error of about 0.7%. The case for the model in this work, however, is quite different as the acoustic pressure fields in the fluid channel were defined by mathematical equations (described in the following subsections), rather than determined by boundary conditions; and (3) the same goes with the acoustic streaming velocities. As plotted in Fig. 3, the size of mesh elements was shown to have little effect on the modelled acoustic streaming velocities. It can be seen that as few as one element in the z-direction of the fluid channel is sufficient to accurately capture the streaming velocity magnitudes (with a mesh-induced numerical error of less than 1%). However, a higher mesh density (three elements in the z-direction, shown in Fig. 2(b)) was used for all the numerical simulations shown below, unless otherwise stated, in order to allow for situations more complex than the simple case explored in the mesh dependency study. For the mesh constitution described in Fig. 2(b), this model resulted in 1425 mesh elements, a peak RAM usage of 1.78 GB, and a running time of about 6 min for solving the sequential steps (1)–(3) described in the following section.

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

461

Fig. 3. (Colour online) Relationship between the normalised acoustic streaming velocities and mesh sizes, where h is the model height, a is the side length of mesh elements, and ASV and MSV are abbreviations of average streaming velocity and maximum streaming velocity magnitudes in the model regime, respectively.

4.2. Numerical procedure The whole numerical procedure used to model the standard and inverse TPS patterns was split into the following sequential steps: 4.2.1. Definition of acoustic fields In this work, we defined the acoustic pressure field from combinations of cavity modes and TW modes. Compared to the linear acoustic model we built in previous work, where a normal distribution of boundary vibration was applied to approximate the 3D first-order acoustic fields in fluid channels [23,25], with this model we are able to obtain more insights into the relative significance of the TW and SW components [26] and the effect of their phase differences on the resulting TPS patterns, which are illustrated in detail below. In our model, based on the structure of our experimental device, we decomposed the 3D first-order acoustic pressure, p1 , established in the fluid channel, into two components, a SW component, p1s , and a TW component propagating along the fluid channel (in the x-direction), p1t ,

p1 = p1s + p1t ,

(10)

p1s = p0s cos (kxs x ) cos (kys y ) sin (kzs z )eiωt ,

(11)

p1t = p0t eikxt x cos (kyt y ) sin (kzt z )ei(ωt+ϕ ) ,

(12)

where the second s and t subscripts indicate the SW and TW components respectively, p0 is the acoustic pressure amplitude, ω is the angular frequency, kx , ky and kz are the wave numbers in the x, y and z directions and ϕ represents the phase difference between the SW and TW components. The reason for defining the TW component only in the x-direction of the fluid channel is that the acoustofluidic channel is made of glass (acoustically hard material), in which acoustic energy travelling down the glass capillary is largely absorbed by the connectors and tubing at the ends.

462

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

Two points need to be clarified here. Firstly, unlike the acoustic step in the reduced-fluid model shown in our previous work [25], where reflecting and radiation boundary conditions were applied to approximate the total acoustic pressure fields, in this model, as the total acoustic pressure field and its SW and TW components are characterised with Eqs. (10)– (12), we do not need an acoustic model (thus no need to apply boundary conditions), but only need to solve Eqs. (10)–(12) in a frequency domain. Then, although we have previously used Eqs. (10)–(12) to demonstrate a TW component is required to generate TPS vortices in thin-layer acoustofluidic devices [26], the effects of phase difference between the SW and TW components of an acoustic field on the variations of TPS patterns, however, have not yet been studied or discussed; and we will show below that, with the help of the inverse four-quadrant TPS patterns observed in our experimental device shown in Section 3, Eqs. (10)–(12) are extremely useful for demonstrating the driving mechanisms of standard and inverse TPS patterns in thin-layer BAW resonators. 4.2.2. Derivation of LV fields In this work, we have applied the LVM [10] based on the perturbation method [35] to model the 3D acoustic streaming fields in the bulk of the fluid channel of the capillary. The LVM was initially proposed by Nyborg [10], who showed that the limit value of the streaming velocity (i.e. the LV just outside the acoustic boundary layer) can be derived from the linear acoustic velocity field and can be effectively applied to solve the outer streaming field, and later modified by Lee and Wang [11] to a more generalised version and by Bach and Bruus [44] for curved boundaries. This method is more computationally efficient than the conventional numerical method by considering the RSF [33,45] as it is not necessary to use extremely fine mesh elements to resolve the acoustic and streaming fields in the thin acoustic boundary layer. We have recently established the viability of this method for solving 3D boundary-driven acoustic streaming fields in thin-layer resonant acoustofluidic manipulation devices [25], vibrating thin-membrane systems [46], and ultrasound-excited phononic crystal systems [47]. The LVs at all driving boundaries were calculated from the first-order acoustic velocity fields. In this model, the main acoustic streaming fields in the bulk of the fluid are driven by LVs on the top and bottom boundaries (z = ±h/2), which take the form: [48]





uL = −

1 d w1 Re qx + u∗1 (2 + i )∇ · u1 − (2 + 3i ) 4ω dz

vL = −

1 d w1 Re qy + v∗1 (2 + i )∇ · u1 − (2 + 3i ) 4ω dz







,

(13)

,

(14)



qx = u1

du∗1 du∗ + v1 1 , dx dy

(15)

qy = u1

dv∗1 d v∗ + v1 1 , dx dy

(16)

where uL and vL are the x- and y-components of LVs, Re{·} represents the real part of the quantity inside, and u1 , v1 and w1 are the x-, y- and z-components of the complex first-order acoustic velocity vector, u1 . The superscript, ∗ , represents the complex conjugate. Based on the acoustic pressure field defined in Eqs. (10)-(12), the acoustic velocity components are calculated following the momentum conservation equation, ∂ u1 /∂ t + ∇ p1 /ρ 0 = 0. In a recent work, Bach and Bruus [44] extended Lee and Wang’s theory [11] and derived more general LV equations including the effect of Stokes drift, which results from the no-slip boundary condition at the wall [49]. An additional term, which was missing in Lee and Wang’s theory, was added to the LV equations. However, this missing term has little effect on the resulting limiting velocities in our resonant system as the quality factor of the system Q is much larger than 1, that is Q  1, which means that the wall velocity is negligible compared to the bulk acoustic velocity |u1 | and thus the missing term is negligible in the LV equations [44]. Therefore, Eqs. (13)–(16) are valid for approximating the outer streaming fields in our transducer-capillary system. 4.2.3. Modelling of acoustic streaming fields The acoustic streaming fields in this thin-layer resonant acoustofluidic device were simulated from a ‘creeping flow’ model, where the LV fields solved in the previous step were applied to as LV boundary conditions. As the acoustic boundary layer is not considered in this work, with the assumption of low velocity, incompressible flow, the term ∇ · ρ1 u1 in the LHS of Eq. (8) is zero [45], and thus the last term in the RHS of Eq. (9), ∇∇ · u2 , is also zero. Then, as discussed by Lighthill [6], the RSF outside the boundary layer can set up hydrostatic stresses, which, however, will not create vortices in the absence of attenuation; hence, the RSF in Eq. (9) is also ignored. Therefore, in this work, the governing Eqs. (8) and (9) for solving the second-order outer streaming velocities, u2 , and the associated pressure fields, p2 , are simplified to [33]

∇ p2 = μ∇ 2 u2 ,

(17)

∇ · u2 = 0.

(18)

Here, the bottom and top surfaces (z = ±h/2) were considered as LV boundary conditions, surfaces x = 0 and y = 0 were symmetric conditions and the remaining were slip boundary conditions (boundary layer is not considered here).

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

463

Fig. 4. (Colour online) Modelling of the inverse four-quadrant transducer-plane streaming pattern: (a) the considered 3D model dimensions (mm), where the dash-dot lines show the symmetry planes; and the total first-order acoustic pressure magnitudes (Pa); (b) the standing wave pressure component (Pa); (c) the travelling wave pressure component (Pa); (d) the active intensity field on the bottom boundary z = −h/2 (W/m2 ); (e) the 3D acoustic streaming velocity field (m/s), where velocity vectors are shown at two heights within the model (z positions of one third and two thirds of the model height); and (f) xy view of the streaming field shown in (e), where the five-pointed star represents the vortex centre. The amplitudes of the standing and travelling wave components used in simulations are 3 and 0.3 MPa, respectively.

4.3. Results and discussion It should be noted that, prior to defining the acoustic field using Eqs. (10)–(12), the linear acoustic model applied in previous work [23,25] was used to predict the total 3D SW mode, which can narrow the searching list of SW and TW combinations for the formation of this inverse four-quadrant TPS pattern. The main principle is based on the fact [23] that the vorticity in LV pattern, which drives the vorticity of TPS patterns in thin-layer acoustofluidic devices, is proportional to the solenoidal part of the active intensity field (which are analysed further below). Thus, to seek an acoustic field for the formation of this inverse four-quadrant TPS pattern, it is only necessary to find the acoustic field which induces an inverse four-quadrant active intensity pattern on the driving boundaries that the limiting velocities are calculated. It would be very difficult to know the TW mode from the experiments alone without this theoretical model. It was found from the linear acoustic model that the total acoustic pressure field shown in Fig. 4(a), which is close to the (1, 2, 1) cavity mode but includes both SW and TW components, gives rise to the inverse four-quadrant active intensity pattern, which is similar to that driving the four-quadrant TPS pattern in another device [25]. To distinguish the contributions of respectively the SW and TW components in this device, we thus define here a (1, 2, 1) SW mode (shown in Fig. 4(b)) according to the modelled total acoustic pressure field and examine the active intensity patterns under different TW modes. Through parametric studies, we found that, for pressure expressions presented in Eqs. (10)–(12), a combination of (1, 2, 1) SW mode and (t, 0, 1) TW mode (see Fig. 4(c)) and a phase difference ϕ = π /2 between them are required to produce the inverse four-quadrant active intensity pattern on the driving boundaries (shown in Fig. 4(d)). The modelled acoustic streaming fields driven by acoustic pressure fields illustrated above are presented in Fig. 4(e)-(f), where a well-formed TPS vortex (vortex centre near y = −w/4) can be seen in the quadrant model, which is the same as that of well-formed fourquadrant TPS pattern [25] but rotates in an opposite direction. It is worth mentioning that the acoustic pressure magnitudes presented in Fig. 4(a)-(c) are exemplified pressure fields (representing p0t  p0s ) which generate the streaming pattern and velocity magnitudes plotted in Fig. 4(e). Until now, as summarised in Fig. 5, in two glass capillaries at different fluid channel sizes (Fig. 5(a)), we have shown the full set of well-formed four-quadrant TPS patterns (Fig. 5(b)), in which the streaming vortex in each quadrant of the vertical (z-direction) pressure nodal plane rotates in opposite directions. Although the first order acoustic pressure fields established in these two devices (shown in Fig. 5(c)) are the same combination of SW and TW components, the phase differences between the SW and TW components in these two devices differ by π (Fig. 5(d)). To obtain more insights into the underlying mechanism, i.e. to understand why a change of phase difference between the SW and TW components can make these two four-quadrant TPS patterns rotate in opposite directions (Fig. 5(b)), we must establish how the phase difference in the driving terms affects the vorticity of the LV vector field in the xy plane (at the boundary where the LVs are calculated, z = ±h/2) as the vorticity of TPS fields is driven by vorticity in the LV patterns themselves.

464

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

Fig. 5. (Colour online) The full set of well-formed four-quadrant transducer-plane streaming that have been observed in two thin-layer acoustofluidic devices: (a) dimensions of fluid channel cross-sections; (b) four-quadrant transducer-plane streaming patterns; (c) acoustic fields driving the corrsponding streaming fields shown in (b); and (d) phase differences between the standing wave (SW) and travelling wave (TW) components for the formation of streaming patterns presented in (b).

According to the pressure distributions expressed in Eqs. (10)–(12), the LV vector field at the driving boundaries can be written as

uL = −

1 1 2 2kxs + 2k2ys Jx + k2xt + k2yt − 2k2zt Ix − Re{qx }, 4ω 2 ρ0 ω 2

(19)

vL = −

1 1 2 2kxs + 2k2ys Jy + k2xt + k2yt − 2k2zt Iy − Re{qy }, 4ω 2 ρ0 ω 2

(20)

where ρ 0 is the fluid density, Ix and Iy are the x- and y-components of the active intensity (the real part of the complex acoustic intensity), and Jx and Jy are the x- and y-components of the reactive intensity (the imaginary part of the complex acoustic intensity). It can be seen that three vector fields contribute to the LV field: (1) reactive intensity, (Jx , Jy ); (2) active intensity, (Ix , Iy ); and (3) (qx , qy ). Firstly, using the relation ∇ × u1 = 0 in a linear acoustic field (thus the spatial derivatives of ∇ × u1 must also be zero), we can derive that the curl of the field Q = (qx , qy , 0) is everywhere zero and hence will not contribute to the vorticity of the LV field and thus four-quadrant TPS patterns. Then, according to Mann [30], only the active intensity, the real part of complex sound intensity, can have a rotational component while reactive intensity is irrotational in an acoustic field. Thus, the rotational component of the streaming field in the xy plane is proportional to the active sound intensity components in Eqs. (19) and (20). In other words, the vorticity of TPS is only related to the active intensity field on the driving boundaries. In the following, we are going to discuss how the phase difference between the SW and the TW influences the active intensity vector field and hence the four-quadrant TPS patterns. As shown in Fig. 5, both the standard and inverse four-quadrant TPS patterns are excited at the same TW mode, (t, 0, 1), so we have kyt = 0 for these two cases. Therefore, for pressure distributions described in Eqs. (10)–(12), the active intensity vector field at the driving boundaries z = ±h/2 can be derived:

Ix = − Iy =

0.5 p1t

ρ0 ω

{kxt p1t + p1s cos (kys y )[kxt cos (kxs x ) cos (kxt x + ϕ ) + kxs sin (kxt x + ϕ ) sin (kxs x )]},

0.5 p1t p1s

ρ0 ω

cos (kxs x ) sin (kxt x + ϕ )[−kys sin (kys y )].

(21) (22)

It can be seen that the phase difference ϕ has effects on both the x- and y-components of the active intensity. It can be easily known from Eq. (22) that, if ϕ is increased/decreased by π , Iy at each position (x, y) in the fluid channel will point to the opposite direction under the effect of the term sin (kxt x + ϕ ) although its magnitude will remain the same. The situation for Ix , however, is a bit more complicated. As described in Eq. (21), the quantity in the curly bracket is a summation of two terms, respectively multiples of p1t and p1s , while only the second term is dependent on the phase difference ϕ . Considering the condition that p1t should be weaker (generally a few times smaller) than p1s for the formation of TPS patterns (otherwise no vortex will be formed in the xy plane) [26], the second term shall have a bigger impact on the magnitudes of the quantity in the curly bracket of Eq. (21) (which is proved from our model). Thus, the sign of Ix at each position (x, y) in the fluid channel will be altered when ϕ differ by π for a sign change of the second term, and the magnitude of Ix at each position (x, y) in the fluid channel will thus also be different, either larger or smaller depending on the relationship between the initial signs of the first and second terms in the curly bracket. Therefore, the mechanism of how the phase difference ϕ between the SW and TW components of 3D acoustic pressure fields excited in thin-layer channels influences the resulting TPS streaming patterns can here be briefly summarised: the

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

465

Fig. 6. (Colour online) Relationship between the phase differences (ϕ) between the standing wave and travelling wave components and the x-positions of centres of the transducer-plane streaming vortices at the vertical acoustic pressure nodal plane (z = 0), where the insets show orientations of streaming vortices in the third quadrant for negative and positive ϕ. The slope of the solid lines is − 1/kxt (kxt is the travelling wave number in the x direction).

active intensity vortex pattern on the top and bottom surfaces of the model (i.e. the driving boundaries where LVs are calculated) can be modified by a change of ϕ and is reversed when ϕ differ by π , and so does the LV field, from which standard and inverse TPS streaming patterns in the bulk of the fluid are thus generated. To further investigate the effect of phase difference between the SW and TW components on the TPS fields, we modelled the acoustic streaming fields at different ϕ ranging from − −π to π in the same fluid channel shown in Fig. 4(a) and found that a change of ϕ can shift the streaming vortex along the fluid channel (in the x direction). The modelled x-positions of the centre of the four-quadrant TPS vortex (indicated with a five-pointed star in Fig. 4(f)) in the third quadrant for phases − π < ϕ < π are presented in Fig. 6, where, for comparison, the calculated values of the vortex centre of active intensity are also plotted. The theoretical line was calculated based on Eq. (22), which indicates that it is only necessary to make sin (kxt x + ϕ ) = 0 if we would like to calculate the x-positions of active intensity vortex, and thus in this third quadrant model we have

x0 = −

ϕ kxt

,

(23)

where x0 is the x-position of active intensity vortex centre on the driving boundaries. It can be seen that the modelled results compare very well with the theoretical predictions. For every change of π on the phase difference between the SW and TW components, an inverse TPS pattern can be formed. Although many of these streaming patterns have not been experimentally observed, the consistence between the experiments and numerical modelling on these two reversed fourquadrant TPS patterns indicates that this is an efficient model for predicting 3D acoustic and streaming fields in thin-layer acoustofluidic devices and the results may provide insights for many acousto-microfluidic applications (discussed in the following sections). In a previous work, we showed direct comparisons between the modelled and measured acoustic streaming velocity magnitudes in a glass capillary device (Fig. 10 in Ref. [25]). Further work, however, are still required to make more precise comparisons between the measured acoustic streaming fields and those obtained from the mathematical model in this work. It was found from the numerical simulations that, for a fixed value of p0s (or p0t ), the magnitudes of acoustic streaming velocities increase with a growth of p0t (or p0s ), making it difficult to directly compare the acoustic streaming velocity magnitudes between the modelling and experiments as it could hardly experimentally separate the SW and TW components in a thin-layer resonant microfluidic channel fabricated on acoustically hard materials, e.g. the one shown in Fig. 1, where the SW and TW components are generated from the same oscillation of channel walls. We hope in future to build an acoustofluidic system that can individually control the SW and TW components and make better comparisons between the numerical simulations and experiments. 5. Revisit of unusual 6 × 6 in-plane streaming patterns Having established the underlying physics of standard and inverse TPS patterns in glass capillaries, it is worth mentioning that the above numerical analyses can also be applied to explain the complex phenomena seen in other acoustofluidic particle manipulation systems. Here, as an example, we revisited the problem of unusual 6 × 6 in-plane acoustic streaming vortex patterns seen in our previous work (Section IV-2 in Ref. [24]). While interpreting the driving mechanism of a 6 × 6

466

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

Fig. 7. (Colour online) Numerical simulations of unusual 6 × 6 in-plane streaming patterns in an acoustofluidic particle manipulation device: (a) – (b) the full device model considered in a previous work [24] and the corresponding dimensions of the fluid channel (unit: mm); (c) the simplified, reduced-fluid quadrant model considered in this work; and (d) – (e) the modelled two 6 × 6 transducer-plane streaming patterns in the central square area of the fluid channel formed under different acoustic conditions, where ϕ is the phase difference between the (0, 3, 0) standing wave mode and the (t, 0, 0) travelling wave mode, and the dashed-dot lines represent symmetric boundary conditions.

in-plane acoustic streaming vortex pattern observed in a silicon-based acoustofluidic particle manipulation device (Fig. 4 in Ref. [22]), through 3D full device (including the PZT transducer, silicon, water and glass layers, shown in Fig. 7(a)) numerical simulations, we previously demonstrated that the relative phases of the total acoustic field in the visualisation plane vary in fluid channels with different dimensions (Section IV-3 in Ref. [24]). And with a slight change on the dimensions of the central square region of the fluid channel, each vortex of the 6 × 6 in-plane acoustic streaming vortex pattern was found to circulate inversely. With the help of the model presented in this paper, we are now in a position to have a deeper understand why the orientation of each vortex of those two 6 × 6 in-plane acoustic streaming patterns in the xy plane are reversed when the fluid channel dimensions are slightly different: it is because the phase differences between the SW and TW components established in these two fluid channels differ by π , which cause different rotations of active sound intensity on the driving boundaries and hence in-plane streaming fields in the bulk of the fluid. As shown in Fig. 7(c), we considered a reduced-fluid model consisting of a quarter of the fluid channel in the full device model (Fig. 7(b)) as the latter is symmetric to planes x = 0 and y = 0, and defined the 3D total acoustic pressure field from combinations of SW and TW components using Eqs. (10)–(12). This is impossible to implement in our previous full device model, which solves the acoustic coupling of different layers of resonant acoustofluidic devices, as both the SW and TW components in the fluid channel were excited from the same harmonic vibration of channel walls. Firstly, with reference to the modelled acoustic pressure field plotted in Fig. 9 (a) in Ref. [24], the SW mode in the fluid channel was determined, which is a (0, 3, 0) mode. Then, through parametric studies, it was found that, as shown in Fig. 7(d)-(e), a combination of a (0, 3, 0) SW mode and a (t, 0, 0) TW mode is required to generate these two well-formed 6 × 6 in-plane vortex patterns. The phase differences of these two cases differ by π , as is the case for the standard and inverse four-quadrant TPS patterns demonstrated in Fig. 5(b). The mechanisms behind this are similar to those discussed in Section 4.3: for acoustic pressure fields whose ϕ differ by π , the circulations of active intensity fields, the sole rotational part in the LV equations, on the driving boundaries are reversed, which thus generate the reversed 6 × 6 in-plane acoustic streaming patterns in the bulk of the fluid channels. 6. Concluding discussion In this paper, we have shown for the first time an inverse four-quadrant TPS pattern from numerical modelling and experimental measurements, and thus have presented the full set of four-quadrant TPS patterns in thin-layer acoustofluidic

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

467

manipulation devices and demonstrated their driving mechanisms. The importance of the phase difference between the SW and TW components, ϕ , on the formation of TPS patterns in thin-layer acoustofluidic particle manipulation devices have been illustrated in detail. In addition to demonstrating the newly observed inverse four-quadrant TPS pattern, we have predicted here that, at various ϕ , a variety of TPS patterns could be excited. Full LV equations have been derived from the first-order acoustic velocity fields containing both the SW and TW components, showing that the vorticity of the LV field that drives the 3D TPS vortices is only related to the rotational active intensity field, which supports our findings from linear acoustic models previously applied for the modelling of TPS patterns [23]. The numerical model for the simulation of 3D acoustic streaming fields presented in this work can accurately capture both the SW and TW features of a 3D acoustic field established in the fluid channel of an acoustofluidic device, which is extremely useful for interpreting acoustofluidic phenomena too complex for analytical solutions. As an example, based on our previous work, through numerical simulations we have presented deeper investigations on the underlying mechanisms of unusual 6 × 6 in-plane streaming patterns observed in a silicon-based acoustofluidic particle manipulation system. As this numerical model solves acoustic streaming fields in acoustic fields containing SW and TW modes that could be generated by any means, we believe that it cannot only be used to simulate TPS patterns in resonant BAW microfluidic particle manipulation devices such as glass capillaries and silicon-based microfluidic channels discussed in this work, but could also be applied to predict those in SAW-actuated particle manipulation devices [50] and those associated with Bessel beams [51]. An important next step is to understand how device structures, ultrasonic excitations and fluid channel dimensions affect the phase difference between the SW and TW components, ϕ , in a 3D acoustic field established in the fluid channel of an acoustofluidic device, which could provide means for the control of ϕ in a single thin-layer acoustofluidic manipulation device. By establishing this, a variety of potential micro-acoustofluidic applications could be realised. For example, (1) in thin-layer resonant acoustofluidic particle manipulation devices working in the MHz region, one would expect the movement of submicron particles exposed in the fluid follows the TPS circulations. As predicted in this work, a change of ϕ could shift the TPS vortices along the fluid channel, which thus would be useful for the transport/pumping of a single particle or a cluster of particles in thin-layer resonant acoustofluidic devices, similar to the fluid and particle pumping using acoustic streaming vortex arrays produced from asymmetrically oscillating sharp edges [52]. (2) Moreover, in a single acoustofluidic manipulation device, sweeping between two reversed TPS streaming patterns could in average dramatically decrease the streaming velocity magnitudes but will not significantly affect the acoustic radiation forces on microparticles, which could open doors for the two-dimensional (y- and z-directions) manipulation of acoustic streaming-dominated motions, such as those of nanometre sized particles, cells, virus, bacteria and so forth. In a thin-layer fluid channel fabricated on sound hard materials such as glass or silicon, it would be very difficult to separate/control the SW and TW components within the device as they both result from the same vibration of channel walls. However, it might be achieved with other designs. In a water tank, for example, it could be possible to control the magnitudes of SW and TW components and their phase relations with isolated transducers or ultrasonic phased arrays. Acknowledgements The authors gratefully acknowledge the financial support for this work received from the National Natural Science Foundation of China (grant number 11804060); the Seed Fund Programme of School of Electromechanical Engineering; and the Youth Hundred-Talents Programme of Guangdong University of Technology (grant number 220413195). References [1] B. Hammarstrom, T. Laurell, J. Nilsson, Seed particle-enabled acoustic trapping of bacteria and nanoparticles in continuous flow systems, Lab. Chip 12 (2012) 4296–4304, doi:10.1039/C2lc40697g. [2] B. Hammarstrom, B. Nilson, T. Laurell, J. Nilsson, S. Ekstrom, Acoustic trapping for bacteria identification in positive blood cultures with MALDI-TOF mS, Anal. Chem. 86 (2014) 10560–10567, doi:10.1021/ac502020f. [3] M. Wiklund, R. Green, M. Ohlin, Acoustofluidics 14: applications of acoustic streaming in microfluidic devices, Lab. Chip 12 (2012) 2438–2451, doi:10. 1039/C2lc40203c. [4] C. Eckart, Vortices and streams caused by sound waves, Phys. Rev. 73 (1948) 68–76, doi:10.1103/PhysRev.73.68. [5] J. Vanneste, O. Buhler, Streaming by leaky surface acoustic waves, Proc. R. Soc. Math. Phys. 467 (2011) 1779–1800, doi:10.1098/rspa.2010.0457. [6] S.J. Lighthill, Acoustic streaming, J. Sound Vib. 61 (1978) 391–418, doi:10.1016/0022- 460X(78)90388- 7. [7] H. Schlichting, Berechnung ebener periodischer grenzschichtstromungen (Calculation of plane periodic boundary layer streaming), Physikalische Zeitschrift 33 (1932) 327–335. [8] L. Rayleigh, On the circulation of air observed in kundt’s Tubes, and on some allied acoustical problems, Philos. Trans. R. Soc. Lond. 175 (1884) 1–21. [9] J.J. Lei, M. Hill, C.P.D. Albarran, P. Glynne-Jones, Effects of micron scale surface profiles on acoustic streaming, Microfluid Nanofluidics 22 (2018) 140, doi:10.1007/S10404- 018- 2161- 2. [10] W.L. Nyborg, Acoustic streaming near a boundary, J. Acoust. Soc. Am. 30 (1958) 329–339, doi:10.1121/1.1909587. [11] C.P. Lee, T.G. Wang, Near-boundary streaming around a small sphere due to two orthogonal standing waves, J. Acoust. Soc. Am. 85 (1989) 1081–1088, doi:10.1121/1.397491. [12] R. Barnkob, P. Augustsson, T. Laurell, H. Bruus, Acoustic radiation- and streaming-induced microparticle velocities determined by microparticle image velocimetry in an ultrasound symmetry plane, Phys. Rev. E 86 (2012) 056307, doi:10.1103/Physreve.86.056307. [13] P.B. Muller, R. Barnkob, M.J.H. Jensen, H. Bruus, A numerical study of microparticle acoustophoresis driven by acoustic radiation forces and streaminginduced drag forces, Lab. Chip 12 (2012) 4617–4627, doi:10.1039/C2LC40612H. [14] P. Hahn, I. Leibacher, T. Baasch, J. Dual, Numerical simulation of acoustofluidic manipulation by radiation forces and acoustic streaming for complex particles, Lab. Chip 15 (2015) 4302–4313, doi:10.1039/c5lc00866b.

468

J. Lei, F. Cheng and Z. Guo / Applied Mathematical Modelling 77 (2020) 456–468

[15] D.A. Gubaidullin, P.P. Ossipov, A.A. Abdyushev, Simulation of aerosol distribution in hyperbolic resonator, Appl. Math. Model 62 (2018) 181–193, doi:10. 1016/j.apm.2018.05.028. [16] F. Li, X.X. Xia, Z.T. Deng, J.J. Lei, Y.X. Shen, Q. Lin, W. Zhou, L. Meng, J.R. Wu, F.Y. Cai, H.R. Zheng, Ultrafast Rayleigh-like streaming in a sub-wavelength slit between two phononic crystal plates, J. Appl. Phys. 125 (2019) 134903, doi:10.1063/1.5058206. [17] V. Daru, I. Reyt, H. Bailliet, C. Weisman, D. Baltean-Carlès, Acoustic and streaming velocity components in a resonant waveguide at high acoustic levels, J. Acoust. Soc. Am. 141 (2017) 563–574, doi:10.1121/1.4974058. [18] V. Daru, C. Weisman, D. Baltean-Carlès, I. Reyt, H. Bailliet, Inertial effects on acoustic Rayleigh streaming flow: transient and established regimes, Wave Motion 74 (2017) 1–17, doi:10.1016/j.wavemoti.2017.06.001. [19] J.T. Karlsen, W. Qiu, P. Augustsson, H. Bruus, Acoustic streaming and its suppression in inhomogeneous fluids, Phys. Rev. Lett. 120 (2018) 054501, doi:10.1103/PhysRevLett.120.054501. [20] W. Qiu, J.T. Karlsen, H. Bruus, P. Augustsson, Experimental characterization of acoustic streaming in gradients of density and compressibility, Phys. Rev. Appl. 11 (2019) 024018, doi:10.1103/PhysRevApplied.11.024018. [21] D. Baltean-Carlès, V. Daru, C. Weisman, S. Tabakova, H. Bailliet, An unexpected balance between outer rayleigh streaming sources, J. Fluid Mech. 867 (2019) 985–1011, doi:10.1017/jfm.2019.111. [22] S.M. Hagsater, T.G. Jensen, H. Bruus, J.P. Kutter, Acoustic resonances in microfluidic chips: full-image micro-PIV experiments and numerical simulations, Lab. Chip 7 (2007) 1336–1344, doi:10.1039/B704864e. [23] J.J. Lei, P. Glynne-Jones, M. Hill, Modal Rayleigh-like streaming in layered acoustofluidic devices, Phys. Fluids 28 (2016) 012004, doi:10.1063/1.4939590. [24] J. Lei, M. Hill, P. Glynne-Jones, Numerical simulation of 3D boundary-driven acoustic streaming in microfluidic devices, Lab. Chip 14 (2014) 532–541, doi:10.1039/c3lc50985k. [25] J. Lei, P. Glynne-Jones, M. Hill, Acoustic streaming in the transducer plane in ultrasonic particle manipulation devices, Lab. Chip 13 (2013) 2133–2143, doi:10.1039/c3lc0 0 010a. [26] J.J. Lei, M. Hill, P. Glynne-Jones, Transducer-Plane streaming patterns in thin-layer acoustofluidic devices, Phys. Rev. Appl. 8 (2017) 014018, doi:10.1103/ PhysRevApplied.8.014018. [27] A.A. Doinikov, P. Thibault, P. Marmottant, Acoustic streaming induced by two orthogonal ultrasound standing waves in a microfluidic channel, Ultrasonics 87 (2018) 7–19, doi:10.1016/j.ultras.2018.02.002. [28] A.A. Doinikov, P. Thibault, P. Marmottant, Acoustic streaming in a microfluidic channel with a reflector: case of a standing wave generated by two counterpropagating leaky surface waves, Phys. Rev. E 96 (2017) 013101, doi:10.1103/PhysRevE.96.013101. [29] A.A. Doinikov, P. Thibault, P. Marmottant, Acoustic streaming generated by two orthogonal standing waves propagating between two rigid walls, J. Acoust. Soc. Am. 141 (2017) 1282, doi:10.1121/1.4976088. [30] J.A. Mann, J. Tichy, A.J. Romano, Instantaneous and time-averaged energy-transfer in acoustic fields, J. Acoust. Soc. Am. 82 (1987) 17–30, doi:10.1121/ 1.395562. [31] H. Bruus, Acoustofluidics 7: the acoustic radiation force on small particles, Lab. Chip 12 (2012) 1014–1021, doi:10.1039/C2lc21068a. [32] W.L. Nyborg, Acoustic streaming, in: W.P. Mason (Ed.), Physical Acoustics, Academic, New York, 1965, pp. 290–295. [33] J.J. Lei, P. Glynne-Jones, M. Hill, Comparing methods for the modelling of boundary-driven streaming in acoustofluidic devices, Microfluid Nanofluidics 21 (2017) 23, doi:10.1007/s10404- 017- 1865- z. [34] H. Bruus, Acoustofluidics 2: perturbation theory and ultrasound resonance modes, Lab. Chip 12 (2012) 20–28, doi:10.1039/C1lc20770a. [35] S.S. Sadhal, Acoustofluidics 13: analysis of acoustic streaming by perturbation methods foreword, Lab. Chip 12 (2012) 2292–2300, doi:10.1039/ C2lc40202e. [36] H. Bruus, Theoretical aspects of microscale acoustofluidics, 2018, arXiv:1802.05597. [37] W.L. Nyborg, Acoustic streaming, in: M.F. Hamilton, D.T. Blackstock (Eds.), Nonlinear Acoustics, Academic, San Diego, 1998. [38] M.W.H. Ley, H. Bruus, Three-Dimensional numerical modeling of acoustic trapping in glass capillaries, Phys. Rev. Appl. 8 (2017) 024020, doi:10.1103/ PhysRevApplied.8.024020. [39] D. Carugo, D.N. Ankrett, P. Glynne-Jones, L. Capretto, R.J. Boltryk, X.L. Zhang, P.A. Townsend, M. Hill, Contrast agent-free sonoporation: the use of an ultrasonic standing wave microfluidic system for the delivery of pharmaceutical agents, Biomicrofluidics 5 (2011) 044108, doi:10.1063/1.3660352. [40] D.N. Ankrett, D. Carugo, J. Lei, P. Glynne-Jones, P.A. Townsend, X. Zhang, M. Hill, The effect of ultrasound-related stimuli on cell viability in microfluidic channels, J Nanobiotechnol 11 (2013) 20, doi:10.1186/1477-3155-11-20. [41] P. Mishra, M. Hill, P. Glynne-Jones, Deformation of red blood cells using acoustic radiation forces, Biomicrofluidics 8 (2014) 034109, doi:10.1063/1. 4882777. [42] S.W. Li, P. Glynne-Jones, O.G. Andriotis, K.Y. Ching, U.S. Jonnalagadda, R.O.C. Oreffo, M. Hill, R.S. Tare, Application of an acoustofluidic perfusion bioreactor for cartilage tissue engineering, Lab. Chip 14 (2014) 4475–4485, doi:10.1039/c4lc00956h. [43] COMSOL multiphysics 5.2, COMSOL multiphysics, http://www.comsol.com/, 2015. [44] J.S. Bach, H. Bruus, Theory of pressure acoustics with viscous boundary layers and streaming in curved elastic cavities, J. Acoust. Soc. Am. 144 (2018) 766, doi:10.1121/1.5049579. [45] M.F. Hamilton, Y.A. Ilinskii, E.A. Zabolotskaya, Acoustic streaming generated by standing waves in two-dimensional channels of arbitrary width, J. Acoust. Soc. Am. 113 (2003) 153–160, doi:10.1121/1.1528928. [46] J. Lei, Formation of inverse chladni patterns in liquids at microscale: roles of acoustic radiation and streaming-induced drag forces, Microfluid Nanofluidics 21 (2017) 50, doi:10.1007/s10404- 017- 1888- 5. [47] F. Li, Y. Xiao, J.J. Lei, X.X. Xia, W. Zhou, L. Meng, L.L. Niu, J.R. Wu, J.Y. Li, F.Y. Cai, H.R. Zheng, Rapid acoustophoretic motion of microparticles manipulated by phononic crystals, Appl. Phys. Lett. 113 (2018) 173503, doi:10.1063/1.5052045. [48] J. Lei, An investigation of boundary-driven streaming in acoustofluidic systems for particle and cell manipulation, Faculty of Engineering and Physical Sciences (Ph.D. thesis), University of Southampton, 2015, p. 165. [49] N. Nama, T.J. Huang, F. Costanzo, Acoustic streaming: an arbitrary lagrangian-eulerian perspective, J. Fluid Mech. 825 (2017) 600–630, doi:10.1017/jfm. 2017.338. [50] A. Ozcelik, J. Rufo, F. Guo, Y.Y. Gu, P. Li, J. Lata, T.J. Huang, Acoustic tweezers for the life sciences, Nat. Methods 15 (2018) 1021–1028, doi:10.1038/ s41592- 018- 0222- 9. [51] Z.Y. Hong, J. Zhang, B.W. Drinkwater, Observation of orbital angular momentum transfer from bessel-shaped acoustic vortices to diphasic liquidmicroparticle mixtures, Phys. Rev. Lett. 114 (2015) 214301, doi:10.1103/Physrevlett.114.214301. [52] P.H. Huang, N. Nama, Z. Mao, P. Li, J. Rufo, Y. Chen, Y. Xie, C.H. Wei, L. Wang, T.J. Huang, A reliable and programmable acoustofluidic pump powered by oscillating sharp-edge structures, Lab. Chip 14 (2014) 4319–4323, doi:10.1039/c4lc00806e.