Strong wave-mean-flow coupling in baroclinic acoustic streaming

TL;DR

This paper develops a coupled wave-mean-flow model for baroclinic acoustic streaming in stratified fluids, revealing strong two-way interactions and potential for enhanced heat transfer, validated by analytical and numerical methods.

Contribution

It introduces a novel coupled system of equations capturing strong wave-mean-flow interactions driven by baroclinic effects, extending classical streaming theory.

Findings

Mean flow arises at first order in wave amplitude, enabling two-way coupling.

Analytical solutions match prior DNS results with high accuracy.

Numerical simulations demonstrate potential for heat transfer enhancement.

Abstract

The interaction of an acoustic wave with a stratified fluid can drive strong streaming flows owing to the baroclinic production of fluctuating vorticity, as recently demonstrated by Chini et al. (J. Fluid Mech., 744, 2014, pp. 329). In the present investigation, a set of wave/mean-flow interaction equations is derived that governs the coupled dynamics of a standing acoustic wave mode of characteristic (small) amplitude {\epsilon} and the streaming flow it drives in a thin channel with walls maintained at differing temperatures. Unlike classical Rayleigh streaming, the resulting mean flow arises at O({\epsilon}) rather than at O({\epsilon^2}). Consequently, fully two-way coupling between the waves and the mean flow is possible: the streaming is sufficiently strong to induce O(1) rearrangements of the imposed background temperature and density fields, which modifies the spatial structure and frequency of the acoustic mode on the streaming time scale. A novel Wentzel-Kramers-Brillouin-Jeffreys analysis is developed to average over the fast wave dynamics, enabling the coupled system to be integrated strictly on the slow time scale of the streaming flow. Analytical solutions of the reduced system are derived for weak wave forcing and are shown to reproduce results from prior direct numerical simulations (DNS) of the compressible Navier Stokes and heat equations with remarkable accuracy. Moreover, numerical simulations of the reduced system are performed in the regime of strong wave mean flow coupling for a fraction of the computational cost of the corresponding DNS. These simulations shed light on the potential for baroclinic acoustic streaming to be used as an effective means to enhance heat transfer.

Tables

Table 1. Table 1: Dimensional variables and parameters

Notation	Definition
$\widetilde{𝐮}=(\tilde{u},\tilde{v})$	Gas velocity
$\tilde{\rho }$	Gas density
$\tilde{p}$	Gas pressure
$\tilde{T}$	Gas temperature
$(\tilde{x},\tilde{y})$	Horizontal, vertical (wall-normal) coordinate
$\tilde{t}$	Time variable
$H_{\ast }$	Channel height
$k_{\ast }$	Acoustic-wave horizontal wavenumber
${\mu }_{\ast }$	Dynamic viscosity
${\kappa }_{\ast }$	Thermal conductivity
$R_{\ast }$	Specific gas constant
$(c_{v_{\ast }},c_{p_{\ast }})$	Constant volume, pressure specific heat coefficient
$a_{\ast }=\sqrt{\left(c_{p_{\ast }}/c_{v_{\ast }}\right)R_{\ast }{T}_{\ast }}$	Background sound speed
$p_{\ast }$	Background pressure field

Table 2. Table 2: Dimensionless variables and parameters

Variable	Scale	Parameter	Definition	Scaling
$x$	$k_{\ast }^{-1}$	Strouhal number $S$	$a_{\ast }/U_{\ast }$	$S=1/ϵ$
$y$	$H_{\ast }$	Aspect ratio $\delta$	$k_{\ast }{H}_{\ast }$	$\delta =\sqrt{ϵ}h$
$t$	${\left(a_{\ast }{k}_{\ast }\right)}^{-1}$	Temperature gradient $\Gamma$	$\Delta {\Theta }_{\ast }/T_{\ast }$	$\Gamma =O\left(1\right)$
$u$	$a_{\ast }$	Reynolds number $Re$	${\rho }_{\ast }{U}_{\ast }/\left(k_{\ast }{\mu }_{\ast }\right)$	$Re=R{e}_s/ϵ$
$v$	$\left(k_{\ast }{H}_{\ast }\right)a_{\ast }$	Péclet number $Pe$	${\rho }_{\ast }{c}_{p_{\ast }}{U}_{\ast }/\left(k_{\ast }{\kappa }_{\ast }\right)$	$Pe=P{e}_s/ϵ$
$\rho$	${\rho }_{\ast }\equiv p_{\ast }/\left(R_{\ast }{T}_{\ast }\right)$	Specific heat ratio $\gamma$	$c_{p_{\ast }}/c_{v_{\ast }}$	$\gamma =O\left(1\right)$
$T$	$T_{\ast }$
$P$	$p_{\ast }$

Equations

$$\tilde{\rho }\left[{\partial }_{\tilde{t}}\tilde{\mathbf{u}}+(\tilde{\mathbf{u}}\cdot \tilde{\nabla })\tilde{\mathbf{u}}\right]=-\tilde{\nabla }\tilde{p}+{\mu }_{\ast }\left[{\tilde{\nabla }}^2\tilde{\mathbf{u}}+\frac{1}{3}\tilde{\nabla }(\tilde{\nabla }\cdot \tilde{\mathbf{u}})\right],$$

$${\partial }_{\tilde{t}}\tilde{\rho }+\tilde{\nabla }\cdot \left(\tilde{\rho }\tilde{\mathbf{u}}\right)=0,$$

$$\tilde{\rho }{c}_v\left[{\partial }_{\tilde{t}}\tilde{T}+(\tilde{\mathbf{u}}\cdot \tilde{\nabla })\tilde{T}\right]=-\tilde{p}\left(\tilde{\nabla }\cdot \tilde{u}\right)+{\kappa }_{\ast }{\tilde{\nabla }}^2\tilde{T},$$

$$\tilde{p}=\tilde{\rho }{R}_{\ast }\tilde{T},$$

$${\tilde{T}}_B=T_{\ast }\left(1+\Gamma \frac{y_{\ast }}{H_{\ast }}\right),{\tilde{p}}_B=p_{\ast },$$

$$(u,v)$$

$$\pi$$

$$\Theta$$

$$\rho$$

$$\varphi (t)=\frac{\Phi (T)}{ϵ},$$

$$\omega (T)=\frac{\mathrm{d}\varphi }{\mathrm{d}t}=\frac{\mathrm{d}\Phi }{\mathrm{d}T},$$

$$\omega ={\omega }_0+ϵ{\omega }_1+O({ϵ}^2).$$

$$\bar{f}(x,y,T)=\frac{1}{2n\pi }{\int }_{\varphi }^{\varphi +2n\pi }f(x,y,s,T)\mathrm{d}s$$

$$f(x,y,\varphi ,T)=\bar{f}(x,y,T)+f^{\prime }(x,y,\varphi ,T),$$

$${\bar{\rho }}_0\left({\partial }_T{\bar{u}}_1+{\bar{u}}_1{\partial }_x{\bar{u}}_1+{\bar{v}}_1{\partial }_y{\bar{u}}_1\right)$$

$${\partial }_y{\bar{\pi }}_2$$

$${\partial }_T{\bar{\rho }}_0+{\partial }_x({\bar{\rho }}_0{\bar{u}}_1)+{\partial }_y({\bar{\rho }}_0{\bar{v}}_1)$$

$${\partial }_T{\bar{\Theta }}_0+{\bar{u}}_1{\partial }_x{\bar{\Theta }}_0+{\bar{v}}_1{\partial }_y\left({\bar{\Theta }}_0+T_B\right)$$

$${\bar{\rho }}_0$$

$${\omega }_0\overline{{\rho }_0}{\partial }_{\varphi }{u}_1^{\prime }+\frac{1}{\gamma }{\partial }_x{\pi }_1^{\prime }=0,{\partial }_y{\pi }_1^{\prime }=0,$$

$${\omega }_0{\partial }_{\varphi }{\rho }_1^{\prime }+{\partial }_x({\bar{\rho }}_0{u}_1^{\prime })+{\partial }_y({\bar{\rho }}_0{v}_1^{\prime })=0,$$

$${\omega }_0{\partial }_{\varphi }{\Theta }_1^{\prime }+u_1^{\prime }{\partial }_x{\bar{\Theta }}_0+v_1^{\prime }{\partial }_y\left({\bar{\Theta }}_0+T_B\right)+(\gamma -1)({\bar{\Theta }}_0+T_B)({\partial }_x{u}_1^{\prime }+{\partial }_y{v}_1^{\prime })=0,$$

$${\pi }_1^{\prime }-{\rho }_1^{\prime }({\bar{\Theta }}_0+T_B)-{\bar{\rho }}_0{\Theta }_1^{\prime }=0.$$

$$f_1^{\prime }(x,y,\varphi ,T)=\frac{A(T)}{2}\left({\hat{f}}_1(x,y,T)e^{i\varphi }+c.c.\right),$$

$${\hat{\pi }}_1=\frac{i\gamma }{{\omega }_0}\left({\partial }_x{\hat{u}}_1+{\partial }_y{\hat{v}}_1\right).$$

$${\partial }_x\left({\partial }_x{\hat{u}}_1+{\partial }_y{\hat{v}}_1\right)+{\omega }_0^2{\bar{\rho }}_0{\hat{u}}_1=0,$$

$${\partial }_y\left({\partial }_x{\hat{u}}_1+{\partial }_y{\hat{v}}_1\right)=0.$$

$${\partial }_y(\text{ref}wave:eq1)+{\partial }_x(\text{ref}wave:eq2)\Rightarrow {\partial }_y\left({\bar{\rho }}_0{\hat{u}}_1\right)=0\Rightarrow {\hat{u}}_1=\frac{q}{{\bar{\rho }}_0},$$

$${\partial }_x{\hat{u}}_1+{\partial }_y{\hat{v}}_1=g,$$

$${\hat{u}}_1=-\frac{g^{\prime }}{{\omega }_0^2{\bar{\rho }}_0},$$

$${\hat{v}}_1=yg+{\partial }_x\left(\frac{g^{\prime }}{{\omega }_0^2}{\int }_0^y\frac{\mathrm{d}y}{{\bar{\rho }}_0}\right).$$

$$\frac{\mathrm{d}}{\mathrm{d}x}\left(g^{\prime }\alpha \right)+{\omega }_0^{2}g=0,$$

$$\alpha (x,T)={\int }_0^1\frac{\mathrm{d}y}{{\bar{\rho }}_0}.$$

$${\int }_0^{2\pi }{g}_A(x)g_B(x)\mathrm{d}x=\frac{1}{{\omega }_A^2-{\omega }_B^2}{\left[\alpha (g_A{g}_B^{\prime }-g_B{g}_A^{\prime })\right]}_0^{2\pi }=0.$$

Full text

1 Introduction

Sound waves can drive Eulerian flows that evolve on a slow time scale compared to the period of the waves. The theoretical study of this phenomenon, called acoustic (or, in other contexts, steady) streaming, can be traced back to Rayleigh in the 19th century (Rayleigh, 1884). Given that ultrasonic power sources are now routinely used in laboratory experiments, acoustic streaming has been widely observed, often as a source of unwanted flow. Nonetheless, streaming also has been recognized as a practical means to enhance transport and mixing and has, for instance, been used to improve the efficiency of chemical reactions occurring near a catalytic solid phase that otherwise would be controlled by molecular diffusion (Bengtsson & Laurell, 2004); to directly mix chemical species (Yaralioglu et al., 2004); and for activated irrigation in medical applications including root-canal procedures (Verhaagen et al., 2014). Heat also can be transported by streaming flows, and acoustic waves therefore can be used to accelerate the cooling of a hot object, as recently reviewed by Legay et al. (2011). Acoustic streaming technologies are of particular interest in the zero-gravity environment, where natural convective flows do not exist and acoustic thermal management systems may provide a reliable, efficient and lightweight alternative to fans.

In a characteristically lucid lecture, Lighthill identified the different regimes of acoustic streaming in a homogeneous medium (Lighthill, 1978). Owing to attenuation mechanisms, acoustic waves generate a Reynolds stress divergence capable of driving a mean flow, which is balanced either by viscous forces (termed “Rayleigh streaming” if, in addition, the sound waves are damped in oscillatory boundary layers) or by inertia (“Stuart streaming”). The former regime occurs for small values of the streaming Reynolds number $Re_{s}=U_{s}L/\nu$, where $U_{s}$ is a characteristic streaming speed, $L$ is a typical dimension of the system and $\nu$ is the kinematic viscosity of the fluid. In Rayleigh streaming, the (laminar) cellular mean flow that is generated is localized within a few wavelengths of any solid boundary (Nyborg, 1958); the streaming can be analytically computed for simple geometries, e.g. in a channel (Rayleigh, 1884; Hamilton et al., 2003) or adjacent to a circular cylinder (Holtsmark et al., 1954). This regime has become important in microfluidics, where the vortices induced by acoustic streaming in microchannels can be used to mix chemicals (see references above). In contrast, for large $Re_{s}$, the streaming flow acquires a jet-like structure and can become turbulent (Stuart, 1966; Lighthill, 1978).

The presence of an inhomogeneous background temperature (or density) field strongly affects the fundamental mechanics and kinematics of acoustic streaming: streaming velocities are significantly enhanced and the flow patterns are substantially altered. These changes are evident in the early experiments of Fand & Kaye (1960) and in subsequent experiments and numerical simulations; see e.g. Loh et al. (2002), Hyun et al. (2005), Lin & Farouk (2008), Nabavi et al. (2008), Atkas & Ozgumus (2010), Dreeben & Chini (2011) and Karlsen et al. (2018). Consequently, the resulting flow and associated transport cannot be computed simply by coupling the corresponding isothermal (e.g. Rayleigh or Stuart) streaming with the heat or other appropriate transport equation. Instead, new physical phenomena occur, which renders this problem both complex and interesting. Experimental challenges arise because natural convection and acoustic streaming may be difficult to disentangle in the laboratory; numerical challenges result from the need to resolve compressible fluid dynamics on temporal scales ranging from the acoustic wave period to the slow time scale over which the streaming flow evolves; while the primary theoretical challenge is to elucidate the novel phenomenology resulting from fully two-way coupling between the sound waves and the mean flow.

This striking change in the character of the streaming has been observed in various contexts in which an agency other than viscosity generates vorticity in the oscillatory flow (e.g. the acoustic waves). Amin (1988) and Riley & Trinh (2001) noted fundamental changes in the steady streaming driven by a non-conservative body force in their study of fluid flow in the presence of $g$-jitter, i.e. a fluctuating gravitational field in an otherwise gravity-free environment. Motivated by the observation that streaming velocities in high-intensity discharge lamps are two orders of magnitude larger than those predicted by Rayleigh streaming theory (Dreeben & Chini, 2011), Chini et al. (2014) derived a theory capable of accounting for both the observed streaming pattern and magnitude. As discussed more fully in § 2, the mechanism underlying this large-amplitude streaming is the baroclinic production of sound-wave vorticity arising from the misalignment of fluctuating isobars and mean isopycnals. Similarly, non-classical streaming phenomena have been observed in microfluidic systems with gradients in density; in particular, Karlsen et al. (2016, 2018) recently obtained a local expression for the acoustic force density driving the streaming flow as a function of the acoustic-wave characteristics.

The primary objective of the present investigation is to systematically extend the recent theory of Chini et al. (2014) to efficiently capture the two-way coupling that can occur in baroclinic acoustic streaming and to quantify the concomitant heat transfer. Accordingly, a novel Wentzel–Kramers–Brillouin–Jeffreys (WKBJ) analysis is performed to enable prediction of the slow evolution of the acoustic wave amplitude, thereby obviating the need to explicitly simulate the fast oscillatory dynamics. We focus on perhaps the most well-documented acoustic streaming configuration: a thin two-dimensional channel with an imposed standing acoustic wave oscillating in the wall-parallel direction. For a homogeneous system, the resulting streaming flow was first described in the pioneering work of Rayleigh (1884) in the limits $Re_{s}\ll 1$ and $\delta_{*}\ll H_{*}\ll k_{*}^{-1}$, where $\delta_{*}=\sqrt{2\nu_{*}/\omega_{*}}$ is the thickness of the oscillatory (Stokes) boundary layers ($\nu_{*}$ is the kinematic viscosity, and $\omega_{*}$ is the wave angular frequency), $H_{*}$ is the channel width and $k_{*}$ is the wavenumber of the acoustic wave. The streaming flow comprises a wall-parallel array of counter-rotating vortices, stacked in the wall-normal direction and having a characteristic velocity $3U_{*}^{2}/(16a_{*})$, where $U_{*}$ is the maximum fluctuating velocity induced by the standing acoustic wave and $a_{*}$ is the speed of sound. Experiments were first performed in a tube and showed quantitative agreement with the predictions of Rayleigh (Andrade, 1931). Hamilton et al. (2003) extended this theoretical study to channels of arbitrary width $H_{*}$, the only restrictions being $Re_{s}\ll 1$ and $\delta_{*}\ll k_{*}^{-1}$. When the upper and lower walls of the channel are maintained at fixed but differing temperatures, both experiments and direct numerical simulations of the compressible Navier–Stokes and heat equations indicate a change in the streaming phenomenology: the stacked vortices merge and their characteristic velocity increases (Loh et al., 2002; Lin & Farouk, 2008). To date, no theory has correctly predicted the resulting streaming-flow pattern and intensity; our study, which focuses on the regime $Re_{s}\gtrsim 1$, fills this gap in the literature.

The remainder of the paper is organized as follows. After formulating the problem for the instantaneous dynamics, we carry out a multiple scale analysis (§ 2) to obtain a reduced but two time-scale system. In § 3, we analyse the wave dynamics to show that this multiscale system can be integrated strictly on the slow time scale. We then consider, in § 4, the limit of weak wave forcing, in which the streaming flow does not produce appreciable feedback on the waves; in particular, we derive an approximate analytical solution and compare it to the streaming flow numerically computed by Lin & Farouk (2008), demonstrating excellent quantitative agreement. In § 5, we perform numerical simulations of our reduced model and characterize the resulting fully-coupled waves and mean flows. We summarize our key findings and suggest possible further extensions in § 6.

2 Two time-scale wave/mean-flow system

2.1 Flow configuration

The problem we consider is similar to that introduced in Chini et al. (2014). Specifically, we analyze the two-dimensional flow of an ideal gas, with specific gas constant $R_{*}$ and constant dynamic viscosity $\mu_{*}$ and thermal conductivity $\kappa_{*}$, in a channel with walls separated by a distance $H_{*}$ 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 $\omega_{*}$, yielding a standing sound wave with spatial wavenumber $k_{*}$. 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_{*}$. All dependent fields are required to satisfy a $2\pi/k_{*}$ 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}$.

In contrast to the study of Chini et al. (2014), the thermal driving is achieved by fixing the temperatures of the lower and upper walls to be $T_{*}$ and $T_{*}+\Delta\Theta_{*}$, respectively, rather than by including a volumetric heat source. For convenience, we take the temperature differential $\Delta\Theta_{*}>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{\rho}$, $\tilde{p}$, $\tilde{T}$ and $\tilde{\mathbf{u}}$, respectively, where $\tilde{\mathbf{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

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where the two-dimensional gradient operator $\tilde{\nabla}=(\partial_{\tilde{x}},\partial_{\tilde{y}})$. Note that dilatational (or ‘bulk’) viscosity has been neglected in (1), and viscous heating has been omitted in (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 et al., 1973). 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{\mathbf{u}}=\mathbf{0}$, are referred to as the background fields (denoted with a subscript ‘$B$’) and are found to be

[TABLE]

where $p_{*}$ is a constant and the dimensionless temperature differential $\Gamma\equiv\Delta\Theta_{*}/T_{*}$. Table 1 summarizes the dimensional fields and parameters used in the following analysis.

2.2 Scaling and non-dimensionalisation

To facilitate the asymptotic analysis, we non-dimensionalise the governing equations by scaling the dependent and independent variables as outlined in table 2. The $\tilde{x}$-velocity component is scaled with $a_{*}$, the sound speed at temperature $T_{*}$. This scaling introduces into the dimensionless governing equations the Strouhal number $S=a_{*}/U_{*}$, where $U_{*}$ (rather than $a_{*}$) 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 $\epsilon\equiv 1/S$ and consider the asymptotic limit $\epsilon\to 0$ with all other dimensionless parameters scaled as appropriate powers of $\epsilon$. Since $\epsilon\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 [$\mathit{O}(1/\epsilon)$] 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_{*}^{-1}$ and $H_{*}$, respectively, so that the gas lies in the domain defined by $x\in[0,2\pi]$ and $y\in[0,1]$. Time $\tilde{t}$ is non-dimensionalised using the inverse reference wave frequency $\omega_{*}^{-1}=(a_{*}k_{*})^{-1}$. The vertical ($\tilde{y}$), or wall-normal, velocity component is scaled by $(k_{*}H_{*})a_{*}$. The domain aspect ratio $\delta\equiv k_{*}H_{*}$ is assumed to be small and, more precisely, is chosen so that $\delta=\sqrt{\epsilon}h$, where $h$ is a dimensionless parameter of order unity. Although the cross-channel heat flux is not expected to be maximized as $\delta\to 0$ (since the streaming flow will be largely horizontal), we follow Chini et al. (2014) 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.

The temperature $T_{*}$ of the lower wall is used to nondimensionalise the temperature field $\tilde{T}$. In the analysis that follows, $\Gamma$ is fixed, i.e. $\mathit{O}(1)$, as $\epsilon\to 0$, although the smallness of $\epsilon$ in acoustic streaming ensures that our multiple scale analysis remains accurate even for $\Gamma\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_{*}$ (in contrast to Rayleigh streaming, where streaming speeds are proportional to $\epsilon U_{*}$). 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 $\delta^{-2}\propto\epsilon^{-1}$. We therefore define the streaming Reynolds and Péclet numbers $Re_{s}=\epsilon Re$ and $Pe_{s}=\epsilon 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 $\epsilon$ 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\epsilon 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:

[TABLE]

where $\pi=(\tilde{p}-\tilde{p}_{B})/p_{*}$ and $\Theta=(\tilde{T}-\tilde{T}_{B})/T_{*}$ 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 (6) follows from the scaling of the Strouhal number ($S=1/\epsilon$). The state equation constrains the temperature and density perturbations to be of the same order. In baroclinic acoustic streaming, the $\mathit{O}(\epsilon)$ streaming flow is sufficiently strong to induce $\mathit{O}(1)$ rearrangements of the background temperature and density fields over an $\mathit{O}(1/\epsilon)$ time period. Thus, in contrast to other studies of acoustic streaming in the presence of inhomogeneous temperature fields (e.g. Červenka & Bednarřik (2017)), it is crucial that the expansions for both $\Theta$ and $\rho$ begin at $\mathit{O}(1)$, as first shown in Chini et al. (2014).

Owing to the large perturbations to the background density field, the natural frequency of an acoustic mode may evolve in time. Here, we extend the analysis of Chini et al. (2014) 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,\phi,T)$, where $\phi$ and $T$ are treated as independent variables. The rapidly-varying phase $\phi$ may be written as

[TABLE]

where $\mathrm{d}\Phi/\mathrm{d}T$ is of order unity. We define the instantaneous angular frequency $\omega(T)$,

[TABLE]

and expand $\Phi=\Phi_{0}+\epsilon\Phi_{1}+O(\epsilon^{2})$, so that

[TABLE]

Finally, in order to distinguish the streaming flow from the acoustic wave, we introduce the fast-time average of a function $f(x,y,\phi,T)$,

[TABLE]

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

[TABLE]

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. (2014) and are reproduced here for completeness. The evolution of the streaming fields is governed by the following equations:

[TABLE]

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 (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 $\epsilon$ 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.

[TABLE]

[TABLE]

[TABLE]

[TABLE]

These equations describe the fast dynamics of approximately linear non-dissipative acoustic waves in a medium whose mean density field $\bar{\rho}_{0}$ evolves slowly in time. Conversely, the evolution of $\bar{\rho}_{0}(x,y,T)$ depends on the acoustic waves, as can be gleaned from inspection of (15)–(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 $\mathit{O}(\epsilon^{2})$ while the streaming flow governed by (15)–(19) is $\mathit{O}(\epsilon)$. 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. (2014), 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. (2016, 2018) also noted this two-way coupling in their computational studies of acoustic streaming in inhomogeneous fluids. Nevertheless, to make analytical progress, Chini et al. (2014) 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. (2014) 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 Averaging over fast wave dynamics

We now characterize the dynamics of the acoustic waves on both the fast and slow time scales with the aim of eliminating the need to explicitly simulate the fast evolution. Inspection of (20)–(23) reveals that the wave dynamics directly depends on a single slowly-varying field: the leading-order mean density $\bar{\rho}_{0}=(\bar{\Theta}_{0}+T_{B})^{-1}$. For purposes of the analysis described in this section, this field is presumed to be given. Consequently, the fluctuation equations (20)–(23) comprise a linear homogeneous system, and we henceforth consider a single eigenvector. A generic fluctuation field $f_{1}^{\prime}$ can be expressed as

[TABLE]

where $f_{1}^{\prime}$ stands for any fluctuation variable ($u_{1}^{\prime}$, $v_{1}^{\prime}$, $\rho_{1}^{\prime}$, $\pi_{1}^{\prime}$, $\Theta_{1}^{\prime}$); $A(T)$ is the slowly-evolving modal amplitude (here taken to be real without loss of generality); $\hat{f}$ is a complex function that describes the spatial structure of the mode; and c.c. denotes the complex conjugate. A normalization condition, specified subsequently, must be imposed on $\hat{f}$ to render this decomposition unique. We next describe the determination of the spatial structure of the mode (as defined by the functions $\hat{u}_{1}$, $\hat{v}_{1}$, etc.) and then derive a novel amplitude equation governing the slow evolution of the generally a priori unknown function $A(T)$.

3.1 Mode structure

Substituting the decomposition (24) into the fluctuation equations (20)–(23) yields a linear but non-separable two-dimensional (partial) differential eigenvalue problem for the spatial structure and frequency $\omega_{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 (24), (21)–(23) can be combined to deduce

[TABLE]

With this expression for the acoustic wave pressure $\hat{\pi}_{1}$, the momentum equations (20) become two coupled partial differential equations for $\hat{u}_{1}$ and $\hat{v}_{1}$:

[TABLE]

[TABLE]

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

[TABLE]

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

[TABLE]

where $g$ is a second unknown function of $x$ and $T$ only. Equations (26), (28) and (29) imply that $q=-g^{\prime}/\omega_{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.

[TABLE]

[TABLE]

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

[TABLE]

where

[TABLE]

To characterize the function $g(x,T)$, we first take, without loss of generality, $\hat{u}_{1}$ to be a real field: (26) and (27) 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 (32) must be solved subject to the boundary conditions $g^{\prime}(0)=g^{\prime}(2\pi)=0$. This requirement leads to an orthogonality condition: let $(g_{A},g_{B})$ be two eigenvectors and $(\omega_{A},\omega_{B})$ their angular eigenfrequencies; then

[TABLE]

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

[TABLE]

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

3.2 Wave amplitude

As explained in the previous subsection, the acoustic mode shape $g(x,T)$ and frequency $\omega_{0}(T)$ can be computed at every time $T$ for a given mean density profile $\bar{\rho}_{0}(x,y,T)$. The advection of hot or cold gas by the streaming flow will cause both $g$ and $\omega_{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.

1. We collect terms in the dimensionless governing equations at $\mathit{O}(\epsilon^{2})$. The resulting equations are reported in Appendix A.

2. We make the ansatz that generic $\mathit{O}(\epsilon^{2})$ fluctuation field $f_{2}^{\prime}$ can be represented as

[TABLE]

and, in direct analogy with the manipulations performed in § 3.1, reduce the $\mathit{O}(\epsilon^{2})$ fluctuation system for the five unknown fields $\hat{u}_{2}$, $\hat{v}_{2}$, $\hat{\pi}_{2}$, $\hat{\rho}_{2}$ and $\hat{\Theta}_{2}$ to a system of two equations for the two fields $\hat{u}_{2}$ and $\hat{v}_{2}$. [Harmonics of the form $e^{i2\phi}$ 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

[TABLE]

[TABLE]

The linear operator acting on the left-hand side of this system is identical to that arising in the leading-order fluctuation equations (26)–(27). 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.

3. We derive a solvability condition for the $\mathit{O}(\epsilon^{2})$ system, which requires determination of the adjoint linear operator, by invoking the Fredholm alternative theorem in the usual manner; see Appendix C.

4. We enforce the solvability condition to obtain an equation for $\mathrm{d}A/\mathrm{d}T$; cf. Appendix D.

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

[TABLE]

where $\int_{D}$ refers to definite integration over the spatial domain. (Note that, since the temperature and velocity fluctuations are out of phase, $\hat{\Theta}_{1}$ is strictly imaginary.) The lack of terms nonlinear in $A$ in (39) 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 (39) 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, 1988). Owing to the occurrence of the mean fields $\bar{u}_{1}$, $\bar{v}_{1}$ and $\bar{\rho}_{0}$ in (39), the amplitude equation is, in fact, nonlinear. A distinguishing feature of (39) is that, unlike amplitude equations derived in numerous other contexts, determination of the coefficients requires evaluation of functional derivatives that capture the $\mathit{O}(1)$ variations in generic eigenfunction $\hat{f}_{1}$ caused by changes in the mean density field $\bar{\rho}_{0}$ that occur on the slow time $T$; see Appendices B and D for details.

The quantity $A\omega_{0}^{-1}$ on the left-hand side of (39) 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

[TABLE]

This expression can be evaluated using (30) to obtain

[TABLE]

The last integral reduces to $\omega_{0}^{2}$ following an integration by parts and ulitisation of the differential equation (32) and the normalization condition (35), yielding

[TABLE]

The amplitude equation therefore can be interpreted as an energy balance for the acoustic wave, with the left-hand side of (39) equalling $(1/\overline{E}_{K})\mathrm{d}\overline{E}_{K}/\mathrm{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 (32), while the evolution of the amplitude can be determined by integrating (39) over the slow time scale. These computations are performed in conjunction with the numerical solution of the streaming equations (15)–(19).

4 One-way coupling

Although the asymptotically-reduced equations constitute a substantial simplification of the full compressible Navier–Stokes equations, they defy analytical solution owing to the occurrence of nonlinearities and two-way coupling between the waves and the streaming flow. For sufficiently weak streaming, however, an approximate steady-state solution can be derived: in this limit, mean advection is weak and the mean density perturbations are small, thereby ameliorating these two difficulties. Here, we derive this approximate analytical solution and demonstrate that it accurately describes results obtained from direct numerical simulations of the compressible Navier–Stokes and heat equations reported in the literature.

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 (5), i.e.

[TABLE]

The coefficient $\alpha$ defined in (33) then does not depend on $x$ and is simply equal to $1+\Gamma/2$. The acoustic-wave eigenfunction $g$ characterizing the shape of the acoustic mode follows from the solution of the second-order differential equation (32) and, with the prescribed boundary conditions $g^{\prime}(0)=g^{\prime}(2\pi)=0$ and the normalization condition (35), is given by

[TABLE]

where the integer $n$ is set to unity for this study. The angular frequency of this mode $\omega_{0}=\sqrt{1+\Gamma/2}$, and the velocity field ($\hat{u}_{1}$, $\hat{v}_{1}$) follows from (30)–(31) and reduces to

[TABLE]

[TABLE]

This velocity field is plotted in figure 2 for two different values of $\Gamma$. A rotational (vortical) component may be discerned, particularly for the $\Gamma=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 (15) and responsible for driving the streaming flow can be explicitly evaluated:

[TABLE]

If $\Gamma=0$, i.e. for a homogeneous fluid, $R(x,y)=\partial_{x}[(A^{2}/4\pi)\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 $\epsilon$ owing to the action of viscous torques within oscillatory boundary layers. If, however, $\Gamma\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 (47) is sufficiently weak that the streaming equations (15)–(19) can be linearised. 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 linearised versions of equations (15)–(16), yielding

[TABLE]

while conservation of mass and internal energy (17)–(18) imply

[TABLE]

Finally, we obtain from equations (47)–(49) a single differential equation for $\bar{\Theta}_{0}$:

[TABLE]

This equation is supplemented with the boundary conditions:

$\bar{\Theta}_{0}(x,y=0)=\bar{\Theta}_{0}(x,y=1)=0$, since the wall temperatures are held constant; 2. 2.

$\partial_{yy}\bar{\Theta}_{0}(x,y=0)=\partial_{yy}\bar{\Theta}_{0}(x,y=1)=0$, to enforce $\bar{v}_{1}=0$ at the walls [see (17)–(19)]; and 3. 3.

$\partial_{yyy}\bar{\Theta}_{0}(x,y=0)=\partial_{yyy}\bar{\Theta}_{0}(x,y=1)=0$, to enforce the no-slip boundary condition at each wall, upon using (49).

With these boundary conditions, a unique solution can be found. For illustration, we focus on the case $\Gamma=1$, for which

[TABLE]

where

[TABLE]

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{\Theta}_{0}\ll 1$, and that the nonlinear terms be negligible, e.g. $\bar{u}_{1}^{2}\ll 1$. From equations (15)–(19) and (51), these constraints imply (for $\Gamma=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 (51), various properties of the streaming flow can be deduced. First, the maximum dimensional baroclinic streaming velocity obtained here

[TABLE]

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

[TABLE]

Clearly, for small values of $\epsilon$, the baroclinic streaming dominates any streaming driven by viscous torques acting in near-wall Stokes layers.

The streaming velocity field ($\bar{u}_{1}$, $\bar{v}_{1}$) computed from (50) for $\Gamma=0.2$ is plotted along with the total temperature $T_{B}+\bar{\Theta}_{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. (2002); Hyun et al. (2005); Lin & Farouk (2008); Nabavi et al. (2008); Atkas & Ozgumus (2010); Dreeben & Chini (2011).

Quantitative comparisons can be made with the results of Lin & Farouk (2008), 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 (2008) are approximately

[TABLE]

The amplitude $A$ of the acoustic waves is not reported, but can reasonably 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$ (resp. $y$) component of the dimensional streaming velocity at $x=3\pi/4$ (resp. $x=\pi/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_{*}=353~{}\mathrm{m}\cdot\mathrm{s}^{-1}$ and $\bar{v}_{1}$ by $\sqrt{\epsilon}ha_{*}$.

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 (2008) 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. (2014) and systematically extended in the present study.

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}_{*}$ to the diffusive flux,

[TABLE]

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

[TABLE]

so that the Nusselt number is given at leading order by

[TABLE]

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

[TABLE]

Moreover, from equations (17)–(19), the following equality can be derived,

[TABLE]

Finally, using the boundary condition $\bar{\Theta}_{0}=0$ at the top and bottom walls, we obtain

[TABLE]

In the limit of one-way coupling and for $\Gamma=1$, we evaluate this quantity with the approximate expression (51) for $\bar{\Theta}_{0}$ derived in § 4.2,

[TABLE]

This expression should be contrasted with that obtained for Rayleigh streaming (i.e. for sufficiently small temperature differences) and computed by Vainshtein et al. (1995):

[TABLE]

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, $\epsilon\ll 10^{-2}$. Equation (63) 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 (39) suggests that an acoustic wave may be amplified via interaction with the solid boundaries. Indeed, the first term on the right-hand side of (39) 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 $\partial_{y}\hat{\Theta}_{1}(x,y=1)-\partial_{y}\hat{\Theta}_{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{\Theta}_{1}$ can be determined from the acoustic wave solution derived in § 4.1 and then used to evaluate

[TABLE]

Since $\gamma>1$, this term is negative and damps the acoustic waves, unless the thermal forcing $\Gamma$ 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.

5 Two-way coupling

The analytical solution derived in § 4 and used to compute the heat flux enhancement is not valid when there is two-way coupling between the acoustic waves and the streaming flow, as occurs for sufficiently large values of the dimensionless parameter groups $A^{2}Re_{s}h^{2}$ and $A^{2}Re_{s}Pe_{s}h^{4}$. To investigate the dynamics in this regime, we perform numerical simulations of the reduced equations derived in § 2 and § 3 using the spectral computing environment Dedalus (Burns et al., 2018).

The spatial domain is discretized on a $512\times 128$ ($x\times y$) grid using a Fourier–Chebyshev pseudospectral scheme. Temporal integration is performed strictly on the slow time scale $T$ with a 2nd-order Runge–Kutta method. At each time step, the shape and angular frequency of the standing acoustic-wave mode are obtained by numerically solving the eigenvalue problem (32) derived in § 3.1 using a Chebyshev spectral method and global QZ algorithm. The amplitude of the acoustic mode $A$ either is taken to be a prescribed time-independent value or evolved according to the amplitude equation with an additional constant input power $P$, resulting in a term $P/\overline{E}_{K}$ on the right-hand side of (39). In the first scenario, the forcing must be imagined to be tuned continuously to maintain the constancy of the wave amplitude, and the amplitude equation can then be utilised to self-consistently evaluate the power $P$ required. In the second scenario, in which the input power $P$ is held constant, (39) must be co-evolved numerically with the streaming equations.

Depending upon the parameter regime, numerical integration of (15)–(19) can result in a grid-scale instability manifesting as high-wavenumber variability in the $x$ direction. To remedy this difficulty, we regularize the streaming dynamics simply by including the next-order terms in $\epsilon$ that, while formally small, introduce a viscous term of the form $\partial_{xx}\bar{u}_{1}$ in (15) and eliminate the numerical instability.

Here, we describe the the numerical simulations performed at constant amplitude $A$. The other parameters are chosen as follows:

[TABLE]

Although not depicted here, for small values of $A$ the resulting steady-state velocity field agrees quantitatively with the one computed in § 4. As $A$ is increased, however, the cellular streamline pattern evident in figure 3 changes into one dominated by narrow jet-like flow structures; see figure 5, which shows the steady-state fields obtained for $A=2.8$. Unlike streaming at high $Re_{s}$ in the absence of background density stratification (Vainshtein et al., 1995), symmetry with respect to the mid-plane $y=1/2$ is broken. Moreover, for this simulation, the angular frequency of the acoustic waves changes by more than $3\%$, providing quantitative evidence of the two-way coupling between the waves and the streaming flow.

The heat-flux enhancement factors achieved by steady streaming states obtained for various values of $A$ (and computed using (58)) are plotted in figure 6. For small values of $A$, the results follow the theoretical prediction derived in § 4, but ultimately deviate as the streaming flow develops localized jets. Even in the regime of strong coupling, the rapid growth of the Nusselt number ($Nu-1\propto A^{4}$) is qualitatively observed, strongly suggesting that significant increases in heat transport may be realized for sufficiently large values of $A$ (see § 6).

6 Discussion

Both experiments and direct numerical simulations confirm that an imposed temperature difference or, more generally, any inhomogeneous background density field strongly affects both the pattern and intensity of acoustic streaming. To correctly predict the streaming dynamics, it is thus necessary to properly account for these density gradients, as in the work of Chini et al. (2014) and in the present study. The striking effect of the inhomogeneous density field on the streaming originates from the baroclinic generation of wave vorticity that occurs even in the absence of viscous torques. This vorticity generation mechanism is readily inferred by taking the curl of the linearised Euler equations describing the leading-order acoustic wave dynamics,

[TABLE]

and from the instantaneous acoustic-wave velocity field shown in figure 2. As already emphasized in Chini et al. (2014), this mechanism enhances the transfer of energy from the waves to the streaming flow; accordingly, this regime is appropriately termed baroclinic acoustic streaming.

A quantitative understanding of baroclinic acoustic streaming is necessary for ascertaining the extent to which acoustics can be used to improve the transport and mixing of heat or of any dense or light solute chemical species. The separation in time scales between the period of the waves and the dynamics of the streaming flow, typically $\gtrsim 10^{3}$, renders multiple scale analysis very attractive if not essential. The analysis reported here has been shown to accurately describe acoustic streaming in a thin channel across which a temperature difference is imposed: the explicit solution derived in the limit of one-way coupling quantitatively fits the data of a previous direct numerical simulation of the instantaneous governing equations (cf. § 4.3). In the same limit, our analysis enables us to readily compute the heat flux enhancement engendered by the streaming and, in particular, to identify a very strong dependence of $Nu$ on the dimensionless channel height $h$: $Nu-1\propto h^{8}$ (see (63)). This result clearly indicates that, to enhance the heat flux in a channel, the frequency of the acoustic wave should be tuned so that the aspect ratio $\delta$ is of order unity.

The description of the $\delta=O(1)$ scaling regime also should be amenable to multiple scale analysis, an extension we are currently pursuing. In this regime, baroclinic streaming continues to dominate viscous streaming, yet it is possible to more naturally account for viscous losses in the oscillatory boundary layers. Indeed, whereas the streaming flow induced by viscous dissipation in the Stokes boundary layers can be shown to be of higher order in the reduced system (Chini et al., 2014), an order of magnitude estimate of the dimensional time scale $\tilde{\tau}_{BL}$ over which the energy dissipation in the boundary layers damps the waves is given by

[TABLE]

where $\delta_{BL}\sim\sqrt{\mu_{*}/(\rho_{*}a_{*}k_{*})}$ is the thickness of the Stokes layers. Therefore, viscous dissipation in the oscillatory boundary layers takes place on a time scale that falls intermediate to the slow and the fast scales, and thus cannot be easily introduced within the framework of the small aspect-ratio analysis. Here, these viscous losses are implicitly presumed to be offset by the wave forcing mechanism. (Note that in classical Rayleigh streaming, i.e. in a homogeneous medium, the dissipation time scale also lies between the wave period and the time required for steady streaming to be established.) Practically, this results in a difference between the heat flux at the upper and lower walls, and this difference is the input power necessary to sustain the out-of equilibrium state (Lin & Farouk, 2008). With the scaling $\delta=O(1)$, however, viscous dissipation in the Stokes layers occurs on the streaming time scale ($\tilde{\tau}_{BL}\omega_{*}\propto 1/\epsilon$) and therefore can be readily incorporated into the analysis.

The present study extends the investigation of Chini et al. (2014), in which a similar multiple scale analysis is performed to predict acoustic streaming in a thin channel with cold top and bottom boundaries and an imposed volumetric heat source, in three significant ways. First, we have derived an analytical solution for baroclinic acoustic streaming in the limit of one-way coupling and $\mathit{O}(1)$ Prandtl number that quantitatively accords with prior direct numerical simulations. Secondly, we performed numerical simulations of the reduced system that captures the feedback from the streaming to the acoustic wave field, i.e. two-way wave/mean-flow coupling, as manifested by the deviation of the resulting heat flux enhancement from that predicted using the analytical solution. In this stronger streaming regime, there is a transition from a rather smoothly-varying cellular streaming flow for which the wave-induced Reynolds stress divergence is balanced by the mean viscous force to a flow exhibiting jet-like structures and thin (streaming) boundary layers, akin to the transition from the Rayleigh to Stuart streaming regimes. The detailed flow pattern, which exhibits only upward jets, differs from that realized at large $Re_{s}$ in the absence of stratification. Finally, perhaps the primary advance of the present study is that the dynamics of a multiple time-scale quasilinear wave/mean-flow system has been reduced to dynamical evolution strictly on the slow time scale: the spatial structure of the wave field is determined at each coarse time step from the solution of a 1D eigenvalue problem (§ 3.1), while the modal amplitude evolves according to (39).

This reduction yields immediate computational and theoretical advantages. Indeed, by numerically integrating the streaming equations (15)–(19) together with the wave amplitude equation (39) and eigenvalue problem (32), the instantaneous dynamics need not be simulated using supercomputing resources (e.g. as in Loh et al. (2002)) nor approximated by alternatively time-advancing the wave dynamics and the streaming flow (Karlsen et al., 2018). Our algorithm thereby enables accurate and inexpensive numerical simulations over several thousand acoustic-wave periods to be performed in a regime where the waves and streaming flow are strongly coupled. For example, in figure 7, we plot the evolution of the Nusselt number $Nu$ and the wave amplitude $A$ for baroclinic acoustic streaming with a small (dimensionless) constant input power $P=10^{-3}$. Evidently, a steady state is not achieved until $T\approx 1000$, a time period corresponding to roughly $T/\epsilon=10^{5}$ acoustic-wave cycles! Even using modern computational capabilities, direct numerical simulation of the instantaneous Navier–Stokes equations for this scenario would be prohibitively expensive. (For comparison, Lin & Farouk (2008) were able to run their case 1C simulations, discussed in § 4.3, for only a few hundered cycles.) This computational challenge highlights the value of the asymptotically-reduced model, including the amplitude equation (39), which obviates the need for directly simulating the fast dynamics. We note that the model of the quasi-biennial oscillation developed by Plumb (1977) has certain similarities with the one derived here for baroclinic acoustic streaming: in the former case, the internal gravity waves are explicitly determined using a WKBJ approximation, ultimately yielding a closed system for the mean (streaming) flow. The analog of the amplitude equation (39), however, has not been derived, with the waves instead being held at constant amplitude. From a theoretical perspective, the derivation of the amplitude equation directly reveals the potential for a thermoacoustic-like instability in which an acoustic wave spontaneously grows in a stably stratified background, thereby increasing the heat flux. It seems reasonable to assert that, without the amplitude equation, this instability would be difficult to anticipate and to mechanistically understand.

To conclude, we emphasize the importance of obtaining a quantitative understanding of the interaction of acoustic waves with a strongly stratified fluid. We are particularly interested in the potential heat-flux enhancement that can be achieved via baroclinic acoustic streaming in the absence of natural convection as a lightweight means of cooling electronics aboard spacecraft. For this purpose, the analysis performed here should be extended to allow for channel aspect ratios of order unity. For Earth-based applications, it also may be necessary to investigate the impact of buoyancy (gravitational) forces, which introduces into the analysis another dimensionless parameter, the Richardson number $Ri=g_{*}/(k_{*}a_{*}^{2})$, where $g_{*}$ is the acceleration of gravity. Using the values $g_{*}=10~{}\mathrm{m}\cdot\mathrm{s}^{-1}$, $a_{*}=333~{}\mathrm{m}\cdot\mathrm{s}^{-1}$ and $k_{*}=100~{}\mathrm{m}^{-1}$, we obtain $Ri=10^{-6}$. Accordingly, if $Ri$ were scaled in proportion to $\epsilon^{3/2}$, the mean-flow dynamics would be modified at leading-order in $\epsilon$; i.e. incorporation of gravity would, indeed, modify the results reported here, implying that quantitative experimental validation of the present study would require micro-gravity or microscale (microfluidic) environments. Another important extension concerns the cooling of a hot object immersed in a fluid (e.g. again in the absence of natural convection) and subjected to acoustic-wave forcing. As for the internal flow configuration studied here, a streaming flow is expected to develop and thereby enhance heat transfer. To date, analysis of this fundamental problem, e.g. for a heated cylinder (first performed by Richardson (1967) and Davidson (1973); see Riley (2001) for additional references), has been implicitly and, indeed, perhaps unknowingly restricted to small temperature anomalies. Again, if the temperature differential between the object and the ambient fluid is sufficiently large, fundamental changes in the streaming dynamics may be anticipated with important ramifications for mixing rates and heat transfer that have yet to be properly explored. Above all, very few experiments of the baroclinic acoustic streaming regime have been performed: well-controlled experiments in micro-gravity or stably stratified systems would be invaluable for quantitative assessment of the theory. Finally, a complementary perspective on acoustic streaming in stratified flows was recently offered by Beisner et al. (2015), who provided evidence that acoustically-induced mixing can suppress combustion in micro-gravity environments. This phenomenon, too, relies on the strong coupling between a density stratified fluid and an acoustic wave.

Acknowledgements

The authors are grateful to Keaton Burns for sharing his expertise on the Dedalus computing environment and for support from the 2017 Woods Hole Summer Program in Geophysical Fluid Dynamics (Woods Hole, MA), where much of this work was completed.

Appendix A

In this appendix, we report the dimensionless set of equations for the waves at second order. (The leading order equations are given in (20)–(23); also see Chini et al. (2014).) In these equations, the linear operator acting on the $\mathit{O}(\epsilon^{2})$ fluctuation fields is identified on the left-hand side and the inhomogeneous forcing arising from the leading-order terms is collected on the right-hand side.

The $x$-component of the momentum equation is

[TABLE] 2. 2.

The $y$-component of the momentum equation is

[TABLE] 3. 3.

Conservation of mass requires

[TABLE] 4. 4.

The internal energy balance is

[TABLE] 5. 5.

The equation of state is

[TABLE]

Appendix B

We can reduce the set of equations given in Appendix A to (37)–(38). This procedure is strictly analogous to the corresponding manipulations performed at leading order. The imaginary parts of $\mathcal{F}$ and $\mathcal{G}$, respectively $\mathcal{F}_{i}$ and $\mathcal{G}_{i}$, are found to be

[TABLE]

and

[TABLE]

where $\delta f/\delta T=(\partial\bar{\rho}_{0}/\partial T)(\delta f/\delta\bar{\rho}_{0})$ is shorthand notation for a differentiation operation involving the functional derivative $\delta f/\delta\bar{\rho}_{0}$. (Strictly, $\delta f/\delta T$ should read $\mathrm{d}f/\mathrm{d}T$, but the $\delta$-notation is used to emphasize that functional differentiation is involved.) $\mathcal{H}_{r}$ is defined by

[TABLE]

Appendix C

In this appendix, we derive the solvability condition for the system (37)–(38). In the vector space $(\mathbb{R}^{2}\rightarrow\mathbb{C})^{2}$, a vector

[TABLE]

The linear operator

[TABLE]

so the set of equations (37)–(38) is of the form $\mathcal{L}\mathbf{V}=\mathbf{F}$. Given two vectors $\mathbf{V}_{A}$ and $\mathbf{V}_{B}$, we define their scalar product

[TABLE]

where $\hat{f}^{*}$ stands for the complex conjugate of $\hat{f}$. Given the $2\pi$-periodicity requirement in $x$ and the kinematic boundary conditions in $y$, the operator $\mathcal{L}$ is self-adjoint ($\mathcal{L}^{\dagger}=\mathcal{L}$); i.e.

[TABLE]

Therefore, the vectors in the kernel of $\mathcal{L}$ consist of the $\mathit{O}(\epsilon)$ acoustic-wave modes already described. Since one of these modes $\mathbf{V}_{1}=(\hat{u}_{1},\hat{v}_{1})$, then

[TABLE]

The imaginary part of (81) requires

[TABLE]

recalling that both $\hat{u}_{1}$ and $\hat{v}_{1}$ are real-valued (so $\hat{u}_{1}=\hat{u}_{1}^{*}$ and $\hat{v}_{1}=\hat{v}_{1}^{*}$).

Appendix D

Evaluation of the solvability condition (82) requires functional forms for $\hat{u}_{1}$ and $\hat{v}_{1}$, as detailed in (30)–(32), and $\mathcal{F}_{i}$ and $\mathcal{G}_{i}$, whose explicit expressions are reported in (B)–(75). Obtaining (39) from (82) is a lengthy but straightforward calculation, except for the treatment of the functional derivatives. Certain functional derivatives can be simplified by exploiting the normalization condition; we illustrate this technical point through the computation of two integrals necessary to derive (39),

[TABLE]

and

[TABLE]

For a real-valued functional $J$ of a function $j$, the functional derivative $\delta J/\delta j$ describes how $J(j)$ evolves when $j\rightarrow j+w$ with $|w|\rightarrow 0$:

[TABLE]

In the present work, the functions we consider ($\hat{u}_{1}$, $\hat{v}_{1}$, etc.) have a functional dependence on $\bar{\rho}_{0}$ and also depend on the real parameters $x$ and $y$. Specifically, we are interested in how the quantity $\hat{f}(x,y;\bar{\rho}_{0})$ evolves at a fixed position $x$ and $y$, while $\bar{\rho}_{0}$ evolves slowly in time. Formally, we therefore should define a functional $\hat{f}_{x,y}$ such that

[TABLE]

and the shorthand $\delta\hat{f}/\delta T$ should therefore be understood as

[TABLE]

The first integral $I_{1}$ is thus given by

[TABLE]

According to the normalization condition (35), each of the integrals on the right-hand side of (88) is equal to unity so that $I_{1}=0$.

Similarly, the second integral is expressed as a function of $g$ with (30),

[TABLE]

The various contributions to (90) can be simplified using the function $\alpha$ introduced in (33) since only $\bar{\rho}_{0}$ depends on $y$:

[TABLE]

The normalization condition (35) with the constitutive relation for $g$ (32) yields

[TABLE]

This result can be used to compute the first integral in (91),

[TABLE]

Finally, we obtain

[TABLE]