Rotating turbulent thermal convection at very large Rayleigh numbers

Abstract We report on turbulent thermal convection experiments in a rotating cylinder with a diameter ($D$) to height ($H$) aspect ratio of $\varGamma =D/H=0.5$. Using nitrogen and pressurised sulphur hexafluoride we cover Rayleigh numbers (Ra) from $8\times 10^{9}$ to $8\times 10^{14}$ at Prandtl numbers $0.72\lesssim {\textit {Pr}}\lesssim 0.94$. For these Ra we measure the global vertical heat flux (i.e. the Nusselt number – Nu), as well as statistical quantities of local temperature measurements, as a function of the rotation rate, i.e. the inverse Rossby number – 1/Ro. In contrast to measurements in fluids with a higher Pr we do not find a heat transport enhancement, but only a decrease of Nu with increasing 1/Ro. When normalised with Nu(0) for the non-rotating system, data for all different Ra collapse and, for sufficiently large 1/Ro, follow a power law ${\textit {Nu}}/{\textit {Nu}}_0\propto (1/{\textit {Ro}})^{-0.43}$. Furthermore, we find that both the heat transport and the temperature field qualitatively change at rotation rates $1/{\textit {Ro}}^*_1=0.8$ and $1/{\textit {Ro}}^*_2=4$. We interpret the first transition at $1/{\textit {Ro}}^*_1$ as change from a large-scale circulation roll to the recently discovered boundary zonal flow (BZF). The second transition at rotation rate $1/{\textit {Ro}}^*_2$ is not associated with a change of the flow morphology, but is rather the rotation rate for which the BZF is at its maximum. For faster rotation the vertical transport of warm and cold fluid near the sidewall within the BZF decreases and hence so does Nu.


Introduction
Thermal convection, the flow driven by a thermal gradient, is one of the most important heat transport mechanisms in many natural and industrial systems and has been studied for many decades in the well-defined Rayleigh-Bénard convection (RBC) system (Bénard 1900;Rayleigh 1916;Kadanoff 2001;Ahlers, Grossmann & Lohse 2009b). There, a horizontal fluid layer of height H is confined by a warm plate at its bottom and a cold plate at its top. Under sufficiently strong thermal driving, the flow in RBC is highly turbulent and as such difficult to predict and calculate analytically. Nevertheless, good progress has been made in recent years, in particular in understanding the functional dependency between thermal driving and global-averaged heat transport (see e.g. Chavanne et al. 1997;Niemela et al. 2000;Grossmann & Lohse 2000, 2001He et al. 2012;Stevens et al. 2013;Whitehead & Wittenberg 2014;Sondak, Smith & Waleffe 2015;Shishkina et al. 2017;Tobasco & Doering 2017).
The most turbulent convection flows in nature occur in planets, stars and other large scale astrophysical systems. In these systems, the flow is strongly influenced by the rotation of their celestial body. For example, the fluid motion in the Earth's outer core (Buffett 2000) is aligned with the Earth's rotation axis, creating a dipolar magnetic field in the same direction. Another example is the Earth's atmosphere, where pressure equilibration in depression areas is suppressed by Coriolis forces, causing hurricanes to develop that can exists for weeks and travel over thousands of kilometres (Holton 2004). It is therefore of crucial importance to study thermal convection in a rotating system.
For the sake of simplicity, one usually studies the RBC system under Oberbeck-Boussinesq (OB) conditions, i.e. the temperature difference between the bottom and top are small enough so that fluid properties are constant everywhere in the fluid (Oberbeck 1879;Boussinesq 1903;Spiegel & Veronis 1960). In the theoretical description of this problem, only the density is assumed to depend linearly on the temperature, whereby the isobaric thermal expansion α is the relevant coefficient. In this case, the non-rotating system is governed by two dimensionless control parameters. These are the Rayleigh number Ra = gα TH 3 νκ , (1.1) and the Prandtl number (1.2) Here, T denotes the temperature difference between the bottom and top plates, g the gravitational acceleration, ν and κ are the kinematic viscosity and the thermal diffusivity, respectively. In simulations and experiments, one is often interested in the heat flux from the bottom to the top, which is expressed as the dimensionless Nusselt number with q being the dimensional heat flux and q cond being the heat flux that would occur purely due to conduction. In the most general case, calculation of q cond demands numerical integration of the heat conductivity λ(T) over the entire cell (Shishkina, Weiss & Bodenschatz 2016). For constant fluid properties (under OB conditions) this simplifies to q cond = λ T/H. Considerable effort has been devoted to understanding how Nu depends on Ra and Pr (Ahlers et al. 2009b), and multiple theoretical models have been proposed that aim to predict these dependencies (see for instance Malkus 1954;Kraichnan 1962;Castaing et al. 1989;Grossmann & Lohse 2000;Stevens et al. 2013;Tobasco & Doering 2017). For a rotating system the rotation rate Ω is an additional control parameter. It is incorporated in the dimensionless convective Rossby number (1.4) and the Ekman number (1.5) We note that, in the literature, other definitions for Ek are often used which differ by a factor of two. However, since one is mostly interested in scaling behaviours, constant coefficients do not influence the resulting power laws. In this paper we will usually consider the inverse Rossby number 1/Ro since it is a dimensionless rotation rate. One effect of rotation is that it increases the critical Rayleigh number Ra c above which a fluid at rest starts convecting. For a laterally infinite fluid layer, Ra c scales as Ek −4/3 (Chandrasekhar 1981). This is the reason why in rotating RBC often a reduced Rayleigh number is used as a control parameter, defined as Ra = RaEk 4/3 .
(1.6) A major question for rotating RBC is to understand how rotation affects the heat transport, i.e. how Nu depends on Ra, Pr and Ro (or Ek). For sufficiently large Ra when the flow is turbulent, one can very roughly distinguish two different regimes. For small rotation rates, the flow in the bulk is still dominated by buoyancy. Then, heat is predominantly transported by thermal plumes that detach from the top and bottom and that self-organise in a large-scale convection (LSC) role. In this regime, the top and bottom boundary layers are the major bottlenecks for the heat transport. Depending on Ra and Pr in this regime, the heat transport can both increase and decrease with increasing rotation rate (increasing 1/Ro) (see e.g. Weiss & Ahlers 2011a;Weiss, Wei & Ahlers 2016). The increase (most significant for large Pr) is attributed to Ekman pumping, which occurs in the vortices that form close to the top and bottom boundaries when plumes emerge from the boundary layers and get twisted by the Coriolis forces. At much larger rotation rates, Nu decreases monotonically with increasing rotation rates (King et al. 2009). This is caused by the suppression of vertical velocity gradients i.e. the Taylor-Proudmann effect (Taylor 1923) and thus a reduced advective transport inside the bulk region. At constant Ra a critical rotation rate is reached when the critical Rayleigh number for the onset of convection Ra c exceeds Ra. Then convection is completely suppressed (Chandrasekhar 1981) and heat transport is solely by conduction, i.e. Nu = 1. For fast rotation rates, geostrophic balance is reached, meaning that pressure gradients are balanced by Coriolis forces (Holton 2004). Numerical simulations have found four different flow morphologies in this regime (Julien et al. 2012b;Stellmach et al. 2014;Plumley et al. 2016). These are (i) cellular convection, (ii) convective Taylor columns, (iii) plumes and (iv) geostrophic turbulence. While cellular convection occurs at relatively small Ra, geostrophic turbulence occurs at large Ra.
In particular in the geo-and astrophysical community, modelling the heat transport under geostrophic conditions via scaling laws of the form Nu ∝ Ra γ Pr β Ek α is of great importance. Finding good exponents γ and α, however, is challenging. While various models that aim to predict γ and α have been proposed in the past (see e.g. Portegies et al. 2008;King et al. 2009;King, Stellmach & Aurnou 2012;Julien et al. 2012a), they cannot easily be verified or falsified due to the rather small range of geostrophic turbulence that can be achieved in simulations and laboratory experiments. To reach geostrophic turbulence one needs strong thermal driving (large Ra) and fast rotation so that Coriolis forces dominate the flow over buoyancy (small Ro). While large Ra can be achieved in experiments easily, the rotation rates are limited by unwanted centrifugal forces that occur in large cylindrical cells (for an overview of experiments and a related discussion see Cheng et al. 2018). Simulations are usually restricted by the strong turbulence that demand high spatial and temporal resolution over many scales, which becomes computationally very expensive at larger Ra. As a result, reliable direct numerical simulation (DNS) data are available close to convection onset (Ra c ), but do not cover sufficiently large Ra ranges in the turbulent regime, in particular when Ekman layers occur at the top and bottom due to no-slip boundaries (see e.g. Stellmach et al. 2014;Plumley et al. 2016).
So far it is also not clear where transitions between different regimes are expected to occur. In particular, at which rotation rates do Coriolis forces dominate over buoyancy? Different propositions have been made for relevant mechanisms that consider either properties of the bulk flow (e.g. Gilman 1977) or the dynamics in the boundary layers (King et al. 2009;Julien et al. 2012b;Cheng et al. 2018;Long et al. 2020). Below, in § 3, where we discuss Nusselt number measurements, we will also compare our measurements with proposed scalings in the geostrophic regime and regime transitions.
While most theoretical models assume laterally infinite or periodic systems, in experiments and many simulations on rotating RBC an upright cylinder is considered of finite aspect ratio Γ = D/H between its diameter D and its height H and with adiabatic no-slip sidewalls. It is known that the onset of convection (i.e. when Ra exceeds Ra c ) in rotating RBC in finite Γ -containers occurs first at the lateral sidewall in the form of travelling waves, so called wall modes (Buell & Catton 1983a,b;Goldstein et al. 1993;Zhong, Ecke & Steinberg 1993;Goldstein et al. 1994;Bajaj, Ahlers & Pesch 2002). These wall modes set in at significantly smaller Ra than the convection in the bulk or in an infinite system. For sufficiently large rotation rates the onset of wall modes was calculated to occur at Ra w ∝ Ek −1 (Goldstein et al. 1993(Goldstein et al. , 1994Herrmann & Busse 1993;Kuo & Cross 1993;Zhang & Liao 2009;Favier & Knobloch 2020).
The lateral sidewall not only plays an important role in convection close to the onset, but also in the turbulent state. There, Stewartson layers form, in which fluid is pumped from the top and bottom towards the midheight of the cell and from there into the bulk (Stewartson 1957(Stewartson , 1966Kunnen et al. 2011). In addition, another flow pattern has been found recently, which develops very close to the lateral sidewall at sufficiently fast rotation rates, termed the boundary zonal flow (BZF) (de Wit et al. 2020;Zhang et al. 2020). A schematic of the BZF is shown in figure 4(a). This flow structure is characterised by a narrow region close to the sidewall with a positive average azimuthal velocity (prograde) that surrounds a core region with negative average azimuthal velocity (retrograde). Inside this narrow region warm fluid rises on one side and cold fluid sinks on the other. The temperature field is periodic with wavenumber Γ (Zhang, Ecke & Shishkina 2021) in the azimuthal direction, and drifts in the retrograde direction. Within the BZF warm and cold fluid is pumped in large areas from the bottom to the top and vice versa, causing the heat flux there to be significantly larger than in the bulk. In fact, in Γ = 1/2 cells, more than 60 % of the heat is transported inside the BZF. The radial width of the BZF decreases with increasing rotation rates. In a recent study by Favier & Knobloch (2020) it was suggested that the BZF is related to the wall modes close to convection onset, and originates out of them via an Eckhaus-like instability.
In this paper we present a comprehensive analysis of global heat flux measurements as well as of measurements of the temperature at different vertical and radial positions for rotating RBC. We measure statistical properties of the temperature field, such as probability density functions (PDFs), time-averaged values, standard deviations or skewness as a function of the radial and vertical coordinates, Ra, and 1/Ro. In particular, we see qualitative differences between these quantities when measured inside the BZF and when measured in the radial bulk.
The paper is structured as follows. In the next section we explain the details of our experimental set-up. Section 3 reports on measurements of the Nusselt number for different Ra and 1/Ro. Section 4 discusses the structure and dynamics of the BZF based on temperature measurements inside the sidewall. In § 5 we report on vertical profiles of the time-averaged temperature, followed by a detailed analysis of the temperature fluctuations in § 6. This analysis considers the entire PDF ( § 6.1) such as its standard deviation or its skewness, since all these quantities help to determine the onset of the BZF. We close the paper with a conclusion in § 7.

Experimental set-up
The experiments were conducted at the high pressure convection facility (HPCF) at the Max-Planck-Institute for Dynamics and Self-Organization in Göttingen. The convection cell has been described in detail before in previous publications (see e.g. Ahlers, Funfschilling & Bodenschatz 2009a;Ahlers et al. 2012b). Here, we use the version HPCF-II, which consists of a cylinder of height H = 2.24 m and a diameter D = 1.12 m, resulting in an aspect ratio of Γ = D/H = 0.50. The sidewall of the cylinder was made of a 9.5 mm thick acrylic. The top and bottom plates were made of high-purity copper with a heat conductivity of 394 W m −1 K −1 . The bottom plate consisted of two such copper plates (35 mm and 25 mm thick), separated by a 5 mm thin Lexan plate in between. The temperature drop across the Lexan plate was used to determine the vertical heat flux.
The bottom plate was heated with an ohmic heater at its bottom side, whereas the top plate was cooled using a circulated water bath which was temperature controlled with a precision of 0.02 K. Thermal shields underneath the bottom plate and around the sidewall ensured minimal heat loss through the bottom and the sides. Additional micro-shields close to the boundary layers at top and bottom aimed to reduce the heat transport through the sidewall close to the vertical boundaries (Ahlers 2000;Stevens, Lohse & Verzicco 2014).
The cell was mounted on a rotating table that could sustain an axial load up to 2800 kg. The table was driven by a direct drive motor (Siemens 1FW6150 SIMOTICS T Torque-motor), which could deliver a torque of up to 1000 Nm even at very low rotation rates (down to 1 rad min −1 ), ensuring very smooth rotation even at such low speeds. Water feed-through and electrical slip rings were used to bring the coolant as well as electrical connections from the laboratory into the rotating frame. The maximal rotation rate in our experiment was 2 rad s −1 to keep the influence of the centrifugal forces small. The maximal Froude number Fr = Ω 2 D/g (with D = HΓ ) in our experiment did not exceed the value of Fr c = 0.5 below which the influence of centrifugal forces is expected to be small compared to Coriolis forces and where their influence on the flow field and the heat transport are expected to be not significant (see Horn & Aurnou 2018. The rotating table and cell were installed into the U-Boot of Göttingen, a 4 m long and up to 4.3 m high pressure vessel (see Ahlers et al. 2009a), which can be filled with nitrogen or sulphur hexafluoride (SF 6 ). The pressure inside the U-Boot could be increased to 19 bar.
Next to the vertical heat flux, we also measured temperatures at various locations inside the sidewall and the fluid. For this, three rows of thermistors were placed into blind holes inside the sidewall at heights H/4, H/2 and 3H/4. Each row consists of eight thermistors distributed azimuthally at equal distance. An additional 62 thermistors were arranged in vertical columns close to the sidewall. An overview of the radial and azimuthal locations of these thermistors is given in figure 1.
During a typical experiment, the rotation rate Ω, as well as T b and T t were held constant while thermistors were read out every 5 s. Experiments lasted at least 12 h, while the first 2 h were discarded during which T b and T t settled to their desired values. We conducted various sets of experiments, with constant Ra and different 1/Ro by changing the rotation rate Ω. Experiments with SF 6 were also conducted for different pressures, causing a change of Pr from 0.78 (at 1 bar) to 0.97 (at 19 bar). While the U-Boot temperature T U was found to have only little impact on the results, we, however, kept it close to the mean temperature T m = (T b + T t )/2. A list of measurements is given in table 1.
Often in this paper, whenever we are investigating how certain quantities change under the influence of rotation at constant Ra, we present data from run E2e with Ra = 4.9 × 10 13 . There is nothing special about this run, except that its Rayleigh number is somewhere in the middle of all Rayleigh numbers.

Nusselt number measurements
We have conducted heat flux measurements for a rather large range of Ra and 1/Ro. The measurements without rotation are in good agreement with previous measurements in a similar set-up (Ahlers et al. 2012b), and follow a power law relation of Nu(1/Ro = 0) ∝ Ra 0.314 (see supplemental material available at https://doi.org/10.1017/jfm.2020.1149). Since we are predominantly interested in the effect of rotation on the heat transport, we show in figure 2 the normalised Nusselt number Nu/Nu 0 = Nu(1/Ro)/Nu(0) as a function of the dimensionless rotation rate, i.e. the inverse Rossby number (1/Ro). As clearly seen in figure 2(a), the data for different Ra collapse rather well over the entire 1/Ro range. Looking at the data, one can distinguish three different 1/Ro regimes, with different monotonic behaviours of Nu/Nu 0 in each of them. (i) For 1/Ro (1/Ro * 1 = 0.8), Nu/Nu 0 is rather constant and close to unity, i.e. rotation has barely any effect on the vertical heat transport. (ii) For 0.8 1/Ro 4.0, Nu/Nu 0 decreases with increasing 1/Ro. (iii) In the range (1/Ro * 2 = 4.0) 1/Ro, Nu/Nu 0 decreases strongly with 1/Ro revealing a fitted exponent of the power law Nu/Nu 0 ∝ 1/Ro β of β = −0.43 ± 0.02. It seems clear that in regime (i) the flow is dominated by the action of buoyancy, while in regime (iii) rotation and the occurring Coriolis forces cause a suppression of vertical fluid motion and hence of vertical convective heat transport (i.e. the Taylor-Proudman effect). While model predictions based on scaling estimates predict power law relationships between Nu/Nu 0 and 1/Ro we note that a logarithmic function appears to fit the data in regime (iii) slightly better (green dash-dotted line in figure 2(a), Nu/Nu 0 ∝ log((0.015 ± 0.002)/Ro)). Furthermore, we note that in regime (iii) the data points for larger Ra are slightly above those for smaller Ra. While the difference is within the margin of uncertainty of the Nu 0 measurements of approximately 1 %-2 % or so, we cannot exclude a Ra-dependency. A possible explanation could be the increasing influence of centrifugal forces. Data points at larger Ra have also larger Fr at the same 1/Ro. However, one then would also observe a slight temperature decrease at the sidewall at midheight, which we did not observe (figure 7).
While we do not see any sign of heat transport enhancement in our measurements, we nevertheless want to estimate at which 1/Ro such an enhancement is expected to occur from previous measurements and how this compares with our transitional rotation rates 1/Ro * 1 and 1/Ro * 2 . The critical inverse Rossby number for the onset of heat transport enhancement (1/Ro c ) (see e.g. Zhong et al. 2009;Weiss et al. 2016) depends on the aspect ratio of the convection cylinder Γ (Weiss et al. 2010). For example, for Γ = 0.5, 1/Ro c = 0.8 was measured for water as the working fluid (Pr = 4.38). Indeed, this value is quite close to our first transition from regime (i) to the intermediate regime (ii). However, we learnt from previous measurements in cylinders of Γ = 1 that 1/Ro c increases with decreasing Pr. In Weiss et al. (2016), an empirical power law fit of 1/Ro c = K 1 Pr α was determined with Nu/Nu 0 1/Ro 7.7 × 10 09 2.1 × 10 10 2.0 × 10 11 8.4 × 10 12 4.9 × 10 13 3.9 × 10 14 8.0 × 10 14 1.71(1/Ro) −0.43 -0.33log(0.015/Ro) 1/Ro 7.7 × 10 9 (THIS-EXP) 6.2 × 10 9 (EN14) 1 × 10 9 (HS15) 1 × 10 9 (ZS20) (a) ( b) Figure 2. (a) Semi-logarithmic plot of Nu/Nu 0 as a function of 1/Ro. Different symbols mark different Ra as given in the legend. The dashed vertical lines mark the approximate points of transition at 1/Ro * 1 = 0.8 and 1/Ro * 2 = 4.0. The dashed red and dot-dashed green lines are power law and logarithmic fits to the data in the range 4 < 1/Ro. (b) Shows our data for the lowest Ra (solid blue circles, error bars mark the estimated 2 % uncertainty) in comparison with experimental and numerical results by others. Filled red squares are results measurements in helium (Pr = 0.7) at Ra = 6.2 × 10 9 (Ecke & Niemela 2014). Open down-pointing black triangles are results from DNS for Ra = 10 9 , Pr = 0.8 (Horn & Shishkina 2015). Orange stars mark most recent results for Ra = 10 9 (Zhang et al. 2021). Dotted vertical lines display the estimated critical Rossby number 1/Ro c,est (from Weiss et al. 2010) and two transition points 1/Ro t and 1/Ro T from Ecke & Niemela (2014).
K 1 = 0.75 and α = −0.41. While it is unclear how α and K 1 depend on Γ , we know that 1/Ro c = a/Γ (1 + b/Γ ) for Pr = 4.38 with a = 0.381 and b = 0.061 (Weiss et al. 2010). Therefore, one could combine the 1/Ro c -dependency of Pr and Γ to 1 This rather rough estimate results for our experiments with Γ = 0.5 and Pr = 0.8 in 1/Ro c,est = 1.74. This value is somewhere in between 1/Ro * 1 and 1/Ro * 2 in the intermediate regime (ii). Neither in the heat transport measurements (figure 2) nor in the local temperature measurements (e.g. figure 8) can any significant changes be observed at this 1/Ro. This suggests that the mechanisms leading to a heat transport enhancement for larger Pr, like Ekman pumping, are absent for Pr < 1 and are not just counteracted by suppression of vertical velocity.
In figure 2(b), we compare our results with other studies that have been conducted in a similar Pr-and Ra-range. As shown, the general trend of our data (blue bullets) agrees well with other results. Quantitative agreement also exists for small and large 1/Ro with the results from cryogenic helium experiments by Ecke & Niemela (2014) (red solid squares) and with DNS at Ra = 10 9 (Horn & Shishkina 2015). For intermediate 1/Ro, however, our heat flux measurements are slightly smaller than those of Ecke & Niemela (2014) and Horn & Shishkina (2015). In fact, Ecke & Niemela (2014) and Horn & Shishkina (2015) observe a very small heat transport enhancement. This might, however, also be within the margin of uncertainty. The enhancement for the black triangles (Horn & Shishkina 2015) is for example only seen in a single data point. Furthermore, no heat transport enhancement is observed in the most recent simulations by Zhang et al. (2021), where finer computational grids were used and the DNS data were collected for much longer than in Horn & Shishkina (2015).
A major goal in research on rotating turbulent convection is to find simple scaling relationships in the regime of geostrophic turbulence, i.e. where Coriolis forces balance horizontal pressure gradients but the flow is still highly turbulent. It is generally difficult to achieve this regime in experiments as both the Rayleigh-and the inverse Rossby number need to be very large.
However, we also note that there are no obvious reasons for the collapse of different Ra-measurements, and various other scaling laws have been proposed in the past for the geostrophic turbulent regime. In the following, we aim to test some of them. In the limit of strong rotation (Ek → 0) the critical Rayleigh number Ra c for the onset of convection (for an infinitely large container) increases according to Ra c ∝ Ek −4/3 (Chandrasekhar 1981). Thus, it is useful to consider RaEk 4/3 as a control parameter. Geostrophic turbulence is expected for RaEk 4/3 1 and it is argued by Julien et al. (2012a) that RaEk 4/3 is the only relevant parameter and therefore the relation Nu ∝ (RaEk 4/3 ) α was proposed.
In figure 3(a), we plot Nu as a function of RaEk 4/3 . The Nu-data in this representation do not collapse. While they are not expected to collapse for large RaEk 4/3 , the data all bend down if RaEk 4/3 decreases and they might collapse for smaller values that were not achievable in our experimental set-up without entering the regime where centrifugal forces start to play a significant role. We plot as solid and dashed lines results from DNS of the full equation and of the asymptotic quasi-geostrophic model by Plumley et al. (2016). While our data are in a completely different range of RaEk 4/3 we see that they do not contradict the data by Plumley et al. (2016). However, it also hints that our measurements might still be far away from the geostrophic turbulent regime.
A good collapse of the Nu/Nu 0 -data was found by Ecke & Niemela (2014) when plotted against Ra/Ek 7/4 . This control parameter was suggested by King et al. (2009King et al. ( , 2012, who proposed that the transition between the buoyancy-dominated and the rotation-dominated regime occurs when the thermal boundary layer and the Ekman boundary layer are of equal thickness. We therefore plot our data points in such a way in figure 3(b). We see in this plot that data taken with N 2 as a working fluid collapse very well for different Ra and data for SF 6 also collapse decently well for larger Ra. However, the N 2 data do not collapse with the SF 6 data. The reason for this mismatch is unclear as both Ra and Pr are different for the two datasets, but also both change within each of the data sets. In particular, Pr ≈ 0.72 for N 2 while 0.8 < Pr < 0.97 for SF 6 , and thus it is unlikely that the different Pr values are causing this mismatch. One of course gets a better collapse when the exponent of Ek is increased to 2, which results in figure 2 since RaEk 2 ∝ Ro 2 . We also want to note that RaEk 7/4 assumes a scaling for the non-rotating case of Nu ∝ Ra 2/7 . In our case the Ra-exponent is larger (Nu ∝ Ra 0.31 ), which would in fact result in a control parameter with an even smaller Ek-exponent (RaEk 1.6 ) and consequently in a worse collapse. Furthermore, the same RaEk 8/5 should collapse the transition from the rotation-dominated regime to the buoyancy-dominated regime according to Julien et al. (2012a), which is when the geostrophic balance breaks down in the thermal boundary layer. We note that in simulations of convection in a spherical geometry using RaEk 8/5 has in fact collapsed the data for different Ra and Ek quite well, however, at significantly smaller Ra (Long et al. 2020).
As just mentioned, RaEk 4/3 was proposed as a control parameter, because the Ra c for the onset of convection under rotation increases as Ek −4/3 . This, however, is only true for a container with an infinite aspect ratio Γ . For any cylinders with finite Γ , convection occurs at much smaller Ra close to the sidewall in the form of periodic wall modes, while no convection occurs in the radial bulk. However, these wall modes already transport heat and one might argue that one should rather compare Ra with its critical value for the onset of wall modes (Ra w ). (In fact this suggestion was made by one of the anonymous referees of this paper, to whom we are grateful.) It has been shown (Herrmann & Busse 1993;Kuo & Cross 1993) that for asymptotically small Ek the critical value scales as Ra w ∝ Ek −1 , i.e. it increases slower with increasing rotation rate than for the laterally infinite case. We therefore plot in figure 3(c) Nu as a function of RaEk. The data in this representation qualitatively look similar to figure 3(a), but the data for sufficiently large rotation rates (filled symbols mark data with 1/Ro ≥ 4) seem to follow a power law better. The best fit (green line in figure 3c) gives Nu ∝ (RaEk) 0.62 . This is quite remarkable and suggests that (i) the influence of the wall is strong, even in the turbulent regime and that mechanisms leading to wall modes at onset prevail also in the highly turbulent regime, perhaps as BZF (see also Favier & Knobloch 2020). It also suggests that (ii) these near wall flow structures (the BZF) play a crucial role for the global heat transport, an observation that has been confirmed by DNS already (Zhang et al. 2020(Zhang et al. , 2021. While the proposed scalings discussed here were based to some extent on physical scaling arguments, it certainly would be interesting to determine the best exponents from a two-dimensional fit of the form Nu ∝ Ra a Ek b to the data. We have done this by using only data with sufficiently large rotation rates (1/Ro ≥ 4). The best fit results in a = 0.54 and b = 0.46. We note that these values are subject to a rather large uncertainty, as there is no very localised narrow minimum in the residuals. Instead, there is a long extended minima valley, such that combinations of a and b with b = 0.85a provide equally good fits (see figure 2 in the supplemental material).
We plot compensated values of Nu in figure 3(d) so that data following this relation lie on a horizontal line. Although the data (solid bullets) are rather close to each other on the vertical scale over two decades on the x-axis, they do show clear deviations for each individual Ra-dataset, which suggests that some other exponents would have been fitted if additional measurements at even larger rotation rates (smaller Ek) would have been available.

The dynamics and strength of the BZF
Most of our measurements here need to be interpreted in light of the previously discovered BZF (de Wit et al. 2020;Zhang et al. 2020). In order to remind the reader, we show in figure 4(a) a schematic based on a projection of simulated temperature data onto the sidewall. In a narrow region close to the sidewall the azimuthal warm fluid rises on one side and sinks on the opposite side. We want to point out that, in fact, the wavenumber of this periodicity is k = 2Γ , so that in a Γ = 1 cell, there are two areas of warm upflow separated by two cold downflow regions (Zhang et al. 2021). While the time-averaged azimuthal velocity in this thin region is positive (prograde), the warm-cold structure itself drifts in the negative (retrograde) direction, as observed in simulations (figure 4(b) and Zhang et al. 2020).
While simulations provide a very detailed view of the temperature structure, in the experiments we can make significantly longer observations for better statistics. For this, we analyse measurements that are taken inside the cylindrical plastic sidewall at heights z = H/4, H/2 and 3H/4 and at 8 equally spaced azimuthal positions θ = (0, π/4, π/2, 3π/4, π, 5π/4, 3π/2, 7π/4). The structure of the BZF can be seen well in figure 4(c), where the temperature (colour code) measured at z = H/2 is plotted as a function of the azimuthal position (x-axis) and time (y-axis). The time span plotted here corresponds to approximately 8 min. Albeit the spatial resolution is significantly smaller than in the simulation (figure 4b), one clearly sees the signature of the BZF, namely a mode k = 1 wave with warm temperature on one side of the cell and cold temperature on the opposite side. In particular, we see here that the temperature structure drifts in the azimuthal direction with negative drift rate ∂θ m /∂t < 0, where θ m is the azimuthal location of the warmest temperature, as defined below in (4.1).
Because the BZF in our Γ = 1/2 cell has an azimuthal wavenumber of k = 1, we can analyse it in the same way as large-scale circulation is usually analysed for the non-rotating (see e.g. Verzicco Figure 4. (a) Schematic of the BZF. Top and bottom show the average azimuthal velocity, the sidewall is colour coded with a snapshot of the simulated fluid temperature. The schematic was created using simulation results that were published already in Zhang et al. (2020). The black circle marks the circumference at midheight, at which temperature was measured and plotted against time in (b). (b) Simulated dimensionless temperature at r = R and z = 0.5H as a function of time for Ra = 10 9 and 1/Ro = 10 (adapted from Zhang et al. 2020).
(c) Azimuthal temperature distribution for Ra = 4.9 × 10 13 , 1/Ro = 9.18, r/R = 1.0, z/H = 0.5 as a function of time. (d) Temperature distribution as a function of the azimuthal angle (blue bullets) and a fit of (4.1) to the data (red solid line). Conditions are as in (c). The fitted parameters θ m and δ m are also marked by a solid vertical line and a down-pointing arrow. Note that for both plots we plot the temperature at θ = 0 also for θ = 2π for a better visual appearance. measurement in time a harmonic function to the data. Here the index i denotes the azimuthal position, and j stands for the vertical position and takes indices j = ('b', 'm', 't'), corresponding to the vertical locations z = (H/4, H/2, 3H/4). The fit parameters are the azimuthally averaged wall temperature T w,j , the amplitude δ j and the orientation θ j . An example of the data and the corresponding fit is shown in figure 4(d). Note that this approach works for both analysing the LSC for small rotation rates as well as analysing the BZF for larger rotation rates. In the following, we will only focus on the dynamics of the azimuthal orientations θ j . The average temperature T w,j will be discussed later in § 5. The amplitudes δ j are well reflected in measurements of d (see § 6.1) and σ (see § 6.2) and therefore will also not be discussed now. Figure 5(a) shows the evolution of θ m (at midheight) as a function of time. Here, we have made θ m continuous by adding or subtracting 2π whenever the difference θ m (t i+1 ) − θ m (t i ) at consecutive time steps was larger than π or smaller than −π, respectively. The data plotted in this way show a nearly linear change of θ m with time, indicating a monotonic drift with fairly constant drift speed ω m = ∂θ m /∂t t .

The azimuthal drift of the BZF
We already see here that the average drift velocity is a non-monotonic function of 1/Ro. Even the direction of rotation changes with increasing 1/Ro. For very small 1/Ro ω m is positive, which means that the observed structure rotates in the same direction (in the rotating frame) as the convection cell, i.e. in the prograde or cyclonic direction. The value of θ m decreases with time for large 1/Ro, which suggest a retrograde or anticyclonic motion of the temperature structure in the region close to the sidewall. This anticyclonic movement of the temperature field corresponds to the reported drift of the BZF in Zhang et al. (2020).
For a more quantitative analysis we calculate the average drift speed ω (b,m,t) as the slope of a linear fit to θ (b,m,t) . The result is plotted in figure 5(b). We see that, for sufficiently small rotation rates (1/Ro < 1/Ro * 1 ), there is a discrepancy between the drift rates measured at heights H/4, H/2 and 3H/4. While ω b and ω t are almost zero, we see a finite drift of the structure in positive (cyclonic) direction at the middle, i.e. ω m > 0. With increasing 1/Ro also ω m increases and reaches a maximum at around 1/Ro ≈ 0.3 or so, before it decreases again and reaches zero at 1/Ro ≈ 1/Ro * 1 . We note that differences in the drift rates at the top and bottom compared to the midheight were also observed in measurements in a Γ = 1/2 cell with water (Pr = 4.38) as the working fluid (Weiss & Ahlers 2011b). There, however, a prograde rotation was observed at all three heights for slow rotation that turned into a retrograde rotation for faster rotation rates.
What is measured for such small 1/Ro is not the BZF, but the structure of the LSC, which is not a coherent structure close to the sidewall. Consider a large-scale circulation roll of elliptical shape that lies diagonally inside the convection cell. Then, the fluid flows at midheight on average towards the cell centre, and the Coriolis force causes an acceleration of the azimuthal velocity, causing a cyclonic drift at midheight but no or only a very slow drift at H/4 and 3H/4. The shape and strength of the LSC is clearly a function or Pr and hence quantitative differences are expected for different Pr.
For 1/Ro > 1/Ro * 1 the drifts ω j are the same for all heights and negative, i.e. in anticyclonic direction, suggesting a vertical coherent temperature structure -the BZF. With increasing rotation rate ω j decreases as well and reaches a minimum at 1/Ro ≈ 6. A further increase in 1/Ro results again in an increase of ω j , so it reaches zero asymptotically.
We will now compare the effect of Ra on this observation. Therefore, we show in figure 5(c) ω m , normalised with the rotation velocity of the convection cylinder Ω for different Ra as a function of 1/Ro. First of all, we see that data for different Ra collapse very well for large 1/Ro, but not so much for small 1/Ro. That is an initial cyclonic drift (ω m /Ω > 0) for small 1/Ro turns into an anticyclonic drift ω m /Ω < 0 for large 1/Ro. The 1/Ro where ω m /Ω = 0 is, albeit close to 1/Ro * 1 , smaller for larger Ra. A similar Ra-dependency has also been observed in previous experiments with water (Pr = 4.38) (Weiss & Ahlers 2011b). For larger 1/Ro, ω m /Ω reaches again a minimum and then increases asymptotically to zero.
The collapse of different Ra-data for sufficiently large 1/Ro occurs when the BZF is present, which suggests that the drift rate of the BZF (ω m ) is for a given 1/Ro proportional to Ω (or 1/Ek), but otherwise independent of Ra. We only cover nearly a decade in 1/Ro, and therefore it is difficult to compare our data with scaling laws observed or predicted by others. Nevertheless, we show in figure 5(d) the same data again (−ω/Ω) for large 1/Ro plotted in a double-log plot. In this logarithmic representation of −ω/Ω, small dependencies of Ra become visible at the largest 1/Ro. We note, however, that also Pr varies slightly and hence it is unclear whether the variations in ω/Ω are caused by Ra or Pr.
In this representation, no single scaling over the entire shown 1/Ro-range is visible. For 2 1/Ro 7 the data seem to follow ∝ (1/Ro) −3/4 or so. For larger 1/Ro data for different Ra seem to diverge slightly and show a larger negative exponent. We note that Zhang et al. (2021) observed a normalised BZF drift rate of ω/Ω ∝ 1/Ro −5/3 , which is not too far from the small-Ra data at larger 1/Ro. We further note that in a recent measurement in slender cylinders (Γ = 0.2) at Pr = 5.2 (de Wit et al. 2020) a relation ω sc H 2 /ν = 6 × Ra 1.16 was measured, however, at constant Ek = 10 −7 . We cannot 10 -1 10 0 10 1 10 -1 10 0 10 1 compare our data with their results, since their dependency on Ek is unknown. We note, however, that a relation ω sc H 2 /ν ∝ Ek 3.32 Ra 1.16 would remove the Ra-dependency when plotted as a function of 1/Ro. This in turn would then lead to ω/Ω ∝ (1/Ro) −2.32 , which is rather different to our data.

4.2.
The azimuthal Fourier modes of the BZF Fitting (4.1) to the thermistor measurements as described above also gives the amplitudes δ j , which we have not discussed so far. By fitting (4.1) to the 8 sidewall thermistors, we gain with δ j basically the energy in the first Fourier mode (E 1 ). However, since we measure the temperature at 8 azimuthal positions, we can also calculate the energy of the second, third and fourth, Fourier modes, namely E 2 , E 3 and E 4 . Figure 6 shows the time average of the first four Fourier modes normalised with the total spectral energy E tot = 4 n=1 E n as a function of 1/Ro for Ra = 2 × 10 11 (a) and Ra = 8 × 10 14 (b). Up to a rotation rate of 1/Ro ≈ 1/Ro * 2 the first mode E 1 is significantly larger than the others. This is due to the existence of the LSC for small 1/Ro and due to the BZF at larger 1/Ro. However, it is interesting that E 1 /E tot reaches a local maximum pretty much at 1/Ro * 1 . After that, it decreases slightly first, before it increases again with increasing 1/Ro. In experiments with water in a cell of Γ = 0.5 we have seen similarly a maximum at 1/Ro ≈ 0.8 of E 1 /E tot (see figure 4a in Weiss & Ahlers 2011b).
We can only speculate, but it seems that E 1 is enhanced right when the BZF starts to appear and the LSC vanishes. The fact that we see the same feature for very different Ra and Pr suggests that the transition of the LSC into the BZF only depends on 1/Ro. For larger 1/Ro, say 1/Ro > 1/Ro * 2 , the first mode E 1 /E tot decreases for small Ra = 2 × 10 11 , while for Ra = 8 × 10 14 it continues to increase. This is in accordance with the behaviour of the standard deviation σ shown in figure 15. This finding is a bit puzzling. While the decreases of E 1 /E tot suggests that the BZF only exists in a finite 1/Ro-range, which depends on Ra, we in fact know from simulation (Zhang et al. 2020) that the size of the BZF layer δ BZF is expected to slightly decrease with Ra at a fixed 1/Ro following δ BZF ∝ Ra −0.08 Ro 0.66 . Since it seems reasonable to believe that the pumping efficiency inside the BZF decreases with decreasing thickness δ BZF one would expect the decrease of E 1 /E tot to appear at smaller 1/Ro for the larger Ra, i.e. the opposite to what is observed in the measurements.

Mean temperature profiles
As shown in figure 1, 62 thermistors were installed inside the cell, predominantly close to the sidewall, to acquire vertical temperature profiles. We now want to analyse these data and see how the vertical temperature profiles as well as fluctuations change as a function of 1/Ro. We note that the thermal and viscous boundary layers close to the top and bottom plates are only of the order of a millimetre or less and are not resolved in these measurements. Data presented here thus are measurements in the (vertical) bulk region.
Before we discuss the effect of rotation on the vertical temperature profiles, we quickly recapitulate what is known for the non-rotating case. It has been shown by Ahlers et al. (2012a) and  that the vertical temperature profile close to the sidewall can be well approximated by a logarithmic function of the form where Θ(z) denotes the normalised temperature with · · · denoting the time average over an entire run. The logarithmic slope of the temperature profiles, i.e. the coefficients a t,b , decrease with Ra following a power law of a t,b ∝ Ra −η , where η depends on the radial position. Our measurements were conducted with the same cell that was used for the comprehensive investigation of the vertical temperature profile ) and our measurements for the non-rotating case are in accordance with those reported in Ahlers et al. (2014). Figure 7 shows Θ as a function of z very close to the sidewall at radial distance r = 0.98R, for Ra = 4.9 × 10 13 and different 1/Ro. While the solid lines in figure 7(a) are just guides to the eyes, we show the same data in figure 7(b) plotted against a logarithmically scaled x-axis, and see that the data indeed follow straight lines for all rotation rates. This indicates that, also under rotation, the vertical temperature profiles are well represented by a logarithmic function ((5.1) and (5.2)) over the measured vertical distance, which is approximately a single order of magnitude. Furthermore, we see that the profiles at the bottom and the top behave rather similarly, at least for the Ra plotted here. The logarithmic temperature gradient increases with increasing rotation, which also suggests that the temperature drop across the thin top and bottom boundary layers is reduced with increasing rotation rate, an observations that is in accordance with previous findings, such as, for instance, those by Julien et al. (2012b). The steeper temperature gradient in the bulk is a result of suppressed vertical fluid motion due to Coriolis forces. We also see that the temperature at midheight is close to the mean temperature T m between top and bottom plates. This is because non-OB effects are rather small in this particular case as the temperature difference between bottom and top was only T = 5 K and the pressure was not too high (10 bar -run E2e). We note that non-OB effects in pressurised gases are expected to lower the midheight temperature compared to T m . The very same effect would occur under the presence of centrifugal forces, which apparently also do not play a significant role here.  Figure 8 shows the logarithmic slope a b as a function of 1/Ro for data taken at r/R = 0.98. We have seen already in figure 7 that, with increasing rotation, the coefficient a t (a b ) increases. In figure 8(a), we see that the change of the logarithmic slope a with 1/Ro is not the same for all 1/Ro. Instead, a for small 1/Ro remains rather constant, but increases for larger 1/Ro. The transition between these two regimes coincides approximately with 1/Ro * 1 . The increase of a t (a b ) is well approximated for large 1/Ro by a power law a t,b ∝ 1/Ro ε , with exponents ε = 0.51 ± 0.02 and 0.50 ± 0.02 for both top and bottom halves of the cell, respectively. The fact that this exponent resembles a pitchfork bifurcation is merely a coincidence for the particular case of Ra = 4.9 × 10 13 . We will see below that the exponent actually depends on Ra.
In figure 8(b) we show how the offset parameters b t and b b change with increasing 1/Ro. From (5.1) and (5.2) we see that b t,b give us the fitted temperature at midheight (z = H/2). Since in the ideal OB case it is θ(H/2) = 0, the b t,b mark whether the fits overshoots or undershoots this value, and hence also contain information on the relative slope of the temperature profile at the top and the bottom compared to the temperature gradient at midheight. For the ideal OB case, we expect b b and b t to have opposite signs. For a better comparison, we thus plot in figure 8(b) b t and −b b as functions of 1/Ro. We see that b t is positive and b b is negative for all 1/Ro, denoting that the logarithmic fits undershoot slightly the temperature at z = H/2. Furthermore, we see that there is an asymmetry between the bottom and the top. For small 1/Ro b t is considerably smaller than b b . This discrepancy might be due to small non-OB effects, or imperfections in the experimental set-up that destroy the up-down symmetry of the system. While b t is rather constant for 1/Ro < 1/Ro * 1 , it increases for 1/Ro > 1/Ro * 1 . This trend is also visible for b b , albeit much weaker.
In general, the fit of (5.1) and (5.2) to the data returns small residuals for not too large rotation rates, i.e. 1/Ro 8. For larger 1/Ro the vertical temperature profile starts to deviate from a logarithmic profile consistently, as shown in figure 8(b). We further note that the residuals are in general larger at the cold top part of the cell than in the warm lower part. However, for sufficiently fast rotation (1/Ro 1/Ro ≈ 8), the residuals at both locations merge and increase consistently.
To compare measurements at different Ra, we show in figure 9(a) the logarithmic slopes normalised with their corresponding values without rotation (i.e. a t (1/Ro)/a t (0)) as a function of the rotation rate 1/Ro. In this way, we can compare the change of the vertical profile for different Ra. We see in this log-log representation that the onset of the a t enhancement occurs at around 1/Ro * 1 for all Ra. The increase of a t for 1/Ro > 1/Ro * 1 is stronger for larger Ra. This comes as a surprise, since the change in Nu seems rather independent of Ra. Moreover, from figure 2, the Nu-reduction with Ra seems to be slightly smaller for larger Ra. We note in this context that the local effective heat conductivity is inversely proportional to the local vertical temperature gradient λ eff = q/( T/ z). In particular, for larger Ra the slope increase can be well approximated using a power law with a Ra-dependent exponent of the form a t (1/Ro)/a t (0) = K a × (1/Ro) η(Ra) , such that all data collapse. The fitted exponents η are plotted as functions of Ra in the inset of figure 9(b) on a semi-logarithmic plot. We note that, in this representation, the data for η do not follow a straight line. In particular, the decrease of η is stronger with decreasing Ra than a logarithmic function would suggest. However, as can be seen in figure 9(a), the normalised slopes deviate from the power law particularly strongly for small Ra and thus the uncertainty of η is also large for these data points. However, a function η = 0.031 ln(Ra) − 0.45 represents the data decently well (black line in the inset in figure 9b) and therefore all data a t /a t (0) collapse onto a single master curve when plotted as a function of (1/Ro) η (Ra) . At this point, these findings are purely empirical, but they hopefully can be used to compare with other measurements, DNS or theoretical models.
From measurements in the non-rotating system  we know that the logarithmic slope of the vertical temperature profiles increases with increasing radial distance from the centre line r/R. This is not surprising, since at least for the non-rotating  case vertical fluid motion and hence vertical fluid transport close to the sidewall are reduced. In figure 10(a) we show how a b changes with 1/Ro at different r/R. For small rotation rates (1/Ro < 1/Ro * 1 ), similar to the non-rotating case, the a b measured closest to the radial centre (r/R = 0.73) have the smallest values and increase with radial distance. It is interesting, however, that a b in this regime is not affected by a change of 1/Ro for r/R = 0.98 and r/R = 0.96, very little for r/R = 0.93 but it changes increasingly for r/R = 0.73. As a result, the a b for different radial distances are getting closer to each other with increasing 1/Ro and finally, somewhere before 1/Ro * 2 , collapse on a single power law curve of the form a b ∝ (1/Ro) 0.63 (for Ra = 4.9 × 10 13 ), so that a b is independent of the radial location.
A similar behaviour is observed for the offset b b shown in figure 10(b). Although the data are more scattered, the values for different r/R differ for 1/Ro < 1/Ro * 1 , but are getting closer for larger 1/Ro. This behaviour suggests that a radial symmetry, which is broken by the LSC in the non-rotating case, is gradually being restored in the rotating case. This is quite surprising. As we will discuss in the next section, close to the wall (r/R = 0.98, 0.96 and 0.93) the temperature and velocity fields are dominated by the BZF and one would expect to see a different average temperature field there and in the radial bulk (r/R = 0.73). It might be that, due to Coriolis forces under rotation, mixing in the horizontal direction is much more effective and thus faster compared to the vertical transport of warm and cold fluid, or in other words, the eddy diffusion time scales in the vertical direction are much smaller than those in horizontal direction. Therefore, vertical temperature gradients increase and horizontal temperature gradients must decrease with increasing rotation rates. Thus, even though the temperature and flow fields are different in the BZF and the radial bulk, the time-averaged temperatures are surprisingly similar.
6. Temperature fluctuations 6.1. The full PDF of the temperature fluctuations So far, we have looked at the time-averaged temperature at different positions as a function of Ra and 1/Ro. In this section, we want to look at temperature fluctuations, as they provide valuable information about the flow field and the specific heat transport mechanisms. Therefore, we first consider the full PDF of the normalised temperatureT = (T − T m )/ T and will discuss some of its quantities later. Figure 11 shows the temperature PDFs p(T) for heights z/H = 0.287 (bottom), 0.493 (midheight) and 0.75 (top) and four different rotation rates, measured close to the sidewall at r/R = 0.98. For the non-rotating case the distribution is symmetric for the midheight (green diamonds) with slightly larger tails than what would be expected for a Gaussian distribution. The distributions for the temperature at the bottom and the top are skewed in the way that at the bottom (top) warmer (colder) temperatures are measured more often than for a Gaussian distribution. The deviation from the Gaussian originates from the warm (cold) plumes that are rising (falling), which are not well mixed with the turbulent fluid. In fact, it has been suggested by Wang, He & Tong (2019) that temperature distributions close to the top and bottom plates can be represented by a superposition of a Gaussian, which reflects fluctuations in the well-mixed turbulent bulk, and an exponential distribution, which represents the contribution from thermal plumes. While we refrain from a detailed analysis of our temperature measurements, we note that the distributions for the non-rotating case are qualitatively consistent with such a superposition.
The temperature distribution clearly changes with increasing rotation rate, not just in quantitative, but also in qualitative terms. At very slow rotation rates, p(T) actually becomes less skewed and follows more a Gaussian distribution, as shown in figure 11(b). This suggests that warm and cold plumes mix faster with the background fluid under very slow rotation. Very warm and very cold plumes lose their heat faster under faster rotation when they move from the bottom and top plates towards the thermistor positions at H/4 and 3H/4. At faster rotation rates, the standard deviation increases as the maximum of the distribution widens. At sufficiently fast rotation, in p(T) the peak splits into two peaks denoting a bimodal temperature distribution (figure 11c). The two peaks move apart from each other with increasing 1/Ro, see figure 11(d). We have already reported on this observation in a previous publication, in which we explain this bimodal PDF as the signature of the BZF that transports warm fluid from the bottom to the top close to the sidewall at one side and cold fluid to the bottom at the opposite side, thus preventing effective mixing between them (Zhang et al. 2020).
For a quantitative analysis we fit a superposition of two Gaussian functions of the form to the data. Equation (6.1) has five independent parameters, the means μ i , the standard deviations σ i (i = 1, 2) and a parameter A ≤ 1, which determines the relative ratio of the amplitudes of the two peaks. From now on we follow the convention that the index 1 refers to the colder mode, such that it is always μ 1 < μ 2 . By fitting (6.1) to the data, we can describe the PDF by the parameters A, μ i and σ i and we can investigate how these parameters change with increasing rotation rates.
We show in figure 12(a) μ 1 and μ 2 as functions of 1/Ro at different vertical positions. We see that for small rotation rates (1/Ro 1/Ro * 1 ) all μ 1,2 are nearly independent of 1/Ro and μ 1 (open symbols) is for any height rather close to μ 2 (closed symbols). In particular, for measurements at midheight, the data points are on top of each other, indicating that Relative distance of the two peaks d as defined in (6.2). Vertical dashed lines mark 1/Ro * 1 and 1/Ro * 2 . Data were acquired at Ra = 4.9 × 10 13 (E2e) and at radial position r/R = 0.98. no bimodal Gaussian is present. As we have mentioned above, the temperature PDFs measured at the top and the bottom half of the cell are skewed. This skewness is detected when fitting (6.2) to the data and hence we get slightly different values for μ 1 and μ 2 . This is not the case for measurements at midheight, where warm and cold plumes are sufficiently mixed. In the regime of small rotation rates, it even seems as if the difference between μ 1 and μ 2 decreases slightly for the top and bottom parts with increasing 1/Ro. When the rotation rates exceed 1/Ro * 1 , for the measurements at midheight, two clearly distinguishable peaks occur, (i.e. a bimodal distribution) resulting in a decreasing μ 1 and an increasing μ 2 . For measurements at the top and bottom, however, only one of the μ 1,2 changes in this regime. These are a decrease of μ 1 for the cold top and an increase of μ 2 for the warm bottom. This observation is quite significant. It suggests that, for example, at the bottom, the pumping of warm fluid from the bottom inside the BZF starts to play a role, while mixing is still small so that the cold flow from the top is not being heated up very strongly. In fact, in this regime, the red open squares (μ 1 at the top) and the blue solid circles (μ 2 at the bottom) approach each other, meaning that the cold areas at the bottom become a bit colder while the warm areas at the top become a bit warmer due to the relative strong pumping inside the BZF.
This changes at 1/Ro * 2 . Now the warm fluid at the top part (μ 2 -blue solid circles) also cools down and the cold fluid at the bottom (μ 1 -red open squares) heats up. While there is still a vertical transport of fluid inside the BZF, the size of the BZF decreases with increasing rotation rate (see Zhang et al. 2020). Therefore the cold down (warm up) flow inside the BZF gets heated (cooled) from the bulk, where the temperature gradient increases. This observation is in accordance with the temperature gradient measurements presented in figure 10(a). There, the gradients in the bulk and close to the sidewall, i.e. in the BZF converge for rotation rates 1/Ro * 1 < 1/Ro < 1/Ro * 2 . It is also insightful to plot the temperature difference between the warm and the cold areas in the BZF. For this, we look at the relative distance between the two peaks as suggested in Holzmann & Vollmer (2008)  We note that d tells us something about how well (small d) or badly (large d) warm and cold fluid is mixed. We plot in figure 12(b) d as a function of 1/Ro for measurements at the top, the middle and the bottom parts of the cell. In this representation, we see, in particular for the midheight of the cell (green diamonds), a bifurcation-like behaviour of d at 1/Ro ≈ 1/Ro * 1 , which we interpret as the onset of the BZF. The PDF consists only of a single peak for smaller 1/Ro -thus d ≈ 0 -and at 1/Ro * 1 the PDF splits into a bimodal distribution, with warm fluid being pumped inside the BZF from the bottom to the top and cold fluid being pumped from the top to the bottom, causing d to sharply increase. The increase of d continues for larger 1/Ro. This implies that, in this regime, the lateral mixing is small compared to the heat transport by vertical pumping from the top and bottom. For the largest 1/Ro, however, we see a decrease of d for just a single point, which we believe is real and not just an outlier. We know (Zhang et al. 2020) that the width of the BZF layer shrinks with ∝ (1/Ro) −0.66 and thus thermal diffusion starts to thermally couple the bulk with the BZF, which is expected to lead to a reduced d at sufficiently large rotation rates.
The behaviour of d for the top (blue circles) and bottom (red squares) is for small 1/Ro qualitatively similar and for larger 1/Ro also quantitatively similar. Just as we discussed above, for small rotation rates (1/Ro < 1/Ro * 1 ) at heights H/4 and 3H/4, d is finite due to the skewness of the PDF. This fact is captured when naively fitting two Gaussian (i.e. (6.1)) to the PDF.
The qualitative changes of p(T) are restricted to a thin region close to the sidewall (the BZF). This is shown in figure 13, where p(T) is plotted as measured at r/R = 0.98 (a) and at r/R = 0.73 (b) for different 1/Ro and at midheight of the cell. We see that the change from a uni-modal to a bimodal distribution with increasing 1/Ro is clearly seen at r/R = 0.98, but not for r/R = 0.73. There, p(T) seems to be unaffected by rotation in the observed 1/Ro-range. Up to large 1/Ro p(T) follows a stretched exponential, which is in accordance with previous findings for the temperature fluctuations in the bulk (e.g. Castaing et al. 1989;Ching 1991;Niemela et al. 2000). An exponential fit of the function e −c|x| β yields β = 1.34 ± 0.03, i.e. the decrease is somewhere in between an exponential (β = 1) and a Gaussian (β = 2).
It is somehow unfortunate that we cannot resolve the BZF sufficiently well in space and therefore cannot plot a detailed radial profile of it, e.g. d(r/R). We nevertheless plot in figure 14 Figure 14. Peak-to-peak distance at midheight d as defined in (6.2) as a function of the radial position. Filled symbols are data measured for Ra = 4.9 × 10 13 . Open symbols show estimates for the width of the BZF δ BZF according to Zhang et al. (2020). Dashed lines are guides to the eyes. and add estimates for the width of the BZF (δ BZF ) from numerical simulations (Zhang et al. 2020(Zhang et al. , 2021 as open symbols. These estimates are δ BZF = 3.55 × Ra −0.08 × (1/Ro) −0.66 . The dashed lines in figure 14 are convenient functions as a guide to the eye to better illustrate the estimated radial profile of d. We note that the δ BZF are in accordance with our observation.
6.2. Amplitude of temperature fluctuations Above, we have already analysed average temperatures ( § 5) and the locations of the peak of the PDF. Now, we want to look at the standard deviations of the temperature, defined as Here, . . . t denotes the average over time. Please note that in this way we have already normalised the temperature fluctuations by the temperature difference T. Before we discuss the rotating case, let us again quickly recapitulate what is known for the non-rotating case.
At a given position and for sufficiently large Ra, σ decreases with increasing Ra, according to σ ∝ Ra −α , as has been observed e.g. by Castaing et al. (1989), Niemela et al. (2000), Daya & Ecke (2002) and Grossmann & Lohse (2004). This decrease is due to stronger mixing as the turbulence intensity increases. The exponent α does slightly depend on the radial position in the flow. In fact we measure α = 0.101 for r/R = 0.73 and α = 0.084 for r/R = 0.98 for measurements at midheight (z = H/2). These values are comparable with other values found in the literature.
Furthermore, the temperature fluctuations depend strongly on the vertical coordinate and are largest when measured close to the top and bottom plate and decrease towards the vertical centre. Similar to the time-averaged temperature, also the temperature fluctuations can be modelled as a logarithmic function of z/H, as reported by He et al. (2014).
The data analysis in this subsection has so far dealt with temperature fluctuations and how they depend on Ra and their radial and vertical locations. We now want to focus on how the fluctuations are influenced by rotation. For this we show in figure 15(a-c) temperature fluctuations at midheight of the cell at two different radial distances and for different Ra. Let us first focus on figure 15(a), where σ is plotted as a function of 1/Ro, for different Ra. One sees that, at very small rotation, σ is larger for smaller Ra, just as discussed above for the non-rotating case. For low rotations rates, σ increases monotonically until it reaches a maximum, beyond which the temperature fluctuations decrease again with increasing 1/Ro. The rotation rate at which σ reaches its maximum depends on Ra and is smallest for the smallest Ra; in our observed Ra range the maximum for the smallest Ra (≈ 7.7 × 10 9 ) is reached at 1/Ro ≈ 4.
For a better comparison of the effect of Ra on σ we normalised for each Ra σ (1/Ro) by σ (0) without rotation and plot the result in figure 15(b). We see that the data normalised in this way collapse rather well for 1/Ro < 2, but start to deviate for larger 1/Ro. In this normalised plot, one also sees that the relative increase in σ is in fact larger for larger Ra. The data plotted here are acquired where the BZF occurs for sufficiently large 1/Ro ( 1/Ro * 1 ). While for the non-rotating case larger Ra leads to a better mixing of warm and cold fluid, this fact is not very relevant inside the BZF, as it is less influenced by the turbulent flow in the radial bulk. In this regard, it is important to note that no sharp transition is observed in these plots. Instead, σ increases rather smoothly until it reaches its maximum. This is all the more surprising since we have seen already in figure 12(b) that the distance between the maxima in the bimodal PDF increases quite abruptly at 1/Ro * 1 . One would guess that the distance between the peaks of the PDF d contributes heavily to the measured temperature fluctuations. While this might nevertheless be the case, the plots here suggests that the PDF widens significantly before it splits into a bimodal distribution, i.e. before the BZF occurs.
For comparison, in figure 15(c) we also plot the normalised fluctuations as measured at r/R = 0.73, outside the BZF in the radial bulk. The difference to the measurements in the BZF are clear, particularly for larger Ra. The fluctuations there are drastically reduced. Instead of a continuous increase with increasing 1/Ro, σ/σ 0 barely changes with 1/Ro up to ≈ 1/Ro * 1 , but then even decreases for the largest Ra with increasing 1/Ro. The Ra-dependency for larger 1/Ro is reversed compared to what was observed at r/R = 0.93. That is the smallest σ/σ 0 is for the largest Ra. For the smallest Ra, we still see a maximum close to 1/Ro * 2 , which even has very similar values of σ/σ 0 as measured for r/R = 0.93. Why are measurements at r/R = 0.73 and r/R = 0.93 so different for large Ra but similar for small Ra? The answer is simple. The size of the BZF also depends on ∝ Ra γ with an exponent of γ ≈ −0.08 according to Zhang et al. (2021). Based on this estimate, at 1/Ro ≈ 1/Ro * 2 the BZF reaches for Ra = 7.7 × 10 9 still up to r/R = 0.71 and the thermistor still measures the BZF. This is also the case for Ra = 2.1 × 10 10 , where 1 − δ BZF ≈ 0.73r/R. For the third smallest Ra, the thermistor at r/R = 0.73 is outside but very close to where the BZF occurs, and one would expect the thermistor to still measure the thermal signature of the BZF, albeit with a smaller amplitude.
We see for Ra = 8.4 × 10 12 and Ra = 4.9 × 10 13 an increase of σ/σ 0 for 8 < 1/Ro, when measured at r/R = 0.73. Further investigations are necessary to interpret this increase; however, it might be related to coherent vertical vortex structures, such as Ekman vortices in which warm and cold fluid is transported from the boundaries deep into the bulk areas.
After we have discussed fluctuations at midheight of the cell, we briefly want to discuss how the vertical temperature fluctuation profile changes with increasing rotation rate. In figure 15(d,e) we therefore plot the temperature fluctuations against a logarithmically scaled z-axis. Measurements were taken inside the BZF at r/R = 0.98. The fluctuations are normalised with the fluctuations at midheight σ 0 to indicate the relative change as a function of the vertical position. Without rotation, fluctuations at the bottom and the top are significantly larger than the midheight fluctuations and, with increasing 1/Ro, the relative temperature fluctuations close to the top and bottom boundaries decrease, in other words the slope ∂(σ (z)/σ (0.5))/∂z decreases. Please recall that this is different to the logarithmic slope of the time-averaged temperature profiles discussed in § 5, which increases for increasing 1/Ro. Here, however, fluctuations at the bottom (top) are just as large as at midheight, since they mainly originate from the warm and cold structures within the BZF.

Skewness γ 1 (T)
In § 6.1 we have seen that the symmetry of p(T) changes with increasing rotation when the measurements are not taken exactly at midheight. The symmetry is in fact a measure of the cause for the temperature fluctuations as shown for example by Wang et al. (2019). When the fluctuations are predominantly caused by the turbulent background a nearly symmetric Gaussian PDF is expected. The presence of warm and cold plumes near the bottom and top causes the PDF to be skewed. Therefore the skewness is a good measure for the existence of thermal plumes, and how well they are mixed with the turbulent interior. For a quantification of these effects, we measure the skewness γ 1 defined as as a function of 1/Ro. Figure 16 shows the change of γ 1 with increasing rotation rate measured at vertical positions z ≈ H/4 (red squares), H/2 (green diamonds) and 3H/4 (blue bullets) at a radial position r/R = 0.98 where the BZF occurs. The trend of reduced skewness with increasing 1/Ro, as we have observed already in figure 11, is clearly visible here. For small rotation rates, the skewness has different signs for the bottom and the top due to the presence of the warm and cold plumes that skew the PDF to warmer and colder temperatures. Only for the four smallest rotation rates is γ 1 unaffected by the rotation of the cylinder. Already for 0.2 1/Ro the absolute values |γ 1 | start to decrease significantly with increasing 1/Ro and reach values close to zero at around 1/Ro ≈ 1/Ro * 2 . The start of the decrease does not correlate with the occurrence of the BZF at 1/Ro * 1 as observed by the sudden increase in d shown in figure 12. There is a small jump visible for the bottom data (red squares) at 1/Ro * 1 , but this might be just an outlier of two points close by. No such step is observed for the data taken at the upper part of the cell (3H/4, blue bullets). While one needs to be careful in interpreting these data, one can clearly say that the signature of the thermal plumes that exist next to a turbulent background is reduced with increasing rotation rates, even before the BZF occurs.
The reduced skewness is in accordance with a better lateral heat exchange, which also leads to a more homogeneous time-averaged temperature (figure 5). That is, because of a slower vertical transport, a relatively better lateral diffusion and, therefore, a reduction of the most extreme temperatures, resulting in a more symmetric PDF.
As discussed above, it seems as though the change in γ 1 does not strongly correlate with the onset of the BZF at 1/Ro * 1 . We therefore also look at measurements at r/R = 0.73, i.e. outside the BZF. These measurements are shown in figure 17. It is unfortunate that we do have a functional thermistor at the upper part of the cell at z/H = 0.75, but not at the corresponding position at z/H = 0.25. We therefore consider a thermistor at z/H = 0.144 which we expect to measure qualitatively similar features, as this is still sufficiently far away from the bottom boundary layer.
We see in figure 17 that, similar to the measurements at r/R = 0.98, above a certain 1/Ro, also here γ 1 decreases (increases) for measurements at z/H = 0.144 (z/H = 0.75) with increasing rotation rate. The explanation for the decrease (increase) of γ 1 here is similar. The lateral mixing is enhanced compared to vertical pumping. While qualitatively similar, there are quantitative differences. Most notably, and somehow surprisingly, we see that |γ 1 | starts to decrease at higher rotation rates, beyond 1/Ro * 1 . Furthermore, γ 1 does not decrease (increase) to nearly zero, but changes signs at around 1/Ro ≈ 2.5 or so. While we cannot explain at this point where this change of sign comes from, we note that a skewness close to zero is expected for a well-mixed system. A non-zero skewness hints that a new pumping mechanism occurs next to the turbulent mixing. This might be related to the formation of Ekman or Taylor vortices inside the bulk.
We also note that the data presented in figure 17 were acquired at Ra = 4.9 × 10 13 . At this Ra we have seen in figure 15(c) that, for the three largest 1/Ro, the fluctuations σ/σ 0 increase again. The same measurements here show a decreased |γ 1 |. At this point, we can only speculate about this change in the monotonic behaviour, but we note that the Froude numbers for the points with the largest rotation rates are Fr = 0.32 (for 1/Ro = 9.8), Fr = 0.4 (for 1/Ro = 11.0) and Fr = 0.5 (for 1/Ro = 12.2), and hence the influence of centrifugal forces on the flow field increases. Whether and in which way centrifugal forces are responsible for the decrease of |γ 1 |, however, demands further investigation.
Since we have seen that the skewness is clearly a function of the vertical location, we plot γ 1 as a function of z/H in figure 18 for different 1/Ro. Without rotation (up triangles in figure 18), γ 1 shows maximal values close to z/H = 0.25 (a) and z/H = 0.75 (b). The reason for this is the role of the large-scale circulation, which carries the plumes along and which has an elliptical shape. Very close to the sidewall at the bottom and top, corner rolls exist (Sun, Xi & Xia 2005) into which the warm and cold plumes usually do not enter. The plumes that detach from the boundary layers and that are carried by the LSC thus 'hit' the sidewall above (or below) these corner rolls. With increasing 1/Ro, |γ 1 | becomes smaller, . The x-axis is scaled logarithmically for better visual appearance. (a,b) Show the same data, but the vertical distance is measured from the bottom in (a) and measured from the top in (b). The dotted vertical line marks the midheight z = H/2. and the maxima and minima close to the top and bottom also disappear. Although in figure 18 we only show seven different 1/Ro, maxima in |γ | only occur for 1/Ro < 1/Ro * 1 , i.e. before the LSC is replaced by the BZF.

Summary and discussion
In this paper we present comprehensive heat flux and temperature measurements in rotating turbulent thermal convection at very high Rayleigh numbers. We use nitrogen and sulphur-hexafluoride at pressures between 1 bar and 19 bar in a two metre tall cylindrical cell of aspect ratio Γ = 0.5, thus covering the Ra-range 7 × 10 9 < Ra < 8 × 10 14 . Due to the slight dependency of the Prandtl number on the pressure, the Prandtl number increases slightly from Pr = 0.69 at 1 bar (nitrogen) to 0.97 at 19 bar (SF 6 ).
We apply slow and moderate rotations of the cylinder around its cylindrical axis, i.e. parallel to gravity and measure characteristic quantities as a function of the rotation rate. In contrast to similar measurements in fluids of larger Pr, we do not find any enhancement of the Nusselt number under rotation. Instead, rotation barely influences the heat transport for 1/Ro 1/Ro * 1 = 0.8, but leads to a strong decrease for rotation rates 1/Ro 1/Ro * 2 = 4. The decrease appears to follow ∝ (1/Ro) −0.43 , but we note that the uncertainty of the exponent is significant as the corresponding fit was done over less then a decade of 1/Ro and a logarithmic function fits the data at least equally well. The reduced Nusselt number Nu/Nu 0 appears to only depend on 1/Ro, but is independent on Ra.
We have also investigated the temperature at different radial and vertical positions as a function of 1/Ro. The resulting measurements need to be interpreted with reference to the recently observed BZF that forms under sufficiently fast rotation in the narrow region close to the lateral sidewall (Zhang et al. 2020). There, a periodic temperature structure occurs in which warm fluid rises along one side, while cold fluid sinks on the opposite side, hence the temperature fluctuations in this region are significantly larger than in the turbulent bulk. Nevertheless, similar to the non-rotating case, also under rotation, the time-averaged temperature at the sidewall can be well represented with a logarithmic dependency of the vertical coordinate, even when the BZF is present. However, the logarithmic slope (a t,b ) increases significantly for 1/Ro > 1/Ro * 1 . 912 A30-29 In addition, while without rotation a depends on the radial position (r/R) of the relevant thermistors, when the BZF is present, a becomes independent of r/R, regardless whether the measurement was conducted inside the BZF or inside the bulk region. That is even more surprising as, in particular, the temperature fluctuations (as measured via the standard deviation in time), show clear differences between the BZF (large fluctuations) and the bulk regions (small fluctuations).
We also analyse the full PDFs of the temperature measurements. As reported before in Zhang et al. (2020), the PDF is qualitatively very different in these states, and thus allows us to clearly distinguish the BZF from the bulk region. While the PDF in the BZF turns from a unimodal PDF at no or slow rotation into a bimodal distribution with two clearly distinguishable maxima for faster rotation; the PDF in the bulk region remains nearly unchanged and unimodal.
The onset of the BZF might not be a sharp one, although the change of certain experimentally accessible properties show a rather sharp bifurcation when 1/Ro increases beyond certain critical values. Although these values depend on which property is considered, it is usually somewhere between 1/Ro * 1 and 1/Ro * 2 . For example, a clear change of the temperature field close to the sidewall can be seen when looking at the distance between the two maxima of the bimodal temperature PDF. This distance starts to increase at approximately 1/Ro * 1 . Similarly, the azimuthal drift of the temperature field close to the sidewall changes direction from being prograde to retrograde at the same 1/Ro. From this study we can state that, at 1/Ro * 1 , a flow state transition takes place. While the flow self-organises into a LSC for smaller 1/Ro, the BZF exists for larger 1/Ro. However, it is not so clear what happens at 1/Ro * 2 . From our measurements, it seems that, at 1/Ro * 2 , the BZF has reached its maximal amplitude, meaning the pumping of warm and cold fluid from the top and bottom boundaries towards the other side causes the largest azimuthal temperature variations. Furthermore, as the BZF also contributes significantly to the vertical heat transport, Nu deceases rather quickly for larger 1/Ro since this pumping also decreases.
To summarise, here we present a very comprehensive dataset about heat flux and temperature measurements in rotating turbulent RBC that is essential to extrapolate the insights from numerical data to larger Rayleigh numbers.