Axiallyhomogeneous RayleighBénard convection in a cylindrical cell
Abstract
Previous numerical studies have shown that the “ultimate regime of thermal convection” can be attained in a RayleighBénard cell when the kinetic and thermal boundary layers are eliminated by replacing the walls with periodic boundary conditions (homogeneous RayleighBénard convection). Then, the heat transfer scales like and turbulence intensity as , where the Rayleigh number indicates the strength of the driving force. However, experiments never operate in unbounded domains and it is important to understand how confinement might alter the approach to this ultimate regime. Here we consider homogeneous RayleighBénard convection in a laterally confined geometry – a small aspectratio vertical cylindrical cell – and show evidence of the ultimate regime as is increased: In spite of the confinement and the resulting kinetic boundary layers, we still find . The system supports exact solutions composed of modes of exponentially growing vertical velocity and temperature fields, with as the critical parameter determining the properties of these modes. Counterintuitively, in the low regime, or for very narrow cylinders, the numerical simulations are susceptible to these solutions which can dominate the dynamics and lead to very high and unsteady heat transfer. As is increased, interaction between modes stabilizes the system, evidenced by the increasing homogeneity and reduced fluctuations in the r.m.s. velocity and temperature fields. We also test that physical results become independent of the periodicity length of the cylinder, a purely numerical parameter, as the aspect ratio is increased.
L. E. Schmidt et al.]LAURA E. SCHMIDT, ENRICO CALZAVARINI, DETLEF LOHSE, FEDERICO TOSCHI and ROBERTO VERZICCO
1 Introduction
There has been longstanding scientific interest in the heat transfer in RayleighBénard convection (RBC) – when a fluid confined between two plates undergoes thermal convection due to a temperature difference between the cold top and hotter bottom plate (Kadanoff (2001); Siggia (1994); Ahlers et al. (2009); Lohse & Xia (2010)). The dynamics of the system depend on the strength of the driving temperature difference (given by the Rayleigh number, ) and the ratio of the kinematic viscosity to the thermal diffusivity (given by the Prandtl number, ). Recently, studies of unconfined thermal convection – starting with the numerical work by Lohse & Toschi (2003) – were prompted by the relevance to natural convection phenomena, as occurs in the Earth’s atmosphere (Celani et al. (2007)) or in the core of stars (Garaud et al. (2010)), and as a way to test theoretical predictions for the scaling of the heat transfer (given by the Nusselt number, ) and turbulence intensity (given by the Reynolds number, ).
In the limit of high it was proposed that the dynamics would reach an ultimate regime of thermal convection, which is dominated by the bulk flow rather than the viscous and thermal boundary layers (Kraichnan (1962)). Within the GrossmannLohse theory of thermal convection, this regime where the dynamics does not depend on thermal plumes is also predicted (Grossmann & Lohse (2000, 2001, 2002)). Under these conditions, the Reynolds and Nusselt numbers are predicted to scale like and , with different scaling depending on the theory (Kraichnan (1962); Spiegel (1971); Grossmann & Lohse (2000, 2001, 2002)). The question of the existence of this regime was addressed in simulations of homogeneous RayleighBénard convection in a triperiodic cell (i.e., without confinement) where evidence was found supporting such a regime when all of the boundaries (and thus boundary layers) were absent (Lohse & Toschi (2003); Calzavarini et al. (2005, 2006)) and confirming (Calzavarini et al. (2005)) the Prandtl number dependence as suggested by Grossmann & Lohse (2000, 2001, 2002). Exponentially growing solutions at low were discovered analytically by Calzavarini et al. (2006) and observed in the simulations, appearing as “elevator modes” composed of strong upwards and downwards jets, and coinciding with an increasing heat transfer until the modes became unstable and disintegrated.
Later, experiments performed in long rectangular channels by Gibert et al. (2006, 2009); Tisserand et al. (2010)) and cylindrical pipes by Cholemari & Arakeri (2009) also found scaling laws consistent with those expected for the ultimate regime, despite the presence of the side walls which cause an additional anisotropy and drag in the flow. They also observed a mean global flow composed of hot upwards and cool downwardsmoving columns. Whereas in standard RBC strong thermal gradients occur at the top and bottom plates, leaving the temperature field within the bulk flow nearly uniform, in the experiments with long cells, there exists a mean, linear temperature gradient throughout the bulk. It is this underlying linear gradient which is used to drive the convection in homogeneous RBC, and in the axiallyperiodic convection cell considered here.
In this work, we numerically and analytically investigate thermal convection in a long vertical cylinder, in both the low and high limits. Though the top and bottom boundary layers are absent, the side walls are still present as in the experimental situations. The objective is to find out whether also in this confined (and presumably more stabilizing) geometry exponentially growing solutions (so called elevatormodes) exist, how the Nusselt and Reynolds number depend on , and how the results compare to recent experiments in confined flow geometries (Perrier et al. (2002); Gibert et al. (2006, 2009); Cholemari & Arakeri (2009)). We first derive exponentially growing mode solutions of the governing Boussinesq equations with laterally confined boundary conditions, and compare them with results from simulations. Then we increase and find that interaction between the modes stabilizes the system. This interaction prevents the growth of the exponential modes, and allows us to define the global heat transfer and the global Reynolds number of the system.
2 The axiallyhomogeneous RayleighBénard system and a class of exact analytic solutions
2.1 Definition of the system
The system under study is axiallyhomogeneous thermal convection in a vertical cylinder, i.e., in a domain with periodic boundary conditions in the axial direction, but with lateral confinement. Standard RayleighBénard convection occurs in a cylinder (or other geometry) with a top and bottom plate kept at fixed temperatures, and this temperature difference drives the flow. In the homogeneous or periodic case, with no boundaries except the radial wall, we drive the flow using a temperature gradient , where is an underlying linear temperature profile in about which fluctuations occur. The absolute local temperature can then be written as . This linear temperature (or density gradient in the case of Cholemari & Arakeri (2009)) is observed in experiments (Gibert et al. (2006, 2009)). The diameter of the cylinder is the only imposed scale, and therefore, we make the hypothesis that is the relevant length scale for the convection (aside from the additional smaller scales as turbulence develops). The periodicity length , which can be combined with the diameter to form the aspect ratio , is much less relevant because it does not change the forcing of the system. We nondimensionalize the governing NavierStokes equations using for lengths, the freefall speed for velocities, and for temperatures. The Prandtl number is and the Rayleigh number is
(2.0) 
Here the material parameters are the coefficient of thermal expansion , the kinematic viscosity , and the thermal diffusivity ; the gravitational acceleration is denoted as . Notice that the same definition of the Rayleigh number was chosen in the experimental studies by Gibert et al. (2006, 2009) and Tisserand et al. (2010). The more usual Rayleigh number based on the height of a closed RayleighBénard cell, , where is the temperature difference between the horizontal plates, is linked to the present one through the relation .
In nondimensional form, the Boussinesq equations read
(2.0) 
where is the velocity field with components in the , , and directions, respectively, and is the dimensionless pressure field.
The velocity field satisfies the noslip condition at the lateral walls located (in nondimensional units) at so that
(2.0) 
For the temperature we consider either a fixed condition at the lateral wall, which in term of the fluctuations means
(2.0) 
or adiabatic conditions at the wall, corresponding to:
(2.0) 
The system defined so far models the physical situation of a vertical pipe connected to two infinite reservoirs of fluids at different temperatures. However, in real experiments homogeneous convective systems are produced by connecting a vertical conduct to two closed and finite chambers filled with fluid at different temperatures (see e.g. Gibert et al. (2006) and Cholemari & Arakeri (2009)). Therefore, in order to mimic such a condition in the following it will be important to take into account the extra boundary condition of no mean vertical momentum and no mean temperature increase in the system
(2.0) 
which shall be satisfied at any point on the vertical axis of the cylinder. Note however that this “no mean value” condition does not prevent from having a positive heat flow, and indeed for the convective regime  as opposed to the conductive one  we expect .
2.2 Exact solutions to governing equations
As in the triperiodic case, we look for exact solutions which are translationally invariant in and have (Calzavarini et al. (2006)). The coupled equations for the vertical velocity and temperature are then:
(2.0) 
An adiabatic side wall boundary condition (Eq. 2.0) is relevant for experiments where tanks or cells are insulated at the sides, and will be implemented throughout the remainder of this paper. Appendix A gives the analytic result for the fixed temperature side wall boundary condition (Eq. 2.0). The mixed (i.e. adiabatic) boundary condition case is more difficult to solve analytically, but solutions can still be written for the special case . This case will be described in detail in the following.
Trying solutions with an exponential growth time dependence (), Eq. 2.0 simplifies and can be uncoupled by introducing two new fields and
(2.0) 
which satisfy
(2.0) 
The azimuthal dependence can be expanded as a combination of and and written in terms of the modes and which satisfy the modified Bessel equation:
(2.0) 
The label for the growth rate is necessary because again there are multiple solutions for a given , , which satisfy the boundary conditions. The radial dependencies of the solutions to Eq. 2.0 are
(2.0) 
where is the modified Bessel function of the first kind. The modified Bessel functions of the second kind, , are also solutions to Eq. 2.0 but are not physically acceptable because they diverge at .
The vertical velocity and temperature fields then have the form
(2.0) 
where the prefactor is unprescribed.
The growth rate of each mode depends on . After applying the boundary conditions at the cylinder wall at (Eqs. 2.0 and 2.0), the solutions correspond to the roots of
(2.0) 
Even though there are multiple solutions for each mode , the maximum growth rate will dominate over the others as time progresses.
Intuitively one might guess that the exponential growth cannot be maintained, and would be arrested by friction on the side walls. However, the solutions to the above system of equations have not ignored viscous drag on the walls. Even though the velocity, and gradients at the walls, are growing exponentially, the frictional force does not overtake the driving buoyancy force.
In the periodic RB cell, the global, timeaveraged dimensionless heat transfer is given as
(2.0) 
where the brackets indicate a volume (surface would be the same) and time average. When and , the instantaneous heat transfer also grows exponentially in time.
0.81  0.022  —  
1  0.52  —  —  
2  0.14  —  —  
0  0.89  0.45  —  
1  0.73  0.12  —  
2  0.52  —  —  
3  0.26  —  —  
0  0.94  0.70  0.28  
1  0.86  0.52  0.0041  
2  0.74  0.31  —  
3  0.60  0.080  —  
4  0.44  —  —  
5  0.25  —  —  
6  0.040  —  — 
2.2.1 Axisymmetric, dipolar modes and critical Rayleigh number
The spatial form of solutions corresponding to the zero mode case () have an axisymmetric profile (see Eq. (2.0)). This means that an increasing amount of fluid is transported through the cylinder and that the mean temperature of the system increases in time. Obviously these zeromode axisymmetric solutions have a nonzero mean vertical flow and nonzero mean temperature fluctuation, therefore they are in conflict with the condition (2.0). Although these solutions are nonphysical (they correspond to an unbounded infinite system) they will be useful to validate our simulations.
From the above considerations we deduce that the first unstable physical mode is , which corresponds to a dipolar flow profile, moving upwards in half of the cylinder and downwards in the other half. Such dipolar modes were also observed in experiments in a cylindrical column (Cholemari & Arakeri (2009)). The critical Rayleigh number corresponding to can be computed from Eq. (2.0), which leads to the value . In Table 1 we list the active modes at different numbers above . For the specific case and , the maximum growth rate is . The first physical dipolar mode has growth rate . The growth rates and their spatial dependencies are confirmed in the following simulations.
3 Numerical Method
Numerical simulations in a cylindrical cell are performed using the code developed by Verzicco & Orlandi (1996) and Verzicco & Camussi (2003), modified to satisfy the periodic conditions at . Furthermore, to satisfy the condition (2.0) the vertical flow and mean temperature fluctuation are subtracted at each time step. A similar procedure was adopted in the earlier triperiodic cell simulations by Lohse & Toschi (2003); Calzavarini et al. (2005, 2006). A finite difference scheme on a staggered grid is used to directly solve the incompressible continuity, momentum equation, and energy equations under the Boussinesq approximation using a fractional step method (Verzicco & Orlandi (1996); Verzicco & Camussi (2003); Kim & Moin (1985)). The grid is nonuniform in the radial direction, but is uniform in due to the periodicity, allowing the pressure solver to take advantage of trigonometric expansions and making the routine more efficient. In the simulations we must fix the axial periodicity length , which is here just a numerical parameter rather than a physical scale for the system dynamics. We expect, and will find, that as the aspect ratio is reduced (by increasing the periodicity length ), the system properties become independent of .
All the simulations have been started with a velocity field at rest and the temperature field consisting only of a random perturbation of maximum amplitude .
The new code implementation with periodic conditions was first checked by direct comparison of the predicted growth rates and radial forms of and from the analysis in the previous section, both for cases where and modes dominated, along with the fixed temperature side wall boundary conditions (Appendix A). For , , and , the axisymmetric and dipolar modes are isolated at early times and are shown in Figs. 1 and 2, comparing well to the theoretical predictions. Here, because the velocity and temperature fields are observed to be independent, we use the maximum vertical velocity field calculated at each time step to represent from Eq. 2.0.
When only a few of the exponentially growing modes are active at low (as for the parameters in Table 1), the velocity and temperature gradients at grow exponentially and thus the resolution unavoidably becomes insufficient at some point. However, this is only an issue for low cells because as more modes are available, their interaction prevents such a situation from occurring. We find that though the resolution can affect the maximum and (and thus, the heat transfer) obtained during the oscillations of the system, the growth rates are unaffected. For example, at , , and where the mode strongly dominates the system, simulations with a resolution of (in , , ) reach a maximum vertical velocity of before stabilizing, compared to for a resolution of (shown in Fig. 3). Even though the early time dynamics depend on the resolution, the growth rate and steadystate that is reached are unaffected. The situation is analogous to what Calzavarini et al. (2006) had found in the unconfined RB system.
The consistency of the code was also verified using the relation between the Nusselt number and the volume and timeaveraged kinetic dissipation rate, . For the periodic cell, which is confirmed within few percents as long as the spatial resolution of the simulation is adequate (see Table 2). It should be noted that for increasing the previous identity is verified less satisfactorily since limitations of the computational resources prevented us from running cases on larger grids.
4 Results
At low we observe individual growing modes which, depending on and , dominate the system dynamics. Already shown in Figs. 1–3 is the exponential growth at early times for simulations with , , and . After the initial growth at , Fig. 5 shows that the system temporarily stabilizes as horizontal velocity field fluctuations increase, which destroys the axiallyuniform modes of Eq. 2.0. This process continues periodically, with the horizontal field growing and decaying. This is visualized in Fig. 4 which shows snapshots of the dimensionless vertical velocity in a vertical cross section through the axis of the cell, just before, during, and after the stabilization process has occurred. Analogous dynamics were observed in the triperiodic homogeneous RB system of Calzavarini et al. (2006), and the process was found to depend sensitively on the numerical resolution and timestep size. Here, we have an additional parameter, the periodicity length which can affect this stabilization process. If the periodicity length is increased so that , instead of observing the longterm periodic behavior, the system settles and stabilizes after the horizontal fluctuations become of the same order as the vertical fluctuations (Fig. 5). This signals that the origin of the physical process which causes the destruction of the exponentially growing modes, is related to the ability of the fluctuations to extend along the direction.
Plotting as a function of time, along with the active , clearly shows convergence of the system to the fastest growing mode for with in Fig. 6.
As is increased, more modes are present and the growth rates become closer together. This leads to interaction between the modes and stabilization of the system by preventing individual modes from dominating the flow. Fig. 7 shows the time variation of the Nusselt number for for and . Though the maximums jump among the lower growthrate modes, there is no prolonged exponential growth as was seen at the lower . There is again a difference between and – the lower aspect ratio simulations do not reach the same maximum growth rates, leading to a more stable system. For reference and comparison to Fig. 4 with , Fig. 8 shows a snapshot of the vertical velocity field in a cell for .
We have already seen that the aspect ratio of the cell, determined by the chosen numerical periodicity length, can affect the stabilization of the system by the breakdown of the exponentially growing modes (Figs. 5 and 7). Because we are interested in understanding the heat transfer and development of turbulence in a long pipe, when going to higher one would hope that the results are unaffected by the choice of the periodicity length. Table 2 shows simulation results at at aspect ratios , , and . The aspect ratio was changed by doubling the periodicity length, and keeping the resolution in the same. It appears that the global system properties indeed start to converge once for . For higher , there is still aspect ratio dependence, but numerical limitations prevent us from increasing the periodicity length further. Since the simulations require half the number of grid points with respect to longer cells with , higher simulations which already require a finer resolution are done with the intermediate aspect ratio .
The spatial variation of the r.m.s. velocity components plotted in Fig. 9 for various shows the horizontal velocity peaks at the axis, and goes to zero at the walls. The vertical velocity fluctuations peak closer to the walls, with a local minimum along the axis, falling to zero at the walls owing to friction. At lower , the horizontal fluctuations do not reach the same size as the vertical fluctuations at the axis. As is increased, increases relative to at the axis, eventually overtaking . The horizontal fluctuations also become more uniform across the cell. The peaks in become sharper and move towards the outside.
Similar behavior was observed in the measured profiles in the buoyancydriven pipe turbulence experiments of Cholemari & Arakeri (2009). In their experiments with and , the profiles resembled those of Fig. 9(a) for the lower and , but the horizontal fluctuations reached only about 60% the magnitude of the vertical fluctuations at the axis. The reason for the resemblance to the lower result here is almost certainly the high value of . To check, we also calculated that the turbulent shear stress was absent across the cell. For all cases, is close to zero compared to the r.m.s. components. As in the experiments, the system remains anisotropic, evidenced by the fact that , especially near the walls where the horizontal fluctuations fall to zero much more rapidly than the vertical fluctuations.
Counterintuitively, the heat transfer at low and smaller, is very large due to single exponentially modes dominating the system in this regime. For these systems, it is inappropriate to define an average or because we have seen that the maximum peaks (as in Fig. 5) depend on the grid resolution. As is increased, the system stabilizes due to mode interaction and the scaling appears consistent with for , as shown in Fig. 10. The values used to determine the scaling law come from the aspect ratio runs (given in Table 2). Errors are estimated by the percentage difference of the two methods to calculate the Nusselt number (one, being the volume average in Eq. 2.0, the other being via the global relation with the kinetic energy dissipation rate). The highest runs appear underresolved but nevertheless the prediction confirms the scaling laws. The mean Reynolds number is calculated as (where is the total dimensional velocity), and also appears consistent with the power law from Fig. 10.
1/4  197  205  3.12  1.02  

1/2  384  428  2.50  0.97  
1/4  271  377  2.23  0.97  

1/2  678  880  –  1.00  
1/4  409  716  1.84  1.00  
1/8  431  730  1.98  1.00  

1/2  4544  4702  2.49  0.832  
1/4  2349  3665  1.80  0.870  
1/8  2196  3664  1.70  0.890  

1/4  2386  4747  1.59  0.83  
1/8  3248  5349  1.76  0.82  
1/4  12191  22909  1.43  0.55 
5 Conclusions
We have presented theoretical results and simulations of thermal convection in a laterally confined cylindrical pipe in the limit of small aspect ratios. At low we find exponentially growing modes consisting of upwards and downwardsflowing columns, analogous to those found in homogeneous RBC with no boundaries and simulations in triperiodic cells (Calzavarini et al. (2005, 2006)), along with experimental observations (Gibert et al. (2006, 2009); Cholemari & Arakeri (2009)). The breakdown of these modes is found to depend on the numerical aspect ratio of the cell (or periodicity length), implying that their destruction is related to the ability of disturbances in the axial direction to grow. At higher , the scaling of the heat transfer and turbulence are consistent with the predicted (Grossmann & Lohse (2000, 2001, 2002); Lohse & Toschi (2003)) scaling of the ultimate, bulkdominated regime of thermal convection, just as in the case of the triperiodic boundary conditions (Lohse & Toschi (2003); Calzavarini et al. (2005)).
It is remarkable that the sidewall boundary conditions and the resulting kinematic boundary layers do not lead to different scaling behavior with a less steep increase of with , as common for the lower regimes, say , in the unifying theory of Grossmann & Lohse (2000, 2001, 2002). The reason must be that these kinetic boundary layers only form at the sidewalls, and not at the top and bottom walls, which do not exist in our setup. We therefore neither expect a less steep increase of with beyond the onset of turbulence in the lateral kinetic boundary layers, which is expected to occur beyond a critical shear Reynolds number (Landau & Lifshitz (1987); Grossmann & Lohse (2000, 2002, 2011)). We note that here we are still considerably below this transition, as even for our largest we only have . It would be worthwhile to push the numerical simulations in this laterally confined geometry to numbers beyond the onset of the shear instability, in order to confirm that also then , without any logarithmic corrections, which seem to be typical for ultimate RB flow (and analogous ultimate TaylorCouette flow) in fully confined (i.e., with boundaries and boundary layers also towards the top and bottom) geometries (Chavanne et al. (1997, 2001); Dubrulle (2001); Grossmann & Lohse (2011); van Gils et al. (2011)).
A Zero temperature fluctuations at side wall solution
Though this boundary condition is not currently relevant for the experiments, it has a simple solution which was used to check the numerical results. Zero temperature fluctuations at the side wall means that physically, the absolute temperature follows the imposed linear profile at the walls. Eq. 2.0 is exactly satisfied by solutions of the form
(A 0) 
where is the th Bessel function of the first kind ( dependence is also acceptable). The growth rates and wave numbers are determined via the boundary conditions, Eqs. 2.0 and 2.0. Here with , is related to the (the th root of ) by
(A 0) 
The dispersion relation determining the growth rate for each mode is given by
(A 0) 
which has the same form as the triperiodic case (up to a constant factor due to the different choices for nondimensionalization), but the allowed values of differ because of the geometry. The critical Rayleigh number for this case, i.e. for , is .
The authors acknowledge the National Computing Facilities (NCF) for awarding us computing time on SARA in Amsterdam. We thank Valentina Lavezzo for careful reading of the manuscript. We also thank the EU COST Action MP0806.
 Ahlers et al. (2009) Ahlers, G., Grossmann, S. & Lohse, D. 2009 Heat transfer and large scale dynamics in turbulent RayleighBénard convection. Rev. Mod. Phys. 81, 503.
 Calzavarini et al. (2006) Calzavarini, E., Doering, C. R., Gibbon, J. D., Lohse, D., Tanabe, A. & Toschi, F. 2006 Exponentially growing solutions of homogeneous RayleighBénard flow. Phys. Rev. E 73, R035301.
 Calzavarini et al. (2005) Calzavarini, E., Lohse, D., Toschi, F. & Tripiccione, R. 2005 Rayleigh and Prandtl number scaling in the bulk of RayleighBénard turbulence. Phys. Fluids 17, 055107.
 Celani et al. (2007) Celani, A., Mazzino, A., Seminara, A. & Tizzi, M. 2007 Droplet condensation in twodimensional bolgiano turbulence. J. Turb. 8, 1–9.
 Chavanne et al. (1997) Chavanne, X., Chilla, F., Castaing, B., Hebral, B., Chabaud, B. & Chaussy, J. 1997 Observation of the ultimate regime in RayleighBénard convection. Phys. Rev. Lett. 79, 3648–3651.
 Chavanne et al. (2001) Chavanne, X., Chilla, F., Chabaud, B., Castaing, B. & Hebral, B. 2001 Turbulent RayleighBénard convection in gaseous and liquid he. Phys. Fluids 13, 1300–1320.
 Cholemari & Arakeri (2009) Cholemari, M. & Arakeri, J. 2009 Axially homogeneous, zero mean flow buoyancydriven turbulence in a vertical pipe. J. Fluid Mech. 621, 69–102.
 Dubrulle (2001) Dubrulle, B. 2001 Momentum transport and torque scaling in TaylorCouette flow from an analogy with turbulent convection. Eur. Phys. J. B 21, 295.
 Garaud et al. (2010) Garaud, P., Ogilvie, G., Miller, N. & Stellmach, S. 2010 A model of the entropy flux and reynolds stress in turbulent convection. Month. Not. Roy. Astro. Soc. 407, 2451–2467.
 Gibert et al. (2006) Gibert, M., Pabiou, H., Chilla, F. & Castaing, B. 2006 HighRayleighnumber convection in a vertical channel. Phys. Rev. Lett. 96, 084501.
 Gibert et al. (2009) Gibert, M., Pabiou, H., Tisserand, J.C., Gertjerenken, B., Castaing, B. & Chillà, F. 2009 Heat convection in a vertical channel: Plumes versus turbulence diffusion. Phys. Fluids 21, 035109.
 van Gils et al. (2011) van Gils, D., Huisman, S. G., Bruggert, G. W., Sun, C. & Lohse, D. 2011 Torque scaling in turbulent TaylorCouette flow with co and counterrotating cylinders. Phys. Rev. Lett. 106, 024502.
 Grossmann & Lohse (2000) Grossmann, S. & Lohse, D. 2000 Scaling in thermal convection: A unifying theory. J. Fluid. Mech. 407, 27–56.
 Grossmann & Lohse (2001) Grossmann, S. & Lohse, D. 2001 Thermal convection for large Prandtl number. Phys. Rev. Lett. 86, 3316–3319.
 Grossmann & Lohse (2002) Grossmann, S. & Lohse, D. 2002 Prandtl and Rayleigh number dependence of the Reynolds number in turbulent thermal convection. Phys. Rev. E 66, 016305.
 Grossmann & Lohse (2011) Grossmann, S. & Lohse, D. 2011 Multiple scaling in the ultimate regime of thermal convection. Phys. Fluids x, y.
 Kadanoff (2001) Kadanoff, L. P. 2001 Turbulent heat flow: Structures and scaling. Phys. Today 54 (8), 34–39.
 Kim & Moin (1985) Kim, J. & Moin, P. 1985 Application of a fractionalstep method to incompressible navierstokes equations. J. Comp. Phys. 59, 308–323.
 Kraichnan (1962) Kraichnan, R. H. 1962 Turbulent thermal convection at arbritrary Prandtl number. Phys. Fluids 5, 1374–1389.
 Landau & Lifshitz (1987) Landau, L. D. & Lifshitz, E. M. 1987 Fluid Mechanics. Oxford: Pergamon Press.
 Lohse & Toschi (2003) Lohse, D. & Toschi, F. 2003 The ultimate state of thermal convection. Phys. Rev. Lett. 90, 034502.
 Lohse & Xia (2010) Lohse, D. & Xia, K.Q. 2010 Smallscale properties of turbulent RayleighBénard convection. Ann. Rev. Fluid Mech. 42, 335–364.
 Perrier et al. (2002) Perrier, F., Morat, P. & LeMouel, J. L. 2002 Dynamics of air avalanches in the access pit of an underground quarry. Phys. Rev. Lett. 89, 134501.
 Siggia (1994) Siggia, E. D. 1994 High Rayleigh number convection. Annu. Rev. Fluid Mech. 26, 137–168.
 Spiegel (1971) Spiegel, E. A. 1971 Convection in stars. Ann. Rev. Astron. Astrophys. 9, 323–352.
 Tisserand et al. (2010) Tisserand, J.C., Creyssels, M., Gibert, M., Castaing, D. & Chillà, F. 2010 Convection in a vertical channel. New J. Physics 12, 075024.
 Verzicco & Camussi (2003) Verzicco, R. & Camussi, R. 2003 Numerical experiments on strongly turbulent thermal convection in a slender cylindrical cell. J. Fluid Mech. 477, 19–49.
 Verzicco & Orlandi (1996) Verzicco, R. & Orlandi, P. 1996 A finitedifference scheme for threedimensional incompressible flow in cylindrical coordinates. J. Comput. Phys. 123, 402–413.