2.1. Governing equation and boundary conditions

The flow satisfies the incompressible Navier–Stokes equations

(2.2)\begin{equation} \rho\left(\frac{\partial\boldsymbol{V}}{\partial T}+\boldsymbol{V}\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{V} \right)={-} \boldsymbol{\nabla} P+\mu\,\nabla^2 \boldsymbol{V}, \quad \boldsymbol{\nabla} \boldsymbol{\cdot}\boldsymbol{V}=0, \end{equation}

where $\boldsymbol {V} (\boldsymbol {X},T)$ is the 3-D fluid velocity, and $P(\boldsymbol {X},T)$ is the pressure. The flow satisfies the no-slip condition ($\boldsymbol {V}=0$) at the cylinder surface and the channel walls.

We non-dimensionalize the equations by choosing $a$ as the characteristic length scale in the $XY$ plane, and $h$ as the characteristic scale in the $Z$ direction. We define dimensionless coordinates (lowercase) according to

(2.3a–d)\begin{equation} X=ax,\quad Y=ay,\quad Z=hz, \quad R=ar. \end{equation}

For later convenience, we write the velocity as $\boldsymbol {V} = \boldsymbol {V}_{\parallel } + V_Z \boldsymbol {e}_Z$, where $\boldsymbol {V}_{\parallel }$ represents velocities in the horizontal ($XY$) plane, and $V_Z$ is the vertical velocity component. The velocity scale $\bar {V}_\infty$ is characteristic of the horizontal velocity $\boldsymbol {V}_{\parallel }$. Continuity then suggests a vertical velocity scale $\bar {V}_\infty h/a$. Substituting these velocity scales into (2.2), and comparing viscous and pressure terms, suggests a pressure scale $\mu \bar {V}_{\infty } a/h^2$. We therefore define the dimensionless velocity components $\boldsymbol {v}_{\parallel }$ and $v_z$, and pressure $p$, according to

(2.4a–c)\begin{equation} \boldsymbol{V}_{{\parallel}}= \bar{V}_\infty \boldsymbol{v}_{{\parallel}},\quad V_Z=\frac{\bar{V}_\infty h}{a}\,v_z, \quad P= \frac{\mu \bar{V}_\infty a}{h^2}\,p. \end{equation}

We also define a rescaled time $t = \omega T$ and introduce the dimensionless parameters

(2.5a–c)\begin{equation} \delta = \sqrt{\frac{2\nu}{h^2 \omega}},\quad \epsilon=\frac{\bar{V}_{\infty}}{a \omega}\quad \mbox{and} \quad \eta=\frac{h}{a}, \end{equation}

where $\nu =\mu /\rho$ is the kinematic viscosity of the fluid. Here, $\epsilon$ is the dimensionless amplitude of the oscillatory flow, and $\delta$ is the ratio of the Stokes layer thickness to the channel height, which characterizes the frequency of the flow, while $\eta$ is the aspect ratio of the cylinder, which characterizes the degree of axial confinement.

The dimensionless 3-D fluid velocity is therefore $\boldsymbol {v} = \boldsymbol {V}/\bar {V}^{\infty }=\boldsymbol {v}_{\parallel } + \eta v_z \boldsymbol {e}_z$. Thus (2.2) rescales (without approximations) as

(2.6a)\begin{gather} \frac{2}{\delta^2} \left(\frac{\partial \boldsymbol{v}_{{\parallel}}}{\partial t}+\epsilon\boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{v}_{{\parallel}}\right)={-}\nabla_{{\parallel}} p+\left(\eta^2\,\nabla^2_{{\parallel}} +\frac{\partial ^2}{\partial z^2}\right) \boldsymbol{v}_{{\parallel}}, \end{gather}

(2.6b)\begin{gather}\frac{2}{\delta^2} \left(\frac{\partial v_z}{\partial t}+\epsilon \boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{\nabla} v_z\right)={-}\frac{1}{\eta^2}\,\frac{\partial p}{\partial z}+\left(\eta^2\, \nabla^2_{{\parallel}}+\frac{\partial ^2}{\partial z^2}\right)v_z, \end{gather}

(2.6c)\begin{gather}\nabla_{{\parallel}}\boldsymbol{\cdot}\boldsymbol{v}_{{\parallel}}+\frac{\partial v_z}{\partial z}=0, \end{gather}

where $\nabla _{\parallel }=\boldsymbol {e}_{x}({\partial }/{\partial x})+\boldsymbol {e}_{y}({\partial }/{\partial y})$ is the gradient in the horizontal plane. Using overbars to indicate the average across the height of the channel (for a dimensionless function $f(z)$, $\bar {f} \equiv (1/2) \int _{-1}^1f(z)\,\mbox {d}z$), the height-averaged velocity at infinity satisfies $\bar {\boldsymbol {v}}(r \rightarrow \infty )=\boldsymbol {e}_{x}\,\mathrm {e}^{{\rm i}t}$. The flow also satisfies no-slip conditions on the cylinder surface and channel walls.

3.1. Small-amplitude approximation and outer flow

We now solve (3.1) and (3.2) to obtain the outer flow in the practically relevant limit of $\epsilon \ll 1$, but with arbitrary $\delta$. As is typical in many streaming calculations (Riley Reference Riley2001), we assume that the oscillation amplitude is smaller than other geometric scales, so $\epsilon \ll \delta$ and $\epsilon \ll \eta \ll 1$. For small $\epsilon$, we develop a solution by applying a perturbation expansion

(3.3)\begin{equation} (p, \boldsymbol{v})=(p_1, \boldsymbol{v}_1)+\epsilon\,(p_2, \boldsymbol{v}_2)+\cdots. \end{equation}

The subscript $1$ identifies the primary flow (oscillatory flow), whereas the subscript $2$ corresponds to the secondary steady flow (streaming). The above series expansion, when substituted into (3.1) and separating orders of $\epsilon$, yields the governing equations for the primary and secondary flow.

3.1.1. Primary flow

The primary flow at $O(\epsilon ^0)$ is governed by

(3.4) \begin{equation} \left.\begin{array}{c@{}} \dfrac{2}{\delta^2}\,\dfrac{\partial \boldsymbol{v}_{1 \parallel}}{\partial t}={-}\nabla_\parallel p_1+\dfrac{\partial^2 \boldsymbol{v}_{1 \parallel}}{\partial z^2}, \quad \dfrac{\partial p_1}{\partial z}=0, \quad \nabla_{{\parallel}}\boldsymbol{\cdot}\boldsymbol{v}_{1 \parallel}+\dfrac{\partial v_{1z}}{\partial z}=0,\\ \mbox{with}\quad\bar{\boldsymbol{v}}_1|_{r\rightarrow\infty} =\boldsymbol{e}_{x}\, \mathrm{e}^{{\rm i}t},\quad \boldsymbol{v}_{1 \parallel}|_{z={\pm} 1}=\boldsymbol{0}, \quad v_{1z}|_{z={\pm} 1}=0,\quad \bar{\boldsymbol{v}}_{1 \parallel }\boldsymbol{\cdot} \boldsymbol{n}|_{r=1}=0. \end{array}\right\} \end{equation}

The system (3.4) is linear, and the only non-homogeneous condition is the one on velocity at infinity. We therefore seek separable solutions proportional to $\mathrm {e}^{\textrm {i}t}$. Observing that pressure is independent of $z$, the general expression for the primary velocity in the outer region is

(3.5)\begin{equation} \boldsymbol{v}_1={-}\frac{ \boldsymbol{\nabla} p_1}{\alpha^2} \left(1-\frac{\cosh{\alpha z}}{\cosh{\alpha}}\right),\quad\mbox{where}\ \alpha=\frac{1+\mathrm{i}}{\delta}.\end{equation}

The continuity equation, along with the normal velocity condition at $r = 1$ (see (3.4)), yields the primary pressure as

(3.6)\begin{equation} p_1= G(\alpha)\left(r+\frac{1}{r}\right)\cos{\theta}\,\mathrm{e}^{{\rm i}t}, \quad\mbox{where}\ G(\alpha) ={-}\frac{\alpha^3}{\alpha-\tanh{(\alpha})},\end{equation}

which agrees with Shih (Reference Shih1970) up to scaling. Note that the result (3.5) is valid in the outer region independent of the details of the primary pressure $p_1$. By contrast, the expression (3.6) explicitly ignores the presence of the inner region. In Appendix A, we will see that a matching condition with the inner region leads to $O(\eta )$ corrections of the outer pressure in (3.4), while leaving (3.5) unchanged. Observe that the primary flow has no vertical velocity component ($v_{1z}=0$), and that the primary pressure varies only in the $r\theta$ plane. The channel-averaged primary (oscillatory) velocity is

(3.7)\begin{equation} \bar{\boldsymbol{v}}_1 = \frac{ \boldsymbol{\nabla} p_1}{G(\alpha)} = \boldsymbol{\nabla} \left[\left(r+\frac{1}{r}\right)\cos{\theta}\right] \mathrm{e}^{{\rm i}t}, \end{equation}

and is used in the calculation of the streaming flow below.

3.1.2. Secondary flow

The convective acceleration term in (3.1a), which is neglected when solving for the primary flow, appears as a body force that drives the secondary flow. It will be understood that all secondary quantities are time-averaged. We assume that $\epsilon \ll \delta$, so that the inertia of the time-averaged secondary flow is negligible Riley (Reference Riley2001), and the governing equations at $O(\epsilon ^1)$ are

(3.8a)\begin{align} \frac{2}{\delta^2}\left<\boldsymbol{v}_{1 \parallel}\boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{v}_{1 \parallel}\right>&={-}\nabla_\parallel p_2+\frac{\partial^2 \boldsymbol{v}_{2 \parallel}}{\partial z^2}, \quad \frac{\partial p_2}{\partial z}=0, \quad \nabla_{{\parallel}}\boldsymbol{\cdot}\boldsymbol{v}_{2 \parallel}+\frac{\partial v_{2z}}{\partial z}=0, \end{align}

(3.8b)\begin{align} \mbox{with}\quad p_2|_{r\rightarrow\infty} &= 0,\quad \boldsymbol{v}_{2 \parallel}|_{z={\pm} 1}=\boldsymbol{0}, \quad v_{2z}|_{z={\pm} 1}=0,\quad\bar{\boldsymbol{v}}_{2 \parallel }\boldsymbol{\cdot} \boldsymbol{n}|_{r=1}=0. \end{align}

Here, $\langle {\cdot } \rangle$ denotes the time average over an oscillation cycle. We note the useful averaging rule $\langle \textrm {Re} [ A \,\mathrm {e}^{\mathrm {i} t}]\,\textrm {Re} [ B \,\mathrm {e}^{\mathrm {i} t}]\rangle = (1/2)\,\textrm {Re}[A^* B]$ for any complex quantities $A$ and $B$, where the asterisk denotes the complex conjugate (Longuet-Higgins Reference Longuet-Higgins1998). We then use (3.5) to evaluate the advective term $\langle \boldsymbol {v}_{1 \parallel }\boldsymbol {\cdot } \boldsymbol {\nabla } \boldsymbol {v}_{1 \parallel }\rangle = (1/2) \boldsymbol {v}_{1 \parallel }^*\boldsymbol {\cdot } \boldsymbol {\nabla } \boldsymbol {v}_{1 \parallel }$. As mentioned before, only the real part is physically relevant. Integrating the horizontal momentum equation in (3.8a) with respect to $z$, applying no-slip conditions and noting that $\alpha ^*=(1-\mathrm {i})/\delta = -\mathrm {i}\alpha$, we obtain

(3.9a)\begin{equation} \boldsymbol{v}_{2 \parallel}=\tfrac{1}{2}(z^2-1)\, \boldsymbol{\nabla} p_2+ g(z,\alpha)\, \boldsymbol{\nabla} | \boldsymbol{\nabla} p_1|^2, \end{equation}

where

(3.9b)\begin{equation} g(z,\alpha)=\frac{\mathrm{i}}{8\alpha^4}\left((z^2-1) \alpha^2+\frac{2\cos{\alpha z}}{\cos \alpha}-\frac{2\cosh\alpha z}{\cosh\alpha}+\frac{\sin\alpha z \sinh\alpha z}{\cos\alpha \cosh\alpha}-\tan\alpha \tanh\alpha\right). \end{equation}

Next, we integrate the continuity equation in (3.8a), and apply boundary conditions at the top and bottom walls, to find that the horizontal flux is divergence-free: ${\boldsymbol {\nabla } \boldsymbol {\cdot }\bar {\boldsymbol {v}}_{2 \parallel }=0}$. We then substitute (3.9) into the above relation to find that the secondary pressure is governed by

(3.10a)\begin{equation} \nabla_{{\parallel}}^2 {p_2} = 3 \bar{g}\, \nabla_{{\parallel}}^2| \boldsymbol{\nabla} p_1|^2, \end{equation}

subject to

(3.10b)\begin{equation} p_2|_{r \to \infty} = 0 \quad \mbox{and} \quad \boldsymbol{\nabla} (p_2 - 3 \bar{g}\,| \boldsymbol{\nabla} p_1|^2) \boldsymbol{\cdot} \boldsymbol{n}|_{r = 1} = 0. \end{equation}

Note that the condition at $r = 1$ corresponds to no penetration of the depth-averaged flow into the cylinder surface; see (3.8b). Solving (3.10), we obtain the secondary pressure as

(3.11a)\begin{align} p_2&=3\bar{g} (| \boldsymbol{\nabla} p_1|^2-|G|^2) \end{align}

(3.11b)\begin{align} &={-}\frac{\mathrm{i}}{16\alpha^5}\,(4 \alpha^3 + 15 (\tanh \alpha - \tan \alpha) +6\alpha\tan\alpha\,\tanh{\alpha})(| \boldsymbol{\nabla} p_1|^2-|G|^2), \end{align}

where we note that $|G|^2 = | \boldsymbol {\nabla } p_1|^2|_{r \to \infty }$ is a constant.

Substituting the expression above into (3.9) yields $\boldsymbol {v}_{2 \parallel }$, and using the continuity equation yields the vertical component of secondary flow $v_{2z}$. Observing that the secondary pressure and velocity both involve $| \boldsymbol {\nabla } p_1|^2=| G(\alpha ) |^2 |\bar {\boldsymbol {v}}_1|^2$ (see (3.7)), we obtain the secondary (time-averaged) velocity

(3.12a)\begin{equation} \boldsymbol{v}_2 = \boldsymbol{v}_{2 \parallel} + \eta v_{2z}\, \boldsymbol{e}_z = F_{{\parallel}}(z,\alpha)\, \boldsymbol{\nabla} (|\bar{\boldsymbol{v}}_1|^2 ) + \eta\,F_z(z, \alpha)\,\nabla^2 (|\bar{\boldsymbol{v}}_1|^2 )\, \boldsymbol{e}_z, \end{equation}

where

(3.12b)\begin{align} F_{{\parallel}}(z,\alpha)&=\frac{\mathrm{i}\, | G(\alpha) |^2}{32\alpha^5}\left(\frac{8\alpha \cos{\alpha z}}{\cos{\alpha}}-\frac{8\alpha \cosh{\alpha z}}{\cosh{\alpha}}+\frac{4\alpha \sin{\alpha z}\,\sinh{\alpha z}}{\cos{\alpha} \cosh{\alpha}}\right.\nonumber\\ &\quad \left. \vphantom{\frac{8\alpha \cosh{\alpha z}}{\cosh{\alpha}}}+ 15(z^2-1) (\tan \alpha -\tanh \alpha)- 2(3 z^2-1) \alpha\tan\alpha \tanh{\alpha}\right), \end{align}

(3.12c)\begin{align} F_z(z,\alpha)&={-}\frac{\mathrm{i}\,| G(\alpha) |^2}{32\alpha^5} \left(\frac{8\sin{\alpha z}}{\cos{\alpha}}- \frac{8\sinh{\alpha z}}{\cosh{\alpha}} +\frac{2(\sin{\alpha z} \cosh{\alpha z}- \cos\alpha z \sinh{\alpha z})}{\cos{\alpha}\,\cosh{\alpha}}\right.\nonumber\\ &\quad \left.\vphantom{\frac{8\alpha \cosh{\alpha z}}{\cosh{\alpha}}}+5z(z^2-3)(\tan \alpha- \tanh \alpha) - 2 z (z^2-1) \alpha \tan \alpha \tanh \alpha\right). \end{align}

The functions $F_{\parallel }(z,\alpha )$ and $F_z(z,\alpha )$ describe horizontal and vertical velocity profiles across the height of the channel, while horizontal variations are captured by $|\bar {\boldsymbol {v}}_1|^2$. We observe that the secondary flow has zero average across the channel everywhere in the flow, not just at the cylinder surface ($\overline {F}_{\parallel } = 0 \implies \bar {\boldsymbol {v}}_2 = \boldsymbol {0}$).

While $\boldsymbol {v}_2$ describes the Eulerian secondary flow, it is the time-averaged motion of material fluid elements (the Lagrangian secondary flow) that is practically relevant and is normally referred to as streaming. The time-averaged Lagrangian velocity is $\boldsymbol {v}_L=\boldsymbol {v}_2+\boldsymbol {v}_d$, where $\boldsymbol {v}_d$ is the Stokes drift defined by $\boldsymbol {v}_d=\langle \int \boldsymbol {v}_1\,\mbox {d}t \boldsymbol {\cdot } \boldsymbol {\nabla } \boldsymbol {v}_1\rangle$ (Riley Reference Riley2001). Using (3.5), the Stokes drift in the outer region has the general form (again using the averaging rules for products of complex oscillating quantities)

(3.13)\begin{equation} \boldsymbol{v}_d = \frac{1}{2}\left|\frac{1}{\alpha^2} \left(1-\frac{\cosh{\alpha z}}{\cosh{\alpha}}\right) \right|^2 (\mathrm{i}\,\nabla_{{\parallel}} p_1^* \boldsymbol{\cdot} \nabla_{{\parallel}} \nabla_{{\parallel}} p_1). \end{equation}

For the $p_1$ given in (3.6), the Stokes drift vanishes identically, so the Lagrangian streaming is given simply by $\boldsymbol {v}_2$. Thus the streaming velocity has zero average across the channel height ($\bar {\boldsymbol {v}}_{L} = 0$) up to leading order in $\eta$. We show in Appendix A that the leading contribution to the channel-averaged streaming occurs at $O(\eta )$ as a result of a correction to the primary pressure $p_1$ to match the flow in the inner region.

3.2. Flow structure

The variation of the streaming velocity across the height of the channel, characterized by the functions $F_{\parallel }$ and $F_z$, is shown in figure 2 for different values of $\delta$. The horizontal flow velocity is characterized by multiple flow reversals across the channel height (changes in the sign of $F_{\parallel }$). Horizontal velocity profiles for $\delta \lesssim 0.1$ have qualitatively similar shapes. Near the top and bottom walls, the horizontal flow increases and decreases rapidly on the scale of the Stokes layer thickness $\delta$ before changing direction (figure 2a). At the outer edge of the Stokes layers the flow achieves a maximal velocity (the maximal negative value of $F_{\parallel }$), which we identify with the time-averaged slip over walls at the edge of thin Stokes layers (Longuet-Higgins Reference Longuet-Higgins1953; Nyborg Reference Nyborg1958). The central core outside the Stokes layers is characterized by a more gradual variation of the flow, which reverses direction another time as one moves towards the central plane ($z=0$) of the channel. For $\delta \gtrsim 0.13$, the second flow reversal near $z = 0$ disappears (see profile for $\delta = 0.2$ in figure 2a). For $\delta$ larger than approximately $0.25$, the shape of the flow remains roughly the same and is characterized by a single flow reversal with a maximum at the central plane where $z=0$. For $\delta \ll 1$, the central core (outside Stokes layers) is described by the parabolic profile with slip, $F_{\parallel } \sim (1/16)(1-3z^2)$, whereas for $\delta \gg 1$, the flow profile is described by $F_{\parallel } \sim (3 \delta ^{-2}/560) (z^2 - 1)(7 z^4 - 28 z^2 +5)$.

Figure 2.

Analytical description of streaming velocity: (a) horizontal and (b) vertical velocity profiles for different values of Stokes layer thickness $\delta$. The flow reverses direction multiple times across the height of the channel depending on the value of $\delta$. (c) The horizontal flow at the central plane $z=0$ reverses direction across $\delta \approx 0.13$. Dashed curves are asymptotic expansions for small and large $\delta$.

The vertical velocity, shown in figure 2(b), exhibits odd symmetry about $z = 0$ but undergoes similar changes in sign as $\delta$ increases. When $\delta$ is small, $F_z$ is mostly negative in the top half of the channel, with a reversal in direction within the Stokes layer. By contrast, for $\delta$ larger than approximately $0.13$, the vertical velocity reverses and $F_z$ is positive throughout the top half of the channel.

Figure 2(c) shows the behaviour of the horizontal velocity at the central plane $z=0$. As $\delta$ increases, the horizontal velocity decreases and then increases rapidly to reach a maximum (negative) value at approximately $\delta =0.35$, with a direction change occurring at $\delta \approx 0.13$. The orange dashed curve represents asymptotic expansions of horizontal central plane velocity function $F_{\parallel }(z = 0) \sim (2-15\delta )/32$ for small $\delta$, whereas the blue dashed curve shows the behaviour of $F_{\parallel }(z = 0) \sim -3/(112\delta ^2)$ for large $\delta$.

In addition to showing velocity distributions, we visualize the 3-D flow using streamline portraits in two-dimensional sections. Streamline plots in different planes – ‘length–width’ (${xy}$), ‘length–depth’ (${xz}$) and ‘width–depth’ (${yz}$) – are shown in figure 3, with the oscillation direction and coordinates marked.

Figure 3.

Theoretical streamline plots for $\delta =1$, with oscillations driven along $x$. (a) Streamlines in the horizontal plane show a vortex structure that reverses direction across $z$ (blue and red streamlines are at $z = 0$ and $z = 0.7$, respectively). The reversal of flow is made clearer by streamlines in the vertical planes (b) and (c).

Figure 3(a) shows streamlines around a cylinder (shaded grey) on the horizontal plane but at different depths for $\delta = 1$. Blue streamlines represent the flow at the axial mid-plane ($z=0$), showing a structure qualitatively similar to two-dimensional streaming. However, streamlines at a depths $z=0.7$, shown in red, indicate flow in the opposite direction. This reversal of flow direction (associated with changes in the sign of $F_{\parallel }$) is a feature that is absent not only from two-dimensional theory (which ignores variations along $z$), but also from previous theory (Shih Reference Shih1970) and experiments (Costalonga et al. Reference Costalonga, Brunet and Peerhossaini2015) in the Hele-Shaw limit. It is interesting to note that previous descriptions of 3-D streaming flows under less extreme vertical confinement (of rigid cylinders and semi-cylindrical bubbles) do not report such a flow reversal either (Lutz et al. Reference Lutz, Chen and Schwartz2005; Marin et al. Reference Marin, Rossi, Rallabandi, Wang, Hilgenfeldt and Kähler2015). In § 4, we confirm that such a flow reversal with $z$ does indeed occur in experiments in Hele-Shaw geometries.

Figures 3(b) and 3(c) indicate flows in vertical planes and also show that the direction of secondary flow reverses across the depth of the channel. Unlike figure 3(b), we see in figure 3(c) that the flow appears to locally penetrate the cylinder surface, even though the depth-averaged velocity normal to the cylinder surface remains identically zero by construction. This is understood by recalling that the present theory is restricted to the ‘outer region’ where $r - 1 \gg \eta$. To accommodate the local velocity normal to the cylinder at $r=1$, it is necessary to acknowledge the presence of an inner region of width $\eta$ around the cylinder that ensures that the streaming velocity vanishes everywhere on the cylinder surface. Though we do not resolve the details of this inner region here, we provide some estimates of its effects in Appendix A.

4.1. Visualization

Matlab and ImageJ are used to process the recorded images of the microparticles as they are transported by the flow. The microparticles are expected to behave as passive tracers of the flow due to their small size ($a_p/a \lesssim 0.005$, $a_p/h \lesssim 0.016$) and inertia (Stokes number $a_p^2 \omega /\nu \lesssim 0.01$). We extract particle traces from the experimental images by using an in-house implementation of the Flowtrace algorithm (Gilpin, Prakash & Prakash Reference Gilpin, Prakash and Prakash2017). Figure 4, with the oscillation direction represented by pink arrows, indicates steady pathlines of tracers around a cylinder of radius $a = 600\ \mathrm {\mu }$m, but with different values of $h$ and driving frequency $f$. The flow structure is robust across our experiments: four steady vortices are created adjacent to the cylinder. Figure 4 shows these vortices for two representative experiments (one pair of vortices per experiment). The flow structure projected in the imaging $xy$ plane resembles two-dimensional streaming around a cylinder in an unbounded or concentrically bounded fluid (Bertelsen et al. Reference Bertelsen, Svardal and Tjøtta1973), and is qualitatively similar to the theoretical streamline portrait of figure 3(a).

Figure 4.

Pathlines of tracer particles around cylinders of radius $a = 600\ \mathrm {\mu }$m for different experimental conditions; double arrows indicate the oscillation direction. Steady pathlines for two experiments with different values of channel height $h$ and frequency $f$. The qualitative flow pattern is robust across our experiments. (a) $h = 240\ {\rm \mu}{\rm m},\ f = 100\ {\rm Hz}\ (\eta = 0.4,\ \delta = 0.23)$; (b) $h = 80\ {\rm \mu}{\rm m},\ f = 300\ {\rm Hz}\ (\eta = 0.13,\ \delta = 0.41)$.

A closer look at the trajectories for small $\eta$, however, reveals a striking difference from two-dimensional streaming. Figure 5(a) shows trajectories of particles in an experiment with $\eta \approx 0.13$ and $\delta \approx 0.41$ (cf. figure 4b), digitally processed such that the brighter part of each trace represents the current position of a particle (head of the trace), and the darker part represents the past position of the same particle (tail of the trace); see Gilpin et al. (Reference Gilpin, Prakash and Prakash2017). This leads to the appearance of thicker heads and thinner tails of the particle tracers (note that this is not an optical effect but a digital one that accounts for temporal information). We observe pairs of apparently nearby tracers that are clearly travelling in opposite directions. We note that this feature is not localized and can be observed everywhere in the flow; see the supplementary movie available at https://doi.org/10.1017/jfm.2023.276. Zooming into one of the marked regions (insets in figure 5(a)) makes the counter-motion of nearby particles clearer. The theory rationalizes this observation as the consequence of the particles being at different $z$ locations in the channel with a reversed flow direction, as shown in figure 5(b); see also figure 3.

Figure 5.

Streaming in a Hele-Shaw cell with $a=600\ \mathrm {\mu }$m, $2h=160 \mathrm {\mu }$m at $f = 300$ Hz ($\eta = 0.13$, $\delta =0.41$). (a) Particle traces, now including temporal information: for each trace, the brighter ‘head’ shows the particle's present position, whereas the darker ‘tail’ represents the particle's past position. The insets detail different parts of the flow where pairs of tracer particles are travelling in opposite directions, but at different depths (yellow arrows indicate the direction of motion to guide the eye). (b) Theoretically predicted streamlines near the central plane ($z=0.19$, blue dashed) and near the boundary wall ($z=0.95$, red solid). (c,d) Experiments with the same parameters as (a) but isolating particles at specific $z$ depths (see main text). Particles at $r\approx 1.6$, $\theta \approx 30^{\circ }$ but at different depths, (c) $z=0.19$ and (d) $z=0.95$, move in opposite directions, consistent with the theoretical streamline portraits in (b).

Figures 5(c) and 5(d) provide direct experimental confirmation of this reversal of flow with $z$. These figures correspond to a single experiment conducted with the same parameters as in figure 5(a) but observed under a higher magnification to obtain a shallower depth of focus, which enables the isolation of particles at different $z$ locations. We also utilize somewhat larger tracer particles ($2 a_p = 10.0\ \mathrm {\mu }$m, density 1.05 g cm$^{-3}$) to aid visualization. The two particles in figures 5(c) and 5(d) are captured at positions $r\approx 1.6$, $\theta \approx 30^{\circ }$ but at different locations along the channel height, $z=0.19$ and $0.95$, respectively. The direction of flow is clearly different at the two $z$ locations, and is consistent with the directions predicted theoretically in figure 5(b) at those $z$ locations. Furthermore, we see experimentally that the flow is faster at $z = 0.19$ (near the channel centre) than at $z = 0.95$ (near the channel wall) by a factor ${\approx }2.8$, which is also consistent with the theory, which predicts a factor ${\approx }2.4$. Thus both the vortex structure and direction of flow are consistent between experiments and theory for narrow channels.

4.2. Velocity decay

Our streaming theory predicts a 3-D velocity field $\boldsymbol {v}_2=\boldsymbol {v}_{2 \parallel }+\eta v_{2z}\, {\boldsymbol {e}}_z$. Substituting $|\bar {\boldsymbol {v}}_1|^2=1+1/r^4-2 \cos {(2 \theta )}/r^2$ from (3.5) into (3.12), we obtain

(4.1)\begin{equation} \boldsymbol{v}_2=(v_{r2}, v_{\theta 2}, \eta v_{2z}) = \left(\frac{4 F_{{\parallel}}}{r^5}\,(r^2\cos{2\theta}-1), \frac{4F_{{\parallel}}}{r^3}\sin{2\theta}, \frac{16 \eta F_z}{r^6}\right). \end{equation}

As is typical for streaming flows, the flow speed decays with the distance from the obstacle. Two-dimensional theory predicts that the velocity decays with the square of the distance to the centre of the cylinder: $|\boldsymbol {v}_2| \propto r^{-2}$ (Bertelsen et al. Reference Bertelsen, Svardal and Tjøtta1973). However, observe from the above solution that horizontal velocities are predicted to have a dominant far-field decay of $|\boldsymbol {v}_2| \propto r^{-3}$.

We analyse experimental streaming velocities using (an adaptation by Blair & Dufresne (Reference Blair and Dufresne2008) of) the particle-tracking algorithms of Crocker & Grier (Reference Crocker and Grier1996). From the experimental videos, we calculate the speeds of microparticles, which we approximate by the magnitude of the horizontal velocity, as the vertical velocity is smaller. The results of our analysis are shown in figure 6. We first identify particles in a narrow sector, concentric with the cylinder, of angular width $2\delta \theta$ about a mean angle $\theta$ (we pick $\theta = {\rm \pi}/4$, $\delta \theta =0.1$). We subdivide this sector into multiple radial sections (typically 30) starting from the cylinder surface. We collect hundreds of particle trajectories in the sector of interest, and calculate their positions and speeds as they pass through the angle $\theta$ (light grey markers in figure 6). We then average the particle positions and speeds within each radial section of the sector, which we identify with the mean speed at the projected averaged locations of the particles $v(r, \theta )$ (coloured markers in figure 6). Deviations from the average thus represent the spread of the individual particle positions and speeds. The maximum (bin-averaged) speed for each experiment is denoted $v_{m}$, and its radial location is denoted $r_{m}$.

Figure 6.

Decay of normalized horizontal speed $v/v_{m}$ with dimensionless distance from the centre of the cylinder $r = R/a$, showing experiments (symbols) and the theoretical prediction (dashed line). The maximal speed for each experiment is $v_{m}$. Light grey symbols are velocity data from individual particles at $\theta = {\rm \pi}/4 \pm 0.1$, whereas coloured symbols are data after a local averaging. (a) Data for $v/v_{m}$ from different experiments decay approximately as $r^{-3}$ but do not overlap due to different aspect ratio $\eta$ and Stokes layer thickness $\delta$. (b) Rescaling velocities as suggested by the present theory causes data from different experiments to collapse onto a universal curve $v/(v_{m} r_{m}^3) \approx 1.1 r^{-3}$. The inset shows the same data on a linear scale.

In figure 6, we plot the ratio of (bin-averaged) speed in each segment and maximal speed among all particles, as a function of the dimensionless distance to the cylinder axis. The experimental data are collected from three characteristic experiments with cylinders confined axially in microchannels: (i) with height $2h=480\ \mathrm {\mu }$m and cylinder radius $a=600\ \mathrm {\mu }$m, operating at frequency $f = 100$ Hz ($\eta =0.4$, $\delta = 0.23$; see panel (a) of figure 4); (ii) $2h=480\ \mathrm {\mu }$m, $a=750\ \mathrm {\mu }$m, driven at $f = 50$ Hz ($\eta =0.32$, $\delta = 0.33$); (iii) $2h=160\ \mathrm {\mu }$m, $a=600\ \mathrm {\mu }$m, operating at $f = 300$ Hz ($\eta =0.13$, $\delta = 0.4$; see panel (b) of figure 4 and figure 5a). As shown in figure 5(a) and discussed earlier, the large depth of field picks up particles at different depths, encompassing a wide range of speeds (including speeds close to zero at depths where flow reverses); our data set includes all of these trajectories.

Figure 6(a) shows the experimental data for $v/v_{m}$ versus $r$. As our theory predicts, the streaming velocity is found to decay as the inverse cube of the distance from the cylinder ($|\boldsymbol {v}_2|\propto r^{-3}$). Parameters $\epsilon$, $\eta$ and $\delta$ vary from experiment to experiment, resulting in different $v_{m}$ and $r_{m}$. We observe that the theory predicts $v_m = r_m = 1$ at $\eta \rightarrow 0$, which is consistent with the experiment for the smallest $\eta$.

Although the theory is formally accurate only for small $\eta$, we find that it nonetheless provides a useful way to organize all of the data onto a universal curve; see figure 6(b). We assume that the maximum binned mean speed from experimental data exists in the outer region, even for moderate $\eta$. The theory predicts a velocity decay behaviour $v\propto r^{-3}$ (with a prefactor that generally depends on $\delta$ and $\eta$), so we anticipate that the maximal velocity satisfies $v_{m}\propto r^{-3}_{m}$, with the same prefactor. Rearranging the equation yields a universal curve $v/(v_{m} r^3_{m}) = 1/r^3$. Applying this relationship, we find that the rescaled experimental data in figure 6(b) closely follow this universal curve across two orders of magnitude, though the agreement is better with an additional prefactor $1.1$.