2.1. Governing hydrodynamical equations

We consider RB convection of a Boussinesq fluid contained in a spherical shell of outer radius $r_{o}$ and inner radius $r_{i}$ . The boundaries are impermeable, no-slip, and at constant temperatures $T_{bot}$ and $T_{top}$ . We adopt a dimensionless formulation using the shell gap $d=r_{o}-r_{i}$ as the reference length scale and the viscous dissipation time $d^{2}/{\it\nu}$ as the reference timescale. Temperature is given in units of ${\rm\Delta}T=T_{top}-T_{bot}$ , the imposed temperature contrast over the shell. Velocity and pressure are expressed in units of ${\it\nu}/d$ and ${\it\rho}_{o}{\it\nu}^{2}/d^{2}$ , respectively. Gravity is non-dimensionalised using its reference value at the outer boundary $g_{o}$ . The dimensionless equations for the velocity $\boldsymbol{u}$ , the pressure $p$ and the temperature $T$ are given by

where $\boldsymbol{e}_{r}$ is the unit vector in the radial direction and $g$ is gravity. Several gravity profiles have been classically considered to model convection in spherical shells. For instance, self-gravitating spherical shells with a constant density correspond to $g\sim r$ (e.g. Tilgner Reference Tilgner1996), while RB convection models with infinite Prandtl number usually assume a constant gravity in the perspective of modelling Earth’s mantle (e.g. Bercovici et al. Reference Bercovici, Schubert, Glatzmaier and Zebib1989). The assumption of a centrally condensed mass has also been frequently assumed when modelling rotating convection (e.g. Gilman & Glatzmaier Reference Gilman and Glatzmaier1981; Jones et al. Reference Jones, Boronski, Brun, Glatzmaier, Gastine, Miesch and Wicht2011) and yields $g\sim 1/r^{2}$ . Finally, the artificial central force field of the microgravity experiments effectively takes the form of $g\sim 1/r^{5}$ (Hart et al. Reference Hart, Glatzmaier and Toomre1986; Feudel et al. Reference Feudel, Bergemann, Tuckerman, Egbers, Futterer, Gellert and Hollerbach2011; Futterer et al. Reference Futterer, Krebs, Plesa, Zaussinger, Hollerbach, Breuer and Egbers2013). To explore the possible impact of these various radial distributions of buoyancy on RB convection in spherical shells, we consider different models with the four following gravity profiles: $g\in [r/r_{o},\,1,\,(r_{o}/r)^{2},\,(r_{o}/r)^{5}]$ . Particular attention will be paid to the cases with $g=(r_{o}/r)^{2}$ , which is the only radial function compatible with an exact analysis of the dissipation rates (see below, § 2.3).

The dimensionless set of (2.1)–(2.3) is governed by the Rayleigh number $\mathit{Ra}$ , the Prandtl number $\mathit{Pr}$ and the radius ratio of the spherical shell ${\it\eta}$ defined by

where ${\it\nu}$ and ${\it\kappa}$ are the viscous and thermal diffusivities and ${\it\alpha}$ is the thermal expansivity.

2.2. Diagnostic parameters

To quantify the impact of the different control parameters on the transport of heat and momentum, we analyse several diagnostic properties. We adopt the following notations regarding different averaging procedures. Overbars $\overline{\cdots \,}$ correspond to a time average

where ${\it\tau}$ is the time averaging interval. Spatial averaging over the whole volume of the spherical shell is denoted by triangular brackets $\langle \cdots \!\,\rangle$ , while $\langle \cdots \!\,\rangle _{s}$ correspond to an average over a spherical surface:

where $V$ is the volume of the spherical shell, $r$ is the radius, ${\it\theta}$ the colatitude and ${\it\phi}$ the longitude.

The convective heat transport is characterised by the Nusselt number $\mathit{Nu}$ , the ratio of the total heat flux to the heat carried by conduction. In spherical shells with isothermal boundaries, the conductive temperature profile $T_{c}$ is the solution of

which yields

For the sake of clarity, we adopt in the following the notation ${\it\vartheta}$ for the time-averaged and horizontally averaged radial temperature profile:

The Nusselt number then reads

The typical root-mean-square (r.m.s.) flow velocity is given by the Reynolds number

while the radial profile for the time and horizontally averaged horizontal velocity is defined by

2.3. Exact dissipation relationships in spherical shells

The mean buoyancy power averaged over the whole volume of a spherical shell is expressed by

Using the Nusselt number definition (2.10) and the conductive temperature profile (2.8) then leads to

The first term in the parentheses becomes identical to the imposed temperature drop ${\rm\Delta}T$ for a gravity $g\sim 1/r^{2}$ :

and thus yields an analytical relation between $P$ and $\mathit{Nu}$ . For any other gravity model, we have to consider the actual spherically symmetric radial temperature profile ${\it\vartheta}(r)$ . Christensen & Aubert (Reference Christensen and Aubert2006) solve this problem by approximating ${\it\vartheta}(r)$ by the diffusive solution (2.8) and obtain an approximate relation between $P$ and $(Ra/\mathit{Pr}^{2})(\mathit{Nu}-1)$ . This motivates our particular focus on the $g=(r_{o}/r)^{2}$ cases which allows us to conduct an exact analysis of the dissipation rates and therefore check the applicability of the GL theory to convection in spherical shells.

Noting that $({\it\eta}/(1-{\it\eta})^{2})\int _{r_{i}}^{r_{o}}g\,\text{d}r=-1/(1-{\it\eta})^{2}$ for $g=(r_{o}/r)^{2}$ , one finally obtains the exact relation for the viscous dissipation rate ${\it\epsilon}_{U}$ :

The thermal dissipation rate can be obtained by multiplying the temperature equation (2.3) by $T$ and integrating over the whole volume of the spherical shell. This yields

These two exact relations (2.16) and (2.17) can be used to validate the spatial resolutions of the numerical models with $g=(r_{o}/r)^{2}$ . To do so, we introduce ${\it\chi}_{{\it\epsilon}_{U}}$ and ${\it\chi}_{{\it\epsilon}_{T}}$ , the ratios of the two sides of (2.16) and (2.17):

2.4. Setting up a parameter study

3.1. Definitions

Several different approaches are traditionally considered to define the thermal boundary layer thickness ${\it\lambda}_{T}$ . They either rely on the horizontally averaged mean radial temperature profile ${\it\vartheta}(r)$ or on the temperature fluctuation ${\it\sigma}$ defined as

Among the possible estimates based on ${\it\vartheta}(r)$ , the slope method (e.g. Verzicco & Camussi Reference Verzicco and Camussi1999; Breuer et al. Reference Breuer, Wessling, Schmalzl and Hansen2004; Liu & Ecke Reference Liu and Ecke2011) defines ${\it\lambda}_{T}$ as the depth where the linear fit to ${\it\vartheta}(r)$ near the boundaries intersects the linear fit to the temperature profile at mid-depth. Alternatively, ${\it\sigma}$ exhibits sharp local maxima close to the walls. The radial distance separating those peaks from the corresponding nearest boundary can be used to define the thermal boundary layer thicknesses (e.g. Tilgner Reference Tilgner1996; King, Stellmach & Buffett Reference King, Stellmach and Buffett2013). Figure 2(a) shows that both definitions of ${\it\lambda}_{T}$ actually yield nearly indistinguishable boundary layer thicknesses. We therefore adopt the slope method to define the thermal boundary layers.

Figure 2. (a) Radial profiles of the time and horizontally averaged mean temperature ${\it\vartheta}(r)$ (solid black line) and the temperature variance ${\it\sigma}$ (dotted black line). The thermal boundary layers ${\it\lambda}_{T}^{i}$ and ${\it\lambda}_{T}^{o}$ are highlighted by the grey shaded area. They are defined as the depths where the linear fit to ${\it\vartheta}(r)$ near the top (bottom) crosses the linear fit to the temperature profile at mid-depth (dashed black lines). (b) Radial profiles of the time and horizontally averaged horizontal velocity $\mathit{Re}_{h}(r)$ . The viscous boundary layers are either defined by the local maxima of $\mathit{Re}_{h}$ (black dotted lines, ${\it\lambda}_{U,m}^{i}$ , ${\it\lambda}_{U,m}^{o}$ ) or by the intersection of the linear fit to $\mathit{Re}_{h}$ near the inner (outer) boundary with the tangent to the local maxima (dashed black lines). This second definition is highlighted by a grey shaded area ( ${\it\lambda}_{U}^{i}$ , ${\it\lambda}_{U}^{o}$ ). These profiles have been obtained from a numerical model with $\mathit{Ra}=10^{7}$ , ${\it\eta}=0.6$ and $g=(r_{o}/r)^{2}$ .

There are also several ways to define the viscous boundary layers. Figure 2(b) shows the vertical profile of the root-mean-square horizontal velocity $\mathit{Re}_{h}$ . This profile exhibits strong increases close to the boundaries that are accompanied by well-defined peaks. Following Tilgner (Reference Tilgner1996) and Kerr & Herring (Reference Kerr and Herring2000), the first way to define the kinematic boundary layer is thus to measure the distance between the walls and these local maxima. This commonly used definition gives ${\it\lambda}_{U,m}^{i}$ ( ${\it\lambda}_{U,m}^{o}$ ) for the inner (outer) spherical boundary. Another possible method to estimate the viscous boundary layer follows a similar strategy as the slope method that we adopted for the thermal boundary layers (Breuer et al. Reference Breuer, Wessling, Schmalzl and Hansen2004). ${\it\lambda}_{U}^{i}$ ( ${\it\lambda}_{U}^{o}$ ) is defined as the distance from the inner (outer) wall where the linear fit to $\mathit{Re}_{h}$ near the inner (outer) boundary intersects the horizontal line passing through the maximum horizontal velocity.

Figure 3. (a) ${\it\vartheta}(r)$ for different radius ratios ${\it\eta}$ with a gravity $g=(r_{o}/r)^{2}$ . These models have approximately the same Nusselt number $12<Nu<14$ . (b) ${\it\vartheta}(r)$ for different gravity profiles with a fixed radius ratio ${\it\eta}=0.6$ . These models have approximately the same Nusselt number $10<Nu<11$ .

Figure 2(b) reveals that these two definitions lead to very distinct viscous boundary layer thicknesses. In particular, the definition based on the position of the local maxima of $\mathit{Re}_{h}$ yields much thicker boundary layers than the tangent intersection method, i.e. ${\it\lambda}_{U,m}^{i,o}>{\it\lambda}_{U}^{i,o}$ . The discrepancies of these two definitions are further discussed in § 4.

3.2. Asymmetric thermal boundary layers and mean bulk temperature

Figure 2 also reveals a pronounced asymmetry in the mean temperature profiles with a much larger temperature drop at the inner boundary than at the outer boundary. As a consequence, the mean temperature of the spherical shell $T_{m}=(1/V)\int _{V}T\,\text{d}V$ is far below ${\rm\Delta}T/2$ . Determining how the mean temperature depends on the radius ratio ${\it\eta}$ has been an ongoing open question in mantle convection studies with infinite Prandtl number (e.g. Bercovici et al. Reference Bercovici, Schubert, Glatzmaier and Zebib1989; Jarvis Reference Jarvis1993; Vangelov & Jarvis Reference Vangelov and Jarvis1994; Jarvis et al. Reference Jarvis, Glatzmaier and Vangelov1995; Sotin & Labrosse Reference Sotin and Labrosse1999; Shahnas et al. Reference Shahnas, Lowman, Jarvis and Bunge2008; Deschamps et al. Reference Deschamps, Tackley and Nakagawa2010; O’Farrell et al. Reference O’Farrell, Lowman and Bunge2013). To analyse this issue in numerical models with $\mathit{Pr}=1$ , we have performed a systematic parameter study varying both the radius ratio of the spherical shell ${\it\eta}$ and the gravity profile $g(r)$ (see table 2). Figure 3 shows some selected radial profiles of the mean temperature ${\it\vartheta}$ for various radius ratios (a) and gravity profiles (b) for cases with similar $\mathit{Nu}$ . For small values of ${\it\eta}$ , the large difference between the inner and the outer surfaces lead to a strong asymmetry in the temperature distribution: nearly 90 % of the total temperature drop occurs at the inner boundary when ${\it\eta}=0.2$ . In thinner spherical shells, the mean temperature gradually approaches a more symmetric distribution to finally reach $T_{m}=0.5$ when ${\it\eta}\rightarrow 1$ (no curvature). Figure 3(b) also reveals that a change in the gravity profile has a direct impact on the mean temperature profile. This shows that both the shell geometry and the radial distribution of buoyancy affect the temperature of the fluid bulk in RB convection in spherical shells.

To analytically access the asymmetries in thickness and temperature drop observed in figure 3, we first assume that heat is purely transported by conduction in the thin thermal boundary layers. Heat flux conservation through spherical surfaces (2.10) then yields

where the thermal boundary layers are assumed to correspond to a linear conduction profile with a temperature drop ${\rm\Delta}T^{i}$ ( ${\rm\Delta}T^{o}$ ) over a thickness ${\it\lambda}_{T}^{i}$ ( ${\it\lambda}_{T}^{o}$ ). As shown in figures 2 and 3, the fluid bulk is isothermal and forms the majority of the fluid by volume. We can thus further assume that the temperature drops occur only in the thin boundary layers, which leads to

Equations (3.2) and (3.3) are nevertheless not sufficient to determine the three unknowns ${\rm\Delta}T^{i}$ , ${\rm\Delta}T^{o}$ , ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ and an additional physical assumption is required.

A hypothesis frequently used in mantle convection models with infinite Prandtl number in a spherical geometry (Jarvis Reference Jarvis1993; Vangelov & Jarvis Reference Vangelov and Jarvis1994) is to further assume that both thermal boundary layers are marginally stable such that the local boundary layer Rayleigh numbers $\mathit{Ra}_{{\it\lambda}}^{i}$ and $\mathit{Ra}_{{\it\lambda}}^{o}$ are equal:

This means that both thermal boundary layers adjust their thickness and temperature drop to yield $\mathit{Ra}_{{\it\lambda}}^{i}\sim \mathit{Ra}_{{\it\lambda}}^{o}\sim \mathit{Ra}_{c}\simeq 1000$ (e.g. Malkus Reference Malkus1954). The temperature drops at both boundaries and the ratio of the thermal boundary layer thicknesses can then be derived using (3.2) and (3.3)

where

is the ratio of the gravitational acceleration between the inner and the outer boundaries. Figure 4(a) reveals that the marginal stability hypothesis is not fulfilled when different radius ratios and gravity profiles are considered. This is particularly obvious for small radius ratios where $\mathit{Ra}_{{\it\lambda}}^{o}$ is more than 10 times larger than $\mathit{Ra}_{{\it\lambda}}^{i}$ . This discrepancy tends to vanish when ${\it\eta}\rightarrow 1$ , where curvature and gravity variations become unimportant. As a consequence, there is a significant mismatch between the predicted mean bulk temperature from (3.5) and the actual values (figure 4 b). Deschamps et al. (Reference Deschamps, Tackley and Nakagawa2010) also reported a similar deviation from (3.5) in their spherical shell models with infinite Prandtl number. They suggest that assuming instead $\mathit{Ra}_{{\it\lambda}}^{o}/\mathit{Ra}_{{\it\lambda}}^{i}\sim {\it\eta}^{2}$ might help to improve the agreement with the data. This however cannot account for the additional dependence on the gravity profile visible in figure 4. We finally note that $\mathit{Ra}_{{\it\lambda}}<400$ for the database of numerical simulations explored here, which suggests that the thermal boundary layers are stable in all our simulations.

Figure 4. (a) Ratio of boundary layer Rayleigh numbers (3.4) for various radius ratios and gravity profiles. The horizontal dashed line corresponds to the hypothetical identity $\mathit{Ra}_{{\it\lambda}}^{i}=\mathit{Ra}_{{\it\lambda}}^{o}$ . (b) Temperature drop at the outer boundary layer. The lines correspond to the theoretical prediction given in (3.5).

Alternatively Wu & Libchaber (Reference Wu and Libchaber1991) and Zhang, Childress & Libchaber (Reference Zhang, Childress and Libchaber1997) proposed that the thermal boundary layers adapt their thicknesses such that the mean hot and cold temperature fluctuations at mid-depth are equal. Their experiments with Helium indeed revealed that the statistical distribution of the temperature at mid-depth was symmetrical. They further assumed that the thermal fluctuations at mid-depth can be identified with the boundary layer temperature scales ${\it\theta}^{i}=({\it\nu}{\it\kappa})/({\it\alpha}g_{i}{{\it\lambda}_{T}^{i}}^{3})$ and ${\it\theta}^{o}=({\it\nu}{\it\kappa})/({\it\alpha}g_{o}{{\it\lambda}_{T}^{o}}^{3})$ , which characterise the temperature scale of the thermal boundary layers in a different way than the relative temperature drops ${\rm\Delta}T^{i}$ and ${\rm\Delta}T^{o}$ . This second hypothesis yields

and the corresponding temperature drops and boundary layer thicknesses ratio

Figure 5(a) shows ${\it\theta}^{o}/{\it\theta}^{i}$ for different radius ratios and gravity profiles, while figure 5(b) shows a comparison between the predicted mean bulk temperature and the actual values. Besides the cases with $g=(r_{o}/r)^{2}$ , which are in relatively good agreement with the predicted scalings, the identity of the boundary layer temperature scale is in general not fulfilled for the other gravity profiles. The actual mean bulk temperature is thus poorly described by (3.8). We note that previous findings by Ahlers et al. (Reference Ahlers, Brown, Fontenele Araujo, Funfschilling, Grossmann and Lohse2006) already reported that the theory by Wu & Libchaber’s does also not hold when the transport properties depend on temperature (i.e. non-Oberbeck–Boussinesq convection).

Figure 5. (a) Ratio of boundary layer temperature scales (3.7) for various radius ratios and gravity profiles. The horizontal dashed line corresponds to the hypothetical identity ${\it\theta}^{i}={\it\theta}^{o}$ . (b) Temperature drop at the outer boundary layer. The lines correspond to the theoretical prediction given in (3.8).

3.3. Conservation of the average plume density in spherical shells

As demonstrated in the previous section, none of the hypotheses classically employed accurately account for the temperature drops and the boundary layer asymmetry observed in spherical shells. We must therefore find a dynamical quantity that could possibly be identified between the two boundary layers.

Figure 6 shows visualisations of the thermal boundary layers for three selected numerical models with different radius ratios and gravity profiles. The isocontours displayed in (a–c) reveal the intricate plume structure. Long and thin sheet-like structures form the main network of plumes. During their migration along the spherical surfaces, these sheet-like plumes can collide and convolute with each other to give rise to mushroom-type plumes (see Zhou & Xia Reference Zhou and Xia2010b ; Chillà & Schumacher Reference Chillà and Schumacher2012). During this morphological evolution, mushroom-type plumes acquire a strong radial vorticity component. These mushroom-type plumes are particularly visible at the connection points of the sheet plumes network at the inner thermal boundary layer (red isosurface in figure 6 a–c). Figure 6(d–f) shows the corresponding equatorial and radial cuts of the temperature fluctuation $T^{\prime }=T-{\it\vartheta}$ . These panels further highlight the plume asymmetry between the inner and the outer thermal boundary layers. For instance, the case with ${\it\eta}=0.3$ and $g=(r_{o}/r)^{5}$ (a,d) features an outer boundary layer approximately 4.5 times thicker than the inner one. Accordingly, the mushroom-type plumes that depart from the outer boundary layer are significantly thicker than the ones emitted from the inner boundary. This discrepancy tends to vanish in the thin shell case ( ${\it\eta}=0.8$ , c,f) in which curvature and gravity variations play a less significant role ( ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}\simeq 1.02$ in that case).

Figure 6. (a–c) Isosurfaces and equatorial cut of the temperature for three numerical models: hot at the inner thermal boundary layer $T(r=r_{i}+{\it\lambda}_{T}^{i})$ in red, cold at the outer thermal boundary layer $T(r=r_{o}-{\it\lambda}_{T}^{o})$ in blue. (d–f) Meridional cuts and surfaces of the temperature fluctuations $T^{\prime }$ . The inner (outer) surface corresponds to the location of the inner (outer) thermal boundary layers. Color levels range from $-$ 0.2 (black) to 0.2 (white). (a,d) Correspond to a model with $\mathit{Ra}=3\times 10^{6},\,{\it\eta}=0.3$ and $g=(r_{o}/r)^{5}$ . (b,e) Correspond to a model with $\mathit{Ra}=10^{8},\,{\it\eta}=0.6$ and $g=(r_{o}/r)^{2}$ . (c,f) Correspond to a model with $\mathit{Ra}=4\times 10^{7},\,{\it\eta}=0.8$ and $g=r/r_{o}$ .

Puthenveettil & Arakeri (Reference Puthenveettil and Arakeri2005) and Zhou & Xia (Reference Zhou and Xia2010b ) performed a statistical analysis of the geometrical properties of thermal plumes in experimental RB convection. By tracking a large number of plumes, their analysis revealed that both the plume separation and the width of the sheet-like plumes follow a log-normal probability density function (PDF).

To further assess how the average plume properties of the inner and outer thermal boundary layers compare with each other in a spherical geometry, we adopt a simpler strategy by only focusing on the statistics of the plume density. The plume density per surface unit at a given radius $r$ is expressed by

where $N$ is the number of plumes, approximated here by the ratio of the spherical surface area to the mean inter-plume area $\bar{\mathscr{S}}$ :

This inter-plume area $\bar{\mathscr{S}}$ can be further related to the average plume separation $\bar{\ell }$ via $\bar{\mathscr{S}}\sim ({\rm\pi}/4)\,\bar{\ell }^{2}$ .

An accurate evaluation of the inter-plume area for each thermal boundary layer however requires the extraction of the plumes from the background fluid. Most of the criteria employed to determine the location of the plume boundaries are based on thresholds of certain physical quantities (see Shishkina & Wagner Reference Shishkina and Wagner2008, for a review of the different plume extraction techniques). This encompasses threshold values on the temperature fluctuations $T^{\prime }$ (Zhou & Xia Reference Zhou and Xia2002), on the vertical velocity $u_{r}$ (Ching et al. Reference Ching, Guo, Shang, Tong and Xia2004) or on the thermal dissipation rate ${\it\epsilon}_{T}$ (Shishkina & Wagner Reference Shishkina and Wagner2005). The choice of the threshold value however remains an open problem. Alternatively, Vipin & Puthenveettil (Reference Vipin and Puthenveettil2013) show that the sign of the horizontal divergence of the velocity $\boldsymbol{{\rm\nabla}}_{H}\boldsymbol{\cdot }\boldsymbol{u}$ might provide a simple and threshold-free criterion to separate the plumes from the background fluid

Fluid regions with $\boldsymbol{{\rm\nabla}}_{H}\boldsymbol{\cdot }\boldsymbol{u}<0$ indeed correspond to local regions of positive vertical acceleration, expected inside the plumes, while the fluid regions with $\boldsymbol{{\rm\nabla}}_{H}\boldsymbol{\cdot }\boldsymbol{u}>0$ characterise the inter-plume area.

To analyse the statistics of $\mathscr{S}$ , we thus consider here several criteria based either on a threshold value of the temperature fluctuations or on the sign of the horizontal divergence. This means that a given inter-plume area at the inner (outer) thermal boundary layer is either defined as an enclosed region surrounded by hot (cold) sheet-like plumes carrying a temperature perturbation $|T^{\prime }|>t$ ; or by an enclosed region with $\boldsymbol{{\rm\nabla}}_{H}\boldsymbol{\cdot }\boldsymbol{u}>0$ . To further estimate the possible impact of the chosen threshold value on $\mathscr{S}$ , we vary $t$ between ${\it\sigma}/4$ and ${\it\sigma}$ . This yields

where the physical criterion $\mathscr{T}_{i}$ ( $\mathscr{T}_{o}$ ) to extract the plume boundaries at the inner (outer) boundary layer is given by

where $r_{{\it\lambda}}^{i}=r_{i}+{\it\lambda}_{T}^{i}$ ( $r_{{\it\lambda}}^{o}=r_{o}-{\it\lambda}_{T}^{o}$ ) for the inner (outer) thermal boundary layer.

Figure 7 shows an example of such a characterisation procedure for the inner thermal boundary layer of a numerical model with $\mathit{Ra}=10^{6}$ , ${\it\eta}=0.6$ , $g=(r_{o}/r)^{2}$ . (b) Illustrates a plume extraction process when using $|T^{\prime }|>{\it\sigma}/2$ to determine the location of the plumes: the black area correspond to the inter-plume spacing while the white area correspond to the complementary plume network location. The fainter emerging sheet-like plumes are filtered out and only the remaining ‘skeleton’ of the plume network is selected by this extraction process. The choice of ${\it\sigma}/2$ is however arbitrary and can influence the evaluation of the number of plumes. The insets displayed in (c–e) illustrate the sensitivity of the plume extraction process on the criterion employed to detect the plumes. In particular, using the threshold based on the largest temperature fluctuations $|T^{\prime }|>{\it\sigma}$ can lead to the fragmentation of the detected plume lanes into several isolated smaller regions. As a consequence, several neighbouring inter-plume areas may be artificially connected when using this criterion. In contrast, using the sign of the horizontal divergence to estimate the plumes location yields much broader sheet-like plumes. As visible in (e), the plume boundaries frequently correspond to local maxima of the thermal dissipation rate ${\it\epsilon}_{T}$ (Shishkina & Wagner Reference Shishkina and Wagner2008).

Figure 7. (a) Temperature fluctuation at the inner thermal boundary layer $T^{\prime }(r=r_{i}+{\it\lambda}_{T}^{i})$ displayed in a Hammer projection for a case with $\mathit{Ra}=10^{6}$ , ${\it\eta}=0.6$ , $g=(r_{o}/r)^{2}$ . (b) Corresponding binarised extraction of the plumes boundaries using (3.12) and $T^{\prime }\leqslant {\it\sigma}/2$ to define the inter-plume area. Zoomed-in contour of the temperature fluctuation $T^{\prime }$ (c), the horizontal divergence $-\boldsymbol{{\rm\nabla}}_{H}\boldsymbol{\cdot }\boldsymbol{u}$ (d) and the thermal dissipation rate ${\it\epsilon}_{T}$ (e). The three black contour lines in (c–e) correspond to several criteria to extract the plume boundaries (3.13).

For each criterion given in (3.13), we then calculate the area of each bounded black surface visible in figure 7(b) to construct the statistical distribution of the inter-plume area for both thermal boundary layers. Figure 8 compares the resulting PDFs obtained by combining several snapshots for a numerical model with $\mathit{Ra}=4\times 10^{7}$ , ${\it\eta}=0.8$ and $g=r/r_{o}$ . Besides the criterion $|T^{\prime }|>{\it\sigma}$ which yields PDFs that are slightly shifted towards smaller inter-plume spacing areas, the statistical distributions are found to be relatively insensitive to the detection criterion (3.13). We therefore restrict the following comparison to the criterion $|T^{\prime }|>{\it\sigma}/2$ only.

Figure 8. Probability density functions (PDFs) of the dimensionless inter-plume area $\mathscr{S}$ at the inner $r=r_{i}+{\it\lambda}_{T}^{i}$ (a) and at the outer $r=r_{o}-{\it\lambda}_{T}^{o}$ (b) thermal boundary layers using different criteria to extract the plumes (3.13) for a model with $\mathit{Ra}=4\times 10^{7}$ , ${\it\eta}=0.8$ and $g=r/r_{o}$ .

Figure 9 shows the PDFs for the three numerical models of figure 6. For the two cases with ${\it\eta}=0.6$ and ${\it\eta}=0.8$ (b,c), the statistical distributions for both thermal boundary layers nearly overlap. This means that the inter-plume area is similar at both spherical shell surfaces. In contrast, for the case with ${\it\eta}=0.3$ (a), the two PDFs are offset relative to each other. However, the peaks of the distributions remain relatively close, meaning that once again the inner and the outer thermal boundary layers share a similar average inter-plume area. Puthenveettil & Arakeri (Reference Puthenveettil and Arakeri2005) and Zhou & Xia (Reference Zhou and Xia2010b ) demonstrated that the thermal plume statistics in turbulent RB convection follow a log-normal distribution (see also Shishkina & Wagner Reference Shishkina and Wagner2008; Puthenveettil et al. Reference Puthenveettil, Gunasegarane, Agrawal, Schmeling, Bosbach and Arakeri2011). The large number of plumes in the cases with ${\it\eta}=0.6$ and ${\it\eta}=0.8$ (figure 6 b,c) would allow a characterisation of the nature of the statistical distributions. However, this would be much more difficult in the ${\it\eta}=0.3$ case (figure 6 a) in which the plume density is significantly weaker. As a consequence, no further attempt has been made to characterise the exact nature of the PDFs visible in figure 9, although the universality of the log-normal statistics reported by Puthenveettil & Arakeri (Reference Puthenveettil and Arakeri2005) and Zhou & Xia (Reference Zhou and Xia2010b ) likely indicates that the same statistical distribution should hold here too.

Figure 9. PDF of the dimensionless inter-plume area $\mathscr{S}$ at the outer (orange upward triangles) and at the inner (blue downward triangles) thermal boundary layers using $|T^{\prime }|>{\it\sigma}/2$ to extract plumes. (a) Corresponds to a numerical model with $\mathit{Ra}=3\times 10^{6},\,{\it\eta}=0.3$ and $g=(r_{o}/r)^{5}$ . (b) Corresponds to a model with $\mathit{Ra}=10^{8},\,{\it\eta}=0.6$ and $g=(r_{o}/r)^{2}$ . (c) Corresponds to a model with $\mathit{Ra}=4\times 10^{7},\,{\it\eta}=0.8$ and $g=r/r_{o}$ . The two vertical lines correspond to the predicted values of the mean inter-plume area $\bar{\mathscr{S}}$ derived from (3.16) for both thermal boundary layers.

The inter-plume area statistics therefore reveal that the inner and the outer thermal boundary layers exhibit a similar average plume density, independently of the spherical shell geometry and the gravity profile. Assuming ${\it\rho}_{p}^{o}\simeq {\it\rho}_{p}^{i}$ would allow us to close the system of (3.2)–(3.3) and thus finally estimate ${\rm\Delta}T^{i}$ , ${\rm\Delta}T^{o}$ and ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ . This however requires us to determine an analytical expression of the average inter-plume area $\bar{\mathscr{S}}$ or equivalently of the mean plume separation $\bar{\ell }$ that depends on the boundary layer thickness and the temperature drop.

Using the boundary layer equations for natural convection (Rotem & Claassen Reference Rotem and Claassen1969), Puthenveettil et al. (Reference Puthenveettil, Gunasegarane, Agrawal, Schmeling, Bosbach and Arakeri2011) demonstrated that the thermal boundary layer thickness follows

where $x$ is the distance along the horizontal direction and $\mathit{Ra}_{x}^{i,o}={\it\alpha}g{\rm\Delta}T^{i,o}x^{3}/{\it\nu}{\it\kappa}$ is a Rayleigh number based on the length scale $x$ and on the boundary layer temperature jumps ${\rm\Delta}T^{i,o}$ . As shown on figure 10, using $x=\bar{\ell }/2$ (Puthenveettil & Arakeri Reference Puthenveettil and Arakeri2005; Puthenveettil et al. Reference Puthenveettil, Gunasegarane, Agrawal, Schmeling, Bosbach and Arakeri2011) then allows the following relation for the average plume spacing to be established

which yields

for both thermal boundary layers. We note that an equivalent expression for the average plume spacing can be derived from a simple mechanical description of the equilibrium between production and coalescence of plumes in each boundary layer (see Parmentier & Sotin Reference Parmentier and Sotin2000; King et al. Reference King, Stellmach and Buffett2013).

Figure 10. Schematic showing two adjacent plumes separated by a distance $\bar{\ell }$ . The thick black arrows indicate the direction of merging of the two plumes.

Equation (3.16) is however expected to be valid only at the scaling level. The vertical lines in figure 9 therefore correspond to the estimated average inter-plume area for both thermal boundary layers using (3.16) and $\bar{\mathscr{S}}_{i,o}=0.3\,\bar{\ell }_{i,o}^{2}$ . The predicted average inter-plume area is in good agreement with the peaks of the statistical distributions for the three cases discussed here. The expression (3.16) therefore provides a reasonable estimate of the average plume separation (Puthenveettil & Arakeri Reference Puthenveettil and Arakeri2005; Puthenveettil et al. Reference Puthenveettil, Gunasegarane, Agrawal, Schmeling, Bosbach and Arakeri2011; Gunasegarane & Puthenveettil Reference Gunasegarane and Puthenveettil2014). The comparable observed plume density at both thermal boundary layers thus yields

Using (3.2) and (3.3) then allows us to finally estimate the temperature jumps and the ratio of the thermal boundary layer thicknesses in our dimensionless units:

Figure 11 shows the ratios $\bar{\ell }_{o}/\bar{\ell }_{i}$ , ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ and the temperature jumps ${\rm\Delta}T^{i}$ and ${\rm\Delta}T^{o}$ . In contrast to the previous criteria, either coming from the marginal stability of the boundary layer ((3.5), figure 4) or from the identity of the temperature fluctuations at mid-shell ((3.18), figure 5), the ratio of the average plume separation $\bar{\ell }_{o}/\bar{\ell }_{i}$ now falls much closer to the unity line. Some deviations are nevertheless still visible for spherical shells with ${\it\eta}\leqslant 0.4$ and $g=r/r_{o}$ (orange circles). The comparable average plume density between both boundary layers allows us to accurately predict the asymmetry of the thermal boundary layers ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ and the corresponding temperature drops for the vast majority of the numerical cases explored here (solid lines in b–d).

Figure 11. (a) Ratio of the thermal plume separation estimated by (3.16) for various radius ratios and gravity profiles. The horizontal dashed line corresponds to the identity of the average plume separation between both thermal boundary layers, i.e.  $\bar{\ell }^{i}=\bar{\ell }^{o}$ . (b) Ratio of thermal boundary layer thicknesses. (c) Temperature drop at the inner boundary layer. (d) Temperature drop at the outer boundary layer. The lines in (b–d) correspond to the theoretical prediction given in (3.18).

As we consider a fluid with $\mathit{Pr}=1$ , the viscous boundary layers should show a comparable degree of asymmetry as the thermal boundary layers. Equation (3.18) thus implies

Figure 12 shows the ratio of the viscous boundary layer thicknesses for the different setups explored in this study. The observed asymmetry between the two spherical shell surfaces is in a good agreement with (3.19) (solid black lines).

Figure 12. Ratio of the viscous boundary layer thicknesses for various aspect ratios and gravity profiles. The lines correspond to the theoretical prediction given in (3.19).

3.4. Thermal boundary layer scalings

Using (3.18) and the definition of the Nusselt number (2.10), we can derive the following scaling relations for the thermal boundary layer thicknesses:

Figure 13(a) demonstrates that the boundary layer thicknesses for the numerical simulations of table 1 ( $g=(r_{o}/r)^{2}$ and ${\it\eta}=0.6$ ) are indeed in close agreement with the theoretical predictions. To further check this scaling for other spherical shell configurations, we introduce the following normalisation of the thermal boundary layer thicknesses

This allows us to derive a unified scaling that does not depend on the choice of the gravity profile or on the spherical shell geometry

Figure 13(b) shows this normalised boundary layer thickness for the different spherical shell configurations explored here. Despite the variety of the physical set-ups, the normalised boundary layer thicknesses are in good agreement with the predicted behaviour.

Figure 13. (a) Thermal boundary layer thicknesses at the outer boundary ( ${\it\lambda}_{T}^{i}$ ) and at the inner boundary ( ${\it\lambda}_{T}^{o}$ ) as a function of the Nusselt number for the cases of table 1. The two lines correspond to the theoretical predictions from (3.20). (b) Normalised boundary layer thicknesses as a function of the Nusselt number for different radius ratios and gravity profiles. For the sake of clarity, the outer boundary layer is only displayed for the cases with $g=(r_{o}/r)^{2}$ . The solid line corresponds to $\tilde{{\it\lambda}}_{T}=0.5\,\mathit{Nu}^{-1}$ (3.22).

Figure 14. Ratio of the thermal boundary layer thicknesses ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ as a function of the Rayleigh number for the numerical models of table 1 with ${\it\eta}=0.6$ and $g=(r_{o}/r)^{2}$ . The horizontal dashed line corresponds to the predicted ratio ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ given in (3.18).

A closer inspection of figure 11(b) and table 1 nevertheless reveals a remaining weak dependence of the ratio ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ on the Rayleigh number. This is illustrated in figure 14 which shows ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ as a function of $\mathit{Ra}$ for the ${\it\eta}=0.6$ , $g=(r_{o}/r)^{2}$ cases of table 1. Although some complex variations are visible, the first-order trend is a very slow increase of ${\it\lambda}_{T}^{o}/{\it\lambda}_{T}^{i}$ with $\mathit{Ra}$ (10 % increase over five decades of $\mathit{Ra}$ ). No evidence of saturation is however visible and further deviations from the predicted ratio (horizontal dashed line) might therefore be expected at larger $\mathit{Ra}$ . This variation might cast some doubts on the validity of the previous derivation for higher Rayleigh numbers. This may imply that either the plume separation is not conserved at higher $\mathit{Ra}$ ; or that the estimate of the average plume spacing is too simplistic to capture the detailed plume physics in turbulent convection.

5.1. Bulk and boundary layer contributions to viscous and thermal dissipation rates

The prerequisite of a laminar boundary layer appears to be fulfilled in our numerical models and we can thus attempt to apply the GL formalism to our dataset. The idea of the GL theory is to separate the viscous and thermal dissipation rates into two contributions, one coming from the fluid bulk (indicated by the superscript $bu$ in the following) and one coming from the boundary layers ( $bl$ ), such that

where the contributions from the bulk are defined by

and the boundary layer contributions are given by

The RB flows are then classified in the $\mathit{Ra}{-}\mathit{Pr}$ parameter space according to the dominant contributions to the viscous and thermal dissipation rates. This defines four regimes depending on $\mathit{Ra}$ and $\mathit{Pr}$ : regime I when ${\it\epsilon}_{U}\simeq {\it\epsilon}_{U}^{bl}$ and ${\it\epsilon}_{T}\simeq {\it\epsilon}_{T}^{bl}$ ; regime II when ${\it\epsilon}_{U}\simeq {\it\epsilon}_{U}^{bu}$ and ${\it\epsilon}_{T}\simeq {\it\epsilon}_{T}^{bl}$ ; regime III when ${\it\epsilon}_{U}\simeq {\it\epsilon}_{U}^{bl}$ and ${\it\epsilon}_{T}\simeq {\it\epsilon}_{T}^{bu}$ ; and finally regime IV when ${\it\epsilon}_{U}\simeq {\it\epsilon}_{U}^{bu}$ and ${\it\epsilon}_{T}\simeq {\it\epsilon}_{T}^{bu}$ .

For a unity Prandtl number, the GL theory predicts that the flows should be dominated by dissipations in the boundary layer regions at low $\mathit{Ra}$ (regime I). A transition to another regime where dissipations in the fluid bulk dominates (regime IV) is expected to happen at approximately $\mathit{Ra}\simeq 10^{8}{-}10^{10}$ (Grossmann & Lohse Reference Grossmann and Lohse2000; Ahlers et al. Reference Ahlers, Grossmann and Lohse2009; Stevens et al. Reference Stevens, van der Poel, Grossmann and Lohse2013).

Figure 18 shows the relative contributions of the bulk and boundary layers to the viscous and thermal dissipation rates. The viscous dissipation ${\it\epsilon}_{U}$ is always dominated by the bulk contribution: starting at approximately 60 % at $\mathit{Ra}=10^{4}$ , it nearly reaches 90 % at $\mathit{Ra}=10^{9}$ . In contrast, the thermal dissipation rate is always dominated by the boundary layer regions, such that ${\it\epsilon}_{T}^{bu}$ slowly increases from 10 % at $\mathit{Ra}=10^{4}$ to 30 % at $\mathit{Ra}=10^{9}$ . According to the GL classification, all of our numerical simulations thus belong to regime II of the $\mathit{Ra}{-}\mathit{Pr}$ parameter space, in which ${\it\epsilon}_{U}^{bu}>{\it\epsilon}_{U}^{bl}$ and ${\it\epsilon}_{T}^{bl}>{\it\epsilon}_{T}^{bu}$ . This seems at odds with the GL theory, which predicts that our DNS should be located either in regime I or in regime IV of the parameter space for the range of $\mathit{Ra}$ explored here ( $10^{3}\leqslant Ra\leqslant 10^{9}$ ).

Figure 18. Measured contributions of the boundary layer (open symbols) and the fluid bulk (filled symbols) to the total viscous dissipation rate ${\it\epsilon}_{U}$ (orange squares) and the total thermal dissipation rate ${\it\epsilon}_{T}$ (blue circles).

Figure 19. (a) Viscous dissipation in the bulk of the fluid as a function of $\mathit{Re}$ . (b) Corresponding compensated ${\it\epsilon}_{U}^{bu}$ scaling. (c) Thermal dissipation in the bulk of the fluid as a function of $\mathit{Re}$ . (d) Corresponding compensated ${\it\epsilon}_{T}^{bu}$ scaling. The solid black lines in the four panels correspond to the least-squares fit to the data for the numerical models with $\mathit{Ra}\geqslant 10^{5}$ .

The dominance of the boundary layer contribution in the thermal dissipation rate has already been reported by Verzicco (Reference Verzicco2003) for the same range of $\mathit{Ra}$ values. This phenomenon may be attributed to the dynamical nature of the plumes which permanently detach from the boundary layers and penetrate into the bulk of the fluid. These thermal plumes have the same typical size as the boundary layer thickness and can thus be thought as ‘detached boundary layers’. Grossmann & Lohse (Reference Grossmann and Lohse2004) have therefore suggested modifying their scaling theory to incorporate these detached boundary layers in the thermal dissipation rate. They propose to decompose ${\it\epsilon}_{T}$ into one contribution coming from the plumes ( ${\it\epsilon}_{T}^{pl}$ ) and one coming from the turbulent background ( ${\it\epsilon}_{T}^{bg}$ )

Such a decomposition is however extremely difficult to conduct in spherical shells in which the very large aspect ratio of the convective layer yields a complex and time-dependent multicellular large scale circulation (LSC) pattern (see for instance Bailon-Cuba, Emran & Schumacher Reference Bailon-Cuba, Emran and Schumacher2010, for the influence of large ${\it\Gamma}$ on the LSC). In the following, we therefore first retain the initial decomposition (5.1) before coming back to the inherent problem of accurately separating the different contributions to the dissipation rate in § 5.3.

5.2. Individual scaling laws for the dissipation rates

Based on the hypothesis of homogeneous and isotropic turbulence, the GL theory assumes that the thermal and viscous dissipation rates in the bulk of the fluid scale according to

in our dimensionless units. Figure 19 shows the bulk dissipation rates as a function of $\mathit{Re}$ for the numerical models of table 1. The best fit to the data (solid lines) for the cases with $\mathit{Ra}\geqslant 10^{5}$ yields ${\it\epsilon}_{U}^{bu}\sim \mathit{Re}^{2.79}$ and ${\it\epsilon}_{T}^{bu}\sim \mathit{Re}^{0.7}$ , are only roughly similar to the prediction (5.5). These deviations from the theoretical exponents are further confirmed by the compensated scalings ${\it\epsilon}_{U}^{bu}\,\mathit{Re}^{-3}$ and ${\it\epsilon}_{T}^{bu}\,\mathit{Re}^{-1}$ shown in (b) and (d), which show a coherent remaining dependence on $\mathit{Re}$ . Even at high Reynolds numbers, there is no evidence of convergence towards the exact expected scalings from the GL theory. This is particularly obvious for ${\it\epsilon}_{T}^{bu}$ which remains in close agreement with ${\it\epsilon}_{T}^{bu}\sim \mathit{Re}^{0.7}$ for the whole range of $\mathit{Re}$ values explored here (solid line in figure 19 d). The dependence of ${\it\epsilon}_{U}^{bu}$ on $\mathit{Re}$ shows a gradual steepening of the slope when $\mathit{Re}$ increases, which implies that ${\it\epsilon}_{U}^{bu}(\mathit{Re})$ cannot be accurately represented by a simple power law. Similar deviations from (5.5) have already been reported in the Taylor–Couette flow experiments of Lathrop, Fineberg & Swinney (Reference Lathrop, Fineberg and Swinney1992) and in the numerical simulations of RB homogeneous turbulence of Calzavarini et al. (Reference Calzavarini, Lohse, Toschi and Tripiccione2005).

A similar procedure can be applied to the dissipation rates in the boundary layers. In spherical shells, the volume fraction occupied by the boundary layers can be approximated by the following Taylor expansion to the first-order

in the limit of thin boundary layers ${\it\lambda}\ll r_{o}-r_{i}$ . The viscous dissipation rate in the boundary layers can then be estimated by

in our dimensionless units. As demonstrated in § 4, the boundary layers are laminar and are in reasonable agreement with the PB boundary layer theory. This implies that ${\it\lambda}_{U}\sim \mathit{Re}^{-1/2}$ and thus yields

Figure 20 shows the viscous dissipation rate in the boundary layers as a function of $\mathit{Re}$ for the numerical models of table 1. The least-squares fit to the data (solid lines) yields ${\it\epsilon}_{U}^{bl}\sim \mathit{Re}^{2.52}$ , in close agreement with the expected theoretical exponent. The compensated scaling displayed in (b) reveals a remaining weak secondary dependence of ${\it\epsilon}_{U}^{bl}$ on $\mathit{Re}$ , which is not accurately captured by the power law. The local slope of ${\it\epsilon}_{U}^{bl}(\mathit{Re})$ initially decreases with $\mathit{Re}$ (when $\mathit{Re}\lesssim 200$ ) before slowly increasing at higher Reynolds numbers ( $\mathit{Re}\gtrsim 200$ ). This behaviour suggests that, although simple power laws are at first glance in very good agreement with the GL theory, they may not account for the detailed variations of ${\it\epsilon}_{U}^{bl}(\mathit{Re})$ .

Figure 20. (a) Viscous dissipation in the boundary layers as a function of $\mathit{Re}$ . (b) Corresponding compensated ${\it\epsilon}_{U}^{bl}$ scaling. The solid black lines in the four panels correspond to the least-squares fit to the data for the numerical models with $\mathit{Ra}\geqslant 10^{5}$ .

The boundary layer contribution to the thermal dissipation rate is estimated in a similar way:

which yields

The laminar nature of the boundary layers also implies ${\it\lambda}_{T}\sim \mathit{Re}^{-1/2}$ and thus

Figure 21 shows ${\it\epsilon}_{T}^{bl}$ as a function of $\mathit{Nu}$ and $\mathit{Re}$ for the cases of table 1. The least-squares fits yield ${\it\epsilon}_{T}^{bl}\sim Nu^{0.95}$ and ${\it\epsilon}_{T}^{bl}\sim \mathit{Re}^{0.57}$ (solid lines), close to the expected exponents. However, the compensated scalings displayed in (b) and (d) reveal that the linear fit to ${\it\epsilon}_{T}^{bl}(\mathit{Nu})$ remains in good agreement with the data, while ${\it\epsilon}_{T}^{bl}(\mathit{Re})$ is not accurately described by such a simple fit. The solutions increasingly deviate from the power law at high Reynolds numbers with the local slope of ${\it\epsilon}_{T}^{bl}(\mathit{Re})$ gradually steepening with $\mathit{Re}$ .

Figure 21. (a) Thermal dissipation in the boundary layers as a function of $\mathit{Nu}$ . (b) Corresponding compensated ${\it\epsilon}_{T}^{bl}$ scaling. (c) Thermal dissipation in the boundary layers as a function of $\mathit{Re}$ . (d) Corresponding compensated ${\it\epsilon}_{T}^{bl}$ scaling. The solid black lines in the four panels correspond to the least-squares fit to the data for the numerical models with $\mathit{Ra}\geqslant 10^{5}$ .

Figure 22. (a) ${\it\epsilon}_{U}$ as a function of $\mathit{Re}$ . (b) ${\it\epsilon}_{U}$ normalised by the predictions coming from (5.13) (orange triangles) and (5.14) (blue circles) as a function of $\mathit{Re}$ . (c) ${\it\epsilon}_{T}$ as a function of $\mathit{Re}$ . (d) ${\it\epsilon}_{T}$ normalised by the predictions coming from (5.13) (orange triangles) and (5.14) (blue circles) as a function of $\mathit{Re}$ . The dashed black lines in (c,d) correspond to the equality between the asymptotic scalings and the data.

5.3. Individual versus global scalings

Despite the overall fair agreement with the GL predictions, a close inspection of the dependence of the four dissipation rates on the Reynolds number reveals some remaining dependence on $\mathit{Re}$ , which cannot be perfectly described by simple power laws. This is particularly obvious in the boundary layer contributions ${\it\epsilon}_{U}^{bl}(\mathit{Re})$ and ${\it\epsilon}_{T}^{bl}(\mathit{Re})$ . In addition, both thermal dissipation rates deviate further from the theoretical exponents than their viscous counterparts (see also Verzicco Reference Verzicco2003). One obvious problem is the inherent difficulty in the separation of the bulk and boundary layer contributions already discussed above. The dynamical plumes constantly departing from the boundary layers obviously complicate matters. To check whether the general idea of a boundary layer and a bulk contribution that both scale with the predicted exponents is at least compatible with the total dissipation rates, we directly fit

This leaves only the four prefactors ( $a,b,c,d$ ) as free fitting parameters and yields

We then compare this direct least-squares fit of the total dissipation rates to the sum of the individual scalings obtained in the previous section

The accuracy of the two scalings (5.13) and (5.14) are compared in figure 22, which shows the total viscous and thermal dissipation rates as a function of $\mathit{Re}$ for the numerical models of table 1. While the two scalings are nearly indistinguishable on (a,c), the corresponding normalised scalings shown in (b,d) reveal some important differences. The scalings based on the sum of the power laws derived in the last (5.14) are in relatively poor agreement with the data (5 %–10 % error for $\mathit{Re}>10^{2}$ ) with no obvious asymptotic behaviour. However, the global scalings (5.13) fall much closer to the actual values for the range $10^{2}<\mathit{Re}<10^{4}$ , and approach an asymptote for $\mathit{Re}>10^{2}$ . The deviations observed for the highest $\mathit{Re}$ cases probably have a numerical origin: the averaging timespan used to estimate the dissipation rates are likely too short in the most demanding cases to perfectly average out all the fluctuations.

The total thermal and viscous dissipation rates in our spherical shell simulations are thus better described by the sum of two power laws that follow the GL theory than by the sum of the asymptotic laws derived from the individual contributions, which suffers from an unclear separation of the boundary layer and bulk dynamics.

This result also sheds a new light on the placement of our numerical simulations in the GL regime diagram (figure 18). Equation (5.13) directly provides the estimated relative contributions of the bulk and boundary layers to the viscous and thermal dissipation rates. Figure 23 shows these different contributions for the numerical models of table 1 and reveals a completely different balance than that shown in figure 18. At low Rayleigh numbers, the estimated boundary layer contributions now dominate both the viscous and thermal dissipation rates. The viscous dissipation rate in the fluid bulk gradually increases with $\mathit{Ra}$ and dominates beyond $\mathit{Ra}>10^{7}$ . The bulk contribution to the thermal dissipation rate exhibits a similar trend, gradually increasing from approximately 5 % at $\mathit{Ra}=10^{4}$ to more than 40 % at $\mathit{Ra}=10^{9}$ . While the thermal dissipation rate in the fluid bulk never dominates in the regime explored here, our scaling predicts that it will do so for $\mathit{Ra}\gtrsim 3\times 10^{9}$ . These values would then locate the numerical models with $\mathit{Ra}\leqslant 10^{7}$ in GL regime I of the $\mathit{Ra}{-}\mathit{Pr}$ parameter space. The cases with $10^{7}<\mathit{Ra}<3\times 10^{9}$ would then belong to regime II. The transition to regime IV, where the bulk contributions dominate the dissipation rates, would then happen around $\mathit{Ra}\simeq 3\times 10^{9}$ .

Figure 23. Estimated relative contributions of the boundary layer (open symbols) and the fluid bulk (filled symbols) to $\widehat{{\it\epsilon}_{U}}$ (orange squares) and $\widehat{{\it\epsilon}_{T}}$ (blue circles) using (5.13).

Figure 24. (a) Nusselt number versus Rayleigh number. (b) Reynolds number versus Rayleigh number. (c) Compensated Nusselt number versus Rayleigh number. (d) Compensated Reynolds number versus Rayleigh number. The power laws given in (a,b) have been derived from a best fit to the cases of table 1 with $\mathit{Ra}\geqslant 10^{5}$ . In (c,d), the dashed lines correspond to the numerical solution of (6.1) and the solid lines to the numerical solution of (6.3).

While it is not clear that the contributions inferred from the total dissipation rates actually reflect the exact bulk and boundary layer contributions to ${\it\epsilon}_{U}$ and ${\it\epsilon}_{T}$ , this separation nevertheless allows us to reconcile the classification of our spherical shell convection models in the $\mathit{Ra}{-}\mathit{Pr}$ parameter space with the prediction of the GL theory for a fluid with $\mathit{Pr}=1$ . Future work that will better characterise and separate the bulk and boundary layer dynamics might help to reconcile the individual scalings (figure 18) with the global ones (figure 23), especially at low $\mathit{Ra}$ . In particular, considering dissipation layers as defined by Petschel et al. (Reference Petschel, Stellmach, Wilczek, Lülff and Hansen2013) instead of classical boundary layers might possibly help to better separate the bulk and boundary layer contributions.