2 Simulation details

We have conducted a set of DNS of 3D rotating RB in a horizontally periodic Cartesian computational box. By using a second-order energy-conserving, finite-difference code with fractional time stepping (Verzicco & Orlandi Reference Verzicco and Orlandi1996), we march in time the Navier–Stokes equations plus an advection–diffusion equation for temperature, with the usual Boussinesq approximations (Chandrasekhar Reference Chandrasekhar1961):

with the incompressibility constraint

Here $\boldsymbol{u}$ is the velocity vector; $t$ is time; $\boldsymbol{e}_{z}$ is the unit vector in the vertical direction; $Ra$ is the Rayleigh number, i.e. the non-dimensional temperature difference, defined as $Ra=g{\it\beta}{\it\Delta}L^{3}/({\it\nu}{\it\kappa})$ , with $L$ the height of the system, ${\it\beta}$ the thermal expansion coefficient of the fluid, $g$ the gravitational acceleration, ${\it\Delta}$ the temperature difference between the bottom and top plates, and ${\it\nu}$ and ${\it\kappa}$ the kinematic viscosity and thermal diffusivity of the fluid, respectively; $Ro$ is the Rossby number, i.e. the inverse rotation rate, defined as $Ro=\sqrt{({\it\beta}g{\it\Delta}/L)}/(2{\it\Omega})$ , with ${\it\Omega}$ the angular rotation rate; $Pr$ is the Prandtl number of the fluid, $Pr={\it\nu}/{\it\kappa}$ ; and ${\it\theta}$ is the non-dimensional temperature. Equations (2.1)–(2.3) are non-dimensionalised by using $L$ , ${\it\Delta}$ and the so-called free-fall velocity scale $U=\sqrt{{\it\beta}g{\it\Delta}L}$ . Note that centrifugal buoyancy is neglected here; this means that we are implicitly making the customary assumption that the Froude number $Fr={\it\Omega}^{2}R/g\ll 1$ (Stevens et al. Reference Stevens, Clercx and Lohse2013a ), where $R$ is the horizontal distance to the rotation axis, i.e. the radius of the cylinder in most experiments and simulations. We also define the Nusselt number, i.e. the non-dimensional heat transfer, as $Nu=(\langle u_{z}{\it\theta}\rangle _{A,t}-{\it\kappa}\partial \langle {\it\theta}\rangle _{A,t})/({\it\kappa}{\it\Delta}L^{-1})$ , with $\langle \cdots \rangle _{A,t}$ representing the averaging operator in time and also spatially over a horizontal plane.

The explored parameter values are given in table 1. We vary the rotation rate $Ro$ at a constant thermal driving $Ra$ for two values of $Ra$ , while fixing the Prandtl number to $Pr=1$ . The Ekman number, defined as $Ek={\it\nu}/(2{\it\Omega}L^{2})=Ro\sqrt{Pr/Ra}$ , is small enough to enter into the geostrophic regime, as can be seen from table 1. For completeness, we also define the Taylor number $Ta=Ek^{-2}$ . The aspect ratio, ${\it\Gamma}=D/L$ , where $D$ is the simulation box periodicity length in the horizontal directions, is set to 10 times the most unstable wavelength for convective instability, $L_{c}$ , i.e. ${\it\Gamma}=10L_{c}$ . Wavelength $L_{c}$ scales asymptotically as $L_{c}=4.82Ek^{1/3}$ , with minor corrections at finite $Ek$ (Chandrasekhar Reference Chandrasekhar1961; Niiler & Bisshopp Reference Niiler and Bisshopp1965; Heard & Veronis Reference Heard and Veronis1971). Here, we just take ${\it\Gamma}=48.2Ek^{1/3}$ . The boundary conditions for temperature are fixed as ${\it\theta}=1$ at the bottom plate and ${\it\theta}=0$ at the top plate. For velocity, we employ both no-slip (NS) boundary conditions, i.e. $\boldsymbol{u}=\mathbf{0}$ at both plates, and stress-free (SF) boundary conditions, i.e. $\partial _{z}u_{x}=\partial _{z}u_{y}=0$ and $u_{z}=0$ at the plates.

Table 1.

Parameter values for the computations. For all runs, $Pr=1$ . Each parameter set has been run with both NS and SF boundary conditions. Included are: Rayleigh number $Ra$ , Ekman number $Ek$ , Rossby number $Ro$ , $\widetilde{Ra}=Ra\,Ek^{4/3}$ (Sprague et al. Reference Sprague, Julien, Knobloch and Werne2006), domain aspect ratio ${\it\Gamma}=D/L$ (ratio of horizontal periodicity length $D$ to domain height $L$ ) and number of grid points $N_{x}\times N_{y}\times N_{z}$ in the periodic directions and the vertical direction, respectively. We also list the resulting Nusselt numbers $Nu_{SF}$ and $Nu_{NS}$ , as well as the number of grid points $N({\it\delta}_{{\it\nu}})$ within the kinetic boundary layer in the NS runs (see also § 4.1). Typical time step sizes are $0.004$ and $0.002$ for $Ra=1\times 10^{10}$ and $5\times 10^{10}$ , respectively, expressed in convective time units $H/U$ .

The grid is uniformly discretised in the horizontal direction. In the vertical direction, a clipped Chebyshev distribution is used to cluster points near the boundary layers. While the resulting boundary-layer thicknesses will be presented in § 4.1, we comment here on the adequacy of the near-wall resolution. The number of grid points found within the kinetic boundary layer is reported in the last column of table 1. The most demanding simulation at $Ra=5\times 10^{10}$ and $Ek=1.34\times 10^{-7}$ has 12 grid points within the kinetic boundary layer; the same is true for the case $Ra=1\times 10^{10}$ and $Ek=4.00\times 10^{-7}$ . Higher Ekman numbers result in larger boundary-layer thicknesses and thus they are better resolved. The grid resolutions used here are in line with recent computations by Stellmach et al. (Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014), and thus we consider them to be adequate.

For statistical convergence, the simulations are run until the temperature gradient at both plates is equal to less than $1\,\%$ , and the temporal convergence error of the Nusselt number measured through a volume average is smaller than $2\,\%$ . In practice, this means running the simulations between 100 and 200 convective time units, which are defined using the system height $L$ and the free-fall velocity $U$ . The total time varies: the lower $Ek$ , the slower the convergence of statistics.

We tested the dependence on ${\it\Gamma}$ of the simulations by running two cases at $Ra=1\times 10^{10}$ and $Ro=0.04$ with a twice larger ${\it\Gamma}$ (thus increasing the computational box size by a factor four, and hence the computational load by at least that amount). The Nusselt number is the same for no-slip plates but shows some difference in the stress-free case. For no-slip plates, we do not expect a strong dependence on ${\it\Gamma}$ as long as it is large enough. However, for stress-free plates, differences may occur that are related to the flow structure. This will be addressed in more detail in § 5.

Figure 1.

Phase diagram of rotating convection (a) in the $(Ek,Ra/Ra_{c})$ parameter space, as suggested by Ecke & Niemela (Reference Ecke and Niemela2014), and (b) in the $(Ta,Ra)$ parameter space (Vorobieff & Ecke Reference Vorobieff and Ecke2002). The blue circles indicate the parameter values for which we have performed simulations using both no-slip (NS) and stress-free (SF) boundary conditions. The circles are either open for $Ra=1\times 10^{10}$ or filled for $Ra=5\times 10^{10}$ . Three different regimes of convection can be discerned: non-rotating convection, rotation-affected convection and (rotation-dominated) geostrophic convection. The lines display various relations suggested in the literature that bound the regimes. Ecke & Niemela (Reference Ecke and Niemela2014) suggested $Ra/Ra_{c}=3$ as a lower bound of the geostrophic regime (dotted black line), $Ra=0.25Ek^{-1.8}$ for the transition between the geostrophic and rotation-affected regimes (solid black line), and $Ro=0.35$ for the transition to the non-rotating regime (dashed black line). The dash-dotted black line is a transition valid for higher $Pr\approx 6$ (Ecke & Niemela Reference Ecke and Niemela2014). Two alternative predictions for the transition to the geostrophic regime are also displayed: upper grey dashed line, $Ra=1.4Ek^{-7/4}$ (King et al. Reference King, Stellmach, Noir, Hansen and Aurnou2009); and lower grey dashed line, $Ra\approx 10Ek^{-3/2}$ (King, Stellmach & Aurnou Reference King, Stellmach and Aurnou2012).

To indicate how the current simulations fit in with the previous work on this topic, we display our current parameter values in the ( $Ek,Ra/Ra_{c}$ ) phase diagram of figure 1, where $Ra_{c}$ is the critical Rayleigh number for the onset of the convective instability with rotation, i.e. $Ra_{c}=8.7Ek^{-4/3}$ (Chandrasekhar Reference Chandrasekhar1961), and in the ( $Ta,Ra$ ) phase diagram. Figure 1(a) is based on figure 4 of Ecke & Niemela (Reference Ecke and Niemela2014), and shows how the parameter values of the simulations in this paper are positioned relative to some of the proposed bounds on the geostrophic regime. Typically, a lower bound for the geostrophic regime is chosen such that the resulting flow is supercritical enough for a turbulent flow to develop. We follow Ecke & Niemela (Reference Ecke and Niemela2014) by choosing $Ra/Ra_{c}=3$ (dotted black line in figure 1). Ecke & Niemela (Reference Ecke and Niemela2014) discern two additional transitions based on their heat-transfer measurements in cryogenic helium with $Pr=0.7$ . When reducing the Rayleigh number at constant $Ek$ , the first transition seen is when rotation starts to reduce the heat transfer. This is well described by $Ro=0.35$ (dashed black line in figure 1). When $Ra$ is reduced even further, a transition to a steeper scaling law relating $Nu$ and $Ra$ was found. Ecke & Niemela (Reference Ecke and Niemela2014) interpreted this as the transition to the geostrophic regime. This transition was best described by the relation $Ra=0.25Ek^{-1.8}$ (solid black line in figure 1) according to their data. For completeness, two other suggested relations for the transfer to the geostrophic regime have also been included in figure 1: $Ra=1.4Ek^{-7/4}$ as suggested by King et al. (Reference King, Stellmach, Noir, Hansen and Aurnou2009) (upper grey dashed line) and $Ra\approx 10Ek^{-3/2}$ from King et al. (Reference King, Stellmach and Aurnou2012) (lower grey dashed line). From the diagram, we expect our simulations to show a transition from the rotation-affected to the geostrophic regime for all criteria, except the one suggested by King et al. (Reference King, Stellmach, Noir, Hansen and Aurnou2009).

3 Heat transfer

In this section, we investigate the convective heat transfer through the fluid layer as a function of the applied control parameters. In the geostrophic regime, no consensus has been reached on the heat-transfer dependence $Nu(Ra,Pr,Ek)$ , in particular because it is challenging to achieve the extreme parameter values for $Ra$ and $Ek$ in both experiments and simulations. We summarise the results from the literature reported earlier in table 2, and indicate the method (experimental, numerical or from theory) and range of parameters considered. It must be emphasised that most of these works are outside of the geostrophic regime; the exceptions are the theories by King et al. (Reference King, Stellmach and Aurnou2012) and Julien et al. (Reference Julien, Knobloch, Rubio and Vasil2012a ), as well as the numerical simulations by the latter authors, which consider the asymptotically reduced equations for rapid rotation.

Table 2.

Proposed and observed scaling laws for the heat transfer in rotation-dominated convection found in the literature, expressed as $Nu(Ra,Ek)$ . The Prandtl number dependence has been omitted, as it is out of the scope for the present paper. We also mention the range of parameters for which the scalings are intended or observed, as well as the method by which the scalings have been obtained (‘exp.’ for experiments, ‘num.’ for numerical simulations and ‘theor.’ for scaling arguments). Finally, the boundary conditions for the numerical studies are mentioned. These are the topic of the current work. The simulations by Julien et al. (Reference Julien, Knobloch, Rubio and Vasil2012a ) are based on the asymptotic equations and employ what amounts to stress-free boundary conditions. Note that Schmitz & Tilgner (Reference Schmitz and Tilgner2009, Reference Schmitz and Tilgner2010) have defined the Ekman number without the factor $2$ , hence the different parameter values in this table compared to their papers.

The two recent experimental investigations that enter into the geostrophic regime show quite different results. Ecke & Niemela (Reference Ecke and Niemela2014) achieved $4\times 10^{9}<Ra<4\times 10^{11}$ and $2\times 10^{-7}<Ek<3\times 10^{-5}$ at $Pr=0.7$ used cryogenic helium gas as the working fluid. Their data could be described as $Nu\sim (Ra/Ra_{c})^{{\it\gamma}}$ , with ${\it\gamma}\approx 1$ using direct measurement or $1.2<{\it\gamma}<1.6$ after rescaling of the original data following existing theoretical arguments. The scaling ranges were not extensive enough to decisively discern between these scalings. On the other hand, Cheng et al. (Reference Cheng, Stellmach, Ribeiro, Grannan, King and Aurnou2015) employed water as the working fluid. In their tall, slender cell, they could achieve $1\times 10^{10}<Ra<1\times 10^{13}$ and $2\times 10^{-8}<Ek<2\times 10^{-6}$ for $3.5<Pr<6.5$ . They also reported scaling as $Nu\sim (Ra/Ra_{c})^{{\it\gamma}}$ , with a monotonically increasing ${\it\gamma}$ from $1.8$ at $Ek=10^{-3}$ to $3.6$ at $Ek=10^{-7}$ , based on a combination of experimental and DNS results.

Figure 2 shows $Nu$ as a function of $Ek$ obtained from the present simulations. By simple observation, it is clear that the boundary conditions (NS or SF) play a decisive role, even in the slope of the graph, i.e. the exponent ${\it\alpha}$ of the local scaling law $Nu\sim Ek^{{\it\alpha}}$ . Both NS and SF boundary conditions display a transition in the scaling law (indicated with arrows in the graph) at $Ek\approx 9\times 10^{-7}$ for $Ra=1\times 10^{10}$ and at $Ek\approx 3\times 10^{-7}$ for $Ra=5\times 10^{10}$ , as evidenced by the slope change. This transition is generally considered the boundary between rotation-affected and rotation-dominated convection (Ecke & Niemela Reference Ecke and Niemela2014). At Ekman numbers below the transition, we observe distinctly different scalings with $Ek$ . The SF exponents ( ${\it\alpha}=2.04$ for $Ra=1\times 10^{10}$ and ${\it\alpha}=2.21$ for $Ra=5\times 10^{10}$ ) match fairly well (especially for the lower $Ra$ case) with the theoretically predicted exponent ${\it\alpha}=2$ of Julien et al. (Reference Julien, Knobloch, Rubio and Vasil2012a ), found to be valid for simulations of the reduced equations, with boundary conditions that can be described as stress-free.

However, the NS runs reveal effective exponents ${\it\alpha}=1.07$ for $Ra=1\times 10^{10}$ and ${\it\alpha}=1.36$ for $Ra=5\times 10^{10}$ , considerably lower than the SF runs and pronouncedly lower than the correlations inferred by King et al. (Reference King, Stellmach and Aurnou2012), which predict ${\it\alpha}=4$ for NS plates. Schmitz & Tilgner (Reference Schmitz and Tilgner2009, Reference Schmitz and Tilgner2010) have reported simulations with both NS and SF boundaries; they reported good agreement with the exponent $1.5$ , in fair agreement with the current NS runs (at least at $Ra=5\times 10^{10}$ ), but somewhat low for SF. The experiments by Ecke & Niemela (Reference Ecke and Niemela2014) have provided approximate scaling exponents between $1.6$ and $2.1$ , depending on the exact plotting convention to attain data collapse.

It is worth noting that there is quite a difference in the applied $Ra$ between the various works. Typically, values of $Ra$ up to $5\times 10^{9}$ are applied in experiments and simulations (King et al. Reference King, Stellmach, Noir, Hansen and Aurnou2009; Schmitz & Tilgner Reference Schmitz and Tilgner2009, Reference Schmitz and Tilgner2010; King et al. Reference King, Stellmach and Aurnou2012; King, Stellmach & Buffett Reference King, Stellmach and Buffett2013), while only recently higher values have been attained – cf. simulations by Stellmach et al. (Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014), experiments by Ecke & Niemela (Reference Ecke and Niemela2014) and Cheng et al. (Reference Cheng, Stellmach, Ribeiro, Grannan, King and Aurnou2015), as well as the current simulations. It is plausible that what we observe is a new scaling regime opening up at such high $Ra$ , strongly affected by rapid rotation (very low $Ek\sim O(10^{-7})$ ) but still vigorously turbulent (highly supercritical, i.e. $Ra/Ra_{c}\gg 1$ ), i.e. the geostrophic regime.

Figure 2.

Heat transfer (Nusselt number $Nu$ ) as a function of the Ekman number $Ek$ : (a) SF plates; (b) NS plates. Open symbols are for $Ra=1\times 10^{10}$ ; filled symbols for $Ra=5\times 10^{10}$ . The dashed lines depict fitted power-law slopes. The arrows indicate the transition points as we inferred them from these graphs.

To further quantify the scaling laws, we show in figure 3 the compensated Nusselt number with the two scaling laws proposed by both King et al. (Reference King, Stellmach and Aurnou2012) and Julien et al. (Reference Julien, Knobloch, Rubio and Vasil2012a ). Again, we can see that the $Nu\sim Ek^{2}$ (Julien et al. Reference Julien, Knobloch, Rubio and Vasil2012a ) captures well the Ekman-number dependence of the free-slip simulations below the transition, but the no-slip simulations present a very different dependence. The $Nu\sim Ek^{4}$ scaling law (King et al. Reference King, Stellmach and Aurnou2012) can be seen to be a clear overestimate of the scaling exponent relating $Nu$ and  $Ek$ .

Figure 3.

Compensated heat transfer (Nusselt number $Nu$ ) as a function of the Ekman number $Ek$ using the proposed scaling laws: (a)  $Nu\sim Ra^{3}Ek^{4}$ (King et al. Reference King, Stellmach and Aurnou2012) and (b)  $Nu\sim Ra^{3/2}Ek^{2}$ (Julien et al. Reference Julien, Knobloch, Rubio and Vasil2012a ). Open symbols are for $Ra=1\times 10^{10}$ ; filled symbols for $Ra=5\times 10^{10}$ . The red squares represent no-slip boundary conditions, and blue circles free-slip boundary conditions.

Another striking feature of this graph is that, at the same $Ra$ , there is a range for which $Nu$ is lower for SF than for NS boundaries (Stellmach et al. Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014). Generally, NS boundaries are expected to reduce the turbulence intensity of the flow by introducing more friction than SF plates. However, the active nature of the Ekman boundaries, present for NS but absent for SF, can affect the dynamics of the entire fluid layer, enhancing the heat transfer instead of reducing it. We will revisit these results in later sections, where further differences between NS and SF runs are revealed and interpreted.

4 Boundary-layer and bulk dissipation

For RB convection, one can derive from the Navier–Stokes equations with the Boussinesq approximation exact relations for the total dissipation of turbulent kinetic energy and thermal variance within the domain (Shraiman & Siggia Reference Shraiman and Siggia1990). The energy equations are obtained by taking the inner product of $\boldsymbol{u}$ with (2.1) and multiplying (2.2) with ${\it\theta}$ , respectively. After applying the boundary conditions, the dissipation relations in dimensional form read

where ${\it\epsilon}_{u}$ is the (time- and volume-averaged) total dissipation of turbulent kinetic energy in the fluid layer and ${\it\epsilon}_{{\it\theta}}$ is the total dissipation of thermal variance in the layer. These relations do not change when rotation is added: rotation only enters in the momentum equation (2.1), where we find for the Coriolis term in the energy equation that $\boldsymbol{u}\boldsymbol{\cdot }(\boldsymbol{e}_{z}\times \boldsymbol{u})=0$ .

The Grossmann–Lohse heat-transfer theory for non-rotating convection (see Ahlers, Grossmann & Lohse (Reference Ahlers, Grossmann and Lohse2009) for an overview) is based on a division of the total dissipations into bulk and boundary-layer (BL) contributions. Several scaling regimes can be found depending on the dominance of dissipation in either bulk or BL regions, for both ${\it\epsilon}_{u}$ and ${\it\epsilon}_{{\it\theta}}$ . The theoretical arguments by Julien et al. (Reference Julien, Knobloch, Rubio and Vasil2012a ) employ such a division for ${\it\epsilon}_{{\it\theta}}$ (no division of ${\it\epsilon}_{u}$ given that in their SF case no kinetic BLs are formed) and show that the bulk limits the overall heat transfer in the geostrophic-turbulence state of the rotation-dominated geostrophic regime.

5 Flow phenomenology

A remarkable phenomenological change upon entering the geostrophic-turbulence state of the rotation-dominated geostrophic regime is the disappearance of convective Taylor columns and plumes as a result of a loss of vertical coherence (Sprague et al. Reference Sprague, Julien, Knobloch and Werne2006; Julien et al. Reference Julien, Rubio, Grooms and Knobloch2012b ; Stellmach et al. Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014). These vortical structures have been frequently reported ever since the first observation in turbulent rotating RB flow by Rossby (Reference Rossby1969). In the geostrophic-turbulence state, however, such coherent structures seem to be absent. The boundary conditions largely determine the flow phenomenology. For SF plates, large barotropic vortices can be formed under the influence of an inverse energy cascade (Favier et al. Reference Favier, Silvers and Proctor2014; Guervilly et al. Reference Guervilly, Hughes and Jones2014; Rubio et al. Reference Rubio, Julien, Knobloch and Weiss2014; Stellmach et al. Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014), eventually growing to the scale of the domain. For NS plates, such a condensate is not formed. We compare the phenomenology in two snapshots shown in figure 6, which depict the spatial distribution of vertical vorticity ${\it\omega}_{z}=\partial _{x}u_{y}-\partial _{y}u_{x}$ in a horizontal cross-section at mid-height. Figure 6(a) clearly reveals the formation of a large cyclonic vortex in the top right, while the bottom-right part and its periodic continuation on opposite sides hint at the formation of an anticyclonic vortex. Note that this flow is still slowly evolving over time (Favier et al. Reference Favier, Silvers and Proctor2014; Guervilly et al. Reference Guervilly, Hughes and Jones2014; Rubio et al. Reference Rubio, Julien, Knobloch and Weiss2014). Figure 6(b) shows no condensate vortices. Instead, a fluctuating state is found without large-scale long-lived coherent structure. Ekman pumping, present only in the case of NS boundary conditions, can thus be a source of small-scale fluctuations that prevent condensation into large-scale vortices.