Similarities between the structure functions of thermal convection and hydrodynamic turbulence

TL;DR

This study demonstrates that the velocity structure functions in turbulent thermal convection exhibit scaling similar to hydrodynamic turbulence, with non-Gaussian velocity increment distributions and differing energy flux and dissipation rates.

Contribution

The paper reveals that thermal convection shares similar structure function scaling with hydrodynamic turbulence, supported by high-resolution simulations and analysis of probability distributions.

Findings

Structure functions scale similarly to hydrodynamic turbulence

Velocity increment distributions are non-Gaussian with wide tails

Energy flux in convection differs from viscous dissipation rate

Abstract

In this paper, we analyze the scaling of velocity structure functions of turbulent thermal convection. Using high-resolution numerical simulations, we show that the structure functions scale similar to those of hydrodynamic turbulence, with the scaling exponents in agreement with She and Leveque's predictions [Phys. Rev. Lett. 72, 336-339 (1994)]. The probability distribution functions of velocity increments are non-Gaussian with wide tails in the dissipative scales and become close to Gaussian in the inertial range. The tails of the probability distribution follow a stretched exponential. We also show that in thermal convection, the energy flux in the inertial range is less than the viscous dissipation rate. This is unlike in hydrodynamic turbulence where the energy flux and the dissipation rate are equal.

Tables

Table 1. Table 1: For the two simulations of RBC: Rayleigh number Ra Ra \mathrm{Ra} , Nusselt number Nu Nu \mathrm{Nu} , kinematic viscosity ν 𝜈 \nu , viscous dissipation rate ϵ u subscript italic-ϵ 𝑢 \epsilon_{u} , and Kolmogorov length scale η 𝜂 \eta .

Case	Ra	$\mathrm{Nu}$	$\nu$	${ϵ}_u$	$\eta$
Free-slip	$1.1\times {10}^{11}$	$582$	$3.02\times {10}^{-6}$	$2.59\times {10}^{-3}$	$3.21\times {10}^{-4}$
No-slip	$1.0\times {10}^8$	$32.8$	$1.00\times {10}^{-4}$	$3.18\times {10}^{-3}$	$4.21\times {10}^{-3}$

Table 2. Table 2: For the free-slip and no-slip simulations of RBC: prefactor 𝒜 𝒜 \mathcal{A} and the scaling exponent ζ q subscript 𝜁 𝑞 \zeta_{q} for the structure functions computed by fitting the relation | S q u ​ ( l ) | = 𝒜 ​ l ζ q superscript subscript 𝑆 𝑞 𝑢 𝑙 𝒜 superscript 𝑙 subscript 𝜁 𝑞 |S_{q}^{u}(l)|=\mathcal{A}l^{\zeta_{q}} to our data.

	Free-slip simulation ($\mathrm{Ra}=1.1\times {10}^{11}$)		No-slip simulation ($\mathrm{Ra}=1.0\times {10}^8$)
$q$	$𝒜$	${\zeta }_q$	$𝒜$	${\zeta }_q$
2	$\left(2.8\pm 0.1\right)\times {10}^{-2}$	$0.70\pm 0.01$	$\left(2.3\pm 0.1\right)\times {10}^{-2}$	$0.71\pm 0.01$
3	$\left(9.3\pm 0.5\right)\times {10}^{-4}$	$0.97\pm 0.01$	$\left(8.5\pm 0.5\right)\times {10}^{-4}$	$0.98\pm 0.02$
4	$\left(2.0\pm 0.1\right)\times {10}^{-3}$	$1.26\pm 0.02$	$\left(1.6\pm 0.1\right)\times {10}^{-3}$	$1.25\pm 0.02$
5	$\left(1.5\pm 0.1\right)\times {10}^{-4}$	$1.45\pm 0.02$	$\left(2.6\pm 0.2\right)\times {10}^{-4}$	$1.60\pm 0.04$
6	$\left(1.8\pm 0.1\right)\times {10}^{-4}$	$1.69\pm 0.02$	$\left(2.6\pm 0.2\right)\times {10}^{-4}$	$1.76\pm 0.03$
7	$\left(2.1\pm 0.1\right)\times {10}^{-5}$	$1.81\pm 0.02$	$\left(7.6\pm 0.8\right)\times {10}^{-5}$	$2.01\pm 0.05$
8	$\left(2.7\pm 0.3\right)\times {10}^{-5}$	$2.09\pm 0.03$	$\left(6.0\pm 0.6\right)\times {10}^{-5}$	$2.16\pm 0.05$
9	$\left(3.9\pm 0.7\right)\times {10}^{-6}$	$2.14\pm 0.05$	$\left(2.6\pm 0.4\right)\times {10}^{-5}$	$2.33\pm 0.07$
10	$\left(3.1\pm 0.5\right)\times {10}^{-6}$	$2.28\pm 0.05$	$\left(2.1\pm 0.3\right)\times {10}^{-5}$	$2.51\pm 0.07$

Table 3. Table 3: For the two simulations of RBC: energy flux Π u subscript Π 𝑢 \Pi_{u} computed using the third-order structure functions, viscous dissipation rate ϵ u subscript italic-ϵ 𝑢 \epsilon_{u} , and the Kolmogorov constant K KO subscript 𝐾 KO K_{\mathrm{KO}} .

Case	${\Pi }_u$	${ϵ}_u$	$K_{\mathrm{KO}}$
Free-slip	$\left(1.29\pm 0.06\right)\times {10}^{-3}$	$2.59\times {10}^{-3}$	$1.59\pm 0.09$
No-slip	$\left(1.09\pm 0.03\right)\times {10}^{-3}$	$3.18\times {10}^{-3}$	$1.53\pm 0.04$

Equations

$$⟨[\{\mathbf{u}(\mathbf{r}+\mathbf{l})-\mathbf{u}(\mathbf{r})\}\cdot \hat{\mathbf{l}}{]}^3⟩$$

$${\Pi }_u={ϵ}_u,$$

$$T(x,y,z)=T_c(z)+\theta (x,y,z),$$

$${\color[rgb]{0,0,0}\rho=-\rho_{0}\alpha\theta,}$$

$$\frac{\partial \mathbf{u}}{\partial t}+(\mathbf{u}\cdot \nabla )\mathbf{u}$$

$$\frac{\partial \theta }{\partial t}+(\mathbf{u}\cdot \nabla )\theta$$

$$\nabla \cdot \mathbf{u}$$

$$\frac{\partial \mathbf{u}}{\partial t}+\mathbf{u}\cdot \nabla \mathbf{u}$$

$$\frac{\partial \theta }{\partial t}+\mathbf{u}\cdot \nabla \theta$$

$$\nabla \cdot \mathbf{u}$$

$$E_u(k)=K_{\mathrm{KO}}({\Pi }_u{)}^{2/3}{k}^{-5/3},$$

$$\frac{\partial }{\partial t}{E}_u(k,t)=-\frac{\partial }{\partial k}{\Pi }_u(k,t)+\hat{\mathcal{F}}(k,t)-{\hat{D}}_u(k,t),$$

$$\frac{d}{dk}{\Pi }_u(k)=\hat{\mathcal{F}}(k)-{\hat{D}}_u(k).$$

$$\frac{d}{dk}{\Pi }_u(k)=0,\Rightarrow {\Pi }_u(k)=\mathrm{constant}.$$

$${\Pi }_u(k_0)-{\Pi }_u(0)={\int }_0^{k_0}\hat{\mathcal{F}}(k)dk-{\int }_0^{k_0}{\hat{D}}_u(k)dk.$$

$${\Pi }_u(k_0)\approx {\int }_0^{k_0}\hat{\mathcal{F}}(k)dk.$$

$${\Pi }_u(\infty )-{\Pi }_u(k_0)={\int }_{k_0}^{\infty }\hat{\mathcal{F}}(k)dk-{\int }_{k_0}^{\infty }{\hat{D}}_u(k)dk.$$

$${\Pi }_u(k_0)={\int }_{k_0}^{\infty }{\hat{D}}_u(k)dk\approx {\int }_0^{\infty }{\hat{D}}_u(k)dk={ϵ}_u.$$

$$\frac{d}{dk}{\Pi }_u(k)={\hat{\mathcal{F}}}_B(k)-{\hat{D}}_u(k).$$

$$⟨\theta (\mathbf{r})u_z(\mathbf{r})⟩>0.$$

$$\frac{d}{dk}{\Pi }_u(k)\approx 0,\Rightarrow {\Pi }_u(k)\approx \mathrm{constant}.$$

$${\Pi }_u(k_0)-{\Pi }_u(0)={\int }_0^{k_0}{\hat{\mathcal{F}}}_B(k)dk-{\int }_0^{k_0}{\hat{D}}_u(k)dk.$$

$${\Pi }_u(k_0)\approx {\int }_0^{k_0}{\hat{\mathcal{F}}}_B(k)dk.$$

$${\Pi }_u(\infty )-{\Pi }_u(k_0)={\int }_{k_0}^{\infty }{\hat{\mathcal{F}}}_B(k)dk-{\int }_{k_0}^{\infty }{\hat{D}}_u(k)dk.$$

$${\Pi }_u(k_0)={\int }_{k_0}^{\infty }{\hat{D}}_u(k)dk-{\int }_{k_0}^{\infty }{\hat{\mathcal{F}}}_B(k)dk.$$

$${\int }_{k_0}^{\infty }{\hat{D}}_u(k)dk\approx {\int }_0^{\infty }{\hat{D}}_u(k)dk={ϵ}_u.$$

$${\Pi }_u(k_0)\approx {ϵ}_u-{\int }_{k_0}^{\infty }{\hat{\mathcal{F}}}_B(k)dk<{ϵ}_u.$$

$$\frac{\partial }{\partial t}\left[\frac{1}{2}⟨u_i(\mathbf{r})u_i(\mathbf{r}+\mathbf{l})⟩\right]=T_u(\mathbf{l})+{\mathcal{F}}_B(\mathbf{l})-D_u(\mathbf{l}),$$

$$T_u(\mathbf{l})$$

$${\mathcal{F}}_B(\mathbf{l})$$

$$D_u(\mathbf{l})$$

$${\mathcal{F}}_B(\mathbf{l})\approx -T_u(\mathbf{l}).$$

$${\mathcal{F}}_B(\mathbf{l})=\sum _{\mathbf{k}}{\hat{\mathcal{F}}}_B(\mathbf{k})\exp (i\mathbf{k}\cdot \mathbf{l}).$$

$${\hat{\mathcal{F}}}_B(\mathbf{k})=\frac{A}{2}({\delta }_{{\mathbf{k},\mathbf{k}}_0}+{\delta }_{{\mathbf{k},-\mathbf{k}}_0})+B{k}^{-5/3}.$$

Full text

I Introduction

Turbulence remains largely an unsolved problem for scientists and engineers even today. The energetics of three-dimensional homogeneous and isotropic turbulence is, however, well understood and was explained by Kolmogorov (1941a, b). Here, the energy supplied at large scales cascades down to intermediate scales and then to dissipative scales. The rate of energy supply equals the energy flux, $\Pi_{u}$, and the viscous dissipation rate $\epsilon_{u}$. Kolmogorov showed that such flows exhibit the following property (Kolmogorov, 1941a, b; Frisch, 1995):

[TABLE]

for $\eta\ll l\ll L$, where $L$ is the length scale at which energy is supplied and is of the order of the domain size, and $\eta$ is the dissipative scale, called Kolmogorov length scale. In Eq. (1), $\langle.\rangle$ represents the ensemble average, and $\mathbf{u(r)}$ and $\mathbf{u(r+l)}$ are the velocity fields at positions $\mathbf{r}$ and $\mathbf{r+l}$ respectively. The left-hand side of Eq. (1), denoted as $S_{3}^{u}(l)$, is the third-order velocity structure function. For any order $q$, one expects, using dimensional analysis, that $S_{q}^{u}(l)=\langle[\{\mathbf{u(r+l)-u(r)}\}\cdot\mathbf{\hat{l}}]^{q}\rangle\sim l^{q/3}$. Using the theory of Obukhov (1949) and Corrsin (1951) on turbulence with passive scalar $\theta$, dimensional analysis yields $S_{q}^{\theta}(l)\sim l^{q/3}$, where $S_{q}^{\theta}(l)=\langle\{\theta(\mathbf{r+l})-\theta(\mathbf{r})\}^{q}\rangle$ is the structure function for the passive scalar. The aforementioned relations for $S_{q}^{u}$ and $S_{q}^{\theta}$ are known as Kolmogorov-Obukhov (KO) scaling in literature. In reality, however, the exponents deviate from $q/3$ (other than for 3) due to intermittency effects. The velocity structure functions scale as $S_{q}^{u}(l)\sim l^{\zeta_{q}}$, where the exponents $\zeta_{q}$ fit well with the model of She and Leveque (1994).

The scaling of structure functions of turbulent convection, however, remains an unsolved problem and hence is the theme of our present paper. We focus on Rayleigh-Bénard Convection (RBC) that deals with a fluid enclosed between two horizontal plates, with the bottom plate kept at a higher temperature than the top plate. In thermal convection, complications arise due to anisotropy introduced by gravity, and also because the temperature $T$ is an active scalar.

For stably stratified turbulence, Bolgiano (1959) and Obukhov (1959) predicted the kinetic energy spectrum $E_{u}(k$) and thermal energy spectrum $E_{T}(k)$ to scale as $k^{-11/5}$ and $k^{-7/5}$ respectively, where $k\sim 1/l$ is the wavenumber. An extension of Bolgiano-Obukov (BO) theory to structure functions gives $S_{q}^{u}(l)\sim l^{3q/5}$ and $S_{q}^{T}(l)\sim l^{q/5}$, where $S_{q}^{T}$ is the temperature structure function. BO scaling occurs above the Bolgiano length scale $l_{B}$, where the buoyancy forces are dominant. Evidences of BO scaling have been observed in recent studies of stably-stratified Kumar, Chatterjee, and Verma (2014); Verma, Kumar, and Pandey (2017) and rotating stratified turbulence Rosenberg et al. (2015). Using theoretical arguments, Procaccia and Zeitak (1989), L’vov (1991), L’vov and Falkovich (1992), and Rubinstein (1994) proposed the applicability of BO scaling to RBC as well. Researchers have attempted to confirm the above theory with the help of experiments and numerical simulations, as well as using theoretical arguments.

Benzi et al. (1994a, b) simulated both 2D and 3D RBC using Lattice Boltzmann method and computed velocity and temperature structure functions up to the sixth order. They could not observe any discernible scaling for the structure functions due to short inertial range. They found them, however, to be self-similar for a wide range of $l$, a phenomenon known as extended self-similarity (ESS) Benzi et al. (1993); Chakraborty, Frisch, and Ray (2010). Further, they claimed BO scaling from the relationship between the velocity and the temperature structure functions. Ching (2000) computed temperature and velocity structure functions of thermal convection using the experimental data of Heslot, Castaing, and Libchaber (1987), and Sano, Wu, and Libchaber (1989), as well as the numerical data of Benzi et al. (1996). Although Ching (2000) observed two distinct scaling regimes separated by the Bolgiano scale, the scaling exponents deviated from BO theory.

Many researchers obtained KO scaling in the bulk and attributed it to the large value of local $l_{B}$, which is of the same order as the box size. Since $l_{B}$ is small near the walls, it is argued that the structure functions in those regions follow BO scaling. Using third-order structure functions calculated using their lattice Boltzmann simulation data, Calzavarini, Toschi, and Tripiccione (2002) claimed BO scaling near the walls and KO scaling at the cell center. High-resolution multipoint measurements of velocity and temperature fields in water were conducted by Sun, Zhou, and Xia (2006). Their exponents of velocity structure functions computed at the cell center fit well with the She-Leveque model, with the lower orders following Kolmogorov scaling. Using refined similarity hypothesis, Ching et al. (2013) derived power-law relations for conditional velocity and temperature structure functions computed at given values of the locally averaged thermal dissipation rate. Ching et al. (2013) further computed the conditional temperature structure functions up to the fourth order using the experimental data of He and Tong (2009). Based on the observed power-law scaling, they concluded BO scaling near walls and KO scaling at the cell center.

Using the experimental data of Castaing et al. (1989) and Shang et al. (2003), Ching (2007) computed the structure functions of plume velocity and found them to scale similar to the temperature structure functions. This is unlike the case of velocity structure functions in BO scaling, where they scale differently from the temperature structure functions. Kunnen et al. (2008) conducted direct numerical simulations of RBC in a grid resolution of $129\times 257\times 257$. The velocity structure functions computed by them follow BO scaling for Rayleigh number $\mathrm{Ra}=10^{8}$ and Kolmogorov scaling for higher $\mathrm{Ra}$. Ching and Cheng (2008) calculated temperature structure functions using shell model of homogeneous RBC and found them to deviate significantly from BO scaling for $q>4$. Kaczorowski and Xia (2013) conducted direct numerical simulations (DNS) of RBC in grids ranging from $64^{3}$ to $770^{3}$, and found that the velocity structure functions computed at cell center approach Kolmogorov scaling for lower orders.

From the conflicting nature of past results, it is clear that the behaviour of the structure functions of turbulent convection has not yet been clearly established. Lohse and Xia (2010) reviewed the experimental, numerical and theoretical results of past works critically and raised doubts on the applicability of BO scaling in RBC. Recently, using phenomenological arguments and numerical simulations, Kumar, Chatterjee, and Verma (2014) and Verma, Kumar, and Pandey (2017) showed Kolmogorov energy spectrum in RBC. Using energetics arguments, they derived that the energy cascade rate in turbulent convection is constant, leading to Kolmogorov scaling. Their predictions are being accepted and acknowledged by several groups as is evident from recent literature Meuel et al. (2018); Alexakis and Biferale (2018); Bruneau et al. (2018); Shestakov, Stepanov, and Frick (2017); Pawar and Arakeri (2016); Bhattacharjee (2015); Schumacher, Götzfried, and Scheel (2015). However, some researchers still believe that BO scaling is applicable to RBC Rosenberg et al. (2015); Shrestha, Mompean, and Calzavarini (2016); Rincon et al. (2017); Rincon and Rieutord (2018). In this paper, using numerical simulations, we reinforce the results of Kumar, Chatterjee, and Verma (2014) and Verma, Kumar, and Pandey (2017) by showing that the velocity structure functions of thermal convection scale similarly as those of 3D hydrodynamic turbulence. We further show that although the energy flux in turbulent convection is constant similar to hydrodynamic turbulence, it differs from viscous dissipation rate. We will discuss the scaling of temperature structure functions in a future work.

The outline of the paper is as follows: In Sec. II, we describe the governing equations of RBC. In Sec. III, we discuss the phenomenology of turbulent convection and derive the scaling of third-order structure functions. In Sec. IV, we briefly discuss the simulation details and the procedure employed to calculate the velocity structure functions. In Sec. V, we present the scaling of the structure functions and discuss the nature of the probability distribution functions of velocity increments. Further, we compare the energy flux and viscous dissipation rate in RBC and show that the flux is less than the dissipation rate. Finally, we conclude in Sec. VI.

II Governing equations

In RBC, under the Boussinesq approximation Chandrasekhar (2013); Chillà and Schumacher (2012), we assume the kinematic viscosity $\nu$, thermal diffusivity $\kappa$, and thermal expansion coefficient $\alpha$ to be constants. Further, the density of the fluid is taken to be constant except for the buoyancy term in the momentum equation. The temperature field $T$ can be split as

[TABLE]

where $T_{c}(z)$ is the conduction temperature profile, and $\theta(x,y,z)$ is the deviation of temperature from the conduction state. Further, the temperature fluctuation $\theta$ is related to the density fluctuation $\rho$ as Chandrasekhar (2013); Verma (2018)

[TABLE]

where $\rho_{0}$ is the mean fluid density. The governing equations of RBC are as follows:

[TABLE]

where $\mathbf{u}$ and $\sigma$ are the velocity and the pressure fields respectively, and $\Delta$ and $d$ are the temperature difference and distance respectively between the top and the bottom plates.

Using $d$ as the length scale, $\sqrt{\alpha g\Delta d}$ as the velocity scale, and $\Delta$ as the temperature scale, we non-dimensionalize Eqs. (3)-(5), which yields

[TABLE]

where $\mathrm{Ra}=\alpha g\Delta d^{3}/(\nu\kappa)$ is the Rayleigh number, and $\mathrm{Pr}=\nu/\kappa$ is the Prandtl number. The Rayleigh and Prandtl numbers are the main governing parameters of RBC.

In the next section we construct a phenomenology for the structure functions of turbulent convection.

III Hydrodynamic turbulence-like phenomenology for turbulent convection

III.1 Energy fluxes and spectra in hydrodynamic turbulence and thermal convection

For 3D hydrodynamic turbulence, the energy cascade rate $\Pi_{u}$ in turbulent flows is constant in the inertial range ($\eta\ll l\ll L$). Dimensional analysis gives the following relation for the energy spectrum $E_{u}(k)$:

[TABLE]

where $K_{\mathrm{KO}}$ is the Kolmogorov constant. The aforementioned $k^{-5/3}$ spectrum is known as Kolmogorov’s spectrum. In this section, we briefly describe the phenomenological arguments of Kumar, Chatterjee, and Verma (2014), Verma, Kumar, and Pandey (2017), and Verma (2018), according to which the energy spectrum in turbulent convection follows Kolmogorov scaling with constant energy flux, contrary to the arguments of L’vov (1991) and L’vov and Falkovich (1992), who propose Bolgiano-Obukhov scaling with variable flux.

In all turbulent flows, the following can be derived using Eq. (3) (see Refs. Frisch (1995); Lesieur (2008); Verma (2018)):

[TABLE]

where $\hat{\mathcal{F}}(k,t)$ is the energy feed due to forcing, and $\hat{D}_{u}(k)$ is the dissipation rate of kinetic energy. For a steady state, we have $\frac{\partial}{\partial t}E_{u}(k,t)\approx 0$ that modifies Eq. (10) to

[TABLE]

Now, we will separately consider hydrodynamic turbulence and RBC and show that the flux is constant for both the cases. However, there is a difference between the two fluxes, as shown below.

III.1.1 Hydrodynamic turbulence

The forcing in hydrodynamic turbulence is supplied at small wavenumbers. In the inertial range, $\hat{\mathcal{F}}(k)=0$ and $\hat{D}_{u}(k)$ is negligible. This results in the following (Kolmogorov, 1941a, b; Lesieur, 2008; Frisch, 1995):

[TABLE]

Note that in hydrodynamic turbulence, the forcing injection $\mathcal{F}(k)$ is modelled numerically in many ways. Refer to Canuto et al. (1988) for details.

Let us consider a small wavenumber $k_{0}$ that lies in the inertial range and is slightly larger than the forcing wavenumber. Integration of Eq. (11) from [math] to $k_{0}$ yields

[TABLE]

Note that $\int_{0}^{k_{0}}\hat{\mathcal{F}}(k)dk$ is the total energy injection rate for hydrodynamic turbulence. Since $\Pi_{u}(0)=0$ and the dissipation at small wavenumbers is negligible, we obtain

[TABLE]

Now, integration of Eq. (11) from $k_{0}$ to $\infty$ yields

[TABLE]

Since $\Pi_{u}(\infty)=0$ and $\hat{\mathcal{F}}(k)=0$ for $k\in[k_{0},\infty)$, we get

[TABLE]

Note that $k_{0}$ is small, and $\hat{D}_{u}(k)$ is small in the forcing band. Therefore the lower limit of the aforementioned integration has been replaced with 0. Thus, using Eqs. (14,16) we deduce that in hydrodynamic turbulence, the energy flux in the inertial range is constant, and is approximately equal to the dissipation rate $\epsilon_{u}$ and the total energy injection rate.

III.1.2 Thermal convection

In turbulent convection, the energy is injected into the system by buoyancy. We denote this energy feed as $\hat{\mathcal{F}}_{B}(k)$. Note that we do not inject energy externally in convection as we do in hydrodynamic turbulence. Further, unlike hydrodynamic turbulence, $\hat{\mathcal{F}}_{B}(k)$ acts at all scales in thermal convection. Replacing $\hat{\mathcal{F}}(k)$ with $\hat{\mathcal{F}}_{B}(k)$, we rewrite Eq. (11) as

[TABLE]

Since hot plumes ascend and the cold plumes descend, $u_{z}$ and $\theta$ are positively correlated, which means that (Kumar, Chatterjee, and Verma, 2014; Verma, Kumar, and Pandey, 2017)

[TABLE]

Using this condition, Kumar, Chatterjee, and Verma (2014) and Verma, Kumar, and Pandey (2017) claimed that $\hat{\mathcal{F}}_{B}(k)>0$, that is, buoyancy feeds energy to the system. Hence, $\frac{d}{dk}\Pi_{u}(k)>0$ in steady state from Eq. (17). It is important to note that in stably-stratified flows, buoyancy depletes energy from the system. Thus, for such flows, $\hat{\mathcal{F}}_{B}<0$, resulting in $\frac{d}{dk}\Pi_{u}(k)<0$. This means that the flux decreases with wavenumber in the inertial range; this is an important ingradient of Bolgiano-Obukhov scaling Bolgiano (1959); Obukhov (1959). Since the flux does not decrease with wavenumber in thermal convection, Bolgiano-Obukhov scaling is ruled out.

Further, in turbulent convection, Pandey and Verma (2016) and Pandey et al. (2016a) showed that buoyancy is strong only at large scales and is weak in the inertial range. Nath et al. (2016) showed that the distribution of velocity field in turbulent convection is nearly isotropic similar to hydrodynamic turbulence, again indicating weak buoyancy.

Based on the above observations, Kumar, Chatterjee, and Verma (2014), and Verma, Kumar, and Pandey (2017) argued that $\hat{\mathcal{F}}_{B}$ does not bring about a noticeable increase in $\Pi_{u}(k)$ (See Fig. 1). Therefore, $\hat{\mathcal{F}}_{B}\approx\hat{D}_{u}\approx 0$, which reduces Eq. (17) to

[TABLE]

Thus, it can be inferred from Eq. (18) that Kolmogorov’s theory of hydrodynamic turbulence is also applicable to thermal convection. Integrating Eq. (17) from [math] to a small wavenumber $k_{0}$ lying in the inertial range yields

[TABLE]

Since $\Pi_{u}(0)=0$ and the dissipation rate is negligible at small wavenumbers, the above equation reduces to

[TABLE]

Since $\hat{\mathcal{F}}_{B}(k)$ is strong at large scales, we deduce from Eq. (20) that a large part of energy is injected by buoyancy at large scales that contributes to the energy flux in the inertial range; this feature is similar to hydrodynamic turbulence.

There is, however, a difference between the energetics of RBC and that of 3D turbulence. Integrating Eq. (17) from $k_{0}$ to $\infty$ yields

[TABLE]

Since $\Pi_{u}(\infty)=0$, the above equation becomes

[TABLE]

Since $k_{0}$ is small compared to the dissipation range wavenumbers, we can write

[TABLE]

Now,$\int_{k_{0}}^{\infty}\hat{\mathcal{F}}_{B}(k)dk$ is the energy injected by buoyancy at small scales. It must be noted that $\int_{k_{0}}^{\infty}\hat{\mathcal{F}}_{B}(k)dk>0$ in RBC, because $\hat{\mathcal{F}}_{B}(k)$, albeit weak, is positive and adds up to a significant amount when integrated over the inertial and dissipation range (see Sec. V.3). Therefore,

[TABLE]

Eq. (23) clearly shows that unlike in hydrodynamic turbulence, the energy flux in the inertial range is smaller than the dissipation rate due to the energy injected by buoyancy at small scales. Recall that in hydrodynamic turbulence, no energy is injected in these regimes. In Sec. V.3, using the results of numerical simulations of turbulent convection, we show that the energy flux is smaller than the dissipation rate by a factor of two to three for our selected cases. Note that this factor likely depends on Ra, Pr, type of boundary conditions, etc. A careful study of the spectra and fluxes of thermal convection for different regimes of Ra and Pr needs to be carried out to ascertain how this factor depends on the aforementioned parameters.

In the next subsection, following the procedure of Kolmogorov (1941a, b), we derive the relation for the third-order velocity structure functions of turbulent convection.

III.2 Velocity structure functions of turbulent convection

Sun, Zhou, and Xia (2006) and Zhou, Sun, and Xia (2008) performed experiments of turbulent thermal convection and observed isotropy in regions away from walls. Using detailed numerical simulations, Nath et al. (2016) computed the modal energy of the inertial-range Fourier modes of turbulent convection as a function of polar angle $\Theta$ (angle between buoyancy direction and the wavenumber), and found it to be approximately independent of $\Theta$. Thus, they showed that turbulent convection is nearly isotropic. In Sec. V.1, we compute the second-order velocity structure functions as functions of $l$ and $\Theta$ ($\Theta$ is the angle between the buoyancy direction and $l$) using our numerical data, and show that they are nearly independent of $\Theta$. This again shows near-isotropy in thermal convection. We believe that isotropy is related to the fact that in turbulent convection, buoyancy “effectively” injects energy at large scales, but it is weak in the inertial range.

Further, at high Rayleigh numbers, the boundary layers are very thin, with the boundary layer thickness $\delta_{u}\ll d$, $d$ being the domain height. Therefore, for simplification, we neglect the effects of boundary layers and consider the system to be homogeneous. In Appendix A we show that in turbulent thermal convection, the planar structure functions and those computed in the entire domain exhibit somewhat similar scaling; this result too validates the assumptions of approximate homogeneity and isotropy for turbulent convection. Using the assumptions of homogeneity, isotropy and steady state, and following similar lines of arguments as Kolmogorov (1941a, b), we derive the relation for third-order structure function for turbulent convection in the following discussion.

For homogeneous and incompressible turbulent flows, the temporal evolution of the second-order velocity correlation function can be written as follows Kolmogorov, 1941a, b; Frisch, 1995:

[TABLE]

where

[TABLE]

Here, $T_{u}(\mathbf{l})$ is the non-linear energy transfer at scale $\mathbf{l}$, $\mathcal{F}_{B}(\mathbf{l})$ is the force correlation at $\mathbf{l}$, and $D_{u}(\mathbf{l})$ is the corresponding dissipation rate. The symbol $\nabla^{\prime 2}$ represents the Laplacian at $\mathbf{r+l}$. Under a steady state, the left-hand side of Eq. (24) disappears. Further, we focus on the inertial range where $D_{u}(\mathbf{l})\approx 0$ that yields

[TABLE]

Now, $\mathcal{F}_{B}(\mathbf{l})$ can be expanded as Fourier series as follows:

[TABLE]

Following Verma, Kumar, and Pandey (2017), we model $\hat{\mathcal{F}}_{B}(\mathbf{k})$ as (Frisch, 1995)

[TABLE]

Substitution of Eq. (27) in Eq. (26) yields

[TABLE]

This is because $\mathbf{k}_{0}\cdot\mathbf{l}\approx 0$ since turbulent convection is essentially forced by large-scale plumes (Verma, Kumar, and Pandey, 2017). Here, $B$ is a small constant. Now, for an isotropic flow, $T_{u}(\mathbf{l})=T_{u}(l)$, and is related to the third-order structure function $S_{3}^{u}(l)$ as (see Frisch (1995))

[TABLE]

Combining Eqs. (25), (28) and (29), we get

[TABLE]

Integrating the above expression twice, and noting that $S_{3}^{u}(0)=0$, we obtain the following relation:

[TABLE]

Now, we assume that the large-scale buoyant energy feed at $k=k_{0}$ equals the energy flux $\Pi_{u}$, and that $B$ is small. Therefore, we have $A\approx\Pi_{u}$, and

[TABLE]

Thus, the scaling of the third-order structure functions of RBC is similar to those of 3D hydrodynamic turbulence, except that $\epsilon_{u}$ of $S_{3}^{u}(l)$ is replaced by $\Pi_{u}$. Note that $\Pi_{u}<\epsilon_{u}$ for RBC. We will verify the above relation in Sec. V using numerical simulations.

It is important to note that for hydrodynamic turbulence, $\hat{\mathcal{F}}(k)$ is provided at small wavenumbers and is equal to the viscous dissipation rate $\epsilon_{u}$. Inverse Fourier transform of $\hat{\mathcal{F}}(k)$ results in a constant value of $\mathcal{F}(l)$ that equals $\epsilon_{u}$. Using the same procedure as shown above, one can derive that $S_{3}^{u}(l)=-(4/5)\epsilon_{u}l$. Note that in RBC, $\epsilon_{u}$ of the above $S_{3}^{u}(l)$ is replaced by $\Pi_{u}$. We also remark that our arguments are consistent with the results of Kunnen and Clercx (2014), who computed the scale-by-scale energy budget in direct numerical simulations of RBC and showed that $S_{3}^{u}(l)\neq-(4/5)\epsilon_{u}l$ for convective turbulence.

Finally, as mentioned previously, it must be noted that Eq. (32) has been derived under the assumption of homogeneity and isotropy, which may not be the case for all regimes of turbulent convection. For example, Nath et al. (2016) has shown that anistropy is stronger for large Prandtl numbers. Thus, we cannot make the assumption of isotropy in this regime.

In the next section, we discuss the numerical techniques involved in the computation of the structure functions.

IV Numerical methods

We use two sets of numerical data to compute the velocity structure functions, each set having different boundary conditions. The first set is the data of Verma, Kumar, and Pandey (2017), who performed direct numerical simulation (DNS) of RBC on a $4096^{3}$ grid. The grid corresponds to a cube of unit dimension. The Rayleigh and Prandtl numbers were chosen as $1.1\times 10^{11}$ and unity respectively. The grid corresponds to a cubical domain of unit dimension. The simulation was performed using a pseudo-spectral code (Verma et al., 2013; Chatterjee et al., 2018). Free-slip and isothermal boundary conditions were employed at the top and bottom plates, and periodic boundary conditions were employed at the side walls. For details, refer to Verma, Kumar, and Pandey (2017).

The second set of data is that of Kumar and Verma (2018). This simulation was performed using a finite volume solver (Jasak et al., 2007) on a non-uniform $256^{3}$ grid that corresponds to a cube of unit dimension. The Rayleigh and Prandtl numbers were chosen as $1\times 10^{8}$ and unity respectively. No-slip boundary conditions were imposed at all the walls; such realistic boundary conditions capture the wall effects. Isothermal boundary conditions were imposed at the top and bottom plates and adiabatic boundary conditions at the side walls. For spatial discretization schemes, time-marching method, and the validation of the code, see Refs. Kumar and Verma (2018); Bhattacharya et al. (2018); Bhattacharya, Samtaney, and Verma (2019). We interpolate the velocity fields to a uniform $256^{3}$ grid.

We compute the velocity structure functions in the entire domain using a combination of shared (OpenMP) and distributed memory (MPI) parallelization (see Pacheco (2011)). The computations involve running six nested for loops: the outer three loops describing the position vector $\mathbf{r}$ and the inner three loops describing $\mathbf{r+l}$. To save computational resources, we condense our free-slip data to $512^{3}$ grid. Note that we are interested only in scales pertaining to the inertial range and not the dissipative scales. After the aforementioned coarsening, we are still able to resolve scales above $6\eta$ and capture the inertial range very well in addition to avoiding unnecessary computational costs. The number of MPI nodes used were equal to the number of grid points in the $x$-direction, while the number of OpenMP threads used were 32.

In the forthcoming section, we will discuss the numerical results.

V Numerical results

In the present section, for turbulent thermal convection, we describe the scaling of the velocity structure functions, the probability distribution functions of velocity increments, and the difference between the energy flux and viscous dissipation rate.

V.1 Structure functions

Before computing the structure functions, we first numerically compute the viscous dissipation rate $\epsilon_{u}$ using the velocity field data of our free-slip and no-slip cases. We use the relation

[TABLE]

to compute the viscous dissipation rate, where $S_{ij}$ is the strain rate tensor, and $\langle.\rangle$ represents the volume average. Further, we compute the Kolmogorov length scale $\eta$ and the Nusselt number Nu using the following relations Frisch (1995); Pope (2000); Lohse and Xia (2010):

[TABLE]

In Table 1, we list the values of Nu, $\nu$, $\epsilon_{u}$ and $\eta$ for both free-slip and no-slip data. Clearly, $\eta$ is larger for the no-slip case because of lower $\mathrm{Ra}$. Further, we remark that the viscous boundary layers are thin for our data, with $\delta_{u}=0$ for the free-slip simulation and $\delta_{u}\approx 2\eta$ for the no-slip simulation Bhattacharya et al. (2018). Thus, most of the flow resides in the bulk.

Next, we validate the assumption of isotropy in turbulent convection. Using both sets of data, we compute the second-order velocity structure functions in the entire domain as functions of $l$ and $\Theta$, where $\Theta$ is the angle between the buoyancy direction and $l$. Figs. 2(a,b) exhibit the polar plots $S_{2}^{u}(l,\Theta)$, with $l$ spanning the inertial-dissipation range ($0<l/\eta<210$ for the free-slip case and $0<l/\eta<40$ for the no-slip case). The figures clearly show that the structure functions are nearly independent of $\Theta$, thereby demonstrating near-isotropy in the inertial-dissipation range.

Now, we compute the magnitude of $S_{q}^{u}$ as a function of $l$ in the entire domain, with $q$ ranging from 2 to 10. Fig. 3 exhibits the plots of structure functions of orders 2, 3, 6, 8 and 10 versus $l/\eta$ for both sets of data. Contrary to the results of Benzi et al. (1994a, b), we observe a discernible scaling range for the third order structure function. The range is found to be $32<l/\eta<200$ for the free-slip data and $19<l/\eta<40$ for the no-slip data. The range is much smaller for the no-slip case because of the higher value of $\eta$. Note that the length scales in the inertial range are much larger than the boundary layer thickness.

We compute the scaling exponents $\zeta_{q}$ and the prefactor $\mathcal{A}$ by fitting the relation $S_{q}^{u}(l)=\mathcal{A}l^{\zeta_{q}}$ to our data within the scaling range. Table 2 lists $\mathcal{A}$ and $\zeta_{q}$ for both sets of data. Note that $\zeta_{3}=0.97$ and $0.98$ for the free-slip and the no-slip cases respectively, which are close to Kolmogorov scaling of $S_{3}^{u}\sim l$. From Table 2 and Figs. 3 and 4, we observe that for lower orders, the scaling exponents $\zeta_{q}$ for free-slip and no-slip boundary conditions are nearly equal, and they are close to $q/3$, which is a generalisation of Kolmogorov’s theory of turbulence. For $q=2$, $\zeta_{2}\approx 2q/3$ that yields $k^{-5/3}$ energy spectrum. These results are consistent with the Kolmogorov energy spectrum in thermal convection observed by Kumar, Chatterjee, and Verma (2014), Verma, Kumar, and Pandey (2017), and Kumar and Verma (2018).

Our results are also consistent with those of Sun, Zhou, and Xia (2006) and Kaczorowski and Xia (2013), who report Kolmogorov scaling of the structure functions of RBC computed at the cell center. On the other hand, our results are contrary to those of Benzi et al. (1994a, b), Calzavarini, Toschi, and Tripiccione (2002), and Kunnen et al. (2008) (for $\mathrm{Ra=10^{8}}$), who deduce Bolgiano-Obukhov scaling based on their simulations. However, it must be noted that Kunnen et al. (2008) could not observe Bolgiano-Obukhov scaling for $\mathrm{Ra>10^{8}}$; rather, they report Kolmogorov scaling, similar to our results. We will discuss more on Bolgiano-Obukhov scaling later in this section.

As illustrated in Table 2 and Fig. 4, higher order $\zeta_{q}$’s for the free-slip data are marginally lower than those for the no-slip data. Also, for higher order structure functions, $\zeta_{q}$ deviates from $q/3$ due to intermittency. To explain intermittency effects in hydrodynamic turbulence, She and Leveque (1994) proposed the following model for $\zeta_{q}$:

[TABLE]

Interestingly, the aforementioned equation describes $\zeta_{q}$ calculated using our RBC data quite well; see Figs. 4(a,b). These results demonstrate similarities between $\zeta_{q}$ scaling in convection and in hydrodynamic turbulence, consistent with earlier results (Kumar, Chatterjee, and Verma, 2014; Verma, Kumar, and Pandey, 2017; Verma, 2018). Our results also match with the experimental work of Sun, Zhou, and Xia (2006), who observed the scaling exponents of structure functions calculated at cell-center to fit with She-Leveque’s model.

In Fig. 5, we plot the logarithms of $S_{2}^{u}$, $S_{6}^{u}$, $S_{8}^{u}$ and $S_{10}^{u}$ versus $\log_{10}|S_{3}^{u}|$ for both free-slip and no-slip cases, and observe the structure functions to be self-similar, that is,

[TABLE]

where $\beta(q,3)=\zeta_{q}/\zeta_{3}$. The computed values of the exponent $\beta(q,3)$ are also shown in the figure. This scaling occurs for $l/\eta$ ranging from $12$ to $530$ for the free-slip case and $9$ to $45$ for the no-slip case. The range of $S_{q}^{u}$ versus $S_{3}^{u}$ plots of Fig. 5 is wider than that of $S_{q}^{u}$ plots of Fig. 3 (In Fig. 5, the range extends well beyond the inertial range to the dissipative scales). This is called extended self-similarity (ESS) Benzi et al. (1993); Chakraborty, Frisch, and Ray (2010). ESS has been observed in previous studies of convection Benzi et al. (1994a, b); Lohse and Xia (2010). Note that ESS was first reported by Benzi et al. (1993) in hydrodynamic turbulence.

According to Pope (2000), the upper limit of the inertial range can be estimated by $l_{EI}^{P}\approx L/6$ and the lower limit $l_{DI}^{P}\approx 60\eta$. Going by this estimate, $l_{EI}^{P}=530\eta$ for our free-slip data. Note that the upper and the lower limits of the power-law range of the structure functions for our free-slip data are of the same order of magnitude as Pope’s estimate. For the no-slip case, because of the large value of $\eta$ and the dissipative nature of OpenFOAM solver, $l_{DI}^{P}(=60\eta)$ is greater than $l_{EI}^{P}(=40\eta)$. Therefore, Pope’s estimate for the lower limit does not hold for the no-slip case; this is expected because Pope’s estimates are expected to work for homogenous and isotropic turbulence, or periodic boundary condition.

An important point to note is that $\zeta_{q}$ curve does not fit with $\zeta_{q}=3q/5$, which is a generalisation of Bolgiano-Obukhov (BO) model. As discussed in Sec. III.1.2, Kumar, Chatterjee, and Verma (2014), and Verma, Kumar, and Pandey (2017) have argued against Bolgiano-Obukhov (BO) model for RBC based on energy flux arguments. This result is contrary to some of the earlier works Kunnen et al. (2008); Benzi et al. (1994a, b); Ching (2000); Procaccia and Zeitak (1989); L’vov (1991); L’vov and Falkovich (1992) that argue in favour of Bolgiano-Obukhov model. Note that Bolgiano length computed using $l_{B}=\mathrm{Nu^{1/2}/(PrRa)^{1/4}}$ are approximately $130\eta$ and $14\eta$ for the free-slip and no-slip boundary conditions respectively. They are marked as vertical lines in Fig. 3. We do not discuss $l_{B}$ in detail because Bolgiano-Obukhov (BO) model has been shown to be inapplicable for RBC (Kumar, Chatterjee, and Verma, 2014; Verma, Kumar, and Pandey, 2017; Verma, 2018)(see Sec. III.1.2).

In Appendix A, we compute the planar structure functions for several horizontal cross sections. We observe that the these structure functions are somewhat similar to those described above, with a difference that planar structure functions exhibit relatively higher fluctuations. This is due to lesser averaging for the planar structure function.

In the next subsection we describe the probability distribution function (PDF) for the velocity difference between two points.

V.2 Probability distribution function for velocity increments

For different values of $l/\eta$, we compute the probability distribution functions (PDFs) of velocity increments, $\delta u=\{\mathbf{u(r+l)-u(r)}\}\cdot\hat{\mathbf{l}}$, using the free-slip and the no-slip data. Fig. 6(a) exhibits the PDFs of $\delta u$ for the free-slip data. For small $l$, the PDFs are non-Gaussian with wide tails. The tails fit with a stretched exponential curve given by $P(\delta u)\sim\exp(-m|\delta u^{*}|^{\alpha})$, where $\delta u^{*}=\delta u/\sqrt{\langle\delta u^{2}\rangle}$. We observe that the stretching exponent $\alpha=0.8$, $1.0$, and $1.8$ for $l/\eta=12$, $62$, and $170$ respectively. Thus, the PDFs become closer to Gaussian (represented by dashed black curve) as $l$ increases. This is expected since the velocities at two largely separated points become independent of each other. Our results are similar to those observed in hydrodynamic turbulence (see Refs Kailasnath, Sreenivasan, and Stolovitzky (1992); Donzis, Yeung, and Sreenivasan (2008)).

Fig. 6(b) exhibits the PDFs of $\delta u$ calculated using the no-slip data. Clearly, the tails are narrower compared to the free-slip case. This is because of the weaker velocity fluctuations owing to the lower Rayleigh number. Moreover, the presence of viscous boundary layers also reduces the fluctuations. Pandey et al. (2016b) show that for the same parameters, the large scale velocity and heat flux are less for convection with no-slip walls than with free-slip walls. Similar to the free-slip case, the tails of the PDFs fit well with a stretched exponential. For $l/\eta=7$, $12$, and $37$, $\alpha$’s are $0.9$, $1.0$, and $1.7$ respectively for the left tail, and $1.0$, $1.2$, and $1.9$ respectively for the right tail. The PDFs become close to Gaussian at large scales, similar to the free-slip case.

V.3 Buoyancy forcing, energy flux and viscous dissipation rate

In this section, we provide a numerical demonstration that the energy flux and the viscous dissipation rate differ in RBC.

Using the third-order velocity structure functions, we calculate the energy flux $\Pi_{u}$ using Eq. (32) as

[TABLE]

We list the values of the energy flux in Table 3. We also compute the Fourier transform of our velocity and temperature field data, and compute the spectral energy flux using the following relation Verma et al. (2013); Chatterjee et al. (2018):

[TABLE]

We plot the flux [computed using Eq. (39)] against $k$ in Fig. 7. We observe the value of the flux to be almost constant in the inertial range and it closely matches with that computed using Eq. (38). In Table 3, we also list the values of $\epsilon_{u}$ computed in Sec.V.1.

From the table, we observe that $\epsilon_{u}\approx 2\Pi_{u}$ for the free-slip case and $\approx 3\Pi_{u}$ for the no-slip case. This is unlike in 3D hydrodynamic turbulence in which flux and viscous dissipation rate are equal. Our results are consistent with our arguments in Sec. III.1.2 where we show that the difference between the flux and the viscous dissipation rate arises due to non-zero buoyancy in the inertial range.

Using the values of $\Pi_{u}(k)$ computed using Eq. (39), we numerically compute $\frac{d}{dk}\Pi_{u}(k)$ using central-difference method. We also compute the energy spectrum $E_{u}(k)$ and obtain the spectrum of viscous dissipation using the relation $\hat{D}_{u}(k)=2\nu k^{2}E_{u}(k)$. Using the values of the dissipation spectrum and $\frac{d}{dk}\Pi_{u}(k)$ and assuming steady state, we compute $\hat{\mathcal{F}}_{B}(k)$ using Eq. (17):

[TABLE]

We plot the values of $\Pi_{u}(k)$, $\hat{\mathcal{F}}_{B}(k)$, and $\int_{k}^{\infty}\hat{\mathcal{F}}_{B}(k^{\prime})dk^{\prime}$ in Fig 7(a) for the free-slip case and in Fig 7(b) for the no-slip case. In each of the plots, we also draw a horizontal line denoting the viscous dissipation rate. As shown in Figs. 7(a,b), in the inertial range,

[TABLE]

and is approximately $\epsilon_{u}/2$ for the free-slip case and $\epsilon_{u}/3$ for the no-slip case. Also, $\hat{\mathcal{F}}_{B}(k^{\prime})$ in the inertial range is weak, consistent with our previous arguments.

In Fig. 8, we plot the cumulative buoyant energy forcing $\int_{0}^{k}\hat{\mathcal{F}}_{B}(k^{\prime})dk^{\prime}$, the cumulative viscous dissipation rate $\int_{0}^{k}\hat{D}_{u}(k^{\prime})dk^{\prime}$, and the energy flux $\Pi_{u}(k)$ against $k$ for both sets of data. The plots clearly show that the cumulative buoyant enegy forcing at small wavenumbers contributes to the energy flux in the inertial range, consistent with our arguments in Sec. III.1.2. For the free-slip data, $\int_{0}^{k}\hat{\mathcal{F}}_{B}(k^{\prime})dk^{\prime}$ remains close to the flux till $k=200$, after which it deviates from $\Pi_{u}(k)$. Similar behavior is also observed for the no-slip data, but with the threshold wavenumber $k=18$. Above these wavenumbers, $\int_{0}^{k}\hat{\mathcal{F}}_{B}(k^{\prime})dk^{\prime}$ increases slowly and merges with the cumulative dissipation rate $\int_{0}^{k}\hat{D}_{u}(k^{\prime})dk^{\prime}$ at dissipation wavenumbers. It is clear that $\int_{0}^{k}\hat{\mathcal{F}}_{B}(k^{\prime})dk^{\prime}$ at small wavenumbers (which contributes to the inertial range energy flux) is respectively $1/2$ and $1/3$ of the total energy injection rate ($\int_{0}^{\infty}\hat{\mathcal{F}}_{B}(k^{\prime})dk^{\prime}$) for the free-slip and the no-slip data.

Lastly, we compute the Kolmogorov constant $K_{\mathrm{KO}}$ by first calculating the constant $C$ using the following relation involving the second-order structure function and the energy flux:

[TABLE]

After this, we compute the Kolmogorov constant using (Pope, 2000)

[TABLE]

We list the values of Kolmogorov constant for both free-slip and no-slip cases in Table 3. Interestingly, $K_{\mathrm{KO}}$ of Table 3 is quite close to that for hydrodynamic turbulence (Frisch, 1995).

VI Conclusions

Using the numerical data of thermal convection, we compute the velocity structure functions $S_{q}^{u}$ for $q=2$ to $10$. The first data set Verma, Kumar, and Pandey (2017) was generated with free-slip boundary conditions for $\mathrm{Ra}=1.1\times 10^{11}$ and $\mathrm{Pr}=1$, and the second set Kumar and Verma (2018) with no-slip boundary conditions with $\mathrm{Ra}=1\times 10^{8}$ and $\mathrm{Pr}=1$. We calculate the scaling exponent $\zeta_{q}$ from $S_{q}^{u}$.

We show that the third-order structure functions, computed using both sets of data, scale according to Kolmogorov’s theory [$S_{3}^{u}=-(4/5)\Pi_{u}l$]. Our results are consistent with Kolmogorov’s energy spectrum observed in turbulent convection. The exponents of the structure functions of thermal convection match well with She-Leveque’s predictions. We demonstrate that the structure functions show extended self-similarity.

We also calculate the probability distribution function (PDF) of velocity increments for different values of the separation distance $l$. We show that for small $l$, the PDFs are non-Gaussian with wide tails. With increasing $l$, the PDFs become closer to Gaussian. The tails of the PDFs follow a stretched exponential, and the stretching exponent increases with $l$. Note that the PDFs of hydrodynamic turbulence show similar behaviour.

We compute the energy flux $\Pi_{u}$ using the third-order structure functions and show that $\Pi_{u}\neq\epsilon_{u}$; instead, it is two to three times less than $\epsilon_{u}$ for our cases. This is unlike in hydrodynamic turbulence where flux equals the dissipation rate. Using phenomenological arguments, we have shown that this difference arises due to non-zero, albeit weak, buoyancy present in the inertial range.

In summary, the scaling behaviour of velocity structure functions of turbulent convection shows similarities with those of 3D hydrodynamic turbulence. We do not analyze the temperature structure functions in this paper. Some of the notable works on temperature structure functions of turbulent convection include those of Ching (2000) and Ching et al. (2013). We will discuss the scaling of temperature structure functions in a future work.

Acknowledgments

We are grateful to A. Kumar and A. Chatterjee for sharing their numerical data with us. We acknowledge R. Samuel and M. Sharma for their contributions in the development of the code to calculate structure functions. We thank S. Chakraborty and S. Vashishtha for useful discussions. Our numerical simulations were performed on Shaheen II at Kaust supercomputing laboratory, Saudi Arabia, under the project k1052. This work was supported by the research grants PLANEX/PHY/2015239 from Indian Space Research Organisation, India, and by the Department of Science and Technology, India (INT/RUS/RSF/P-03) for the Indo-Russian project.

Appendix A Extent of homogeneity in turbulent convection

For very high Rayleigh number RBC, the boundary layers are quite thin. Hence, the flow, mostly residing in the bulk, is nearly homogeneous. However, for relatively smaller Ra (around $10^{8}$), there can be some inhomogeneity due to plumes and large-scale structures. To test the extent of inhomogeneity, we compute the third-order velocity structure functions for three horizontal slices of the free-slip and no-slip flow profiles detailed in the main text. The three slices are at $z=0.25,0.5,0.75$. Note that the $z=0.5$ corresponds to the mid plane.

Figure 9 exhibits the plots of $|S_{3}^{u}(l)|$ vs. $l/\eta$ for the three planes. For the free-slip data with higher Ra [Fig. 9(a)], $|S_{3}^{u}(l)|\sim l^{\zeta_{3}}$, where $\zeta_{3}\approx 1$ for $z=0.25$ and $0.5$. However, for the $z=0.75$ plane, $\zeta_{3}=0.86$, which is slightly below unity. The upper and the lower limits of the scaling range are nearly same as those for the structure functions computed in the entire domain [see Fig. 3(a)]. Thus, it is reasonable to assume the free-slip data to be homogeneous. However we observe stronger spatial inhomogeneities for the no-slip data, as shown in Fig. 9(b). Here, the scaling regime of $|S_{3}^{u}(l)|\sim l^{\zeta_{3}}$ is observed for all the three planes, with $\zeta_{3}$ ranging from 1.02 to 1.15. However, the range of the scaling regime differs for the three planes. Note that the spatial inhomogeneities are stronger for no-slip boundary condition due to the relatively stronger plumes for the no-slip boundaries.

As mentioned earlier, the observed inhomogeneity, which is more prominent for no-slip data, can be attributed to localized plumes. Thus, the structure functions are required to be averaged over more points to cancel out the effects of the plumes. That is why bulk structure functions are smoother than those for the planes, and they are closer to the predictions of She-Leveque She and Leveque (1994).