Hostname: page-component-77c89778f8-swr86 Total loading time: 0 Render date: 2024-07-16T13:25:06.450Z Has data issue: false hasContentIssue false
  • English
  • Français

$Nu\sim Ra^{1/2}$ scaling enabled by multiscale wall roughness in Rayleigh–Bénard turbulence

Published online by Cambridge University Press:  23 April 2019

Xiaojue Zhu Open the ORCID record for Xiaojue Zhu [Opens in a new window] ,
Richard J. A. M. Stevens ,
Olga Shishkina Open the ORCID record for Olga Shishkina [Opens in a new window] ,
Roberto Verzicco Open the ORCID record for Roberto Verzicco [Opens in a new window]  and
Detlef Lohse Open the ORCID record for Detlef Lohse [Opens in a new window]
Xiaojue Zhu*
Affiliation:
Physics of Fluids Group and Max Planck Center Twente for Complex Fluid Dynamics, MESA+ Institute and J. M. Burgers Centre for Fluid Dynamics, University of Twente, P.O. Box 217, 7500AE Enschede, The Netherlands Center of Mathematical Sciences and Applications, and School of Engineering and Applied Sciences, Harvard University, Cambridge, MA 02138, USA
Richard J. A. M. Stevens
Affiliation:
Physics of Fluids Group and Max Planck Center Twente for Complex Fluid Dynamics, MESA+ Institute and J. M. Burgers Centre for Fluid Dynamics, University of Twente, P.O. Box 217, 7500AE Enschede, The Netherlands
Olga Shishkina
Affiliation:
Max Planck Institute for Dynamics and Self-Organization, 37077 Göttingen, Germany
Roberto Verzicco
Affiliation:
Physics of Fluids Group and Max Planck Center Twente for Complex Fluid Dynamics, MESA+ Institute and J. M. Burgers Centre for Fluid Dynamics, University of Twente, P.O. Box 217, 7500AE Enschede, The Netherlands Department of Industrial Engineering, University of Rome ‘Tor Vergata’, Via del Politecnico 1, Roma 00133, Italy
Detlef Lohse
Affiliation:
Physics of Fluids Group and Max Planck Center Twente for Complex Fluid Dynamics, MESA+ Institute and J. M. Burgers Centre for Fluid Dynamics, University of Twente, P.O. Box 217, 7500AE Enschede, The Netherlands Max Planck Institute for Dynamics and Self-Organization, 37077 Göttingen, Germany
*
†Email address for correspondence: xjzhu@g.harvard.edu
Rights & Permissions [Opens in a new window]

Abstract

In turbulent Rayleigh–Bénard (RB) convection with regular, mono-scale, surface roughness, the scaling exponent $\unicode[STIX]{x1D6FD}$ in the relationship between the Nusselt number $Nu$ and the Rayleigh number $Ra$, $Nu\sim Ra^{\unicode[STIX]{x1D6FD}}$ can be ${\approx}1/2$ locally, provided that $Ra$ is large enough to ensure that the thermal boundary layer thickness $\unicode[STIX]{x1D706}_{\unicode[STIX]{x1D703}}$ is comparable to the roughness height. However, at even larger $Ra$, $\unicode[STIX]{x1D706}_{\unicode[STIX]{x1D703}}$ becomes thin enough to follow the irregular surface and $\unicode[STIX]{x1D6FD}$ saturates back to the value for smooth walls (Zhu et al., Phys. Rev. Lett., vol. 119, 2017, 154501). In this paper, we prevent this saturation by employing multiscale roughness. We perform direct numerical simulations of two-dimensional RB convection using an immersed boundary method to capture the rough plates. We find that, for rough boundaries that contain three distinct length scales, a scaling exponent of $\unicode[STIX]{x1D6FD}=0.49\pm 0.02$ can be sustained for at least three decades of $Ra$. The physical reason is that the threshold $Ra$ at which the scaling exponent $\unicode[STIX]{x1D6FD}$ saturates back to the smooth wall value is pushed to larger $Ra$, when the smaller roughness elements fully protrude through the thermal boundary layer. The multiscale roughness employed here may better resemble the irregular surfaces that are encountered in geophysical flows and in some industrial applications.

JFM classification

Turbulent Flows: Turbulent convection
Type
JFM Rapids
Information
Journal of Fluid Mechanics , Volume 869 , 25 June 2019 , R4
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution, and reproduction in any medium, provided the original work is properly cited.
Copyright
© 2019 Cambridge University Press

1 Introduction

Rayleigh–Bénard (RB) convection (Ahlers, Grossmann & Lohse Reference Ahlers, Grossmann and Lohse2009; Lohse & Xia Reference Lohse and Xia2010; Chillà & Schumacher Reference Chillà and Schumacher2012; Xia Reference Xia2013), a flow in a container heated from below and cooled from above, is a paradigmatic system in thermally driven turbulence. The key control parameters are the Rayleigh number and Prandtl number, which are respectively defined as $Ra=\unicode[STIX]{x1D6FC}g\unicode[STIX]{x1D6E5}L^{3}/(\unicode[STIX]{x1D708}\unicode[STIX]{x1D705})$ and $Pr=\unicode[STIX]{x1D708}/\unicode[STIX]{x1D705}$ , where $\unicode[STIX]{x1D6FC}$ is the thermal expansion coefficient, $g$ the gravitational acceleration, $\unicode[STIX]{x1D6E5}$ the temperature drop across the container, $L$ the height of the fluid domain, $\unicode[STIX]{x1D708}$ the kinematic viscosity, and $\unicode[STIX]{x1D705}$ the thermal diffusivity of the fluid. The most relevant response of the system is the heat transfer, which in dimensionless form is expressed as the Nusselt number $Nu$ .

Over the years, much attention has been paid to the scaling relation between $Nu$ and $Ra$ , i.e.  $Nu\sim Ra^{\unicode[STIX]{x1D6FD}}$ . Two of the early attempts were made by Malkus (Reference Malkus1954) and Priestley (Reference Priestley1954), both of whom independently proposed $\unicode[STIX]{x1D6FD}=1/3$ , which reflects their assumption that the heat flux is independent of the distance between the two plates and controlled only by the boundary layers (BLs). Grossmann & Lohse (Reference Grossmann and Lohse2000, Reference Grossmann and Lohse2001), based on an analysis and decomposition of the kinetic and thermal energy dissipation rates into bulk and BL contributions, proposed that there are no pure scaling laws but rather a superposition of various ones. For extremely large $Ra$ , Kraichnan (Reference Kraichnan1962) predicted a so-called ultimate regime with turbulent shear BLs, which led to the relation $Nu\sim Ra^{1/2}(\ln Ra)^{-3/2}$ , where the logarithmic correction becomes negligible with increasing $Ra$ (Spiegel Reference Spiegel1963). Yet there are still debates on the various claims of evidence for this regime. With a low-temperature helium RB experiment, Chavanne et al. (Reference Chavanne, Chilla, Castaing, Hebral, Chabaud and Chaussy1997, Reference Chavanne, Chilla, Chabaud, Castaing and Hebral2001) found that $\unicode[STIX]{x1D6FD}$ increases to 0.38 for $Ra=(2\times 10^{11},10^{14})$ . Taking into account the effects of turbulent BLs, Grossmann & Lohse (Reference Grossmann and Lohse2011) derived a scaling law with a different logarithmic correction as compared to Kraichnan (Reference Kraichnan1962) and formulated this relation as an effective power law with a locally determined effective scaling exponent $\unicode[STIX]{x1D6FD}$ . In particular, they derived that $\unicode[STIX]{x1D6FD}$ should be approximately 0.38 when $Ra$ is approximately $10^{14}$ . This was demonstrated experimentally by He et al. (Reference He, Funfschilling, Bodenschatz and Ahlers2012a ,Reference He, Funfschilling, Nobach, Bodenschatz and Ahlers b ). For more information on general aspects of RB convection, we refer the readers to the reviews by Ahlers et al. (Reference Ahlers, Grossmann and Lohse2009), Lohse & Xia (Reference Lohse and Xia2010), Chillà & Schumacher (Reference Chillà and Schumacher2012) and Xia (Reference Xia2013).

To avoid the influence of the BLs, and therefore to avoid the logarithmic corrections, several successful model systems have been proposed throughout the years. In numerical experiments with periodic boundary conditions in all directions, Lohse & Toschi (Reference Lohse and Toschi2003) and Calzavarini et al. (Reference Calzavarini, Lohse, Toschi and Tripiccione2005) proposed ‘homogeneous’ RB turbulence; Gibert et al. (Reference Gibert, Pabiou, Chilla and Castaing2006) and Pawar & Arakeri (Reference Pawar and Arakeri2018) performed corresponding RB experiments in a ‘cavity’; Lepot, Aumaître & Gallet (Reference Lepot, Aumaître and Gallet2018) proposed radiative heating convection, in which heat is input directly inside an absorption layer. When this absorption length is thicker than the BLs, radiative heating is allowed to bypass the BLs and heat up the bulk turbulent flow directly. In all these cases a scaling exponent of $1/2$ was achieved, because the BLs no longer played a role. We call this regime the ‘asymptotic ultimate regime’.

For conventional RB convection, in which BLs close to the bottom and top plate are formed, wall roughness has been introduced in an attempt to trigger an earlier onset of a turbulent BL; see the reviews Ahlers et al. (Reference Ahlers, Grossmann and Lohse2009), Chillà & Schumacher (Reference Chillà and Schumacher2012) and Xia (Reference Xia2013) for detailed discussions. The results for three-dimensional (3-D) simulations and experiments can be divided into two main categories. First, there are studies which show that roughness can increase the scaling exponent $\unicode[STIX]{x1D6FD}$ from a value slightly below $1/3$ to a value between $1/3$ and $1/2$ (Roche et al. Reference Roche, Castaing, Chabaud and Hebral2001; Qiu, Xia & Tong Reference Qiu, Xia and Tong2005; Stringano & Verzicco Reference Stringano and Verzicco2006; Tisserand et al. Reference Tisserand, Creyssels, Gasteuil, Pabiou, Gibert, Castaing and Chilla2011; Salort et al. Reference Salort, Liot, Rusaouen, Seychelles, Tisserand, Creyssels, Castaing and Chillá2014; Wei et al. Reference Wei, Chan, Ni, Zhao and Xia2014; Xie & Xia Reference Xie and Xia2017). Shen, Tong & Xia (Reference Shen, Tong and Xia1996), Du & Tong (Reference Du and Tong2000), Wei et al. (Reference Wei, Chan, Ni, Zhao and Xia2014) and Xie & Xia (Reference Xie and Xia2017) found that the scaling exponent $\unicode[STIX]{x1D6FD}$ remains roughly the same when roughness is introduced. Whether an increase in the scaling exponent is observed or not depends on the roughness configuration and the explored $Ra$ and $Pr$ regime. Roche et al. (Reference Roche, Castaing, Chabaud and Hebral2001) designed the first experiment to possibly reach the ultimate regime without the logarithmic correction (asymptotic ultimate regime) using a cylindrical cell with grooved roughness in both plates and in the sidewall. They observed a scaling exponent $\unicode[STIX]{x1D6FD}$ of 0.51 in the Rayleigh number range $Ra=(2\times 10^{12},5\times 10^{13})$ . Wagner & Shishkina (Reference Wagner and Shishkina2015) showed in direct numerical simulations (DNS) for rectangular roughness that the scaling exponent $\unicode[STIX]{x1D6FD}$ increases compared to the smooth case for low $Ra$ before it saturates back to the smooth wall value at large $Ra$ . For RB with pyramid-shaped roughness, Xie & Xia (Reference Xie and Xia2017) varied the roughness aspect ratio $\unicode[STIX]{x1D706}$ , which they define as the height of a roughness element over its base width, from $0.5$ to $4.0$ . With increasing $Ra$ they identified three regimes. The transition between regime I and II occurs when the thermal BL becomes thinner than the roughness height and the transition between regime II and III occurs when the viscous BL thickness becomes smaller than the roughness height. They found that in regime II the scaling exponent $\unicode[STIX]{x1D6FD}$ increases from $0.36$ to $0.59$ when $\unicode[STIX]{x1D706}$ is increased from $0.5$ to $4.0$ . In regime III they found that these scaling exponents saturate to $0.30$ to $0.50$ , respectively, with increasing $\unicode[STIX]{x1D706}$ . Rusaouën et al. (Reference Rusaouën, Liot, Castaing, Salort and Chillà2018) performed RB experiments with water in cylindrical containers for $Ra$ up to $10^{12}$ . They performed a set of measurements using smooth and rough plates with cubic roughness elements in a square lattice. In these experiments several regimes were identified for the rough case. With increasing $Ra$ they first observed a regime in which the heat transfer is similar to the smooth case, followed by a regime in which the heat transfer is enhanced by a modification of the $Nu$ versus $Ra$ number scaling, before a third regime is obtained in which the heat transfer scaling is similar to the smooth case, but with a larger prefactor.

For all the rough wall RB studies that we mentioned above, a 3-D geometry of the cell has been adopted. Recently, using DNS of two-dimensional (2-D) RB convection with roughness of varying heights and wavelengths for $Pr=1$ , Toppaladoddi, Succi & Wettlaufer (Reference Toppaladoddi, Succi and Wettlaufer2017) observed the existence of $\unicode[STIX]{x1D6FD}=0.483$ by fitting the data in the range $Ra=(4.6\times 10^{6},3\times 10^{9})$ and interpreted this exponent as an achievement of the ultimate regime. In contrast, Zhu et al. (Reference Zhu, Stevens, Verzicco and Lohse2017) showed that: (i) there is no pure scaling exponent in that $Ra$ range; (ii) although $\unicode[STIX]{x1D6FD}$ can locally reach $1/2$ in the range $Ra=(10^{8},3\times 10^{9})$ , this should not be interpreted as the attainment of the ultimate regime, because a further increase of $Ra$ leads to another regime where a thin thermal BL uniformly follows the rough surfaces, and thus the classical BL-controlled regime is recovered, causing the scaling to saturate to the classical effective $Nu$ versus $Ra$ scaling exponent close to $1/3$ .

The main question we want to address in this paper is: can the range of $Ra$ where the effective $1/2$ scaling exponent manifests be extended? We note that in all the studies mentioned above, uniform roughness of a single length scale was adopted. For this situation, the $1/2$ effective exponent can be observed when roughness starts to perturb the thermal BL, as mentioned before. If, with increasing $Ra$ , smaller and smaller roughness length scales are introduced, the different size roughness elements will protrude through the thermal BL one by one. Therefore, the flow can be maintained in a transition state and the $1/2$ effective exponent can be sustained. In this manuscript, we will demonstrate this conjecture by means of multiscale wall roughness. In fact, two decades ago Villermaux (Reference Villermaux1998) theoretically pioneered the research of RB convection with multiscale cubic roughness, with power-law-distributed asperity heights. He formulated a new scaling relation and found that the heat transfer scaling exponent can be significantly enhanced. Later, Ciliberto & Laroche (Reference Ciliberto and Laroche1999) experimentally explored multiscale roughness by gluing glass spheres of controlled diameter on the bottom copper plate, and found that $\unicode[STIX]{x1D6FD}$ increases to 0.45. In contrast, for a periodic roughness case, they found that the scaling exponent is similar to that in the smooth case.

Another motivation for this study is that for real-world applications and geophysical flows, the situation is far more complex, with surface roughness often containing different length scales. For example, in cities, there is huge difference among the heights of the buildings, and also natural terrains contain multiscale structure. Assuming roughness is multiscale provides a practically useful simplification (Rodriguez-Iturbe et al. Reference Rodriguez-Iturbe, Marani, Rigon and Rinaldo1994).

The paper is organized as follows: in § 2 we describe the numerical method and the parameter set-up used in the simulations. In § 3 we show how multiscale roughness alters RB turbulence. In § 4 we briefly summarize the results and give an outlook to potential future work.

2 Numerical details

We solve the Boussinesq equations with the second-order staggered finite-difference code AFiD (Verzicco & Orlandi Reference Verzicco and Orlandi1996; van der Poel et al. Reference van der Poel, Ostilla-Mónico, Donners and Verzicco2015; Zhu et al. Reference Zhu, Phillips, Spandan, Donners, Ruetsch, Romero, Ostilla-Mónico, Yang, Lohse, Verzicco, Massimiliano and Stevens2018b ) in 2-D. The reason why we resorted to 2-D simulations is that they are much less expensive than the 3-D case and thus we can cover a much wider range of $Ra$ . The details of the numerical methods, the parallelization and the different versions (CPU and GPU) can be found in Verzicco & Orlandi (Reference Verzicco and Orlandi1996), van der Poel et al. (Reference van der Poel, Ostilla-Mónico, Donners and Verzicco2015) and Zhu et al. (Reference Zhu, Phillips, Spandan, Donners, Ruetsch, Romero, Ostilla-Mónico, Yang, Lohse, Verzicco, Massimiliano and Stevens2018b ). The code has been extensively validated and used under various conditions (Zhu et al. Reference Zhu, Stevens, Verzicco and Lohse2017, Reference Zhu, Mathai, Stevens, Verzicco and Lohse2018a ,Reference Zhu, Phillips, Spandan, Donners, Ruetsch, Romero, Ostilla-Mónico, Yang, Lohse, Verzicco, Massimiliano and Stevens b ). The governing equations in the dimensionless form read:

(2.1) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\boldsymbol{u}}{\unicode[STIX]{x2202}t}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\boldsymbol{u}=-\unicode[STIX]{x1D735}p+\sqrt{\frac{Pr}{Ra}}\unicode[STIX]{x1D6FB}^{2}\boldsymbol{u}+\unicode[STIX]{x1D703}\hat{\boldsymbol{z}}, & \displaystyle\end{eqnarray}$$
(2.2) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}=0, & \displaystyle\end{eqnarray}$$
(2.3) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D703}}{\unicode[STIX]{x2202}t}+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}=\frac{1}{\sqrt{RaPr}}\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D703}, & \displaystyle\end{eqnarray}$$

where $\hat{\boldsymbol{z}}$ is the unit vector pointing in the direction opposite to gravity, $\boldsymbol{u}$ the velocity vector normalized by the free-fall velocity $\sqrt{g\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D6E5}L}$ , $t$ the dimensionless time normalized by $\sqrt{L/(g\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D6E5})}$ , $\unicode[STIX]{x1D703}$ the temperature normalized by $\unicode[STIX]{x1D6E5}$ , and $p$ the pressure normalized by $g\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D6E5}/L$ . As shown in the above equations, the control parameters of the system are $Ra$ and $Pr$ . The boundary conditions on the top and bottom plates are no-slip for the velocity and constant for the temperature. Periodic conditions are applied to the horizontal boundaries. In all our simulations, $Pr$ is fixed to 1. The aspect ratio is chosen as $\unicode[STIX]{x1D6E4}\equiv W/L=2$ , where $W$ is the width of the computational domain.

For the rough cases, the characteristic length scale that we use to express $Ra$ is the equivalent smooth wall height $L^{\prime }$ . This height is defined by determining the height of a domain with smooth boundaries that would have the same fluid volume. $Nu$ is calculated from $Nu=\sqrt{RaPr}\langle u_{z}\unicode[STIX]{x1D703}\rangle _{A}-\langle \unicode[STIX]{x2202}_{z}\unicode[STIX]{x1D703}\rangle _{A}$ , where $u_{z}$ denotes the instantaneous vertical velocity and $\langle \cdot \rangle _{A}$ the average over any horizontal plane between the rough plates.

Figure 1. (a) A sketch of the computational domain and the roughness elements. (b) Roughness element $R_{1}$ is the base element and the length scale is $R_{1}=0.1$ . The structure is multiscale as $R_{n+1}=2^{-n}R_{n}$ , $n=1,2,3$ .

An immersed boundary method (IBM) has been implemented to cope with rough surfaces (Fadlun et al. Reference Fadlun, Verzicco, Orlandi and Mohd-Yusof2000). The basic idea of the IBM is that a body force term in the Navier–Stokes equation can mimic the effects of the boundaries. For more information on IBM, we refer to the reviews by Peskin (Reference Peskin2002) and Mittal & Iaccarino (Reference Mittal and Iaccarino2005).

Figure 2. The instantaneous temperature fields at (a) $Ra=10^{8}$ , (b) $Ra=10^{9}$ , (c) $Ra=10^{10}$ and (d) $Ra=10^{11}$ . It can be seen that with $Ra$ increasing, plumes are ejected also from smaller and smaller roughness elements.

We now describe the multiscale roughness pattern: we choose a series of wall-mounted sinusoidal elements distributed on both the top and bottom plates. The sinusoidal elements all have the same aspect ratio 1. The multiscale roughness implementation is similar to that in Yang & Meneveau (Reference Yang and Meneveau2017), although in that study square roughness elements were adopted and positioned randomly. The size of the largest roughness element is used as the reference scale, $R_{1}=0.1$ . At the second and third generation, we have the rough elements size as $R_{n+1}=2^{-n}R_{n}$ . No roughness elements of intermediate sizes are included, and the roughness height spectrum is thus discrete (Yang & Meneveau Reference Yang and Meneveau2017). Figure 1 gives an overview of the computational domain and the roughness elements. Adequate resolution was ensured for all cases, i.e. the mesh is stretched in the wall normal direction with the finest grid implemented around the tips of the biggest roughness elements. There are at least 12 points inside the BL. The statistics were averaged over 200 free-fall time units. In the rough case for $Ra=10^{11}$ , $10\,240\times 5120$ grid points, in the horizontal and vertical direction, respectively, were used. Further details about the simulation parameters can be found in appendix A.

3 Results

We first compare the flow structures for increasing $Ra$ . Figure 2 shows the instantaneous temperature snapshots for four $Ra$ , ranging from $10^{8}$ to $10^{11}$ . At the lowest $Ra=10^{8}$ , within the cavity regions, the flow is viscosity dominated. Interestingly, figure 3 shows that the heat transfer for the case with multiscale roughness is approximately $15\,\%$ lower than for the case with smooth walls. The same phenomenon of heat transfer reduction due to roughness was observed by Shishkina & Wagner (Reference Shishkina and Wagner2011) and Zhang et al. (Reference Zhang, Sun, Bao and Zhou2018). The physical reason for the heat transfer reduction is that the hot/cold fluid is trapped between the roughness elements, which thus leads to a thicker thermal BL and therefore to a lower overall heat transport. For larger $Ra$ , plumes start to develop from the tips of the roughness elements and eventually, at the largest $Ra=10^{11}$ , plumes are formed even in the sloping surfaces of the smallest rough elements. This observation indicates the impact of multiscale roughness on the flow structure and heat transfer for increasing $Ra$ .

Figure 3. (a) $Nu$ as a function $Ra$ for the smooth case and the multiscale rough case. For the smooth case, the scaling exponent is $\unicode[STIX]{x1D6FD}=0.29\pm 0.01$ . For the multiscale rough case, the scaling exponent is $\unicode[STIX]{x1D6FD}=0.49\pm 0.02$ . As a reference, the results for mono-scale roughness are also included (Zhu et al. Reference Zhu, Stevens, Verzicco and Lohse2017), which clearly show two scaling regimes. Note that $Ra$ is defined based on the equivalent smooth wall height. In the mono-scale roughness case, 20 sinusoidal roughness elements of the same height (0.1) were adopted. For the multiscale roughness cases considered here, 10 of these large roughness elements are replaced by one $R_{2}$ and two $R_{3}$ generation roughness elements. Therefore, the total number of roughness elements for the multiscale roughness geometry is $40$ . $Nu$ is smaller for the multiscale roughness case than for the mono-scale roughness case, because the latter has larger roughness elements. (b) Same as in (a) but in a compensated way for the multiscale rough case. Note that we use only one specific aspect ratio for the roughness elements. If the aspect ratio changes, the scaling exponent will also change.

Figure 4. Sketches on why regular periodic roughness with the same height leads to scaling saturation and why multiscale roughness increases the exponent in a wider range of $Ra$ . Orange parts are the regions where the thermal BL are. (a) At lower $Ra$ , the roughness is below the thermal BL and has little impact on scaling relations. (b,c) At intermediate $Ra$ , the roughness starts to protrude through the thermal BL, but not to the valley of the roughness elements. For multiscale roughness, it is easy to imagine that the range of $Ra$ is wider in this stage, as only with increasing $Ra$ will the smaller and smaller roughness elements protrude through the thermal BL. (d) When Ra is large enough, a thin thermal BL is uniformly distributed along the rough surfaces and the scaling exponent will saturate back to the value close to the smooth case. This case is not reached in this study.

Next, we check how the scaling relation evolves with increasing $Ra$ , comparing the smooth and the multiscale case. Figure 3 shows the $Nu$ scaling behaviour as a function of $Ra$ , in a log–log plot and in a compensated plot. As was shown before in 2-D RB (DeLuca et al. Reference Deluca, Werne, Rosner and Cattaneo1990; Johnston & Doering Reference Johnston and Doering2009; Zhu et al. Reference Zhu, Stevens, Verzicco and Lohse2017), the smooth case has an effective scaling exponent $\unicode[STIX]{x1D6FD}=0.29\pm 0.01$ , extending over four decades, from $Ra=10^{8}$ to $Ra=10^{12}$ . For the mono-scale roughness case, two distinct effective scaling exponent can be observed, i.e.  $\unicode[STIX]{x1D6FD}=0.50\pm 0.02$ , for one and half decade; then $\unicode[STIX]{x1D6FD}$ saturates back to 0.33. With the introduction of multiscale roughness, the heat transfer is greatly enhanced. Within 95 % of the confidence bound, we get the fit of $Nu\sim 0.00257Ra^{0.49\pm 0.02}$ for three decades of $Ra$ , from $Ra=10^{8}$ to $Ra=10^{11}$ . A root mean square error 2.89 is found for the fit. To our knowledge, this is the first realization of such a large scaling exponent in such a wide range of $Ra$ in RB. It is remarkable that it is realized in spite of the relatively low $Ra$ numbers of the simulations. Obviously, as for the smooth RB, an asymptotic $1/2$ exponent is expected only when $Ra$ approaches infinity.

We will now explain the physical mechanism which leads to this considerable enhancement of the exponent for a wide range of $Ra$ . Let us first recall why in the case of periodic roughness with one single height, the regime with an effective scaling exponent close to $1/2$ survives for a limited range of $Ra$ and then saturates back to a value close to the smooth case. In the former regime, the roughness elements start to protrude into the thermal BL. Only weak secondary vortices are generated in the cavities and the resulting mixing is not efficient. Therefore, the flow there is still dominated by viscosity. In the latter regime, secondary vortices are strong enough to mix all the fluids inside the cavities. Thus the roughness elements are covered by a thin thermal BL which is uniformly distributed along the rough surfaces, effectively mimicking an increased wet surface area. Therefore, the effect of BL is restored and the classical BL-controlled regime is retrieved (Zhu et al. Reference Zhu, Stevens, Verzicco and Lohse2017).

Multiscale roughness essentially extends this effective $1/2$ scaling regime further. Figure 4 shows sketches of thermal BLs for increasing $Ra$ . As $Ra$ increases, the thermal BL becomes thinner, the smaller roughness elements start to perturb the thermal BL, and this process continues until the smallest roughness elements perturb the thermal BL. Therefore, the system stays in the transitioning state and the $1/2$ exponent is observed over a wider range of $Ra$ , compared to the case of periodic roughness with the same height. To give further evidence for this explanation, in figure 5, we show the averaged mean temperature profiles for the valley points in the cavity regions of the roughness elements close to the wall. From the temperature profile for $Ra=10^{8}$ we can detect the influence of the largest roughness $R_{1}$ . At $Ra=10^{9}$ also the effect of the $R_{2}$ roughness can be identified. Last, but not least, at $Ra=10^{10}$ , the influence of smallest $R_{3}$ roughness starts to manifest.

Figure 5. Mean temperature profile as a function of the wall normal coordinate averaged at the $x$ -locations of the valley points in the cavity regions, i.e. the locations where no roughness is added on top of the plate. As explained in the caption of figure 1 there are 40 roughness elements, and these profiles are averaged over all 40 corresponding valley locations.

4 Summary, discussion, and outlook

In this manuscript, we use 2-D DNS to study the effects of multiscale roughness on RB turbulence. In our case, the multiscale roughness is composed of three different roughness length scales of sinusoidal shapes. We show that with this implementation the $Nu$ versus $Ra$ scaling relation $Nu\sim Ra^{0.49\pm 0.02}$ can be observed for at least three decades of $Ra$ , while for mono-scale roughness the scaling could be observed only for only one and half decades in $Ra$ (Zhu et al. Reference Zhu, Stevens, Verzicco and Lohse2017). However, there are still several open issues, for example:

  1. (i) We stress that even though we found $Nu\sim Ra^{0.5}$ over an extended $Ra$ range, this is probably still a transitional regime and not the (asymptotic) ultimate regime in which the BLs are fully turbulent. This means that it is very likely that for the considered roughness geometry the heat transfer scaling exponent saturates back to the value of the smooth case for larger $Ra$ . The situation may also change once the BL become turbulent. Recently, MacDonald et al. (Reference MacDonald, Hutchins, Lohse and Chung2019) found an effective exponent of $\unicode[STIX]{x1D6FD}\approx 0.42$ in the large- $Ra$ regime for forced convection in channel flow under the assumption that the BL profile become logarithmic.

  2. (ii) In this work, we modelled three different roughness length scales. We expect that with more length scales, the $Ra$ range in which the $1/2$ scaling exponent manifests might be more extended. Simulations for RB with more roughness length scales are needed to settle this question.

  3. (iii) The current simulations are for 2-D only and it would be interesting to see results in a fully 3-D case. Such simulations would be much more computationally intensive, but very interesting comparisons to the experimental results by, for example, Roche et al. (Reference Roche, Castaing, Chabaud and Hebral2001) would be possible. In that study a scaling exponent of approximately $1/2$ up to approximately $Ra=5\times 10^{13}$ was observed.

  4. (iv) Ciliberto & Laroche (Reference Ciliberto and Laroche1999) observed $\unicode[STIX]{x1D6FD}\approx 0.45$ for RB with multiscale glass spheres glued on a copper plate. Understanding the effect of the different heat conductivities for the two roughness elements (glue and glass, both with smaller heat conductivity than copper) will be very helpful.

  5. (v) Here we consider only the situation for $Pr=1$ . Recent experiments by Xie & Xia (Reference Xie and Xia2017) suggest that also the $Pr$ number might play an important role on the effect of roughness on the overall heat transport. Simulations to investigate this effect would be very interesting.

  6. (vi) Finally we note that for turbulent Taylor–Couette flow with mono-scale roughness that aligns in the azimuthal direction, simulations (Zhu et al. Reference Zhu, Ostilla-Monico, Verzicco and Lohse2016) yielded an intermediate regime where ‘Nusselt number’ $Nu_{\unicode[STIX]{x1D714}}\sim Ta^{0.5}$ (angular velocity transport versus Taylor number). When $Ta$ is large enough, the exponent saturates back to the smooth case value. The question is whether the $Ta$ range where the $1/2$ exponent shows up can also be extended with multiscale roughness. For turbulent Taylor–Couette flow with mono-scale roughness that aligns in the axial direction, the asymptotic ultimate regime $1/2$ scaling has already been achieved (Cadot et al. Reference Cadot, Couder, Daerr, Douady and Tsinober1997; van den Berg et al. Reference van den Berg, Doering, Lohse and Lathrop2003; Zhu et al. Reference Zhu, Verschoof, Bakhuis, Huisman, Verzicco, Sun and Lohse2018c ), and pressure drag has been identified as the origin thereof (Zhu et al. Reference Zhu, Verschoof, Bakhuis, Huisman, Verzicco, Sun and Lohse2018c ).

Acknowledgements

This work was financially supported by the European Research Council (ERC) under the European Unions Horizon 2020 research and innovation programme (grant agreement no. 740479) and the Netherlands Center for Multiscale Catalytic Energy Conversion (MCEC), an NWO Gravitation programme funded by the Ministry of Education, Culture and Science of the Government of the Netherlands. We acknowledge support by the Priority Program SPP 1881 ‘Turbulent Superstructures’ of the Deutsche Forschungsgemeinschaft (DFG). This work was partially carried out on the Dutch national e-infrastructure with the support of SURF Cooperative. We also acknowledge PRACE for awarding us access to MareNostrum based in Italy at the Barcelona Supercomputing Center (BSC) and JUWELS at the Jülich Supercomputing Centre under PRACE project number 2017174146.

Appendix. Numerical details

The numerical details are given in table 1.

Table 1. $Ra$ , resolution in the horizontal ( $n_{x}$ ) and wall normal ( $n_{z}$ ) directions, and $Nu$ number for the multiscale roughness cases considered in this study. For all cases the domain aspect ratio is 2 and $Pr=1$ . The uncertainties in $Nu$ is smaller than $1\,\%$ for all cases. Corresponding information for the mono-scale roughness cases has been reported in Zhu et al. (Reference Zhu, Stevens, Verzicco and Lohse2017) and for the smooth case in Zhu et al. (Reference Zhu, Mathai, Stevens, Verzicco and Lohse2018a ).

References

Ahlers, G., Grossmann, S. & Lohse, D. 2009 Heat transfer and large scale dynamics in turbulent Rayleigh–Bénard convection. Rev. Mod. Phys. 81, 503–538.Google Scholar
van den Berg, T. H., Doering, C. R., Lohse, D. & Lathrop, D. 2003 Smooth and rough boundaries in turbulent Taylor–Couette flow. Phys. Rev. E 68, 036307.Google Scholar
Cadot, O., Couder, Y., Daerr, A., Douady, S. & Tsinober, A. 1997 Energy injection in closed turbulent flows: stirring through boundary layers versus inertial stirring. Phys. Rev.  E 56, 427–433.Google Scholar
Calzavarini, E., Lohse, D., Toschi, F. & Tripiccione, R. 2005 Rayleigh and Prandtl number scaling in the bulk of Rayleigh–Bénard turbulence. Phys. Fluids 17, 055107.Google Scholar
Chavanne, X., Chilla, F., Castaing, B., Hebral, B., Chabaud, B. & Chaussy, J. 1997 Observation of the ultimate regime in Rayleigh–Bénard convection. Phys. Rev. Lett. 79, 3648–3651.Google Scholar
Chavanne, X., Chilla, F., Chabaud, B., Castaing, B. & Hebral, B. 2001 Turbulent Rayleigh–Bénard convection in gaseous and liquid He. Phys. Fluids 13, 1300–1320.Google Scholar
Chillà, F. & Schumacher, J. 2012 New perspectives in turbulent Rayleigh–Bénard convection. Eur. Phys. J. E 35, 58.Google Scholar
Ciliberto, S. & Laroche, C. 1999 Random roughness of boundary increases the turbulent convection scaling exponent. Phys. Rev. Lett. 82, 3998–4001.Google Scholar
Deluca, E. E., Werne, J., Rosner, R. & Cattaneo, F. 1990 Numerical simulations of soft and hard turbulence: preliminary results for two-dimensional convection. Phys. Rev. Lett. 64 (20), 2370.Google Scholar
Du, Y. B. & Tong, P. 2000 Turbulent thermal convection in a cell with ordered rough boundaries. J. Fluid Mech. 407, 57–84.Google Scholar
Fadlun, E. A., Verzicco, R., Orlandi, P. & Mohd-Yusof, J. 2000 Combined immersed-boundary finite-difference methods for three-dimensional complex flow simulations. J. Comput. Phys. 161, 35–60.Google Scholar
Gibert, M., Pabiou, H., Chilla, F. & Castaing, B. 2006 High-Rayleigh-number convection in a vertical channel. Phys. Rev. Lett. 96, 084501.Google Scholar
Grossmann, S. & Lohse, D. 2000 Scaling in thermal convection: a unifying view. J. Fluid. Mech. 407, 27–56.Google Scholar
Grossmann, S. & Lohse, D. 2001 Thermal convection for large Prandtl number. Phys. Rev. Lett. 86, 3316–3319.Google Scholar
Grossmann, S. & Lohse, D. 2011 Multiple scaling in the ultimate regime of thermal convection. Phys. Fluids 23, 045108.Google Scholar
He, X., Funfschilling, D., Bodenschatz, E. & Ahlers, G. 2012a Heat transport by turbulent Rayleigh–Bénard convection for Pr = 0. 8 and 4 × 1011 < Ra < 2 × 1014 : ultimate-state transition for aspect ratio 𝛤 = 1. 00. New J. Phys. 14, 063030.Google Scholar
He, X., Funfschilling, D., Nobach, H., Bodenschatz, E. & Ahlers, G. 2012b Transition to the ultimate state of turbulent Rayleigh–Bénard convection. Phys. Rev. Lett. 108, 024502.Google Scholar
Johnston, H. & Doering, C. R. 2009 Comparison of turbulent thermal convection between conditions of constant temperature and constant flux. Phys. Rev. Lett. 102, 064501.Google Scholar
Kraichnan, R. H. 1962 Turbulent thermal convection at arbritrary Prandtl number. Phys. Fluids 5, 1374–1389.Google Scholar
Lepot, S., Aumaître, S. & Gallet, B. 2018 Radiative heating achieves the ultimate regime of thermal convection. Proc. Natl Acad. Sci. USA 115 (36), 8937–8941.Google Scholar
Lohse, D. & Toschi, F. 2003 The ultimate state of thermal convection. Phys. Rev. Lett. 90, 034502.Google Scholar
Lohse, D. & Xia, K.-Q. 2010 Small-scale properties of turbulent Rayleigh–Bénard convection. Annu. Rev. Fluid Mech. 42, 335–364.Google Scholar
MacDonald, M., Hutchins, N., Lohse, D. & Chung, D.2019 Heat transfer in fully-rough-wall-bounded turbulent flow in the ultimate regime. Phys. Rev. Fluids (submitted).Google Scholar
Malkus, M. V. R. 1954 The heat transport and spectrum of thermal turbulence. Proc. R. Soc. Lond. A 225, 196–212.Google Scholar
Mittal, R. & Iaccarino, G. 2005 Immersed boundary methods. Annu. Rev. Fluid Mech. 37, 239–261.Google Scholar
Pawar, S. S. & Arakeri, J. H. 2018 Two regimes of flux scaling in axially homogeneous turbulent convection in vertical tube. Phys. Rev. Fluids 1 (4), 042401(R).Google Scholar
Peskin, C. S. 2002 The immersed boundary method. Acta Numer. 11, 479–517.Google Scholar
van der Poel, E. P., Ostilla-Mónico, R., Donners, J. & Verzicco, R. 2015 A pencil distributed finite difference code for strongly turbulent wall–bounded flows. Comput. Fluids 116, 10–16.Google Scholar
Priestley, C. H. B. 1954 Convection from a large horizontal surface. Aust. J. Phys. 7, 176–201.Google Scholar
Qiu, X. L., Xia, K.-Q. & Tong, P. 2005 Experimental study of velocity boundary layer near a rough conducting surface in turbulent natural convection. J. Turbul. 6, 1–13.Google Scholar
Roche, P. E., Castaing, B., Chabaud, B. & Hebral, B. 2001 Observation of the 1/2 power law in Rayleigh–Bénard convection. Phys. Rev. E 63, 045303.Google Scholar
Rodriguez-Iturbe, I., Marani, M., Rigon, R. & Rinaldo, A. 1994 Self-organized river basin landscapes: fractal and multifractal characteristics. Water Resour. Res. 30 (12), 3531–3539.Google Scholar
Rusaouën, E., Liot, O., Castaing, B., Salort, J. & Chillà, F. 2018 Thermal transfer in Rayleigh–Bénard cell with smooth or rough boundaries. J. Fluid Mech. 837, 443–460.Google Scholar
Salort, J., Liot, O., Rusaouen, E., Seychelles, F., Tisserand, J.-C., Creyssels, M., Castaing, B. & Chillá, F. 2014 Thermal boundary layer near roughnesses in turbulent Rayleigh–Bénard convection: flow structure and multistability. Phys. Fluids 26, 015112.Google Scholar
Shen, Y., Tong, P. & Xia, K.-Q. 1996 Turbulent convection over rough surfaces. Phys. Rev. Lett. 76, 908–911.Google Scholar
Shishkina, O. & Wagner, C. 2011 Modelling the influence of wall roughness on heat transfer in thermal convection. J. Fluid Mech. 686, 568–582.Google Scholar
Spiegel, E. A. 1963 A generalization of the mixing-length theory of turbulent convection. Astrophys. J. 138, 216–225.Google Scholar
Stringano, G. & Verzicco, R. 2006 Mean flow structure in thermal convection in a cylindrical cell of aspect-ratio one half. J. Fluid Mech. 548, 1–16.Google Scholar
Tisserand, J. C., Creyssels, M., Gasteuil, Y., Pabiou, H., Gibert, M., Castaing, B. & Chilla, F. 2011 Comparison between rough and smooth plates within the same Rayleigh–Bénard cell. Phys. Fluids 23 (1), 015105.Google Scholar
Toppaladoddi, S., Succi, S. & Wettlaufer, J. S. 2017 Roughness as a route to the ultimate regime of thermal convection. Phys. Rev. Lett. 118, 074503.Google Scholar
Verzicco, R. & Orlandi, P. 1996 A finite-difference scheme for three-dimensional incompressible flow in cylindrical coordinates. J. Comput. Phys. 123, 402–413.Google Scholar
Villermaux, E. 1998 Transfer at rough sheared interfaces. Phys. Rev. Lett. 81 (22), 4859–4862.Google Scholar
Wagner, S. & Shishkina, O. 2015 Heat flux enhancement by regular surface roughness in turbulent thermal convection. J. Fluid Mech. 763, 109–135.Google Scholar
Wei, P., Chan, T.-S., Ni, R., Zhao, X.-Z. & Xia, K.-Q. 2014 Heat transport properties of plates with smooth and rough surfaces in turbulent thermal convection. J. Fluid Mech. 740, 28–46.Google Scholar
Xia, K.-Q. 2013 Current trends and future directions in turbulent thermal convection. Theor. Appl. Mech. Lett. 3 (5), 052001.Google Scholar
Xie, Y.-C. & Xia, K.-Q. 2017 Turbulent thermal convection over rough plates with varying roughness geometries. J. Fluid Mech. 825, 573–599.Google Scholar
Yang, X. I. A. & Meneveau, C. 2017 Modelling turbulent boundary layer flow over fractal-like multiscale terrain using large-eddy simulations and analytical tools. Phil. Trans. R. Soc. Lond. A 375 (2091), 20160098.Google Scholar
Zhang, Y.-Z., Sun, C., Bao, Y. & Zhou, Q. 2018 How surface roughness reduces heat transport for small roughness heights in turbulent Rayleigh–Bénard convection. J. Fluid Mech. 836, R2.Google Scholar
Zhu, X., Mathai, V., Stevens, R. J. A. M., Verzicco, R. & Lohse, D. 2018a Transition to the ultimate regime in two-dimensional Rayleigh–Bénard convection. Phys. Rev. Lett. 120 (14), 144502.Google Scholar
Zhu, X., Ostilla-Monico, R., Verzicco, R. & Lohse, D. 2016 Direct numerical simulation of Taylor–Couette flow with grooved walls: torque scaling and flow structure. J. Fluid Mech. 794, 746–774.Google Scholar
Zhu, X., Phillips, E., Spandan, V., Donners, J., Ruetsch, G., Romero, J., Ostilla-Mónico, R., Yang, Y., Lohse, D., Verzicco, R., Massimiliano, F. & Stevens, R. J. A. M. 2018b AFiD-GPU: a versatile Navier–Stokes solver for wall-bounded turbulent flows on GPU clusters. Comput. Phys. Commun. 229, 199–210.Google Scholar
Zhu, X., Stevens, R. J. A. M., Verzicco, R. & Lohse, D. 2017 Roughness-facilitated local 1/2 scaling does not imply the onset of the ultimate regime of thermal convection. Phys. Rev. Lett. 119 (15), 154501.Google Scholar
Zhu, X., Verschoof, R. A., Bakhuis, D., Huisman, S. G., Verzicco, R., Sun, C. & Lohse, D. 2018c Wall roughness induces asymptotic ultimate turbulence. Nat. Phys. 14 (4), 417–423.Google Scholar
Figure 0

Figure 1. (a) A sketch of the computational domain and the roughness elements. (b) Roughness element $R_{1}$ is the base element and the length scale is $R_{1}=0.1$. The structure is multiscale as $R_{n+1}=2^{-n}R_{n}$, $n=1,2,3$.

Figure 1

Figure 2. The instantaneous temperature fields at (a) $Ra=10^{8}$, (b) $Ra=10^{9}$, (c) $Ra=10^{10}$ and (d) $Ra=10^{11}$. It can be seen that with $Ra$ increasing, plumes are ejected also from smaller and smaller roughness elements.

Figure 2

Figure 3. (a) $Nu$ as a function $Ra$ for the smooth case and the multiscale rough case. For the smooth case, the scaling exponent is $\unicode[STIX]{x1D6FD}=0.29\pm 0.01$. For the multiscale rough case, the scaling exponent is $\unicode[STIX]{x1D6FD}=0.49\pm 0.02$. As a reference, the results for mono-scale roughness are also included (Zhu et al.2017), which clearly show two scaling regimes. Note that $Ra$ is defined based on the equivalent smooth wall height. In the mono-scale roughness case, 20 sinusoidal roughness elements of the same height (0.1) were adopted. For the multiscale roughness cases considered here, 10 of these large roughness elements are replaced by one $R_{2}$ and two $R_{3}$ generation roughness elements. Therefore, the total number of roughness elements for the multiscale roughness geometry is $40$. $Nu$ is smaller for the multiscale roughness case than for the mono-scale roughness case, because the latter has larger roughness elements. (b) Same as in (a) but in a compensated way for the multiscale rough case. Note that we use only one specific aspect ratio for the roughness elements. If the aspect ratio changes, the scaling exponent will also change.

Figure 3

Figure 4. Sketches on why regular periodic roughness with the same height leads to scaling saturation and why multiscale roughness increases the exponent in a wider range of $Ra$. Orange parts are the regions where the thermal BL are. (a) At lower $Ra$, the roughness is below the thermal BL and has little impact on scaling relations. (b,c) At intermediate $Ra$, the roughness starts to protrude through the thermal BL, but not to the valley of the roughness elements. For multiscale roughness, it is easy to imagine that the range of $Ra$ is wider in this stage, as only with increasing $Ra$ will the smaller and smaller roughness elements protrude through the thermal BL. (d) When Ra is large enough, a thin thermal BL is uniformly distributed along the rough surfaces and the scaling exponent will saturate back to the value close to the smooth case. This case is not reached in this study.

Figure 4

Figure 5. Mean temperature profile as a function of the wall normal coordinate averaged at the $x$-locations of the valley points in the cavity regions, i.e. the locations where no roughness is added on top of the plate. As explained in the caption of figure 1 there are 40 roughness elements, and these profiles are averaged over all 40 corresponding valley locations.

Figure 5

Table 1. $Ra$, resolution in the horizontal ($n_{x}$) and wall normal ($n_{z}$) directions, and $Nu$ number for the multiscale roughness cases considered in this study. For all cases the domain aspect ratio is 2 and $Pr=1$. The uncertainties in $Nu$ is smaller than $1\,\%$ for all cases. Corresponding information for the mono-scale roughness cases has been reported in Zhu et al. (2017) and for the smooth case in Zhu et al. (2018a).