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):

(2.1) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial \boldsymbol{u}}{\partial t}+(\boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}})\boldsymbol{u}+\frac{1}{Ro}\boldsymbol{e}_{z}\times \boldsymbol{u}=-\boldsymbol{{\rm\nabla}}p+\sqrt{\frac{Pr}{Ra}}{\rm\nabla}^{2}\boldsymbol{u}+{\it\theta}\boldsymbol{e}_{z}, & \displaystyle\end{eqnarray}$$

(2.2) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial {\it\theta}}{\partial t}+(\boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}}){\it\theta}=\frac{1}{\sqrt{Ra\,Pr}}{\rm\nabla}^{2}{\it\theta}, & \displaystyle\end{eqnarray}$$

with the incompressibility constraint

(2.3) $$\begin{eqnarray}\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\boldsymbol{u}=0.\end{eqnarray}$$

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

(4.1a,b ) $$\begin{eqnarray}{\it\epsilon}_{u}=\frac{{\it\nu}^{3}}{L^{4}}(Nu-1)Ra\,Pr^{-2},\quad {\it\epsilon}_{{\it\theta}}={\it\kappa}\frac{{\it\Delta}^{2}}{L^{2}}Nu,\end{eqnarray}$$

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.

4.1 Boundary-layer scales

In this section we want to compare the distribution of dissipation for both NS and SF plates. To this end, we first need to discern between BL and bulk. Several BL scales have already been introduced in the RB literature. The thickness of the thermal BLs in the non-rotating case is well described by assuming that the bulk is isothermal, and that the temperature drop is fully accommodated by the BLs. This leads to the definition ${\it\delta}_{{\it\theta},Nu}/L=1/(2Nu)$ . This relation is not appropriate for rotating RB convection, given that a mean temperature drop across the bulk is sustained (Julien et al. Reference Julien, Legg, McWilliams and Werne1996). We therefore rely on the common definition of BL thicknesses in turbulence that uses the position of the peak value of the root mean square of temperature fluctuations, denoted by ${\it\delta}_{{\it\theta}}$ . Julien et al. (Reference Julien, Rubio, Grooms and Knobloch2012b ) found this definition to be the most appropriate one. For the kinetic (velocity) BLs we use the positions of the peak root mean square of horizontal velocities, marked  ${\it\delta}_{{\it\nu}}$ .

Figure 4. Boundary-layer thicknesses and their dependence on $Ek$ : (a) SF plates, thermal BL thickness ( ${\it\delta}_{{\it\theta}}$ ); (b) NS plates, thermal BL thickness ( ${\it\delta}_{{\it\theta}}$ , red squares) and kinetic BL thickness ( ${\it\delta}_{{\it\nu}}$ , black triangles). The arrows indicate the transition points as inferred from these graphs, when the steeply decreasing ${\it\delta}_{{\it\theta}}$ begins to flatten out. A power-law fit ${\it\delta}_{{\it\nu}}/L=4.0Ek^{0.51}$ is also included (black dashed line). In both panels, open symbols are for $Ra=1\times 10^{10}$ , and filled symbols for $Ra=5\times 10^{10}$ .

A comparison of these BL scales is presented in figure 4. Starting from the kinetic BLs (black symbols), it is clear that they follow a single scaling, independent of $Ra$ , i.e. their thickness is exclusively determined by $Ek$ . A power law yields the relation ${\it\delta}_{{\it\nu}}/L=4.0Ek^{0.51}$ . Within error, the slope is consistent with the prediction ${\it\delta}_{{\it\nu}}\sim Ek^{1/2}$ for linear Ekman BLs (Greenspan Reference Greenspan1968). It is worth noting that this scaling also matches nicely with the BL scaling laws reported by Kunnen et al. (Reference Kunnen, Geurts and Clercx2010b ), even though the geometry (cylinder instead of periodic cube) and the Prandtl number ( $Pr=6.4$ instead of $Pr=1$ ) are completely different.

On the other hand, the thermal BL thicknesses do reveal some variation with $Ek$ . Before the transition, the thermal BL thickness steeply decreases when $Ek$ is increased, in contrast to the kinetic BLs. In that $Ek$ range, the local scaling laws relating $Nu$ and $Ek$ are steeper for SF than for NS. King et al. (Reference King, Stellmach, Noir, Hansen and Aurnou2009) and King et al. (Reference King, Stellmach and Aurnou2012) proposed that the transition to the rotation-dominated heat-flux scaling is described by the crossing of the kinetic and thermal BL thicknesses, which happens around $Ek=7\times 10^{-7}$ for NS and $Ra=5\times 10^{10}$ in this case. However, the slope change in the scaling laws (figure 2) is found at lower values of $Ek$ . For $Ra=1\times 10^{10}$ a similar mismatch is observed.

A more natural definition of the transition, based on the BL scales as plotted in figure 4, would be the evident slope change of ${\it\delta}_{{\it\theta}}$ at $Ek=8\times 10^{-7}$ for $Ra=1\times 10^{10}$ and at $Ek=3\times 10^{-7}$ for $Ra=5\times 10^{10}$ , which remarkably occurs around the same $Ek$ value for both NS and SF. This transition value matches better with the slope change in the heat-transfer statistics, even if they do not exactly coincide. The transition to the geostrophic regime thus appears to be gradual; different statistics display a change in behaviour at different values of  $Ek$ .

Figure 5. Distribution of dissipations ${\it\epsilon}_{{\it\theta}}$ and ${\it\epsilon}_{u}$ between bulk and BL: (a) SF plates, ${\it\epsilon}_{{\it\theta}}$ ; (b) NS plates, ${\it\epsilon}_{{\it\theta}}$ (red squares) and ${\it\epsilon}_{u}$ (black triangles). The fraction of total dissipation located in the BL region is displayed. Open symbols are for $Ra=1\times 10^{10}$ ; filled symbols for $Ra=5\times 10^{10}$ . The arrows indicate the transition points as we inferred them from these graphs.

4.2 Distribution of dissipation

Using the BL scales of § 4.1, we can now assess how the total dissipation is distributed between BL and bulk regions. This is shown in figure 5, which confirms the picture that, under rapid rotation, the dissipation is mostly concentrated in the bulk. This is the case even more for ${\it\epsilon}_{u}$ than for ${\it\epsilon}_{{\it\theta}}$ . However, we also note that the fraction of ${\it\epsilon}_{{\it\theta}}$ in the BLs appears to start growing when $Ek$ is reduced below ${\sim}\,5\times 10^{-7}$ ( ${\sim}\,2\times 10^{-7}$ ) for $Ra=1\times 10^{10}$ ( $5\times 10^{10}$ ), with an earlier growth appearing for SF than for NS. Looking back to figure 4, it is obvious that the thermal BLs are expanding when $Ek$ is reduced. The larger part of the volume inside the thermal BLs along with a persistent input of thermal fluctuations from the Ekman BLs enhances the BL fraction of ${\it\epsilon}_{{\it\theta}}$ at the lowest considered $Ek$ under NS conditions. For SF plates the thermal BLs are growing even more as $Ek$ is reduced; the increased volume of the BLs appears to be enough for a higher fraction of ${\it\epsilon}_{{\it\theta}}$ there.

Furthermore, the contribution of the kinetic BLs to the total ${\it\epsilon}_{u}$ is remarkable: for $Ra=5\times 10^{10}$ , between $Ek=1.3\times 10^{-7}$ and $Ek=5\times 10^{-7}$ the BL thickness changes by a factor $2$ , yet the fraction of ${\it\epsilon}_{u}$ in the BL remains roughly constant. This confirms that the Ekman BLs, first thought to become passive at low enough $Ek$ (Niiler & Bisshopp Reference Niiler and Bisshopp1965; Heard & Veronis Reference Heard and Veronis1971; Julien & Knobloch Reference Julien and Knobloch1998), are still significantly affecting the flow dynamics (Stellmach et al. Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014).

Figure 6. Snapshot at mid-height ( $z=0.5$ ) of the vertical vorticity from runs at $Ra=5\times 10^{10}$ and $Ek=1.34\times 10^{-7}$ : (a) SF plates; (b) NS plates. Red is positive (cyclonic) vorticity, while blue is negative (anticyclonic) vorticity. Both plots have the same colour scale.

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.

5.1 Relation with spectral energy transfer

In previous works (Favier et al. Reference Favier, Silvers and Proctor2014; Rubio et al. Reference Rubio, Julien, Knobloch and Weiss2014), the energy transfer as a function of wavenumber has been considered to indicate the presence of an inverse energy cascade. Given that figure 6 reveals such a significant difference in flow phenomenology between SF and NS, we anticipate that the spectral energy transfer must also be considerably different. Following Favier et al. (Reference Favier, Silvers and Proctor2014), we define the spectral energy equation as

(5.1) $$\begin{eqnarray}\frac{\text{d}E(K)}{\text{d}t}=\mathop{\sum }_{Q}T(Q,K)-D(K)+F(K),\end{eqnarray}$$

which gives the temporal evolution of the energy $E(K)$ as a function of horizontal wavenumber $K$ in terms of the energy transfer $T(Q,K)$ from wavenumber shell $Q$ to shell $K$ , the dissipation $D(K)$ and the buoyant forcing $F(K)$ . We are particularly interested in the transfer term $T(Q,K)$ , which we evaluate and average vertically over the entire computational domain minus the BLs. We note that this quantity is antisymmetric by definition, i.e. $T(K_{1},K_{2})=-T(K_{2},K_{1})$ and $T(K_{1},K_{1})=0$ for any two wavenumbers $K_{1}$ and $K_{2}$ .

Figure 7. Spectral energy transfer $T(Q,K)$ between wavenumber shells $Q$ and $K$ : (a) for SF and (b) for NS. Positive values indicate that energy is taken from shell $Q$ and transferred to shell $K$ ; negative values imply that $Q$ receives energy from $K$ . Both plots have the same colour scale.

The spectral transfer at $Ra=5\times 10^{10}$ and $Ek=1.34\times 10^{-7}$ is depicted in figure 7, for both SF and NS conditions. We note that we are only showing the first 30 modes, of a total of 256 in the simulation, so the dispersive errors of finite-difference schemes are negligible in this range. The SF picture compares favourably to the earlier studies by Favier et al. (Reference Favier, Silvers and Proctor2014) and Rubio et al. (Reference Rubio, Julien, Knobloch and Weiss2014): a strong localised exchange of energy between neighbouring wavenumbers across the diagonal, i.e. a certain mode $P$ interacts predominantly with its neighbouring modes $P-1$ and $P+1$ . The sign of the transfer reveals that this part of the cascade is direct: energy is transferred to larger wavenumbers. At low $K\lesssim 7$ and higher $Q\gtrsim 7$ , a range of inverse transfer is found where small- $K$ modes receive energy from modes $Q>K$ . The input of energy into the large-scale vortex of figure 6(a) can be recognised as an interaction of the smallest wavenumbers receiving energy from a broad range of higher-order wavenumbers. Such a non-local inverse cascade has been observed before in rotating convection by Favier et al. (Reference Favier, Silvers and Proctor2014), who could display a clearer effect due to a considerably larger domain. It is also similar to the inverse energy cascade in rotating turbulence (Mininni, Alexakis & Pouquet Reference Mininni, Alexakis and Pouquet2009). The non-locality indicates that the largest flow scale directly receives energy from the smaller scales without passing through intermediate scales. In contrast, in the direct cascade, energy is transferred from large to small while passing through all intermediate scales.

For the NS case (figure 7 b), the picture changes. In particular, the ‘staircase’ of energy transfer along the diagonal, which is prominent for SF, is not as strongly present for the larger wavenumbers. Instead, the interactions between modes are less localised, meaning that the interactions are spread more and mode combinations farther from the diagonal are transferring energy. A direct cascade is formed along the diagonal for $K,Q\gtrsim 5$ . Curiously, some signs of an inverse cascade remain: (i) for $K\leqslant 5$ , and (ii) transfers from modes around $Q\approx 5$ to $K\approx 10$ , which are non-localised. The absence of a large-scale structure can be explained by the fact that the lowest wavenumber is receiving considerably less energy and from a narrower range of scales than in SF.

5.2 Relation with mean temperature gradient

The characteristic flow phenomenology of rotating RB convection has been related to the occurrence (and strength) of the persistent mean temperature gradient across the bulk (Julien et al. Reference Julien, Legg, McWilliams and Werne1996), unlike the statistically isothermal bulk of non-rotating RB. The origin of this temperature gradient has been proposed to be increased by lateral mixing, induced by interactions of like-signed vortical plumes. As in the geostrophic-turbulence state the flow phenomenology is altered (cf. figure 6), we can expect that this also affects the strength of the mean temperature gradient. Figure 8 shows the mean temperature gradient as a function of $Ek$ . In the rotation-affected RB, lateral mixing being stronger than vertical mixing leads to mean temperature gradients as large as $-0.5$ for NS and $-0.4$ for SF. However, upon entering the geostrophic-turbulence state by further reducing $Ek$ , the magnitude of the gradient is gradually diminished. We thus expect that in the geostrophic-turbulence state the mixing can become slightly more 3D again.

This behaviour is consistent with previous simulations of the asymptotic equations (Julien et al. Reference Julien, Rubio, Grooms and Knobloch2012b ), which predict that the geostrophic-turbulence state indeed still has a mean temperature gradient, but less steep than when coherent vortical plumes are present. Finally, even though the behaviour of the mean temperature gradient is qualitatively similar for SF and NS plates, the location of the minimum temperature gradient, which could be taken as an additional indicator for the transition, certainly does not coincide between the two cases.

Figure 8. Mean temperature gradient at midplane: (a) for SF and (b) for NS. Open symbols are for $Ra=1\times 10^{10}$ ; filled symbols for $Ra=5\times 10^{10}$ . The arrows indicate the transition points as inferred from these graphs.