Roughness as a Route to the Ultimate Regime of Thermal Convection

TL;DR

This study uses high-resolution simulations to identify optimal wall roughness parameters that enhance heat transfer in turbulent convection, approaching theoretical upper bounds and shedding light on achieving the ultimate convection regime.

Contribution

It demonstrates how wall roughness wavelength influences heat transfer scaling, revealing an optimal wavelength that maximizes the Nusselt number and approaches theoretical growth limits.

Findings

Optimal roughness wavelength enhances heat flux.

Maximum Nusselt number scales as Ra^{0.444}.

Roughness can facilitate the ultimate convection regime.

Abstract

We use highly resolved numerical simulations to study turbulent Rayleigh-B\'enard convection in a cell with sinusoidally rough upper and lower surfaces in two dimensions for $Pr = 1$ and $Ra = \left[4 \times 10^6, 3 \times 10^9\right]$. By varying the wavelength $\lambda$ at a fixed amplitude, we find an optimal wavelength $\lambda_{\text{opt}}$ for which the Nusselt-Rayleigh scaling relation is $\left(Nu-1 \propto Ra^{0.483}\right)$ maximizing the heat flux. This is consistent with the upper bound of Goluskin and Doering \cite{Goluskin:2016} who prove that $Nu$ can grow no faster than ${\cal O} (Ra^{1/2})$ as $Ra \rightarrow \infty$, and thus the concept that roughness facilitates the attainment of the so-called ultimate regime. Our data nearly achieve the largest growth rate permitted by the bound. When $\lambda \ll \lambda_{\text{opt}}$ and $\lambda \gg \lambda_{\text{opt}}$, the planar case is recovered, demonstrating how controlling the wall geometry manipulates the interaction between the boundary layers and the core flow. Finally, for each $Ra$ we choose the maximum $Nu$ among all $\lambda$, and thus optimizing over all $\lambda$, to find $Nu_{\text{opt}} - 1 = 0.01 \times Ra^{0.444}$.

Figures

Equations

$$Nu=\frac{{\left[-\kappa \frac{\partial \overline{T}}{\partial z}+\overline{wT}\right]}_{z=z_e}}{\kappa \Delta T/H},$$

Full text

The ubiquity and importance of thermal convection in many natural and man-made settings is well known Kadanoff (2001); Worster (2000); Wettlaufer (2011). The simplest scenario that has been used to study the fundamental aspects of thermal convection is the Rayleigh-Bénard system Chandrasekhar (2013). The flow in this system is governed by three non-dimensional parameters: (1) the Rayleigh number $Ra=g\alpha\Delta TH^{3}/\nu\kappa$, which is the ratio of buoyancy to viscous forces, where $g$ is the acceleration due to gravity, $\alpha$ the thermal expansion coefficient of the fluid, $\Delta T$ the temperature difference across a layer of fluid of depth $H$, $\nu$ the kinematic viscosity (or momentum diffusivity) and $\kappa$ the thermal diffusivity; (2) the Prandtl number, $Pr=\nu/\kappa$; and (3) the aspect ratio of the cell, $\Gamma$, defined as the ratio of its width to height.

The primary aim of the corpus of studies of turbulent Rayleigh-Bénard convection has been to determine the Nusselt number, $Nu$, defined as the ratio of total heat flux to conductive heat flux (Eq. 1), as a function of the three governing parameters, viz., $Nu=Nu(Ra,Pr,\Gamma)$. For $Ra\gg 1$ and fixed $Pr$ and $\Gamma$, this relation is usually sought in the form of a power law: $Nu=A(Pr,\Gamma)Ra^{\beta}$, where $\beta$ has a fundamental significance for the mechanisms underlying the transport of heat.

The classical theory of Priestley Priestley (1954), Malkus Malkus (1954) and Howard Howard (1966) is based on the argument that as $Ra\rightarrow\infty$ the dimensional heat flux should become independent of the depth of the cell, resulting in $\beta=1/3$. A consequence of this scaling is that the conductive boundary layers (BLs) at the upper and lower surfaces, which are separated by a well mixed interior, do not interact.

However, Kraichnan Kraichnan (1962) reasoned that for extremely large $Ra$ the BLs undergo a transition leading to the generation of smaller scales near the boundaries that increase the system’s efficiency in transporting the heat, predicting that $Nu\sim\left[Ra/\left(\ln Ra\right)^{3}\right]^{1/2}$. In this, “Kraichnan-Spiegel” or “ultimate regime” ($\beta=1/2$), it is argued that the heat flux becomes independent of the molecular properties of the fluid (e.g., Spiegel, 1971; Grossmann and Lohse, 2000). Experimental Niemela et al. (2000); Niemela and Sreenivasan (2006); Urban et al. (2011, 2012) and numerical Verzicco and Camussi (2003); Stevens et al. (2011); Gayen et al. (2013) studies have found $\beta\approx 1/3$. Chavanne et al. Chavanne et al. (1997) and He et al. He et al. (2012) have reported observing transitions to $\beta=0.39$ and $\beta=0.38$ in their respective experiments and these findings continue to stimulate discussion Skrbek and Urban (2015); He et al. (2016). Motivated by studies of shear flow, Borue and Orszag Borue and Orszag (1997) used pseudo-spectral methods at three resolutions (643, 1283, 2563 and hence values of $Ra$) to study “homogeneous” convection, in which the BL’s are effectively removed. Whilst the highest resolution was not numerically converged, the other two resolutions led to a range of $\beta=0.40\pm 0.05$. This idea was later used in Lattice Boltzmann simulations for $Ra=\left[8.64\times 10^{5},1.38\times 10^{7}\right]$, to find $\beta=0.51\pm 0.06$ Lohse and Toschi (2003), ascribing this to the ultimate regime.

Recently, Waleffe et al. Waleffe et al. (2015) and Sondak et al. Sondak et al. (2015) numerically computed the steady solutions to the Oberbeck-Boussinesq equations for $Ra\leq 10^{9}$ and $1\leq Pr\leq 100$ in two dimensions. By fixing $Ra$ and $Pr$, steady solutions for different horizontal wavenumbers, $\alpha$, were computed. The solution that maximized heat transport, $Nu\equiv Nu_{opt}$, was called optimal, for which $\alpha\equiv\alpha_{opt}$ and $Nu_{opt}-1=0.115\times Ra^{0.31}$, which is in agreement with experiments Niemela et al. (2000). Although they found that $\beta$ was independent of $Pr$, the Prandtl number did have considerable effect on the geometry of the coherent structures that transported heat. For $Pr>7$, the scaling for the optimal wavenumber was found to be $\alpha_{\text{opt}}=0.257\times Ra^{0.256}$. The horizontally averaged optimal temperature profiles had the following features: (a) The BLs were always unstably stratified. (b) The core region was either stably ($Pr\leq 7$) or unstably ($Pr>7$) stratified. (c) The transition regions between the core and BLs were always stably stratified. Thus, with small departures, these profiles correspond to the marginally stable profile of Malkus Malkus (1954), with $\beta=1/3$.

An important aspect emerging from the study of planar Rayleigh-Bénard convection in two dimensions for $Pr\geq 1$ is that the flow field Schmalzl et al. (2004) and the $Nu$-$Ra$ scaling relations Johnston and Doering (2009); Waleffe et al. (2015); Sondak et al. (2015) are similar to those in three dimensions. Thus, this correspondence permits one to understand the processes driving the heat transport using well resolved two-dimensional simulations.

It is clear that the value $\beta$ takes in the limit $Ra\rightarrow\infty$ depends on the interaction between the BLs and the core flow. To understand the role of BLs in thermal convection, Shen et al. Shen et al. (1996) introduced rough upper and lower surfaces made of pyramidal elements in a cylindrical cell. They found that these elements enhanced the production of plumes, which were directly injected into the core flow, leading to an increase in $Nu$. The increase in $Nu$ was due to an increase in the pre-factor in the $Nu$–$Ra$ scaling relation. Whereas subsequent experiments found no effect of periodic roughness on $\beta$ Du and Tong (1998, 2000); Ciliberto and Laroche (1999), later studies confirmed that the changes in the flow field brought about by surface roughness do increase the value of $\beta$ from the planar value Roche et al. (2001); Qiu et al. (2005); Stringano et al. (2006); Tisserand et al. (2011); Salort et al. (2014); Wei et al. (2014); Wagner and Shishkina (2015). In our previous study, we used roughness to break the top/bottom boundary layer symmetry, and found that a periodic upper surface with $\lambda_{\text{opt}}=0.154$ maximized the heat transport with $\beta=0.359$ for a smooth lower surface in high resolution numerical simulations Toppaladoddi et al. (2015a). As is the case with the present geometry, when $\lambda\ll\lambda_{\text{opt}}$ and $\lambda\gg\lambda_{\text{opt}}$, the planar results are recovered. For each $Ra$ we determined the maximum $Nu$ among all $\lambda$, thereby optimizing over all $\lambda$, to find $Nu_{\text{opt}}-1=0.058\times Ra^{0.334}$.

The first experimental attempt to use roughness to reach the ultimate regime at $Ra$ accessible in the laboratory was made by Roche et al. Roche et al. (2001), who used V-shaped grooves to cover the entire interior of their cylindrical cell of $\Gamma=0.5$. They observed a transition in $Nu(Ra)$ at $Ra\approx 2\times 10^{12}$, and that the data beyond the point of transition could be fit with a power law with $\beta=0.51$. A similar transition was observed at $Ra=7\times 10^{9}$ in the simulations of Stringano et al. Stringano et al. (2006), who used a cylindrical geometry with V-shaped grooves at the upper and lower surfaces and imposed axisymmetry on the flow. This artificial symmetry had two important effects on the flow field: (1) The production and release of the plumes from the roughness elements was in tandem, resulting in larger plumes; and (2) The plumes traversed the vertical distance without encountering a large scale circulation in the interior region. Both these effects resulted in an increase in the efficiency of the heat transfer. As summarized by Ahlers et al., Ahlers et al. (2009), it was first noted by Niemela & Sreenivasan Niemela and Sreenivasan (2006) that the results of Roche et al. Roche et al. (2001) can be understood as a transition between when the groove depth is less than the BL thickness to a regime where the groove depth is larger than the BL thickness. Ahlers et al., Ahlers et al. (2009) state “More work is needed to resolve this issue.” Here we present results from well resolved numerical simulations of Rayleigh-Bénard convection in a cell with rough upper and lower surfaces in two dimensions. The roughness profiles chosen are sinusoidal. By keeping the amplitude fixed and varying the wavelength of the rough surfaces, we study their effects on the heat transport.

The geometry and the dimensionless equations of motion studied here are shown in figure 1. The aspect ratio of the cell, $\Gamma\equiv L_{x}/L_{z}$, is fixed at 2. The rough surfaces have a wavelength $\lambda\equiv\lambda^{*}/L_{z}$ and an amplitude $h\equiv h^{*}/L_{z}$. The equations of motion for thermal convection are the Oberbeck-Boussinesq (O-B) equations Chandrasekhar (2013), and are non-dimensionalized by choosing $H=L_{z}-2h^{*}$ as the length scale and $U_{0}=\sqrt{g\alpha\Delta TH}$ as the velocity scale. Hence, the time scale is $t_{0}=H/U_{0}$. Here, $\boldsymbol{u}(\boldsymbol{x},t)=\left(u(\boldsymbol{x},t),w(\boldsymbol{x},t)\right)$ is the velocity field, $T(\boldsymbol{x},t)$ is the temperature field, $\boldsymbol{k}$ is the unit vector along the vertical, and $p(\boldsymbol{x},t)$ is the pressure field. No-slip and Dirichlet conditions for $\boldsymbol{u}$ and $T$ are imposed on the rough surfaces, and periodic conditions are used in the horizontal.

The O-B equations were solved using the Lattice Boltzmann method with separate distributions for the momentum and temperature fields Benzi et al. (1992); Massaioli et al. (1993); Chen and Doolen (1998); Shan (1997); Guo et al. (2002). Our code has been extensively tested against results from numerical simulations for a wide range of different flows, and the details of the validation can be found in Toppaladoddi et al. (2015b, a).

For each of ten $\lambda$’s (see Fig. 2) we simulated over the range $Ra=\left[4\times 10^{6},3\times 10^{9}\right]$. The planar wall case is $\lambda=0$, the amplitude of the roughness is fixed at $h=0.1$ and $Pr=1$ for all simulations. We ran the simulations for at least $143\,t_{0}$, where $t_{0}$ is the turnover time, and statistics were collected only after $100\,t_{0}$. The Nusselt number was computed as

[TABLE]

where the overbar represents horizontal and temporal average. We should note here that this definition of $Nu$, in general, does not reduce to unity in the static case for arbitrary roughness geometries Goluskin and Doering (2016); however, for the sinusoidal geometries used here this choice gives $Nu\approx 1$ when $Ra=0$. To give an example of the spatial resolutions in the simulations, for $\lambda=1$ and $Ra=2\times 10^{9}$ the number of grid points used are $N_{x}=2800$ and $N_{z}=1400$. Grid independence was ascertained from simulations at $Ra=2\times 10^{9}$ for $\lambda=0.03$ and $0.2$ using two grids: (a) $N_{x}=2400$, $N_{z}=1200$ and (b) $N_{x}=2000$, $N_{z}=1000$. The difference between $Nu$ computed at $z_{e}=L_{z}/2$ for the two grids was less than $1.2\%$. As an additional check, $Nu$ was computed at three different depths $z_{e}=L_{z}/4,L_{z}/2$, and $3\,L_{z}/4$; and the difference between $Nu$ at any two depths was less than $0.5\%$. More simulation details are provided in the Supplementary Material.

For each $\lambda$, we obtained $\beta$ from a linear least squares fit to the $Nu(Ra)$ simulation data. Figure 2 shows $\beta$ in the scaling relation $Nu-1=A\times Ra^{\beta}$ as a function of $\lambda$. At the optimal wavelength $\lambda_{\text{opt}}=0.1$, $\beta$ attains a maximum value of $0.483$, which indicates that the influence of BLs on heat transport has been minimized. It is clear that in the limits $\lambda\ll\lambda_{\text{opt}}$ and $\lambda\gg\lambda_{\text{opt}}$, the planar case is approached.

The $Nu$-$Ra$ scaling relations for different $\lambda$ are shown in figure 3. The linear least-squares fit for $\lambda_{\text{opt}}=0.1$ giving $Nu-1=0.0042\times Ra^{0.483}$ is shown in figure 3. The roughness elements are ‘submerged’ inside the thermal BLs for $Ra<10^{8}$ (not shown), and hence, as seen in figure 3, the values of $Nu$ for these $Ra$ are close to those for larger $\lambda$. The increase in $\beta$ for $\lambda=0.1$ relative to other $\lambda$ is clear from figure 3(b). Figure 3 also shows the fit obtained for $Nu_{\text{opt}}$($Ra$), which is obtained in the following manner: for each $Ra$ we choose the maximum $Nu$ among all $\lambda$, effectively optimizing over all $\lambda$. This data is described by $Nu_{\text{opt}}-1=0.01\times Ra^{0.444}$.

The flow field for the case of $\lambda_{\text{opt}}=0.1$ and $Ra=2\times 10^{9}$ is shown in Fig. 4, where the following features are apparent:

Two large convection rolls in the cell interior. 2. 2.

The ‘unstable’ BLs at the upper and lower surfaces. 3. 3.

The production of plumes from the fluid moving along the rough surfaces and their ejection from the tips of the roughness elements.

By varying $\lambda$, we have achieved a state in which the interaction between the core flow and the BLs over the roughness elements has been enhanced. This results in an unstable state for the BLs, which then leads to the generation and ejection of plumes from the roughness tips. As noted above, in the case of a single rough wall, the maximum value of $\beta$ was found to be $\beta\approx 0.36$ Toppaladoddi et al. (2015a) but at a slightly larger $\lambda$. This highlights the role played by the second rough wall in further decreasing the role of the BLs in transporting heat. We should note here that in spite of the differences in geometry, our results have a correspondence with those of Waleffe et al. Waleffe et al. (2015) and Sondak et al. Sondak et al. (2015) in that there is a length scale in each setting ($\lambda_{\text{opt}}$ in ours and $\alpha_{\text{opt}}$ in theirs) that optimizes heat transport. The optimization occurs through the manipulation of the coherent structures that transport heat, though in detail it is accomplished in different ways.

Our results are consistent with those of Goluskin & Doering Goluskin and Doering (2016), who used the background method to compute upper bounds 111A detailed discussion of upper-bound studies can be found in Kerswell Kerswell (1998) and Hassanzadeh et al. Hassanzadeh et al. (2014). on $Nu$ for R-B convection in a domain with rough upper and lower surfaces that have square-integrable gradients. They prove that $Nu\leq CRa^{1/2}$, where $C$ depends on the geometry of roughness. Our results show that for the optimal wavelength the heat transport is $Nu-1=0.0042\times Ra^{0.483}$, with the value of $C$ being four orders of magnitude larger than ours, but an exponent approaching their result. Importantly, their approach provides a key framework for exploring a range of amplitudes and wavelengths using our methodology. Finally, our findings demonstrate that the scaling of the ultimate regime is nearly achieved in two dimensions using rough walls. Roche et al. Roche et al. (2001) interpreted their observation of $\beta=1/2$ as being due to a laminar to turbulent transition of the BLs. Here, this state is achieved by the enhanced BL–core flow interaction driven by the roughness, which generates a larger number of intense plumes.

In summary, we have studied convection in a rectangular cell of $\Gamma=2$ with rough upper and lower surfaces. At a fixed roughness amplitude, varying the wavelength $\lambda$ results in a spectrum of exponents in the $Nu$-$Ra$ scaling relation. At $\lambda_{\text{opt}}$ the maximum exponent $\beta_{\text{max}}=0.483$ is achieved, and in the limits $\lambda\ll\lambda_{\text{opt}}$ and $\lambda\gg\lambda_{\text{opt}}$, the planar value of $\beta$ is recovered, which may underlie why certain experiments found no effect of periodic roughness on $\beta$ Du and Tong (1998, 2000); Ciliberto and Laroche (1999). The observation of $\beta_{\text{max}}\approx 0.5$ here has been facilitated by the use of very large amplitude roughness relative to existing studies Roche et al. (2001); Stringano et al. (2006); Wei et al. (2014), indicating the promise of examining this state experimentally for more moderate values of $Ra$ than have been previously necessary. Indeed, by varying both amplitude and wavelength over a significant range, the systematic effects of the BLs, and thus the molecular properties of the fluid, may be realized, comparing and contrasting the concept of a laminar-to-turbulent BL transition, with the enhanced forcing associated with unstable BL’s triggered by the roughness as seen here.

Appendix 1: Optimizing heat transport over wavelength

In figure 5, we show the compensated plot for $Nu_{\text{opt}}$, and it is apparent that the exponent for the $Nu_{\text{opt}}$($Ra$) scaling law is indeed $0.444$ and that the prefactor is $0.01$. Whence, it provides a different means for reaching the same conclusion as described in the manuscript.

Figure 6 shows the variation of $\lambda_{\text{opt}}$ with $Ra$. As can be seen, $\lambda_{\text{opt}}$ decreases from $0.67$ to $0.286$ and finally saturates to $\lambda_{\text{opt}}=0.1$, implying that the wavelength for which $Nu$ is maximum for $Ra\geq 10^{9}$ is $0.1$. This is again consistent with figures 2 and 3b in the manuscript that show that the exponent attains a maximum value for $\lambda=0.1$.

Appendix 2: Simulation Details

The details of all the simulations are provided here. The roughness wavelength is $\lambda$; $Ra$ is the Rayleigh number; $N_{x}$ and $N_{z}$ are the number of grid points along the horizontal and vertical, respectively; $T_{s}$ is the total run time in terms of the turn-over time $t_{0}$; and $Nu$ is the Nusselt number.

$\lambda=0.03$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $1600$ $800$ $287.2$ $10.50$

$10^{7}$ $1600$ $800$ $347.0$ $12.48$

$2\times 10^{7}$ $1600$ $800$ $351.8$ $14.34$

$4\times 10^{7}$ $1600$ $800$ $349.3$ $16.65$

$10^{8}$ $2000$ $1000$ $328.1$ $20.92$

$2\times 10^{8}$ $2000$ $1000$ $327.9$ $26.59$

$5\times 10^{8}$ $2000$ $1000$ $324.7$ $40.18$

$10^{9}$ $2400$ $1200$ $230.7$ $59.18$

$2\times 10^{9}$ $2400$ $1200$ $231.4$ $90.64$

$3\times 10^{9}$ $2800$ $1400$ $193.3$ $113.44$ 2. 2.

$\lambda=0.05$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $1600$ $800$ $347.2$ $9.64$

$10^{7}$ $1600$ $800$ $345.8$ $11.87$

$2\times 10^{7}$ $1600$ $800$ $348.5$ $13.86$

$4\times 10^{7}$ $1600$ $800$ $347.5$ $16.36$

$10^{8}$ $2000$ $1000$ $328.5$ $21.99$

$2\times 10^{8}$ $2000$ $1000$ $329.0$ $30.68$

$5\times 10^{8}$ $2000$ $1000$ $325.9$ $55.32$

$10^{9}$ $2400$ $1200$ $231.2$ $82.36$

$2\times 10^{9}$ $2400$ $1200$ $234.2$ $120.26$

$3\times 10^{9}$ $2800$ $1400$ $194.3$ $127.05$ 3. 3.

$\lambda=0.1$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $615.5$ $9.02$

$10^{7}$ $800$ $400$ $612.6$ $11.52$

$2\times 10^{7}$ $800$ $400$ $613.0$ $13.75$

$4\times 10^{7}$ $800$ $400$ $614.1$ $17.02$

$10^{8}$ $1600$ $800$ $345.8$ $29.67$

$2\times 10^{8}$ $1600$ $800$ $351.6$ $45.02$

$5\times 10^{8}$ $2000$ $1000$ $315.1$ $72.27$

$10^{9}$ $2400$ $1200$ $235.7$ $100.81$

$2\times 10^{9}$ $2400$ $1200$ $233.8$ $135.81$

$3\times 10^{9}$ $2800$ $1400$ $201.1$ $164.62$ 4. 4.

$\lambda=0.154$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $628.9$ $9.17$

$10^{7}$ $800$ $400$ $630.2$ $11.87$

$2\times 10^{7}$ $800$ $400$ $617.7$ $14.25$

$4\times 10^{7}$ $800$ $400$ $627.1$ $18.99$

$10^{8}$ $1600$ $800$ $354.1$ $36.63$

$2\times 10^{8}$ $1600$ $800$ $355.9$ $49.60$

$5\times 10^{8}$ $2000$ $1000$ $357.6$ $76.02$

$10^{9}$ $2400$ $1200$ $233.1$ $99.08$

$2\times 10^{9}$ $2400$ $1200$ $236.1$ $127.52$

$3\times 10^{9}$ $2800$ $1400$ $143.1$ $139.88$ 5. 5.

$\lambda=0.2$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $311.7$ $9.27$

$10^{7}$ $800$ $400$ $311.9$ $12.06$

$2\times 10^{7}$ $800$ $400$ $312.6$ $15.04$

$4\times 10^{7}$ $800$ $400$ $311.3$ $20.22$

$10^{8}$ $1200$ $600$ $545.4$ $33.38$

$2\times 10^{8}$ $1200$ $600$ $549.8$ $46.99$

$5\times 10^{8}$ $2000$ $1000$ $350.2$ $68.30$

$10^{9}$ $2400$ $1200$ $182.3$ $89.55$

$2\times 10^{9}$ $2400$ $1200$ $235.5$ $111.15$

$3\times 10^{9}$ $2800$ $1400$ $201.5$ $127.29$ 6. 6.

$\lambda=0.286$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $624.1$ $9.45$

$10^{7}$ $800$ $400$ $614.0$ $12.34$

$2\times 10^{7}$ $800$ $400$ $605.8$ $16.40$

$4\times 10^{7}$ $800$ $400$ $628.1$ $22.67$

$10^{8}$ $1600$ $800$ $356.3$ $33.59$

$2\times 10^{8}$ $1600$ $800$ $359.5$ $43.71$

$5\times 10^{8}$ $2000$ $1000$ $328.7$ $60.06$

$10^{9}$ $2400$ $1200$ $232.2$ $76.82$

$2\times 10^{9}$ $2400$ $1200$ $230.9$ $94.98$

$3\times 10^{9}$ $2800$ $1400$ $188.5$ $106.92$ 7. 7.

$\lambda=0.4$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $621.3$ $9.81$

$10^{7}$ $800$ $400$ $611.2$ $13.28$

$2\times 10^{7}$ $800$ $400$ $621.1$ $16.94$

$4\times 10^{7}$ $800$ $400$ $623.1$ $21.94$

$10^{8}$ $1600$ $800$ $364.4$ $30.77$

$2\times 10^{8}$ $1600$ $800$ $358.0$ $38.59$

$5\times 10^{8}$ $2000$ $1000$ $325.7$ $51.05$

$10^{9}$ $2400$ $1200$ $218.9$ $63.94$

$2\times 10^{9}$ $2400$ $1200$ $222.2$ $78.21$

$3\times 10^{9}$ $2800$ $1400$ $179.5$ $85.82$ 8. 8.

$\lambda=0.5$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $620.4$ $10.38$

$10^{7}$ $800$ $400$ $625.6$ $13.57$

$2\times 10^{7}$ $800$ $400$ $626.1$ $16.80$

$4\times 10^{7}$ $800$ $400$ $612.2$ $21.40$

$10^{8}$ $1600$ $800$ $361.5$ $28.65$

$2\times 10^{8}$ $1600$ $800$ $357.9$ $36.00$

$5\times 10^{8}$ $2000$ $1000$ $325.8$ $45.93$

$10^{9}$ $2400$ $1200$ $228.4$ $55.98$

$2\times 10^{9}$ $2400$ $1200$ $215.5$ $70.26$

$3\times 10^{9}$ $2800$ $1400$ $185.5$ $77.39$ 9. 9.

$\lambda=0.67$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $630.3$ $10.60$

$10^{7}$ $800$ $400$ $628.0$ $13.29$

$2\times 10^{7}$ $800$ $400$ $606.3$ $16.34$

$4\times 10^{7}$ $800$ $400$ $620.3$ $19.99$

$10^{8}$ $1600$ $800$ $357.5$ $26.20$

$2\times 10^{8}$ $1600$ $800$ $363.0$ $30.32$

$5\times 10^{8}$ $2000$ $1000$ $304.4$ $40.36$

$10^{9}$ $2400$ $1200$ $217.7$ $50.01$

$2\times 10^{9}$ $2400$ $1200$ $213.1$ $62.21$

$3\times 10^{9}$ $2400$ $1200$ $213.2$ $68.66$ 10. 10.

$\lambda=1.0$

 

$Ra$ $N_{x}$ $N_{z}$ $T_{s}$ $Nu$

$4\times 10^{6}$ $800$ $400$ $508.1$ $10.29$

$10^{7}$ $800$ $400$ $506.8$ $12.54$

$2\times 10^{7}$ $800$ $400$ $509.9$ $14.46$

$4\times 10^{7}$ $800$ $400$ $502.8$ $17.58$

$10^{8}$ $1600$ $800$ $351.6$ $22.93$

$2\times 10^{8}$ $1600$ $800$ $368.9$ $28.10$

$5\times 10^{8}$ $2000$ $1000$ $297.6$ $37.51$

$10^{9}$ $2400$ $1200$ $207.6$ $44.34$

$2\times 10^{9}$ $2400$ $1200$ $203.3$ $57.66$

$3\times 10^{9}$ $2800$ $1400$ $173.4$ $65.18$

Appendix 3: Grid Independence Tests

The following tests were performed to ascertain the grid independence of the results:

$\lambda=0.1$ and $Ra=2\times 10^{8}$

Two grids were used: (1) $N_{x}=2400$ and $N_{z}=1200$ and (2) $N_{x}=1600$ and $N_{z}=800$. The difference in $Nu$ from these two runs was $0.3\%$. 2. 2.

$\lambda=0.2$ and $Ra=10^{9}$

Two grids were used: (1) $N_{x}=2400$, $N_{z}=1200$ and (2) $N_{x}=1600$, $N_{z}=800$. The difference in $Nu$ between these two runs was $3.7\%$. 3. 3.

$\lambda=0.03$ and $Ra=2\times 10^{9}$

Two grids were used: (1) $N_{x}=2400$ and $N_{z}=1200$ and (2) $N_{x}=2000$ and $N_{z}=1000$. The difference in $Nu$ from these two runs was $0.9\%$. 4. 4.

$\lambda=0.2$ and $Ra=2\times 10^{9}$

Two grids were used: (1) $N_{x}=2400$, $N_{z}=1200$ and (2) $N_{x}=2000$, $N_{z}=1000$. The difference in $Nu$ between these two runs was $0.1\%$.

We note here that the smaller grid used in test run 2 for $Ra=10^{9}$ was mainly to check the robustness of the code. Such a large $Ra$ in general requires more number of grid points to resolve the flow in the roughness region. The difference in $Nu$ of $3.7\%$ with the higher resolution run demonstrates that the numerical method employed is adequately robust.

Additionally, the code has been thoroughly validated against results from spectral codes for both Rayleigh-Bénard convection and transitional flows in two-dimensional channels Toppaladoddi et al. (2015a, b).