2.1 Flow configuration

Figure 1.

Schematic of the flow configuration. A thermally stratified ideal gas is confined between the plane parallel walls of a long thin channel of height $H_{\ast }$ . A standing acoustic wave of wavelength $2\unicode[STIX]{x03C0}/k_{\ast }$ interacts with the thermal stratification to drive a time mean, or streaming, flow that is sufficiently strong to modify the wave dynamics. The thermal driving is imposed by fixing the temperature $\tilde{T}$ of the lower wall to be $T_{\ast }$ and that of the upper wall to be $T_{\ast }+\unicode[STIX]{x0394}\unicode[STIX]{x1D6E9}_{\ast }$ .

The problem we consider is similar to that introduced in Chini et al. (Reference Chini, Malecha and Dreeben2014). Specifically, we analyse the two-dimensional flow of an ideal gas, with specific gas constant $R_{\ast }$ and constant dynamic viscosity $\unicode[STIX]{x1D707}_{\ast }$ and thermal conductivity $\unicode[STIX]{x1D705}_{\ast }$ , in a channel with walls separated by a distance $H_{\ast }$ in the ${\tilde{y}}$ coordinate direction (see figure 1). Here and throughout, tildes refer to dimensional variables, while asterisks are used to denote dimensional parameters. Subsequently, overbars will be used to designate dimensionless time-averaged fields, while primes will be reserved for dimensionless oscillatory fields. The gas is presumed to be driven in an approximately time-periodic fashion, with frequency $\unicode[STIX]{x1D714}_{\ast }$ , yielding a standing sound wave with spatial wavenumber $k_{\ast }$ . The velocity field is required to satisfy no-slip and zero normal-flow boundary conditions along the channel walls located at ${\tilde{y}}=0$ and ${\tilde{y}}=H_{\ast }$ . All dependent fields are required to satisfy a $2\unicode[STIX]{x03C0}/k_{\ast }$ periodicity condition in the horizontal ( $\tilde{x}$ ) coordinate. In addition, to fix the spatial phase of the sound wave, we impose a symmetry condition along (or, equivalently, a zero-mass-exchange condition across) $\tilde{x}=0$ ; i.e. $\tilde{u} (0,{\tilde{y}},\tilde{t})=0$ , where $\tilde{u}$ is the $\tilde{x}$ -velocity component and $\tilde{t}$ is the time variable. This additional boundary condition holds for any flow developing from initial conditions that are symmetric with respect to $\tilde{x}=0$ , e.g. for the quiescent diffusive state, since then both the initial conditions and the governing equations are invariant with respect to the transformation $\tilde{x}\rightarrow -\tilde{x}$ , $\tilde{u} \rightarrow -\tilde{u}$ and $\tilde{v}\rightarrow \tilde{v}$ .

Table 1.

Dimensional variables and parameters.

In contrast to the study of Chini et al. (Reference Chini, Malecha and Dreeben2014), the thermal driving is achieved by fixing the temperatures of the lower and upper walls to be $T_{\ast }$ and $T_{\ast }+\unicode[STIX]{x0394}\unicode[STIX]{x1D6E9}_{\ast }$ , respectively, rather than by including a volumetric heat source. For convenience, we take the temperature differential $\unicode[STIX]{x0394}\unicode[STIX]{x1D6E9}_{\ast }>0$ , but note that this restriction is not dynamically significant since we do not consider the influence of buoyancy (gravity) in this investigation. Denoting the density, pressure, temperature and velocity fields by $\tilde{\unicode[STIX]{x1D70C}}$ , $\tilde{p}$ , $\tilde{T}$ and $\tilde{\boldsymbol{u}}$ , respectively, where $\tilde{\boldsymbol{u}}=(\tilde{u} ,\tilde{v})$ and $\tilde{v}$ is the ${\tilde{y}}$ -velocity component, the governing (compressible) Navier–Stokes, continuity and energy equations and the ideal gas equation of state can be written as

(2.1) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{\unicode[STIX]{x1D70C}}[\unicode[STIX]{x2202}_{\tilde{t}}\tilde{\boldsymbol{u}}+(\tilde{\boldsymbol{u}}\boldsymbol{\cdot }\tilde{\unicode[STIX]{x1D735}})\tilde{\boldsymbol{u}}]=-\tilde{\unicode[STIX]{x1D735}}\tilde{p}+\unicode[STIX]{x1D707}_{\ast }[\tilde{\unicode[STIX]{x1D6FB}}^{2}\tilde{\boldsymbol{u}}+{\textstyle \frac{1}{3}}\tilde{\unicode[STIX]{x1D735}}(\tilde{\unicode[STIX]{x1D735}}\boldsymbol{\cdot }\tilde{\boldsymbol{u}})], & \displaystyle\end{eqnarray}$$

(2.2) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{\tilde{t}}\tilde{\unicode[STIX]{x1D70C}}+\tilde{\unicode[STIX]{x1D735}}\boldsymbol{\cdot }(\tilde{\unicode[STIX]{x1D70C}}\tilde{\boldsymbol{u}})=0, & \displaystyle\end{eqnarray}$$

(2.3) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{\unicode[STIX]{x1D70C}}c_{v}[\unicode[STIX]{x2202}_{\tilde{t}}\tilde{T}+(\tilde{\boldsymbol{u}}\boldsymbol{\cdot }\tilde{\unicode[STIX]{x1D735}})\tilde{T}]=-\tilde{p}(\tilde{\unicode[STIX]{x1D735}}\boldsymbol{\cdot }\tilde{u} )+\unicode[STIX]{x1D705}_{\ast }\tilde{\unicode[STIX]{x1D6FB}}^{2}\tilde{T}, & \displaystyle\end{eqnarray}$$

(2.4) $$\begin{eqnarray}\displaystyle & \displaystyle \tilde{p}=\tilde{\unicode[STIX]{x1D70C}}R_{\ast }\tilde{T}, & \displaystyle\end{eqnarray}$$

where the two-dimensional gradient operator $\tilde{\unicode[STIX]{x1D735}}=(\unicode[STIX]{x2202}_{\tilde{x}},\unicode[STIX]{x2202}_{{\tilde{y}}})$ . Note that dilatational (or ‘bulk’) viscosity has been neglected in (2.1), and viscous heating has been omitted in (2.3). In practice, the bulk viscosity vanishes for a monatomic gas and, according to early experiments, is smaller than the dynamic (or shear) viscosity for the specific diatomic ideal gas (i.e. nitrogen) studied here (Prangsma, Alberga & Beenakker Reference Prangsma, Alberga and Beenakker1973). More importantly, although significant variations in the dynamic viscosity may be expected owing to the temperature dependence of this coefficient, these variations are neglected to facilitate the analysis. The steady-state pressure and temperature fields in the absence of acoustic waves and streaming flow, that is for $\tilde{\boldsymbol{u}}=\mathbf{0}$ , are referred to as the background fields (denoted with a subscript ‘ $B$ ’) and are found to be

(2.5a,b ) $$\begin{eqnarray}\displaystyle \tilde{T}_{B}=T_{\ast }\left(1+\unicode[STIX]{x1D6E4}\frac{y_{\ast }}{H_{\ast }}\right),\quad \tilde{p}_{B}=p_{\ast }, & & \displaystyle\end{eqnarray}$$

where $p_{\ast }$ is a constant and the dimensionless temperature differential $\unicode[STIX]{x1D6E4}\equiv \unicode[STIX]{x0394}\unicode[STIX]{x1D6E9}_{\ast }/T_{\ast }$ . Table 1 summarizes the dimensional fields and parameters used in the following analysis.

2.2 Scaling and non-dimensionalization

To facilitate the asymptotic analysis, we non-dimensionalize the governing equations by scaling the dependent and independent variables as outlined in table 2. The $\tilde{x}$ -velocity component is scaled with $a_{\ast }$ , the sound speed at temperature $T_{\ast }$ . This scaling introduces into the dimensionless governing equations the Strouhal number $S=a_{\ast }/U_{\ast }$ , where $U_{\ast }$ (rather than $a_{\ast }$ ) is a characteristic oscillatory velocity induced by the standing acoustic wave. The Strouhal number is large ( $10^{3}$ or larger) in many applications, and therefore we introduce $\unicode[STIX]{x1D716}\equiv 1/S$ and consider the asymptotic limit $\unicode[STIX]{x1D716}\rightarrow 0$ with all other dimensionless parameters scaled as appropriate powers of $\unicode[STIX]{x1D716}$ . Since $\unicode[STIX]{x1D716}\ll 1$ , the leading-order acoustic-wave dynamics is linear. Nevertheless, weak wave–wave nonlinearities are crucial for acoustic streaming, as their cumulative effect can be significant over sufficiently many [ $O(1/\unicode[STIX]{x1D716})$ ] acoustic-wave periods (i.e. over the slow time scale). The implied temporal scale separation between the wave and streaming dynamics is readily achieved in both laboratory experiments and streaming-enabled technologies.

The $\tilde{x}$ and ${\tilde{y}}$ coordinates are scaled with $k_{\ast }^{-1}$ and $H_{\ast }$ , respectively, so that the gas lies in the domain defined by $x\in [0,2\unicode[STIX]{x03C0}]$ and $y\in [0,1]$ . Time $\tilde{t}$ is non-dimensionalized using the inverse reference wave frequency $\unicode[STIX]{x1D714}_{\ast }^{-1}=(a_{\ast }k_{\ast })^{-1}$ . The vertical ( ${\tilde{y}}$ ), or wall-normal, velocity component is scaled by $(k_{\ast }H_{\ast })a_{\ast }$ . The domain aspect ratio $\unicode[STIX]{x1D6FF}\equiv k_{\ast }H_{\ast }$ is assumed to be small and, more precisely, is chosen so that $\unicode[STIX]{x1D6FF}=\sqrt{\unicode[STIX]{x1D716}}h$ , where $h$ is a dimensionless parameter of order unity. Although the cross-channel heat flux is not expected to be maximized as $\unicode[STIX]{x1D6FF}\rightarrow 0$ (since the streaming flow will be largely horizontal), we follow Chini et al. (Reference Chini, Malecha and Dreeben2014) and continue to focus on the small aspect-ratio regime for the following reasons. First, as noted in the introduction, most theoretical and computational studies have been performed in this regime, so meaningful comparisons to prior investigations can be made. Secondly, the analysis of the acoustic wave is simplified in a domain that is thin relative to the wavelength of the sound wave. Indeed, the acoustic wave then is dynamically constrained to maintain its first-mode wall-normal structure. Finally, in the small aspect-ratio regime, the leading-order fluctuating pressure gradient is orthogonal to the imposed background density gradient, resulting in a crucial baroclinic contribution to the production of fluctuating vorticity.

Table 2.

Dimensionless variables and parameters.

The temperature $T_{\ast }$ of the lower wall is used to non-dimensionalize the temperature field $\tilde{T}$ . In the analysis that follows, $\unicode[STIX]{x1D6E4}$ is fixed, i.e. $O(1)$ , as $\unicode[STIX]{x1D716}\rightarrow 0$ , although the smallness of $\unicode[STIX]{x1D716}$ in acoustic streaming ensures that our multiple scale analysis remains accurate even for $\unicode[STIX]{x1D6E4}\approx 0.1$ , as will be evident in § 4, where we compare our theoretical predictions with the results of direct numerical simulations.

The Reynolds and Péclet numbers characterizing the acoustic waves are denoted $Re$ and $Pe$ , respectively. Since these parameters are very large compared to unity, both momentum and thermal diffusion can be neglected in the leading-order wave dynamics, at least in the domain interior. Note that in baroclinic acoustic streaming, typical streaming velocities also are of size $U_{\ast }$ (in contrast to Rayleigh streaming, where streaming speeds are proportional to $\unicode[STIX]{x1D716}U_{\ast }$ ). Consequently, $Re$ and $Pe$ also would appear to characterize the streaming flow. Because we consider the limit of small aspect ratio, however, diffusion in the wall-normal ( $y$ ) direction is enhanced by a factor $\unicode[STIX]{x1D6FF}^{-2}\propto \unicode[STIX]{x1D716}^{-1}$ . We therefore define the streaming Reynolds and Péclet numbers $Re_{s}=\unicode[STIX]{x1D716}Re$ and $Pe_{s}=\unicode[STIX]{x1D716}Pe$ as in traditional Rayleigh streaming as the proper measure of mean inertia to the dominant mean diffusive effects.

2.3 Asymptotic analysis

Using the scalings described above (and summarized in table 2), the governing equations and boundary conditions can be recast in dimensionless form. The occurrence of the small parameter $\unicode[STIX]{x1D716}$ in the dimensionless system prompts a multiple scale asymptotic analysis in which the single time variable $t$ characterizing the fast dynamics of the acoustic waves is augmented with a slow-time variable $T\equiv \unicode[STIX]{x1D716}t$ to capture the cumulative effect of weakly nonlinear wave dynamics that ultimately drives streaming. Furthermore, we posit the following asymptotic expansions for the various fields:

(2.6) $$\begin{eqnarray}\displaystyle & \displaystyle (u,v)=\unicode[STIX]{x1D716}(u_{1},v_{1})+\unicode[STIX]{x1D716}^{2}(u_{2},v_{2})+O(\unicode[STIX]{x1D716}^{3}), & \displaystyle\end{eqnarray}$$

(2.7) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x03C0}=\unicode[STIX]{x1D716}\unicode[STIX]{x03C0}_{1}+\unicode[STIX]{x1D716}^{2}\unicode[STIX]{x03C0}_{2}+O(\unicode[STIX]{x1D716}^{3}), & \displaystyle\end{eqnarray}$$

(2.8) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E9}=\unicode[STIX]{x1D6E9}_{0}+\unicode[STIX]{x1D716}\unicode[STIX]{x1D6E9}_{1}+O(\unicode[STIX]{x1D716}^{2}), & \displaystyle\end{eqnarray}$$

(2.9) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70C}=\unicode[STIX]{x1D70C}_{0}+\unicode[STIX]{x1D716}\unicode[STIX]{x1D70C}_{1}+O(\unicode[STIX]{x1D716}^{2}), & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x03C0}=(\,\tilde{p}-\tilde{P}_{B})/p_{\ast }$ and $\unicode[STIX]{x1D6E9}=(\tilde{T}-\tilde{T}_{B})/T_{\ast }$ are, respectively, the dimensionless perturbation pressure and temperature fields. Note that, at each order in these expansions, the field variables can have both fluctuating and mean components, which subsequently will be disentangled via the introduction of a fast-time averaging operation. Expansion (2.6) follows from the scaling of the Strouhal number ( $S=1/\unicode[STIX]{x1D716}$ ). The state equation constrains the temperature and density perturbations to be of the same order. In baroclinic acoustic streaming, the $O(\unicode[STIX]{x1D716})$ streaming flow is sufficiently strong to induce $O(1)$ rearrangements of the background temperature and density fields over an $O(1/\unicode[STIX]{x1D716})$ time period. Thus, in contrast to other studies of acoustic streaming in the presence of inhomogeneous temperature fields (e.g. Červenka & Bednarřik Reference Červenka and Bednarřik2017), it is crucial that the expansions for both $\unicode[STIX]{x1D6E9}$ and $\unicode[STIX]{x1D70C}$ begin at $O(1)$ , as first shown in Chini et al. (Reference Chini, Malecha and Dreeben2014).

Owing to the large perturbations to the background density field, the natural frequency of the acoustic mode may evolve in time. Here, we extend the analysis of Chini et al. (Reference Chini, Malecha and Dreeben2014) by employing a WKBJ approximation to properly capture this slow temporal evolution. Specifically, a generic dependent field $f(x,y,t)$ is re-expressed as $f(x,y,\unicode[STIX]{x1D719},T)$ , where $\unicode[STIX]{x1D719}$ and $T$ are treated as independent variables. The rapidly varying phase $\unicode[STIX]{x1D719}$ may be written as

(2.10) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D719}(t)=\frac{\unicode[STIX]{x1D6F7}(T)}{\unicode[STIX]{x1D716}}, & & \displaystyle\end{eqnarray}$$

where $\text{d}\unicode[STIX]{x1D6F7}/\text{d}T$ is of order unity. We define the instantaneous angular frequency $\unicode[STIX]{x1D714}(T)$ ,

(2.11) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D714}(T)=\frac{\text{d}\unicode[STIX]{x1D719}}{\text{d}t}=\frac{\text{d}\unicode[STIX]{x1D6F7}}{\text{d}T}, & & \displaystyle\end{eqnarray}$$

and expand $\unicode[STIX]{x1D6F7}=\unicode[STIX]{x1D6F7}_{0}+\unicode[STIX]{x1D716}\unicode[STIX]{x1D6F7}_{1}+O(\unicode[STIX]{x1D716}^{2})$ , so that

(2.12) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D714}=\unicode[STIX]{x1D714}_{0}+\unicode[STIX]{x1D716}\unicode[STIX]{x1D714}_{1}+O(\unicode[STIX]{x1D716}^{2}). & & \displaystyle\end{eqnarray}$$

Finally, in order to distinguish the streaming flow from the acoustic wave, we introduce the fast-time average of a function $f(x,y,\unicode[STIX]{x1D719},T)$ ,

(2.13) $$\begin{eqnarray}\displaystyle \bar{f}(x,y,T)=\frac{1}{2n\unicode[STIX]{x03C0}}\int _{\unicode[STIX]{x1D719}}^{\unicode[STIX]{x1D719}+2n\unicode[STIX]{x03C0}}f(x,y,s,T)\,\text{d}s & & \displaystyle\end{eqnarray}$$

for sufficiently large positive integer $n$ , so that any function can be decomposed such that

(2.14) $$\begin{eqnarray}\displaystyle f(x,y,\unicode[STIX]{x1D719},T)=\bar{f}(x,y,T)+f^{\prime }(x,y,\unicode[STIX]{x1D719},T), & & \displaystyle\end{eqnarray}$$

where $\overline{f^{\prime }}=0$ . Thus, $f^{\prime }$ represents the acoustic wave and $\bar{f}$ the streaming flow.

2.4 Leading-order multiscale wave/mean-flow interaction equations

The multiple time scale governing equations for the coupled acoustic-wave/streaming-flow system were first derived in Chini et al. (Reference Chini, Malecha and Dreeben2014) and are reproduced here for completeness. The evolution of the streaming fields is governed by the following equations:

(2.15) $$\begin{eqnarray}\displaystyle & \displaystyle \bar{\unicode[STIX]{x1D70C}}_{0}(\unicode[STIX]{x2202}_{T}\bar{u}_{1}+\bar{u}_{1}\unicode[STIX]{x2202}_{x}\bar{u}_{1}+\bar{v}_{1}\unicode[STIX]{x2202}_{y}\bar{u}_{1})=-\frac{\unicode[STIX]{x2202}_{x}\bar{\unicode[STIX]{x03C0}}_{2}}{\unicode[STIX]{x1D6FE}}-\unicode[STIX]{x2202}_{x}(\bar{\unicode[STIX]{x1D70C}}_{0}\overline{u_{1}^{\prime 2}})-\unicode[STIX]{x2202}_{y}(\bar{\unicode[STIX]{x1D70C}}_{0}\overline{u_{1}^{\prime }v_{1}^{\prime }})+\frac{\unicode[STIX]{x2202}_{yy}\bar{u}_{1}}{Re_{s}h^{2}}, & \displaystyle\end{eqnarray}$$

(2.16) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{y}\bar{\unicode[STIX]{x03C0}}_{2}=0, & \displaystyle\end{eqnarray}$$

(2.17) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{T}\bar{\unicode[STIX]{x1D70C}}_{0}+\unicode[STIX]{x2202}_{x}(\bar{\unicode[STIX]{x1D70C}}_{0}\bar{u}_{1})+\unicode[STIX]{x2202}_{y}(\bar{\unicode[STIX]{x1D70C}}_{0}\bar{v}_{1})=0, & \displaystyle\end{eqnarray}$$

(2.18) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{T}\bar{\unicode[STIX]{x1D6E9}}_{0}+\bar{u}_{1}\unicode[STIX]{x2202}_{x}\bar{\unicode[STIX]{x1D6E9}}_{0}+\bar{v}_{1}\unicode[STIX]{x2202}_{y}(\bar{\unicode[STIX]{x1D6E9}}_{0}+T_{B})=(1-\unicode[STIX]{x1D6FE})(\bar{\unicode[STIX]{x1D6E9}}_{0}+T_{B})(\unicode[STIX]{x2202}_{x}\bar{u}_{1}+\unicode[STIX]{x2202}_{y}\bar{v}_{1})+\frac{\unicode[STIX]{x1D6FE}\unicode[STIX]{x2202}_{yy}\bar{\unicode[STIX]{x1D6E9}}_{0}}{Pe_{s}h^{2}\bar{\unicode[STIX]{x1D70C}}_{0}},\qquad & \displaystyle\end{eqnarray}$$

(2.19) $$\begin{eqnarray}\displaystyle & \displaystyle \bar{\unicode[STIX]{x1D70C}}_{0}=\frac{1}{\bar{\unicode[STIX]{x1D6E9}}_{0}+T_{B}}. & \displaystyle\end{eqnarray}$$

In this slow-time system, the sole – but crucial – effect of the waves arises from the Reynolds stress divergence in the mean $x$ -momentum equation (2.15) that drives the streaming flow. For a homogeneous fluid, this wave-induced force can be offset by a pressure gradient in the bulk: a Rayleigh streaming flow is driven at next order in $\unicode[STIX]{x1D716}$ by the Reynolds stress divergence acting in the oscillatory boundary layers that arise near the no-slip channel walls. In contrast, for baroclinic acoustic streaming, the mean temperature gradient causes this wave-induced forcing to be rotational even within the bulk interior of the domain. Consequently, the bulk force that is created cannot be balanced by a mere adjustment of the mean pressure gradient and instead induces a (strong) mean flow even in the absence of diffusive boundary layers. (Further discussion of this distinction is given in § 6.)

To evaluate this force, we employ the equations governing the leading-order acoustic-wave dynamics, viz.

(2.20a,b ) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D714}_{0}\bar{\unicode[STIX]{x1D70C}_{0}}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D719}}u_{1}^{\prime }+\frac{1}{\unicode[STIX]{x1D6FE}}\unicode[STIX]{x2202}_{x}\unicode[STIX]{x03C0}_{1}^{\prime }=0,\quad \unicode[STIX]{x2202}_{y}\unicode[STIX]{x03C0}_{1}^{\prime }=0, & & \displaystyle\end{eqnarray}$$

(2.21) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D714}_{0}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D70C}_{1}^{\prime }+\unicode[STIX]{x2202}_{x}(\bar{\unicode[STIX]{x1D70C}}_{0}u_{1}^{\prime })+\unicode[STIX]{x2202}_{y}(\bar{\unicode[STIX]{x1D70C}}_{0}v_{1}^{\prime })=0, & \displaystyle\end{eqnarray}$$

(2.22) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D714}_{0}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D6E9}_{1}^{\prime }+u_{1}^{\prime }\unicode[STIX]{x2202}_{x}\bar{\unicode[STIX]{x1D6E9}}_{0}+v_{1}^{\prime }\unicode[STIX]{x2202}_{y}(\bar{\unicode[STIX]{x1D6E9}}_{0}+T_{B})+(\unicode[STIX]{x1D6FE}-1)(\bar{\unicode[STIX]{x1D6E9}}_{0}+T_{B})(\unicode[STIX]{x2202}_{x}u_{1}^{\prime }+\unicode[STIX]{x2202}_{y}v_{1}^{\prime })=0, & \displaystyle\end{eqnarray}$$

(2.23) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x03C0}_{1}^{\prime }-\unicode[STIX]{x1D70C}_{1}^{\prime }(\bar{\unicode[STIX]{x1D6E9}}_{0}+T_{B})-\bar{\unicode[STIX]{x1D70C}}_{0}\unicode[STIX]{x1D6E9}_{1}^{\prime }=0. & \displaystyle\end{eqnarray}$$

These equations describe the fast dynamics of approximately linear non-dissipative acoustic waves in a medium whose mean density field $\bar{\unicode[STIX]{x1D70C}}_{0}$ evolves slowly in time. Conversely, the evolution of $\bar{\unicode[STIX]{x1D70C}}_{0}(x,y,T)$ depends on the acoustic waves, as can be gleaned from inspection of (2.15)–(2.19). Note further that, here, unlike in Rayleigh streaming, the oscillatory Stokes layers are dynamically passive because the streaming induced by near-wall fluctuating viscous torques is $O(\unicode[STIX]{x1D716}^{2})$ while the streaming flow governed by (2.15)–(2.19) is $O(\unicode[STIX]{x1D716})$ . Therefore, the leading-order wave field is required to satisfy only a zero normal-flow boundary condition at each wall, and the details of the oscillatory flow within the Stokes (boundary) layers do not have to be determined at this order (see § 6).

Taken together, these two sets of equations form a closed but two time scale system. As emphasized in Chini et al. (Reference Chini, Malecha and Dreeben2014), the fully two-way coupling between the waves and mean flow captured by this multiscale system renders it fundamentally distinct from classical Rayleigh streaming theory, in which the acoustic-wave field can be computed first and then the response of the streaming flow to the acoustic-wave forcing self-consistently determined (i.e. one-way coupling). Subsequently, Karlsen et al. (Reference Karlsen, Augustsson and Bruus2016, Reference Karlsen, Qiu, Augustsson and Bruus2018) also noted this two-way coupling in their computational studies of acoustic streaming in inhomogeneous fluids. Nevertheless, to make analytical progress, Chini et al. (Reference Chini, Malecha and Dreeben2014) considered a small Prandtl number limit in which the coupling is effectively one way. A primary contribution of the present investigation is to extend the analysis of Chini et al. (Reference Chini, Malecha and Dreeben2014) to systematically treat fully two-way wave/mean-flow interactions. This extension requires the derivation of a novel amplitude equation governing the slow evolution of the acoustic waves, which can only be determined by carrying the asymptotic analysis to next order and imposing an appropriate solvability condition. Heuristically, the slow evolution of the waves is controlled by higher-order terms that e.g. account for energy exchanges with the streaming flow or with the solid boundaries. In the next section, we show how these effects can be self-consistently incorporated.

3.1 Mode structure

Substituting the decomposition (3.1) into the fluctuation equations (2.20)–(2.23) yields a linear but non-separable two-dimensional (partial) differential eigenvalue problem for the spatial structure and frequency $\unicode[STIX]{x1D714}_{0}$ of the leading-order fluctuation fields. By continuing to exploit the small aspect-ratio limit, we nevertheless are able to reduce the required computation to the solution of a one-dimensional eigenvalue problem – a crucial simplification.

To proceed, we note that, using (3.1), (2.21)–(2.23) can be combined to deduce

(3.2) $$\begin{eqnarray}\displaystyle \hat{\unicode[STIX]{x03C0}}_{1}=\frac{\text{i}\unicode[STIX]{x1D6FE}}{\unicode[STIX]{x1D714}_{0}}(\unicode[STIX]{x2202}_{x}\hat{u} _{1}+\unicode[STIX]{x2202}_{y}\hat{v}_{1}). & & \displaystyle\end{eqnarray}$$

With this expression for the acoustic-wave pressure $\hat{\unicode[STIX]{x03C0}}_{1}$ , the momentum equations (2.20) become two coupled partial differential equations for $\hat{u} _{1}$ and $\hat{v}_{1}$ :

(3.3) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{x}(\unicode[STIX]{x2202}_{x}\hat{u} _{1}+\unicode[STIX]{x2202}_{y}\hat{v}_{1})+\unicode[STIX]{x1D714}_{0}^{2}\bar{\unicode[STIX]{x1D70C}}_{0}\hat{u} _{1}=0, & \displaystyle\end{eqnarray}$$

(3.4) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{y}(\unicode[STIX]{x2202}_{x}\hat{u} _{1}+\unicode[STIX]{x2202}_{y}\hat{v}_{1})=0. & \displaystyle\end{eqnarray}$$

A further reduction to a single ordinary differential equation is possible by formally integrating a linear combination of these equations,

(3.5)

where $q$ is an unknown function of $x$ and $T$ only. Integration of (3.4) gives

(3.6) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x2202}_{x}\hat{u} _{1}+\unicode[STIX]{x2202}_{y}\hat{v}_{1}=g, & & \displaystyle\end{eqnarray}$$

where $g$ is a second unknown function of $x$ and $T$ only. Equations (3.3), (3.5) and (3.6) imply that $q=-g^{\prime }/\unicode[STIX]{x1D714}_{0}^{2}$ , where a prime is used to denote partial differentiation of the function $g$ with respect to $x$ , since the $T$ dependence is parametric. The general solution of this system of equations can be obtained using the kinematic boundary condition $\hat{v}_{1}(x,y=0,T)=0$ , viz.

(3.7) $$\begin{eqnarray}\displaystyle & \displaystyle \hat{u} _{1}=-\frac{g^{\prime }}{\unicode[STIX]{x1D714}_{0}^{2}\bar{\unicode[STIX]{x1D70C}}_{0}}, & \displaystyle\end{eqnarray}$$

(3.8) $$\begin{eqnarray}\displaystyle & \displaystyle \hat{v}_{1}=gy+\unicode[STIX]{x2202}_{x}\left(\frac{g^{\prime }}{\unicode[STIX]{x1D714}_{0}^{2}}\int _{0}^{y}\frac{\text{d}y}{\bar{\unicode[STIX]{x1D70C}}_{0}}\right). & \displaystyle\end{eqnarray}$$

Finally, the upper boundary condition $\hat{v}_{1}(x,y=1,T)=0$ provides a constraint on $g$ and $\unicode[STIX]{x1D714}_{0}$ in the form of the ordinary differential eigenvalue problem

(3.9) $$\begin{eqnarray}\displaystyle \frac{\text{d}}{\text{d}x}(g^{\prime }\unicode[STIX]{x1D6FC})+\unicode[STIX]{x1D714}_{0}^{2}g=0, & & \displaystyle\end{eqnarray}$$

where

(3.10) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6FC}(x,T)=\int _{0}^{1}\frac{\text{d}y}{\bar{\unicode[STIX]{x1D70C}}_{0}}. & & \displaystyle\end{eqnarray}$$

To characterize the function $g(x,T)$ , we first take, without loss of generality, $\hat{u} _{1}$ to be a real field: (3.3) and (3.4) then imply that $\hat{v}_{1}$ and, thence, $g$ are also real-valued fields. Moreover, in order to ensure $\hat{u} _{1}(x=0,y,T)=0$ , the ordinary differential equation (3.9) must be solved subject to the boundary conditions $g^{\prime }(0)=g^{\prime }(2\unicode[STIX]{x03C0})=0$ . This requirement leads to an orthogonality condition: let $(g_{A},g_{B})$ be two eigenvectors and $(\unicode[STIX]{x1D714}_{A},\unicode[STIX]{x1D714}_{B})$ their angular eigenfrequencies; then

(3.11) $$\begin{eqnarray}\displaystyle \int _{0}^{2\unicode[STIX]{x03C0}}g_{A}(x)g_{B}(x)\,\text{d}x=\frac{1}{\unicode[STIX]{x1D714}_{A}^{2}-\unicode[STIX]{x1D714}_{B}^{2}}[\unicode[STIX]{x1D6FC}(g_{A}g_{B}^{\prime }-g_{B}g_{A}^{\prime })]_{0}^{2\unicode[STIX]{x03C0}}=0. & & \displaystyle\end{eqnarray}$$

Equation (3.11) provides a convenient scalar product on eigenvectors, and therefore we normalize them according to

(3.12) $$\begin{eqnarray}\displaystyle \int _{0}^{2\unicode[STIX]{x03C0}}g(x)^{2}\,\text{d}x=1. & & \displaystyle\end{eqnarray}$$

This normalization condition resolves the ambiguity in the definition of $A(T)$ and $\hat{f}_{1}$ in (3.1).

3.2 Wave amplitude

As explained in the previous subsection, the acoustic mode shape $g(x,T)$ and frequency $\unicode[STIX]{x1D714}_{0}(T)$ can be computed at every time $T$ for a given mean density profile $\bar{\unicode[STIX]{x1D70C}}_{0}(x,y,T)$ . The advection of hot or cold gas by the streaming flow will cause both $g$ and $\unicode[STIX]{x1D714}_{0}$ to evolve on the slow time scale and induce a two-way coupling between the waves and the streaming flow. The amplitude $A(T)$ of the acoustic mode is also expected to evolve on this slow time scale owing to dissipation by viscosity, energy exchanges with the streaming flow and walls and/or external forcing. To obtain an evolution equation for $A(T)$ , we proceed as follows, relegating details to the appendices.

(i) We collect terms in the dimensionless governing equations at $O(\unicode[STIX]{x1D716}^{2})$ . The resulting equations are reported in appendix A.

(ii) We make the ansatz that generic $O(\unicode[STIX]{x1D716}^{2})$ fluctuation field $f_{2}^{\prime }$ can be represented as

(3.13) $$\begin{eqnarray}\displaystyle f_{2}^{\prime }(x,y,T,\unicode[STIX]{x1D719})=\frac{B(T)}{2}\left(\,\hat{f}_{2}(x,y,T)\text{e}^{\text{i}\unicode[STIX]{x1D719}}+\text{c.c.}\right) & & \displaystyle\end{eqnarray}$$

and, in direct analogy with the manipulations performed in § 3.1, reduce the $O(\unicode[STIX]{x1D716}^{2})$ fluctuation system for the five unknown fields $\hat{u} _{2}$ , $\hat{v}_{2}$ , $\hat{\unicode[STIX]{x03C0}}_{2}$ , $\hat{\unicode[STIX]{x1D70C}}_{2}$ and $\hat{\unicode[STIX]{x1D6E9}}_{2}$ to a system of two equations for the two fields $\hat{u} _{2}$ and $\hat{v}_{2}$ . (Harmonics of the form $\text{e}^{\text{i}2\unicode[STIX]{x1D719}}$ also exist at this order, but are non-resonant and thus do not contribute to the slow-time dynamics; accordingly, these harmonics need not be explicitly computed.) The resulting system of equations has the form

(3.14) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{x}(\unicode[STIX]{x2202}_{x}\hat{u} _{2}+\unicode[STIX]{x2202}_{y}\hat{v}_{2})+\unicode[STIX]{x1D714}_{0}^{2}\bar{\unicode[STIX]{x1D70C}}_{0}\hat{u} _{2}={\mathcal{F}}, & \displaystyle\end{eqnarray}$$

(3.15) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{y}(\unicode[STIX]{x2202}_{x}\hat{u} _{2}+\unicode[STIX]{x2202}_{y}\hat{v}_{2})={\mathcal{G}}. & \displaystyle\end{eqnarray}$$

The linear operator acting on the left-hand side of this system is identical to that arising in the leading-order fluctuation equations (3.3)–(3.4). The right-hand side functions ${\mathcal{F}}$ and ${\mathcal{G}}$ include resonant forcing terms involving the leading-order fluctuations fields ( $f_{1}^{\prime }$ ). Analytical expressions for the imaginary parts of ${\mathcal{F}}$ and ${\mathcal{G}}$ , which are needed for the derivation of the (real) amplitude equation, are given in appendix B.

(iii) We derive a solvability condition for the $O(\unicode[STIX]{x1D716}^{2})$ system, which requires determination of the adjoint linear operator, by invoking the Fredholm alternative theorem in the usual manner; see appendix B.

(iv) We enforce the solvability condition to obtain an equation for $\text{d}A/\text{d}T$ ; cf. appendix D.

Employing this procedure, we derive the following novel amplitude equation:

(3.16) $$\begin{eqnarray}\displaystyle \frac{2}{A\unicode[STIX]{x1D714}_{0}^{-1}}\frac{\text{d}(A\unicode[STIX]{x1D714}_{0}^{-1})}{\text{d}T} & = & \displaystyle -\frac{\text{i}\unicode[STIX]{x1D714}_{0}}{Pe_{s}h^{2}}\int _{D}\,\text{d}x\,\text{d}yg\unicode[STIX]{x2202}_{yy}\hat{\unicode[STIX]{x1D6E9}}_{1}-\frac{1}{2\unicode[STIX]{x1D714}_{0}^{2}}\int _{D}\,\text{d}x\,\text{d}yg^{\prime 2}(\bar{u}_{1}\unicode[STIX]{x2202}_{x}\bar{\unicode[STIX]{x1D70C}}_{0}^{-1}+\bar{v}_{1}\unicode[STIX]{x2202}_{y}\bar{\unicode[STIX]{x1D70C}}_{0}^{-1})\nonumber\\ \displaystyle & & \displaystyle +\,\int _{D}\,\text{d}x\,\text{d}y(\unicode[STIX]{x2202}_{x}\bar{u}_{1}+\unicode[STIX]{x2202}_{y}\bar{v}_{1})\left[(1-\unicode[STIX]{x1D6FE})g^{2}+\frac{g^{\prime 2}}{\unicode[STIX]{x1D714}_{0}^{2}\bar{\unicode[STIX]{x1D70C}}_{0}}\left(\frac{Pe_{s}}{Re_{s}}-\frac{1}{2}\right)\right],\end{eqnarray}$$

where $\int _{D}$ refers to definite integration over the spatial domain. (Note that, since the temperature and velocity fluctuations are out of phase, $\hat{\unicode[STIX]{x1D6E9}}_{1}$ is strictly imaginary.) The lack of terms nonlinear in $A$ in (3.16) confirms that phenomena such as shock-wave formation or harmonic generation are sub-dominant dynamical processes in the given parameter regime relative to, for example, heat exchange with the boundaries or energy exchange with the evolving stratified environment. In particular, the first term on the right-hand side of (3.16) accounts for the time-mean heat transfer between the waves and the boundaries. In § 4.5, we demonstrate that this term can be positive: in this scenario, the waves are driven by a process that is loosely akin to the classical thermoacoustic instability in which acoustic waves in a channel can be excited when a temperature gradient is imposed along the channel walls (Swift Reference Swift1988). Owing to the occurrence of the mean fields $\bar{u}_{1}$ , $\bar{v}_{1}$ and $\bar{\unicode[STIX]{x1D70C}}_{0}$ in (3.16), the amplitude equation is, in fact, nonlinear. A distinguishing feature of (3.16) is that, unlike amplitude equations derived in numerous other contexts, determination of the coefficients requires evaluation of functional derivatives that capture the $O(1)$ variations in generic eigenfunction $\hat{f}_{1}$ caused by changes in the mean density field $\bar{\unicode[STIX]{x1D70C}}_{0}$ that occur on the slow time $T$ ; see appendices B and D for details.

The quantity $A\unicode[STIX]{x1D714}_{0}^{-1}$ on the left-hand side of (3.16) is proportional to the square root of the dimensionless energy of the acoustic wave. Indeed, the leading-order dimensionless kinetic energy $E_{K}$ of the acoustic wave averaged over the fast time scale is given by

(3.17) $$\begin{eqnarray}\displaystyle \overline{E}_{K}\equiv \frac{1}{2}\int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}x\int _{0}^{1}\,\text{d}y\bar{\unicode[STIX]{x1D70C}}_{0}\overline{u_{1}^{\prime 2}}=\frac{A(T)^{2}}{4}\int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}x\int _{0}^{1}\,\text{d}y\bar{\unicode[STIX]{x1D70C}}_{0}\hat{u} _{1}^{2}. & & \displaystyle\end{eqnarray}$$

This expression can be evaluated using (3.7) to obtain

(3.18) $$\begin{eqnarray}\displaystyle \overline{E}_{K}=\frac{A(T)^{2}}{4\unicode[STIX]{x1D714}_{0}^{4}}\int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}x\int _{0}^{1}\,\text{d}y\frac{g^{\prime 2}}{\bar{\unicode[STIX]{x1D70C}}_{0}}=\frac{A(T)^{2}}{4\unicode[STIX]{x1D714}_{0}^{4}}\int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}xg^{\prime 2}\unicode[STIX]{x1D6FC}. & & \displaystyle\end{eqnarray}$$

The last integral reduces to $\unicode[STIX]{x1D714}_{0}^{2}$ following an integration by parts and utilization of the differential equation (3.9) and the normalization condition (3.12), yielding

(3.19) $$\begin{eqnarray}\displaystyle \overline{E}_{K}=\left(\frac{A(T)}{2\unicode[STIX]{x1D714}_{0}}\right)^{2}. & & \displaystyle\end{eqnarray}$$

The amplitude equation therefore can be interpreted as an energy balance for the acoustic wave, with the left-hand side of (3.16) equalling $(1/\overline{E}_{K})\text{d}\overline{E}_{K}/\text{d}T$ .

An important outcome of this study is that even allowing for two-way wave/mean-flow coupling, the WKBJ analysis enables the streaming flow to be computed without evolving the sound waves over the fast time scale. (Of course, in the absence of two-way coupling, as in classical Rayleigh streaming, this averaging is trivial.) Instead, the evolving spatial structure of the waves can be computed at every coarse time step (in a numerical simulation) by solving the one-dimensional eigenvalue problem (3.9), while the evolution of the amplitude can be determined by integrating (3.16) over the slow time scale. These computations are performed in conjunction with the numerical solution of the streaming equations (2.15)–(2.19).

4.1 Acoustic waves

In this section, we assume that the mean density profile varies little and thus can be accurately approximated by the diffusive solution (2.5), i.e.

(4.1) $$\begin{eqnarray}\displaystyle \bar{\unicode[STIX]{x1D70C}}_{0}=\frac{1}{T_{B}}=\frac{1}{1+\unicode[STIX]{x1D6E4}y}. & & \displaystyle\end{eqnarray}$$

The coefficient $\unicode[STIX]{x1D6FC}$ defined in (3.10) then does not depend on $x$ and is simply equal to $1+\unicode[STIX]{x1D6E4}/2$ . The acoustic-wave eigenfunction $g$ characterizing the shape of the acoustic mode follows from the solution of the second-order differential equation (3.9) and, with the prescribed boundary conditions $g^{\prime }(2\unicode[STIX]{x03C0})=g^{\prime }(1)=0$ and the normalization condition (3.12), is given by

(4.2) $$\begin{eqnarray}\displaystyle g(x)=\frac{\cos (nx)}{\sqrt{\unicode[STIX]{x03C0}}}, & & \displaystyle\end{eqnarray}$$

where the integer $n$ is set to unity for this study. The angular frequency of this mode $\unicode[STIX]{x1D714}_{0}=\sqrt{1+\unicode[STIX]{x1D6E4}/2}$ , and the velocity field ( $\hat{u} _{1}$ , $\hat{v}_{1}$ ) follows from (3.7)–(3.8) and reduces to

(4.3) $$\begin{eqnarray}\displaystyle & \displaystyle \hat{u} _{1}=\frac{(1+\unicode[STIX]{x1D6E4}y)\sin (x)}{(1+\unicode[STIX]{x1D6E4}/2)\sqrt{\unicode[STIX]{x03C0}}}, & \displaystyle\end{eqnarray}$$

(4.4) $$\begin{eqnarray}\displaystyle & \displaystyle \hat{v}_{1}=\frac{\unicode[STIX]{x1D6E4}y(1-y)\cos (x)}{2(1+\unicode[STIX]{x1D6E4}/2)\sqrt{\unicode[STIX]{x03C0}}}. & \displaystyle\end{eqnarray}$$

Figure 2.

Standing acoustic-wave velocity field corresponding to the first (i.e. $n=1$ ) eigenmode for $\unicode[STIX]{x1D6E4}=1$ (upper) and $\unicode[STIX]{x1D6E4}=0.2$ (lower), with $\unicode[STIX]{x1D6FE}=1.4$ . Particularly in (a), it is evident that the wall-normal ( $y$ ) velocity component is non-zero owing to the background density stratification and, crucially, the wave is rotational. The lower plot is included to illustrate the dependence of the wave field on the stratification and because $\unicode[STIX]{x1D6E4}=0.2$ is the value used in the direct numerical simulation performed by Lin & Farouk (Reference Lin and Farouk2008) to which we compare our theoretical predictions.

This velocity field is plotted in figure 2 for two different values of $\unicode[STIX]{x1D6E4}$ . A rotational (vortical) component may be discerned, particularly for the $\unicode[STIX]{x1D6E4}=1$ scenario, even though viscous torques are absent since the dynamics in oscillatory boundary layers has been self-consistently omitted. Crucially, the Reynolds stress divergence terms, here denoted $R(x,y)$ , arising in (2.15) and responsible for driving the streaming flow can be explicitly evaluated:

(4.5) $$\begin{eqnarray}\displaystyle R(x,y)=-\unicode[STIX]{x2202}_{x}(\bar{\unicode[STIX]{x1D70C}}_{0}\overline{u_{1}^{\prime 2}})-\unicode[STIX]{x2202}_{y}(\bar{\unicode[STIX]{x1D70C}}_{0}\overline{u_{1}^{\prime }v_{1}^{\prime }})=-\frac{A^{2}}{2\unicode[STIX]{x03C0}(1+\unicode[STIX]{x1D6E4}/2)^{2}}\left(1+\frac{\unicode[STIX]{x1D6E4}}{4}+\frac{\unicode[STIX]{x1D6E4}y}{2}\right)\sin (2x). & & \displaystyle\end{eqnarray}$$

If $\unicode[STIX]{x1D6E4}=0$ , i.e. for a homogeneous fluid, $R(x,y)=\unicode[STIX]{x2202}_{x}[(A^{2}/4\unicode[STIX]{x03C0})\cos (2x)]$ . Clearly, in that case, the wave-induced Reynolds stress divergence can be balanced by a mean pressure gradient, so that streaming is not directly driven; instead, the associated Rayleigh streaming flow arises at next order in $\unicode[STIX]{x1D716}$ owing to the action of viscous torques within oscillatory boundary layers. If, however, $\unicode[STIX]{x1D6E4}\neq 0$ , then $R(x,y)$ can no longer be reduced to gradient form and, consequently, directly drives a streaming flow.

4.2 Streaming flow

We assume that the steady streaming driven by $R(x,y)$ given in (4.5) is sufficiently weak that the streaming equations (2.15)–(2.19) can be linearized. Note that, unlike Rayleigh streaming, the mean flow is compressible even though the streaming Mach number is negligible.

In a steady state, the pressure field can be eliminated from the linearized versions of equations (2.15)–(2.16), yielding

(4.6) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x2202}_{xy}R=-\frac{\unicode[STIX]{x2202}_{yyyx}\bar{u}_{1}}{Re_{s}h^{2}}, & & \displaystyle\end{eqnarray}$$

while conservation of mass and internal energy (2.17)–(2.18) imply

(4.7) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x2202}_{x}\bar{u}_{1}=-\frac{(1+\unicode[STIX]{x1D6E4}y)}{\unicode[STIX]{x1D6E4}Pe_{s}h^{2}}\unicode[STIX]{x2202}_{yyy}\bar{\unicode[STIX]{x1D6E9}}_{0}. & & \displaystyle\end{eqnarray}$$

Finally, we obtain from equations (4.5)–(4.7) a single differential equation for $\bar{\unicode[STIX]{x1D6E9}}_{0}$ :

(4.8) $$\begin{eqnarray}\displaystyle (1+\unicode[STIX]{x1D6E4}y)\unicode[STIX]{x2202}_{yyyyyy}\bar{\unicode[STIX]{x1D6E9}}_{0}+3\unicode[STIX]{x1D6E4}\unicode[STIX]{x2202}_{yyyyy}\bar{\unicode[STIX]{x1D6E9}}_{0}=-\frac{2A^{2}\unicode[STIX]{x1D6E4}^{2}Re_{s}Pe_{s}h^{4}}{\unicode[STIX]{x03C0}(2+\unicode[STIX]{x1D6E4})^{2}}\cos (2x). & & \displaystyle\end{eqnarray}$$

This equation is supplemented with the boundary conditions:

1. (i) $\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=0)=\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=1)=0$ , since the wall temperatures are held constant;

2. (ii) $\unicode[STIX]{x2202}_{yy}\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=0)=\unicode[STIX]{x2202}_{yy}\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=1)=0$ , to enforce $\bar{v}_{1}=0$ at the walls (see (2.17)–(2.19)); and

3. (iii) $\unicode[STIX]{x2202}_{yyy}\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=0)=\unicode[STIX]{x2202}_{yyy}\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=1)=0$ , to enforce the no-slip boundary condition at each wall, upon using (4.7).

With these boundary conditions, a unique solution can be found. For illustration, we focus on the case $\unicode[STIX]{x1D6E4}=1$ , for which

(4.9) $$\begin{eqnarray}\displaystyle \bar{\unicode[STIX]{x1D6E9}}_{0}(x,y)=-\frac{2A^{2}Re_{s}Pe_{s}h^{4}}{9\unicode[STIX]{x03C0}}G(y)\cos (2x), & & \displaystyle\end{eqnarray}$$

where

(4.10) $$\begin{eqnarray}\displaystyle G(y) & = & \displaystyle \frac{1}{1080(-3+\log (16))}\left[\vphantom{x_{y}^{2^{o}}}60(1+y)^{2}\log (1+y)\right.\nonumber\\ \displaystyle & & \displaystyle +\,y\left.\left(94-222\log (2)-90y-20y^{2}-5y^{3}(-5+\log (64))+3y^{4}(-3+\log (16))\right)\!\vphantom{x_{y}^{2^{o}}}\right].\nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

For the self-consistency of this approximation, the assumptions of one-way coupling and linear streaming dynamics require that the mean density change little during the evolution, i.e. $\bar{\unicode[STIX]{x1D6E9}}_{0}\ll 1$ , and that the nonlinear terms be negligible, e.g. $\bar{u}_{1}^{2}\ll 1$ . From equations (2.15)–(2.19) and (4.9), these constraints imply (for $\unicode[STIX]{x1D6E4}=1$ ) upper bounds on the dimensionless parameter combinations $A^{2}Re_{s}h^{2}$ and $A^{2}Re_{s}Pe_{s}h^{4}$ .

4.3 Comparison with previous work

Using the analytical solution (4.9), various properties of the streaming flow can be deduced. First, the maximum dimensional baroclinic streaming velocity obtained here

(4.11) $$\begin{eqnarray}\displaystyle \tilde{u} _{B}\propto (\max |G|)A^{2}Re_{s}h^{2}U_{\ast } & & \displaystyle\end{eqnarray}$$

can be compared to the corresponding velocity $\tilde{u} _{R}$ in Rayleigh streaming resulting from dissipation in the Stokes boundary layers, i.e.

(4.12) $$\begin{eqnarray}\displaystyle \tilde{u} _{R}\propto \unicode[STIX]{x1D716}A^{2}U_{\ast }. & & \displaystyle\end{eqnarray}$$

Clearly, for small values of $\unicode[STIX]{x1D716}$ , the baroclinic streaming dominates any streaming driven by viscous torques acting in near-wall Stokes layers.

Figure 3.

Analytically predicted baroclinic streaming flow obtained in the limit of one-way coupling (i.e. for sufficiently weak sound waves, which implies upper bounds on $A^{2}Re_{s}h^{2}$ and $A^{2}Re_{s}Pe_{s}h^{4}$ ). Colour is used to indicate the total temperature (blue is cold). Recall that the lower boundary is colder than the upper boundary, and that gravitational forces have been neglected. The parameter values are the same as those used to compute the sound wave shown in figure 2 (lower plot); additionally, $A=6$ , $h=2.3$ , $Re_{s}=5.7$ and $Pe_{s}=4.1$ . Collectively, these parameters correspond closely to the values employed by Lin & Farouk (Reference Lin and Farouk2008) in their full numerical simulations of stratified acoustic streaming. In contrast to Rayleigh streaming, the cells span the channel.

The streaming velocity field ( $\bar{u}_{1}$ , $\bar{v}_{1}$ ) computed from (4.8) for $\unicode[STIX]{x1D6E4}=0.2$ is plotted along with the total temperature $T_{B}+\bar{\unicode[STIX]{x1D6E9}}_{0}$ in figure 3. Clearly, baroclinic streaming differs from Rayleigh streaming not only in intensity but also in spatial structure; specifically, here the streaming cells span the channel while in Rayleigh streaming the cells are stacked in the wall-normal direction. These various distinguishing properties accord with prior experiments and numerical simulations; e.g. see Loh et al. (Reference Loh, Hyun, Ro and Kleinstreuer2002), Hyun et al. (Reference Hyun, Lee and Loh2005), Lin & Farouk (Reference Lin and Farouk2008), Nabavi et al. (Reference Nabavi, Siddiqui and Dargahi2008), Atkas & Ozgumus (Reference Atkas and Ozgumus2010) and Dreeben & Chini (Reference Dreeben and Chini2011).

Quantitative comparisons can be made with the results of Lin & Farouk (Reference Lin and Farouk2008), who performed direct numerical simulations of the compressible Navier–Stokes and heat equations specifically to investigate the impact of acoustic streaming in a thin channel on cross-channel heat transport. Thus, the system they considered is very similar to that studied here. In the absence of thermal driving, the streaming flow is accurately predicted by using Rayleigh’s formulation and, accordingly, exhibits a pattern of counter-rotating cells stacked in the $y$ direction. When a temperature difference is imposed, however, the stacked cells merge, resulting in counter-rotating cells that span the channel. In particular, for their case 1C, corresponding to the largest imposed temperature differential, the dimensionless parameters used by Lin & Farouk (Reference Lin and Farouk2008) are approximately

(4.13a-f ) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D716}=10^{-2},\quad \unicode[STIX]{x1D6FE}=1.4,\quad \unicode[STIX]{x1D6E4}=0.2,\quad h=2.3,\quad Re_{s}=5.7,\quad Pe_{s}=4.1.\quad & & \displaystyle\end{eqnarray}$$

The amplitude $A$ of the acoustic waves is not reported, but reasonably can be assumed to be similar to the value arising in the absence of thermal driving, for which $A\approx 6$ . The authors report the values of the $x$ (respectively $y$ ) component of the dimensional streaming velocity at $x=3\unicode[STIX]{x03C0}/4$ (respectively $x=\unicode[STIX]{x03C0}/2$ ). We compute $\bar{u}_{1}$ and $\bar{v}_{1}$ for the same parameters in this one-way coupling limit, and then obtain the corresponding dimensional velocities by multiplying $\bar{u}_{1}$ by the sound speed $a_{\ast }=353~\text{m}~\text{s}^{-1}$ and $\bar{v}_{1}$ by $\sqrt{\unicode[STIX]{x1D716}}ha_{\ast }$ . Figure 4 shows the resulting comparison. The evident quantitative agreement – with no adjustable parameters and despite the fact that the numerical simulations of Lin & Farouk (Reference Lin and Farouk2008) include several physical effects (oscillatory boundary-layer dynamics, viscous heating, inertia and temperature-dependent viscosity and diffusivity coefficients) not incorporated in our analysis – provides strong confirmation of the baroclinic streaming theory developed by Chini et al. (Reference Chini, Malecha and Dreeben2014) and systematically extended in the present study.

Figure 4.

Comparison of the $x$ and $y$ components of the streaming velocity field extracted from the numerical simulations of the instantaneous (compressible) Navier–Stokes and heat equations performed by Lin & Farouk (Reference Lin and Farouk2008) (their case 1C) with the corresponding streaming velocities predicted from the present theory in the limit of one-way wave/mean-flow coupling and linear streaming dynamics. The parameters correspond to: $A=6$ , $\unicode[STIX]{x1D6E4}=0.2$ , $\unicode[STIX]{x1D6FE}=1.4$ , $h=2.3$ , $Re_{s}=5.7$ and $Pe_{s}=4.1$ . Excellent quantitative agreement is observed without the use of any fitting parameters.

4.4 Heat flux enhancement

The streaming flow enhances cross-channel heat transport, as we now demonstrate using the leading-order solutions derived in this section. For this purpose, we introduce the steady-state Nusselt number $Nu$ as the ratio of the (dimensional) total heat flux $\dot{Q}_{\ast }$ to the diffusive flux,

(4.14) $$\begin{eqnarray}\displaystyle Nu=\frac{\dot{Q}_{\ast }}{2\unicode[STIX]{x03C0}\unicode[STIX]{x1D705}_{\ast }\unicode[STIX]{x0394}\unicode[STIX]{x1D6E9}_{\ast }(k_{\ast }H_{\ast })^{-1}}, & & \displaystyle\end{eqnarray}$$

and compute $Nu-1$ . The total heat flux $\dot{Q}_{\ast }$ is evaluated at the top boundary,

(4.15) $$\begin{eqnarray}\displaystyle \dot{Q}_{\ast }=\!\int _{0}^{2\unicode[STIX]{x03C0}k_{\ast }^{-1}}\!\!\unicode[STIX]{x1D705}_{\ast }\unicode[STIX]{x2202}_{{\tilde{y}}}\tilde{T}(\tilde{x},{\tilde{y}}=H_{\ast })\,\text{d}\tilde{x}=\unicode[STIX]{x1D705}_{\ast }(k_{\ast }H_{\ast })^{-1}T_{\ast }\!\left[\int _{0}^{2\unicode[STIX]{x03C0}}\!\!\unicode[STIX]{x2202}_{y}(T_{B}+\bar{\unicode[STIX]{x1D6E9}}_{0})(x,y=1\!)\,\text{d}x+O(\unicode[STIX]{x1D716})\!\right]\!\!, & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

so that the Nusselt number is given at leading order by

(4.16) $$\begin{eqnarray}\displaystyle Nu=1+\frac{1}{2\unicode[STIX]{x03C0}\unicode[STIX]{x1D6E4}}\int _{0}^{2\unicode[STIX]{x03C0}}\unicode[STIX]{x2202}_{y}\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=1)\,\text{d}x+O(\unicode[STIX]{x1D716}). & & \displaystyle\end{eqnarray}$$

Since viscous heating is neglected, the top and bottom heat fluxes are equal in a steady state; hence, the following result holds:

(4.17) $$\begin{eqnarray}\displaystyle \int _{0}^{2\unicode[STIX]{x03C0}}[T\unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D6E9}]_{0}^{1}\,\text{d}x & = & \displaystyle \int _{0}^{2\unicode[STIX]{x03C0}}[(1+\unicode[STIX]{x1D6E4})\unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D6E9}(x,y=1)-\unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D6E9}(x,y=0)]\,\text{d}x,\nonumber\\ \displaystyle & = & \displaystyle \unicode[STIX]{x1D6E4}\int _{0}^{2\unicode[STIX]{x03C0}}\unicode[STIX]{x2202}_{y}\bar{\unicode[STIX]{x1D6E9}}_{0}(x,y=1)\,\text{d}x+O(\unicode[STIX]{x1D716}).\end{eqnarray}$$

Moreover, from equations (2.17)–(2.19), the following equality can be derived,

(4.18) $$\begin{eqnarray}\displaystyle \int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}x\int _{0}^{1}\,\text{d}yT\unicode[STIX]{x2202}_{yy}\bar{\unicode[STIX]{x1D6E9}}_{0}=0\Rightarrow \int _{0}^{2\unicode[STIX]{x03C0}}[T\unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D6E9}]_{0}^{1}\,\text{d}x=\int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}x\int _{0}^{1}\,\text{d}y(\unicode[STIX]{x2202}_{y}T)(\unicode[STIX]{x2202}_{y}\bar{\unicode[STIX]{x1D6E9}}_{0})+O(\unicode[STIX]{x1D716}). & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

Finally, using the boundary condition $\bar{\unicode[STIX]{x1D6E9}}_{0}=0$ at the top and bottom walls, we obtain

(4.19) $$\begin{eqnarray}\displaystyle Nu-1 & = & \displaystyle \frac{1}{2\unicode[STIX]{x03C0}\unicode[STIX]{x1D6E4}^{2}}\int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}x\int _{0}^{1}\,\text{d}y(\unicode[STIX]{x2202}_{y}T)(\unicode[STIX]{x2202}_{y}\bar{\unicode[STIX]{x1D6E9}}_{0})+O(\unicode[STIX]{x1D716})\nonumber\\ \displaystyle & = & \displaystyle \frac{1}{2\unicode[STIX]{x03C0}\unicode[STIX]{x1D6E4}^{2}}\int _{0}^{2\unicode[STIX]{x03C0}}\,\text{d}x\int _{0}^{1}\,\text{d}y(\unicode[STIX]{x2202}_{y}\bar{\unicode[STIX]{x1D6E9}}_{0})^{2}+O(\unicode[STIX]{x1D716}).\end{eqnarray}$$

In the limit of one-way coupling and for $\unicode[STIX]{x1D6E4}=1$ , we evaluate this quantity with the approximate expression (4.9) for $\bar{\unicode[STIX]{x1D6E9}}_{0}$ derived in § 4.2,

(4.20) $$\begin{eqnarray}\displaystyle Nu-1\simeq 3.2\times 10^{-10}(A^{2}Re_{s}Pe_{s}h^{4})^{2}. & & \displaystyle\end{eqnarray}$$

This expression should be contrasted with that obtained for Rayleigh streaming (i.e. for sufficiently small temperature differences) and computed by Vainshtein, Fichman & Gutfinger (Reference Vainshtein, Fichman and Gutfinger1995):

(4.21) $$\begin{eqnarray}\displaystyle (Nu-1)_{R}=6.2\times 10^{-5}(\unicode[STIX]{x1D716}A^{2}Pe_{s}h^{2})^{2}. & & \displaystyle\end{eqnarray}$$

This comparison indicates that, although the Nusselt number derived here is very small (as noted in § 4.2, for the validity of the one-way coupling assumption the dimensionless parameter combination $A^{2}Re_{s}Pe_{s}h^{4}$ cannot be too large), it nevertheless is orders of magnitude larger than $Nu$ resulting from boundary-layer-driven acoustic streaming provided that, numerically, $\unicode[STIX]{x1D716}\ll 10^{-2}$ . Equation (4.20) also suggests that significant heat transport enhancement may be achieved in the limit $A^{2}Re_{s}Pe_{s}h^{4}\gg 1$ , i.e. when there is strong two-way coupling between the acoustic waves and the streaming flow. This strong coupling scenario is investigated in § 5.

4.5 Stability of the quiescent background state

Intriguingly, the amplitude equation (3.16) suggests that an acoustic wave may be amplified via interaction with the solid boundaries. Indeed, the first term on the right-hand side of (3.16) does not depend on the streaming flow, and quantifies how the divergence of the acoustic-wave velocity field (i.e. $g(x)$ ) and the net heat flux to/from the solid boundaries [proportional to $\unicode[STIX]{x2202}_{y}\hat{\unicode[STIX]{x1D6E9}}_{1}(x,y=1)-\unicode[STIX]{x2202}_{y}\hat{\unicode[STIX]{x1D6E9}}_{1}(x,y=0)$ ] are coupled. In a certain parameter regime, this system may act as a thermoacoustic engine and spontaneously generate an acoustic wave that, in turn, would drive a streaming flow and increase the cross-channel heat transport.

To see this, note that $\hat{\unicode[STIX]{x1D6E9}}_{1}$ can be determined from the acoustic-wave solution derived in § 4.1 and then used to evaluate

(4.22) $$\begin{eqnarray}\displaystyle -\frac{\text{i}\unicode[STIX]{x1D714}_{0}}{Pe_{s}h^{2}}\int _{D}\,\text{d}x\,\text{d}yg\unicode[STIX]{x2202}_{yy}\hat{\unicode[STIX]{x1D6E9}}_{1}=\frac{-2\unicode[STIX]{x1D6E4}[1+(2-\unicode[STIX]{x1D6E4})(\unicode[STIX]{x1D6FE}-1)]}{Pe_{s}h^{2}(1+\unicode[STIX]{x1D6E4}/2)}. & & \displaystyle\end{eqnarray}$$

Since $\unicode[STIX]{x1D6FE}>1$ , this term is negative and damps the acoustic waves, unless the thermal forcing $\unicode[STIX]{x1D6E4}$ is sufficiently strong. In that case, a proper discussion of the stability of the quiescent background state would require the evaluation of the acoustic-wave energy dissipation and temperature profile in the oscillatory boundary layers, which has been neglected here (see § 6). Nevertheless, this calculation suggests that such an instability, in which an acoustic wave spontaneously grows by a thermoacoustic effect even in the absence of an imposed wall-parallel ( $x$ ) temperature gradient, is possible.