Passive scalar transport in Couette flow

Abstract

A scaling theory for the passive scalar transport in Couette flow, i.e. the flow between two parallel plates moving with different velocities, is proposed. This flow is determined by the bulk Reynolds number $Re_b$ and the Prandtl number $Pr$. In the turbulent regime, for moderate shear Reynolds number $Re_{\tau }$ and moderate $Pr$, we derive that the passive scalar transport characterised by the Nusselt number $Nu$ scales as $Nu \sim Pr^{1/2}Re_{\tau }^{2}Re_b^{-1}$. We then use the well-established scaling for the friction coefficient $C_f \sim Re_b^{-1/4}$ (corresponding to a shear Reynolds number $Re_{\tau } \sim Re_b^{7/8}$) which holds reasonably well within the range $3\times 10^{3} \leqslant Re_b \leqslant 10^{5}$, to obtain $Nu \sim Pr^{1/2}Re_b^{3/4}$ for the Nusselt number scaling. The theoretical results are tested against direct numerical simulations of Couette flows for the parameter ranges $81 \leqslant Re_b \leqslant 22361$ and $0.1 \leqslant Pr \leqslant 10$, finding good agreement. Analyses of the numerically obtained turbulent flow fields confirm logarithmic mean wall-parallel profiles of the velocity and the passive scalar in the inertial sublayer.

1. Introduction

The transport of passive scalars in shear-driven turbulent flows (Warhaft 2000) is highly relevant to various natural phenomena, for example, in oceanic mixed layers (Large, McWilliams & Doney 1994; Kantha & Clayson 2004) and weakly stable stratified atmospheric boundary layers (Deusebio, Caulfield & Taylor 2015). Understanding turbulent mixing in such phenomena is of major importance for environmental and meteorological processes.

Early experimental work on turbulent passive transport of temperature by Subramanian & Antonia (1981) and Nagano & Tagawa (1988) has shown the existence of characteristic logarithmic profiles of temperature and streamwise velocity for flow in weakly heated channels. These papers also studied the profiles of higher-order statistics. Moreover, Nagano & Tagawa (1988) showed that the probability distribution function of the temperature and streamwise velocity fluctuations in the log-layer are well approximated by Gaussian probability distribution functions. The empirically found log-law of the passive temperature by Kader (1981), based on various experimental data, has proven to be quite successful. Kays & Crawford (1993) not only provide an extensive, in-depth overview of heat and momentum transport in boundary layers, channel flows and pipe flows, but also propose a conduction model for the turbulent Prandtl number with further improvements proposed by Weigand, Ferguson & Crawford (1997) whereas Jischa & Rieke (1979a,b) propose a model based on transport equations for the turbulent kinetic energy and for the turbulent heat flux. Obviously, various efforts to describe the kinetic and thermal boundary layers in shear flows and their (relative) thicknesses have been made from early on, dating back to Prandtl (1910, 1925, 1932) and von Kármán (1921, 1934). Excellent summaries of both the experimental data and the theoretical work are given in the famous textbook on boundary layers by Schlichting & Gersten (2016).

More recently there have been various numerical studies of passive scalar transport in wall-bounded flows. Direct numerical simulations (DNS) of channel flows have revealed the existence of streaky structures in the log-layer (Kim & Moin 1989; Kawamura, Abe & Matsuo 1999, 2004; Debusschere & Rutland 2004; Schwertfirm & Manhart 2007) and large streamwise vortical structures in the bulk which dominate the heat transport (Debusschere & Rutland 2004) and lead to the nonlinear effects (Kawamura et al. 2004). Large eddy simulations using coarse grids and filtered flow fields for the passive scalar have shown that the small scales hardly influence the macroscopic turbulence (Robert & Tiselj 2006). It has been observed that the mean profiles of velocity and passive scalar can be accurately predicted even if the flow is not fully resolved down to the Kolmogorov and Batchelor scales (Bergant & Tiselj 2007). This suggests that the scalar flux is also dominated by the large-scale circulation and not the smaller-scale fluctuations. Similar observations made by Na & Hanratty (2000) show that the high wavenumber fluctuations are damped close to the wall.

It has also been observed that the fluctuations of streamwise velocity and temperature are highly correlated in the boundary layer (Kim & Moin 1989; Na & Hanratty 2000; Tiselj et al. 2001; Liu 2003; Schwertfirm & Manhart 2007; Antonia, Abe & Kawamura 2009; Pirozzoli, Bernardini & Orlandi 2014, 2016). The energy budgets for turbulent production, dissipation and Reynolds stresses have also been well studied (Kim & Moin 1989; Lyons, Hanratty & McLaughlin 1991; Papavassiliou & Hanratty 1997; Wikström & Johansson 1998; Kawamura et al. 1999; Debusschere & Rutland 2004; Schwertfirm & Manhart 2007). For channel flows, Pirozzoli et al. (2016) noted that in the log-layer the production of turbulent kinetic energy is larger than the dissipation. In contrast, in the bulk, the dissipation is larger than the production (Liu 2003), which leads to a transport of turbulent kinetic energy out of the log-layer into the bulk. With the availability of more computational resources, an exploration of channel flows at larger shear Reynolds numbers has revealed that the local dominance of production over dissipation leads to large-scale eddies and an additional term in the log-layer of the scalar field which scales linearly with outer wall distance (Pirozzoli et al. 2016). A comparison between the transport of passive scalars in channel flow and in Couette flow reveals similarities in aspects of flow organisation and flow structures in the bulk as well as in the log-layer (Kawamura, Abe & Shingai 2000; Liu 2003; Debusschere & Rutland 2004). Indeed, the mechanism of scalar transport in both flows is driven by large streamwise vortices (Debusschere & Rutland 2004) and the scalar profiles are very similar (Kawamura et al. 2000). However, in contrast to the channel flow where the turbulent shear stress vanishes at the channel centre, the production term of turbulent kinetic energy is non-zero at the channel centre for Couette flow (Liu 2003) due to the non-zero turbulent shear stress at the mid-height. Therefore, the scalar transport in Couette flows is up to 20 % greater than the scalar transport of the channel flow (Debusschere & Rutland 2004). The value of the von Kármán constant for the scalar profiles is also case dependent and, a priori, cannot be assumed to be universal.

In this paper, we will focus on obtaining the dependence between the transport of passive scalar and momentum transport in Couette flow. The control parameters for the system are the bulk Reynolds number

(1.1)\begin{equation} Re_b \equiv {U_b H}/{\nu}, \end{equation}

characterising the strength of the kinetic driving due to the bottom plate (at $z=0$) moving with speed $(-U_b)$ and the top plate (at distance $z=H$) with speed $U_b$, and the Prandtl number

(1.2)\begin{equation} Pr \equiv {\nu}/{\kappa}, \end{equation}

as the ratio between the kinematic viscosity $\nu$ and the scalar diffusivity $\kappa$ (i.e. a material property of the fluid)

The crucial global response parameters of that system are the Nusselt number ($Nu$), which is the non-dimensionalised scalar flux from the bottom plate to the top plate,

(1.3)\begin{equation} Nu \equiv \frac{ Q H }{\kappa\varDelta }, \end{equation}

with $Q$ being the scalar flux, and the shear Reynolds number

(1.4)\begin{equation} Re_{\tau} \equiv \frac{u_{\tau} H}{2\nu}, \end{equation}

where $u_{\tau } \equiv \sqrt {{\tau _w}/{\rho }} = \sqrt {\nu \left \langle \left.\partial _z u_x\right |_{z=0}\right \rangle _{A,t}}$ is the friction velocity, with $\tau _w$ being the mean wall shear stress, $\rho$ being the density, $u_x$ being the streamwise component of the velocity, $z$ being the wall-normal distance from the bottom wall and $\left \langle \dots \right \rangle _{A,t}$ indicating the mean over time and a wall-parallel plane $A$ at a distance $z$ from the bottom wall (here, $z = 0$). Alternatively, the same response information that is expressed in the shear Reynolds number $Re_\tau$, can also be expressed in terms of the friction coefficient

(1.5)\begin{equation} C_f \equiv \frac{2\tau_w}{\rho U_b^{2}} = 8 \,\frac{ Re_\tau^{2}}{Re_b^{2}} . \end{equation}

The key question is: How do the global response parameters $Nu$ and $Re_\tau$ depend on the control parameters $Re_b$ and $Pr$? That is, we wish to understand the dependencies

(1.6a,b)\begin{equation} Nu(Re_b, Pr) \quad \hbox{and} \quad Re_\tau (Re_b, Pr) . \end{equation}

Na, Dimitrios & Hanratty (1999) found an effective scaling of $Nu \sim Pr^{0.464}$ at $Re_{\tau } \approx 150$. For large Prandtl numbers $Pr > 700$, the effective scaling changes to $Nu \sim Pr^{0.297}$ (Shaw & Hanratty 1977). However, these scaling relations for $Pr$ are obtained empirically through the analysis of numerical and experimental data, respectively, and lack a sound theoretical explanation. In this current work, we derive the Nusselt number dependence through scaling arguments and show its validity through numerical simulations. The ranges of $Pr$ and $Re_{\tau}$ for the present DNS are shown in figure 1.

Figure 1. Parameter space of experiments and simulations on scalar transport taken from various literature sources. The current simulations shown with the black star markers are performed for the range of bulk Reynolds numbers $81 \leqslant Re_b \leqslant 22\,361$ corresponding to shear Reynolds numbers of $7.07 < Re_{\tau } < 546$ using Prandtl numbers $Pr = 0.1$, 0.3, 1.0, 3.0 and 10.0. The other colours represent other experiments and numerical studies from literature: [1] Subramanian & Antonia (1981), [2] Kader (1981), [3] Nagano & Tagawa (1988), [4] Kim & Moin (1989), [5] Lyons et al. (1991), [6] Papavassiliou & Hanratty (1997), [7] Wikström & Johansson (1998), [8] Kawamura et al. (1999), [9] Na et al. (1999), [10] Kawamura et al. (2000), [11] Na & Hanratty (2000), [12] Tiselj et al. (2001), [13] Liu (2003), [14] Schmitt (2003), [15] Kawamura et al. (2004), [16] Debusschere & Rutland (2004), [17] Robert & Tiselj (2006), [18] Schwertfirm & Manhart (2007), [19] van den Berg et al. (2007), [20] Antonia et al. (2009) and [21] Pirozzoli et al. (2016). The star markers correspond to studies of Couette flow, the plus markers correspond to studies of channel flow and the filled circle markers correspond to the studies of pipe flow. A black circular outline around the marker represents an experimental study.

The manuscript is organised in the following way. In § 2 we present the underlying dynamical equations and develop the scaling for passive scalar transport. Section 3 briefly presents the numerical method. In § 4, we check that the numerical results agree with the assumptions made in § 2. Section 5 consists of the comparison between the numerical results and the scaling laws obtained from the theoretical predictions in § 2. Section 6 focuses on the mean velocity and temperature profiles. Section 7 demonstrates that $Nu$ and $Re_{\tau }$ for case of the sheared Rayleigh–Bénard system in the limit of the Rayleigh number going to zero converge to the $Nu$ and $Re_{\tau }$ scaling laws of passive scalar transport. The paper ends with conclusions in § 8.

2. Underlying equations and theory

2.1. Underlying equations

Passive scalar transport in incompressible Couette flow is governed by the Navier–Stokes equation for the velocity field $\boldsymbol {u}({{\boldsymbol {x}}},t) \equiv (u_x,u_y,u_z)$, the continuity equation and the transport equation for the scalar $\theta ({{\boldsymbol {x}}},t)$. In Cartesian coordinates $\boldsymbol {x}\equiv (x,y,z)\equiv (x_1,x_2,x_3)$, they read

(2.1)$$\begin{gather} {\partial_t u_i} + u_j {\partial_j u_i} ={-} {\partial_i p} + \nu {\partial^{2}_j u_i}, \end{gather}$$

(2.2)$$\begin{gather}{\partial_i u_i} = 0, \end{gather}$$

(2.3)$$\begin{gather}{\partial_t \theta} + u_j {\partial_j \theta} = \kappa {\partial^{2}_j \theta}, \end{gather}$$

where $p({{\boldsymbol {x}}},t)$ is the kinematic pressure and $\theta ({{\boldsymbol {x}}},t)$ the passive scalar with the arithmetic mean of the top and bottom wall values subtracted. We consider $x,y$ to be the wall-parallel directions and $z$ to be the wall-normal direction. The boundary conditions for the equations are

(2.4a–d)$$\begin{gather} \left. \theta \right|_{z=0} = \varDelta/2, \quad \left. \theta \right|_{z=H} ={-}\varDelta/2, \quad \left. u_x \right|_{z=0} ={-}U_b, \quad \left. u_x \right|_{z=H} = U_b, \end{gather}$$

(2.5a,b)$$\begin{gather}\left. u_y \right|_{z=0} = \left. u_y \right|_{z=H} = 0, \quad \left. u_z \right|_{z=0} = \left. u_z \right|_{z=H} = 0. \end{gather}$$

The scalar flux $Q$ from the bottom plate to the top plate is

(2.6)\begin{equation} Q= \left\langle u_z \theta \right\rangle_{A,t} - \kappa\left\langle \partial_z \theta \right\rangle_{A,t}. \end{equation}

2.2. Laminar case

In the laminar regime, steady uniform gradients of the streamwise velocity and temperature develop in the flow. Therefore, we have the profiles of the streamwise velocity and scalar given by

(2.7a,b)\begin{equation} u_x = U_b \left(\frac{2z}{H}-1\right), \quad \theta = \frac{\varDelta}{2} \left(1-\frac{2z}{H}\right), \end{equation}

with $z$ as the wall-normal coordinate which lies between $z=0$, corresponding to the bottom plate, and $z=H$, corresponding to the top plate. This gives the well-known and trivial laminar results ${\partial _z u_x} = 2U_b/H$, $\tau _w = 2\rho \nu {U_b}/{H}$, or in dimensionless form

(2.8a,b)\begin{equation} Re_{\tau} = \frac{1}{\sqrt{2} } Re_b^{1/2}, \quad C_f = 4 Re_b^{{-}1}, \end{equation}

and ${\partial \theta }/{\partial z} = -\varDelta /{H}$, $Q = \kappa \varDelta /{H}$, or in dimensionless form

(2.9)\begin{equation} Nu = 1. \end{equation}

2.3. Turbulent regime

From (2.1), we use the Reynolds (1895) averaging to split the velocity and scalar into the mean and fluctuating values as follows:

(2.10a,b)\begin{equation} u_i = \left\langle u_i\right\rangle_{t} + u_i^{\prime}, \quad \theta = \left\langle\theta\right\rangle_{t} + \theta^{\prime}, \end{equation}

with $\left \langle \dots \right \rangle _{t}$ indicating the Reynolds averaging. Assuming that $x$ is the streamwise direction and using (2.1) averaged in time and over the wall-parallel directions, we obtain

(2.11)\begin{equation} \partial_z \left(-\left\langle u_x^{\prime} u_z^{\prime} \right\rangle_{A,t} + \nu \partial_z \left\langle u_x\right\rangle_{A,t}\right) = 0, \end{equation}

which simply reflects momentum conservation in the $z$ direction. Integrating this relation $z$ between $z=0$ and $z$ gives the well-known Reynolds-averaged Navier–Stokes equation (Pope 2000)

(2.12)\begin{equation} -\left\langle u_x^{\prime} u_z^{\prime} \right\rangle_{A,t} + \nu \partial_z \left\langle u_x\right\rangle_{A,t} = \left.\nu\partial_z \left\langle u_x\right\rangle_{A,t}\right|_{z=0} \equiv u_{\tau}^{2}. \end{equation}

Here, $-\left \langle u_x^{\prime } u_z^{\prime } \right \rangle _{A,t}$ represents the Reynolds shear stress and $\nu \partial _z \left \langle u_x\right \rangle _{A,t}$ represents the viscous shear stress. These quantities can be viewed as the convective and diffusive momentum fluxes. This equation means that the total momentum flux given by the sum of convective and diffusive fluxes at any height equals the diffusive momentum flux at the wall in the case of statistically steady flow, i.e. the diffusive flux at the wall is the bottleneck for the total flux through the system. At the wall, we introduce a characteristic length scale in terms of the kinetic boundary layer thickness $\lambda _u$ (Prandtl 1904; Blasius 1908) such that

(2.13)\begin{equation} \left.\nu\partial_z \left\langle u_x\right\rangle_{A,t}\right|_{z=0} = u_{\tau}^{2} = \nu U_b/\lambda_u. \end{equation}

The diffusive time scale of the momentum flux at the wall is then $\lambda _u^{2} / \nu$. Since the convective heat transfer through turbulent motions in the bulk is much quicker than the diffusive heat transfer at the wall, this diffusive time scale represents the bottleneck for the for momentum transport through the system.

Similarly, the Reynolds-averaged form of (2.3) yields

(2.14)\begin{equation} \partial_z \left(-\left\langle \theta^{\prime} u_z^{\prime} \right\rangle_{A,t} + \kappa \partial_z\left\langle\theta\right\rangle_{A,t}\right) = 0, \end{equation}

which when integrated between the limits $z=0$ and $z$ gives

(2.15)\begin{equation} -\left\langle \theta^{\prime} u_z^{\prime} \right\rangle_{A,t} + \kappa \partial_z \left\langle \theta\right\rangle_{A,t} = \left.\kappa\partial_z \left\langle \theta\right\rangle_{A,t}\right|_{z=0} \equiv {\theta_Q}^{2}, \end{equation}

with ${\theta _Q}$ defined analogously to $u_{\tau }$ (Schlichting & Gersten 2016). This equation again reflects that the total scalar flux through the system is determined by the scalar flux at the wall. Therefore, we introduce another characteristic length scale in terms of the thermal boundary layer thickness $\lambda _{\theta }$ such that

(2.16)\begin{equation} {\theta_Q}^{2} = \kappa\varDelta/\lambda_{\theta} = (2\kappa\varDelta/H)Nu. \end{equation}

The diffusive time scale of heat flux at the wall is then $\lambda _{\theta }^{2}/\kappa$, which must be the bottleneck time scale of the system, as the diffusive scalar flux at the wall is the bottleneck for the total scalar flux through the system.

However, since the scalar is a passive quantity, the time scale of the system should be solely governed by the diffusive time scale of the momentum flux at the wall, in which case

(2.17)\begin{equation} \lambda_u^{2}/\nu \approx \lambda_{\theta}^{2}/\kappa. \end{equation}

Postulating this gives,

(2.18)\begin{equation} \lambda_u/\lambda_{\theta} \approx Pr^{1/2}. \end{equation}

Using (1.4), (1.5), (2.13) and (2.16) results in the relation

(2.19)\begin{equation} Nu \approx 2Pr^{1/2} Re_{\tau}^{2} Re_b^{{-}1} = \frac{C_f}{4} Pr^{1/2} Re_b. \end{equation}

Note that the relation (2.19) is similar to the Reynolds (1874) analogy for turbulent boundary layer over a heated flat plate (Schlichting & Gersten 2016) and turbulent flow in heated pipes (Kays 1994; McEligot & Taylor 1996). According to the Reynolds analogy, the shear stress and the heat flux are analogous if the average value of the turbulent Prandtl number $\left \langle Pr_t\right \rangle _L \approx 1$ in the log-layer. Here, the value of the $Pr_t$ is defined at any given point in the flow by

(2.20)\begin{equation} Pr_t = \left. \frac{\left\langle u_x^{\prime} u_z^{\prime} \right\rangle_{t}}{\partial_z \left\langle u_x\right\rangle_{t}} \right/ \frac{\left\langle u_z^{\prime} \theta^{\prime} \right\rangle_{t}}{\partial_z \left\langle \theta\right\rangle_{t}} , \end{equation}

and $\left \langle \dots \right \rangle _L$ indicates the average value in the log-layer. For alternate definitions used to compute $Pr_t$ from the numerical simulations, we refer the reader to §4. Therefore, from the Reynolds analogy,

(2.21)\begin{equation} Nu \approx f(Pr) C_f Re_b. \end{equation}

The comparison of (2.19) and (2.21) thus confirms that the Reynolds analogy is also applicable to the case of passive transport in Couette flow with the value of $f(Pr) \approx Pr^{1/2}/8$.

For the scaling of the shear Reynolds number with the bulk Reynolds number in the turbulent regime, we use the empirical relation $C_f\sim Re_b^{-1/4}$. This relation was first suggested by Blasius (1913) for pipe flows in an intermediate range of $3 \times 10^{3} \leqslant Re_b \leqslant 10^{5}$ and is widely accepted as a good representation of the friction coefficient of pipe flows for the aforementioned range of bulk Reynolds numbers (McKeon, Zargola & Smits 2005). This scaling relation has since been well established through experiments of flows in non-circular pipes (Nikuradse 1930, 1950), rectangular ducts (Dean 1978) and channel flows (Schultz & Flack 2013). Also, the DNS results by Bernardini, Pirozzoli & Orlandi (2014) and by Orlandi, Bernardini & Pirozzoli (2015) show that this power law is a good description of the friction coefficient in turbulent Poiseuille and Couette flows for the considered $Re_b$ range.

Using the scaling relation $C_f\sim Re_b^{-1/4}$ with relation (1.5), we get

(2.22)\begin{equation} Re_{\tau} \sim Re_b^{7/8}. \end{equation}

Since $C_f$ is independent of $Pr$, combining relations (2.22) and (2.19) gives

(2.23)\begin{equation} Nu \sim Pr^{1/2} Re_b^{3/4}. \end{equation}

The scaling relation $C_f\sim Re_b^{-1/4}$ (i.e. $Re_{\tau } \sim Re_b^{7/8}$) holds well for the range of $Re_b$ considered here but it is not valid at high bulk Reynolds numbers. The value of the friction coefficient deviates more and more from the empirical power law for $Re \geqslant 10^{5}$ (McKeon et al. 2005; Scheel, Emran & Schumacher 2013). We note that, for extremely large $Re_b$, it might be more appropriate to use the logarithmic law to obtain the friction coefficient (similar to the Prandtl (1932) turbulent friction law). However, for the range of $Re_b$ for which we claim the validity of the scaling laws derived in this work, the Blasius law and the Prandtl turbulent friction law show minimal differences and therefore the use of the Blasius scaling is sufficient. At higher $Re_b$, we still expect the scaling given by (2.19) to hold, however with a different scaling for $C_f$.

For extremely large $Re_b$, we expect that the scaling law given by relation (2.23) may not be valid and further studies are required to determine the scaling relations for $Nu$ at very high $Re_b$. Here, we draw an analogy from the ultimate regime of turbulent thermal convection, where the scaling relations of $Nu$ change for extremely large thermal forcing (i.e. Rayleigh Number) (Kraichnan 1962; Spiegel 1971; Chavanne et al. 1997; He et al. 2012; Zhu et al. 2018b).

3. Numerical method

We performed numerical simulations in a domain of dimensions $8H \times 4H \times H$ in the streamwise (along the $x$ axis), spanwise (along the $y$ axis) and wall-normal (along the $z$ axis) directions, respectively. We impose periodic boundary conditions in the wall-parallel directions and no-slip boundary conditions at the top and bottom plates.

The non-dimensional form of the incompressible Navier–Stokes equations (2.1) and (2.3) are integrated numerically using the AFiD GPU package (Zhu et al. 2018a) which uses a second-order finite-difference scheme (van der Poel et al. 2015). The code has been validated and verified many times (Verzicco & Orlandi 1996; Verzicco & Camussi 1997, 2003; Stevens, Verzicco & Lohse 2010; Stevens, Lohse & Verzicco 2011; Ostilla-Mónico et al. 2014; Kooij et al. 2018). We use a uniform discretisation in the wall-parallel periodic directions and a non-uniform grid, with a clipped Chebyshev-like clustering of nodes in the wall-normal direction.

The grid for all cases except for ($Re_b=14\,142, Pr=1$) and ($Re_b=22\,360, Pr=1$) consists of 1536, 768 and 256 nodes in the streamwise, spanwise and wall-normal directions, respectively, whereas the grid for ($Re_b=14\,142, Pr=1$) and ($Re_b=22\,360, Pr=1$) consists of 2048, 1024 and 384 nodes in the streamwise, spanwise and wall-normal directions, respectively, to ensure sufficient resolution. We note that the computational domain used in the current study is smaller than those typically used in studies of Couette flow (Avsarkisov et al. 2014; Pirozzoli et al. 2014; Lee & Moser 2018). However, here, we focus on the global quantities $Nu$ and $Re_{\tau }$ (and thus also $C_f$), for which a smaller domain is sufficient. The domain used is, however, larger than the minimal size mentioned in Sekimoto, Atkinson & Soria (2018) to ensure ‘healthy’ turbulence in the log-layer region. To further confirm that the domain size is sufficient, we have verified that these integral properties obtained from an $8H \times 4H \times H$ domain agree to within 1 % of the corresponding values obtained from test domain of size $48H \times 24H \times H$.

4. Boundary layer thickness and turbulent Prandtl number

Figure 2(a) shows that the ratio of the kinetic and scalar boundary layer thicknesses in the turbulent flow simulations indeed scales as $Pr^{1/2}$ while figures 2(b) and 2(c) show the wall-normal location of the kinetic and scalar boundary layers. It is interesting to note that, for $Pr > 0.1$, the scalar boundary layer lies completely inside the buffer region whereas, for $Pr = 0.1$, the scalar boundary layer overlaps with the log-layer. This results in the minor disagreement with the $Pr$ dependence that can be seen in 2(a) for $Re_b = 10\,000$ and $Pr = 0.1$, suggesting that the predictions from the theory are not accurate for $Pr \ll 1$.

Figure 2. (a) The ratio of the kinetic and scalar boundary layer thicknesses normalised with $Pr^{1/2}$ plotted against bulk Reynolds number; (b) the thickness of the kinetic boundary layer in wall units; (c) the thickness of the scalar boundary layer in wall units.

Figure 3 shows $Pr_t$ obtained from our simulations. The average value of $Pr_t$ at a given distance from the wall can be computed from the numerical simulations using either

Figure 3. Turbulent Prandtl number $Pr_t^{(1)}$ computed using (4.1) and averaged over wall-parallel directions for (a) $Pr=0.1$, (b) $Pr=0.3$, (c) $Pr=1.0$, (d) $Pr=3.0$, (e) $Pr=10.0$. ( f) Turbulent Prandtl number $Pr_t^{(1)}$ and $Pr_t^{(2)}$ computed using (4.1) and (4.2), respectively, and averaged over wall-parallel directions in the log-layer i.e. for $z^{+} \geqslant 30$.

(4.1)\begin{equation} Pr_t^{(1)} = \left. \frac{\left\langle u_x^{\prime} u_z^{\prime} \right\rangle_{A,t}}{\partial_z \left\langle u_x\right\rangle_{A,t}} \right/ \frac{\left\langle u_z^{\prime} \theta^{\prime} \right\rangle_{A,t}}{\partial_z \left\langle \theta\right\rangle_{A,t}} , \end{equation}

or

(4.2)\begin{equation} Pr_t^{(2)} = \left. \left\langle \frac{\left\langle u_x^{\prime} u_z^{\prime} \right\rangle_{t}}{\partial_z \left\langle u_x\right\rangle_{t}}\right\rangle_{A} \right/ \left\langle \frac{\left\langle u_z^{\prime} \theta^{\prime} \right\rangle_{t}}{\partial_z \left\langle \theta\right\rangle_{t}}\right\rangle_{A} . \end{equation}

To obtain the profile for $Pr_t^{(1)}$, we first average the time-averaged profiles of $u_x^{\prime } u_z^{\prime }, \partial _z u_x, u_z^{\prime } \theta ^{\prime }$ and $\partial _z \theta$ over the wall-parallel directions and then compute the $Pr_t$ using (4.1). On the other hand, to obtain the profile for $Pr_t^{(2)}$, we use time-averaged three-dimensional flow fields to compute the ratios $\left \langle u_x^{\prime } u_z^{\prime } \right \rangle _{t} / \partial _z \left \langle u_x\right \rangle _{t}$ and $\left \langle u_z^{\prime } \theta ^{\prime } \right \rangle _{t}/\partial _z \left \langle \theta \right \rangle _{t}$ at each grid point, average them over the wall-parallel directions and then compute the $Pr_t$ using (4.2). For a statistically steady flow, the difference between the values computed using the equations (4.1) and (4.2) should be very small.

The wall-normal variation of $Pr_t^{(1)}$ is shown in figures 3(a)–3(e). From figure 3( f), we see that the mean value of the turbulent Prandtl number in the log-layer ($z^{+} \geqslant 30$) is quite sensitive to the definition used in the computation. The range of $z^{+}$ over which $Pr_t(z^{+})$ is averaged also affects the value of the mean $Pr_t$. The average $Pr_t^{(1)}$ obtained by using (4.1) approaches unity for larger $Re_b$ while for larger viscous Prandtl number ($Pr=3.0, Pr=10.0$), the average $Pr_t^{(2)}$ obtained by using (4.2) drops to lower than unity and approaches 0.85, as suggested in the literature (Kader 1981; Malhotra & Kang 1984; Kays 1994; McEligot & Taylor 1996). Once again, it is noteworthy that the data for $Pr = 0.1$ show a slower convergence of $Pr_t$ than the cases for larger $Pr$, suggesting that the predictions of the theory would become inaccurate for $Pr \ll 0.1$. However, for the control parameters of this paper, and for a large range of wall-normal distances $z^{+}$, we can say that, indeed, $Pr_t \approx 1$ is a reasonable assumption consistent with the numerical data.

5. Global scalar transport and wall shear

Now we investigate the variation of $Nu$ and $Re_{\tau }$ (and thus also $C_f$, or vice versa) with increasing $Re_b$. From figures 4(a) and 4(b) in the turbulent regime, we indeed find the Blasius scaling $C_f \sim Re_b^{-1/4}$ (Blasius 1913; Nikuradse 1930, 1950; Orlandi et al. 2015) (indicated with the black dashed line) from the present numerical simulations as well as from the data taken from Pirozzoli et al. (2014), Orlandi et al. (2015), Avsarkisov et al. (2014) and Lee & Moser (2018), closely approximating the Prandtl (1932) turbulent friction law (von Kármán 1934)

(5.1)\begin{equation} \sqrt{\frac{2}{C_f}} = \frac{1}{k} \ln \left(Re_b \sqrt{\frac{C_f}{8}}\right)+B. \end{equation}

It is shown as black dash dotted line where $k = 0.41$ (Pirozzoli et al. 2014) is the von Kármán (1934) constant and $B = 5$.

Figure 4. (a) Value of $C_f$ vs $Re_b$ for various $Pr$. The $C_f \sim Re_b^{-1/4}$ scaling in the turbulent regime (represented with black dashed line) is consistent with the Prandtl (1932) turbulent friction law (shown as the purple dash dot line). (b) Value of $Re_{\tau }$ vs $Re_b$.

In figures 4(c) and 4(d), the scaling for $Re_{\tau }$ changes sharply with the transition from the laminar to the turbulent regime. In the laminar regime, the effect of steady linear gradients can be seen from the $Re_{\tau } \sim Re_b^{1/2}$ scaling. In the turbulent flow regime, we see a good agreement with the suggested scaling law given by (2.22).

A similar sharp transition is also observed for $Nu$ in figure 5(c). Figure 5(d) shows that $Nu$ varies as $Pr^{1/2}$ for the range of $Pr$ considered in the numerical simulations, which reflects the analogy between turbulent transport of the passive scalar in Couette flow and turbulent transport of heat in Rayleigh–Bénard flow in the limit of zero thermal driving. It can be observed that the values of $Nu$ obtained from the numerical simulations follow the scaling laws given, respectively, by (2.9) and (2.23) quite well. We note that Kays & Crawford (1993) give an empirical relation $Nu = 0.022 Pr^{0.5} Re^{0.8}$ for turbulent pipe flows with gases ($0.5< Pr<1$) which is very close to our fit of $Nu = 0.015 Pr^{0.5} Re^{0.75}$.

Figure 5. (a) Value of $Nu$ vs $Re_b$ for different $Pr$. (b) Value of $Nu$ compensated with $Pr^{1/2}$ vs $Re_b$.

The sharpness of the transition from the laminar regime to the turbulent one in the Couette-type shear flow under consideration here is in vast contrast to the very smooth transition from the laminar regime to the turbulent one in Rayleigh–Bénard flow (Ahlers, Grossmann & Lohse 2009; Stevens et al. 2013). The reason for this difference lies in the different type of flow instability: whereas in Rayleigh–Bénard flow linear instabilities are very crucial (Landau & Lifshitz 1987), in Couette-type shear flow, the onset of turbulence is of nonlinear non-normal type (Trefethen et al. 1993; Grossmann 2000; Barkley 2016; Lemoult et al. 2016; Shi, Avila & Hof 2013).

6. Mean velocity and scalar profiles

We now come to the local flow properties. The existence of a logarithmic inertial sublayer in turbulent wall-bounded flows with passive scalar transport has been studied extensively (Yaglom 1979; Jischa & Rieke 1979a; Kader 1981; Pirozzoli et al. 2016). The log-laws for the mean streamwise velocity and the mean temperature may be written as

(6.1a,b)\begin{equation} u^{+} = \frac{1}{k} \ln (z^{+})+B, \quad \theta^{+} = \frac{Pr_t}{k} \ln (z^{+})+D(Pr,Pr_t). \end{equation}

The values of $z^{+}$, $u^{+}$ and $\theta ^{+}$ are defined as usual as

(6.2a,b)$$\begin{gather} z^{+} = {zu_{\tau}}/{\nu}, \quad u^{+} = u/u_{\tau}, \end{gather}$$

(6.3a,b)$$\begin{gather}\theta^{+} = {\left|\theta_w-\theta\right|}/{\theta_{\tau}}, \quad \theta_{\tau} = Q/u_{\tau}. \end{gather}$$

Here, $\theta _{\tau }$ is a quantity analogous to the friction velocity, $\theta _w$ is the mean value of $\theta$ at the wall and $D(Pr,Pr_t)$ is an a priori unknown offset temperature, which in general, depends on $Pr$ and $Pr_t$.

The value of $D(Pr,Pr_t)$ as a function of $Pr$ and $Pr_t$ as obtained from the numerical simulations is shown in table 1. Kader (1981) proposed the following empirical fitting relation:

(6.4)\begin{equation} D = \left( 3.85Pr^{1/3}-1.3 \right)^{2} + \frac{Pr_t}{\kappa} \ln (Pr), \end{equation}

along with $Pr_t = 0.85$ as a modelling parameter.

Table 1. Comparison between the value of the intercept $D$ for the log law obtained from the numerical simulations and the value computed from the empirical relation from Kader (1981).

Figure 6(a) shows the mean streamwise velocity profile normalised with the friction velocity vs height in wall units for the various $Re_b$ obtained from the numerical simulations. The laminar cases (shown by the blue lines) follow the linear profile whereas the turbulent cases (green–red lines) display the log-layer beyond the viscous sublayer. This is also seen when plotting the diagnostic function in figure 6(b) where we see an appreciably well developed log-layer for the turbulent cases. Figure 6(c) shows the mean scalar profiles where we see a similar family of lines, each associated with a particular value of $Pr$.

Figure 6. (a) Mean streamwise velocity normalised with the friction velocity $(u_{\tau })$ vs height in wall units. (b) The diagnostic function for the mean streamwise velocity. (c) Mean scalar value normalised with $\theta _{\tau }$ vs height in wall units. The value of the intercept of the log-law is obtained from the empirical relation (6.4) given by Kader (1981). (d) The diagnostic function for the mean temperature.

The dependence of the turbulent Prandtl number $Pr_t$ on $Pr$ has been well studied for pipe flows (Malhotra & Kang 1984; McEligot & Taylor 1996), channel flows (Kim & Moin 1989) as well as for stratified Couette flows (Zhou, Taylor & Caulfield 2017; Glazunov et al. 2019) with models suggested by Jischa & Rieke (1979a,b). It is considered that, for $Pr<1$, $Pr_t>1$ and for $Pr>1$, $Pr_t<1$ with $Pr_t \to 0.85$ for large $Pr$ and $Re_b$ (Kays 1994). This observation is consistent with the numerical results shown in figure 3. We note that a constant value of $Pr_t = 0.85$ suggested by Kader (1981) seems to be a good fit for the mean temperature profiles of all the values of $Pr$ considered in this study, as shown by the diagnostic function in figure 6(d).

7. Passive scalar transport in Couette flow as limiting case of heat transport in sheared Rayleigh–Bénard

For weak temperature fluctuations and small thermal driving (i.e. small or even zero Rayleigh number), the transport of heat may be considered as passive transport of temperature (Subramanian & Antonia 1981). Therefore, we may view the passive scalar transport problem as a limiting case of the heat transport in the sheared Rayleigh–Bénard system (Ahlers et al. 2009; Lohse & Xia 2010; Chilla & Schumacher 2012). The sheared Rayleigh–Bénard system consists of the standard Rayleigh–Bénard set-up with an additional Couette-type shear forcing. The strength of the thermal driving in the sheared Rayleigh–Bénard system, due to the temperature difference $\varDelta$ between the hot bottom plate and the cold top plate, is given by the Rayleigh number as

(7.1)\begin{equation} Ra \equiv \frac{\beta g H^{3} \varDelta}{\nu \kappa}. \end{equation}

Here, $\beta$ is the thermal expansion coefficient and $g$ is the gravitational acceleration.

The relative strength of thermal and shear driving is quantified by the Richardson number given by

(7.2)\begin{equation} Ri \equiv \frac{Ra}{Pr Re_b^{2}}. \end{equation}

Viewing the passive transport in Couette flow as the limiting case of the sheared Rayleigh–Bénard system with $Ra \to 0$ and $Ri \to 0$, the numerical simulations and results shown in this paper correspond to

(7.3a,b)\begin{equation} Nu(Re_b, Pr, Ra=0) \quad \mbox{and} \quad Re_\tau (Re_b, Pr, Ra=0) . \end{equation}

For the passive scalar transport in Couette flow, given that there is a temperature difference $\varDelta \ne 0$ between the bottom and top plates, the limit $Ra \to 0$ is achieved by setting the thermal expansion coefficient $\beta \to 0$, or alternatively the gravitational acceleration $g \to 0$. The standard Rayleigh–Bénard case of purely thermally driven convective flow (Ahlers et al. 2009; Lohse & Xia 2010; Chilla & Schumacher 2012) is another limiting case, namely

(7.4a,b)\begin{equation} Nu(Re_b = 0 , Pr, Ra) \quad \mbox{and} \quad Re_\tau (Re_b=0, Pr, Ra) . \end{equation}

For that case, Grossmann & Lohse (2000, 2001, 2002, 2004), Stevens et al. (2013) and Shishkina et al. (2017) have developed a unifying theory, based on the decomposition of the kinetic and thermal dissipation rates into boundary layer and bulk contribution, which very successfully describes the experimental and numerical data for the control parameter dependencies of the Nusselt and wind Reynolds numbers (Ahlers et al. 2009; Stevens et al. 2013). Note that some other limiting cases of problem (1.6a,b) have already been analysed before, namely

(7.5a,b)\begin{equation} Nu(Re_b , Pr = 1, Ra) \quad \mbox{and} \quad Re_\tau (Re_b, Pr=1, Ra), \end{equation}

by Blass et al. (2020) and

(7.6a,b)\begin{equation} Nu(Re_b , Pr, Ra= 10^{6}) \quad \mbox{and} \quad Re_\tau (Re_b, Pr, Ra = 10^{6}), \end{equation}

by Blass et al. (2021). In these two papers, the focus was on the interplay and the competition between thermal driving and shear driving, identifying the transitions from dominance of one to dominance of the other, in order to understand the dependencies (7.5a,b) and (7.6a,b).

With increasing shear forcing as compared with thermal forcing, i.e. $Pr Re_b^{2} \gg Ra$, one approaches $Ri \to 0$, which is the limiting case of passive scalar transport in Couette flows. For such large shear forcing (compared with thermal forcing), we expect the $Nu$ and $Re_{\tau }$ to follow the scaling dependencies derived in this work. From figure 7, it can indeed be observed that $Nu(Re_b , Pr, Ra)$ and $Re_\tau (Re_b, Pr, Ra)$ for the sheared Rayleigh–Bénard system as obtained by Blass et al. (2020, 2021) seem to, in the limit of large shear forcing, asymptotically converge to the $Nu(Re_b , Pr, Ra=0)$ and $Re_\tau (Re_b, Pr, Ra=0)$ scaling laws discussed here in the context of passive scalar transport in Couette flows.

Figure 7. (a) Shear Reynolds number ($Re_{\tau }$) plotted against the bulk Reynolds number $Re_b$. (b) Nusselt number ($Nu$) plotted against the bulk Reynolds number ($Re_b$). The data are taken from by Blass et al. (2020, 2021). In the limiting case of large shear, the scaling relations assumed ($Re_{\tau } \sim Re_b^{7/8}$) or derived ($Nu \sim Re_b^{3/4}$) in this paper are recovered.

8. Conclusions

In conclusion, building on the Blasius's empirical scaling relation $C_f \sim Re_b^{-1/4}$ for the friction coefficient, which holds for intermediate $Re_b$, we have studied the scaling relations for passive scalar transport $Nu(Re_b, Pr)$ in turbulent Couette flow. We then performed DNS for the large control parameter ranges $81\leqslant Re_b \leqslant 22\,360$ and $0.1 \leqslant Pr \leqslant 10$, and identified the laminar and turbulent regimes and the transition between them. These numerical results are in good agreement with the derived scaling laws for the whole range of considered control parameters. Some disagreement is observed in the transitional regime where it is very difficult to obtain well-converged results. We verify the validity of the assumptions made to derive the scaling of $Nu(Re_b, Pr)$ and show that this scaling arises due to the presence of strong streamwise fluctuations close to the wall. The scaling for $Re_{\tau }(Re_b, Pr)$ is taken from the empirical Blasius scaling for the friction coefficient $C_f \sim Re_b^{-1/4}$ for pipe flows, which not only shows minimal deviation from the Prandtl (1932) turbulent friction law but also from experimental results for non-circular channels (Nikuradse 1930, 1950), rectangular ducts (Dean 1978) as well as from numerical studies of channel flows (Bernardini et al. 2014) and Couette flows (Orlandi et al. 2015) for the bulk Reynolds number range $3\times 10^{3} \leqslant Re_b \leqslant 10^{5}$.

Next, the mean profiles of streamwise velocity and temperature in the turbulent regime show the existence of the viscous sublayer close to the wall and the inertial sublayer which follows the log-law (Yaglom 1979; Kader 1981; Jischa & Rieke 1979a), as shown by the respective diagnostic functions.

Obviously, it is desirable to further extend the regime of control parameters to even larger shear Reynolds numbers $Re_b$ and even smaller and larger $Pr$, to observe the possible emergence of new regimes. Such new regimes can also emerge for rough walls, which will enhance the scalar transfer, analogous to the enhanced heat transfer in rough-wall Rayleigh–Bénard convection (Ciliberto & Laroche 1999; Emran & Shishkina 2020; Xie & Xia 2017; Zhu et al. 2017; Jiang et al. 2018; Zhu et al. 2019).