Jump rope vortex flow in liquid metal Rayleigh–Bénard convection in a cuboid container of aspect ratio $\varGamma = 5$

Abstract

We study the topology and the temporal dynamics of turbulent Rayleigh–Bénard convection in a liquid metal with a Prandtl number of 0.03 located inside a box with a square base area and an aspect ratio of $\varGamma = 5$. Experiments and numerical simulations are focused on the Rayleigh number range $6.7 \times 10^{4} \leqslant Ra \leqslant 3.5 \times 10^{5}$, where a new cellular flow regime has been reported previously (Akashi et al., Phys. Rev. Fluids, vol. 4, 2019, 033501). This flow structure shows symmetries with respect to the vertical planes crossing at the centre of the container. The dynamic behaviour is dominated by strong three-dimensional oscillations with a period length that corresponds to the turnover time. Our analysis reveals that the flow structure in the $\varGamma = 5$ box corresponds in key features to the jump rope vortex structure, which has recently been discovered in a $\varGamma = 2$ cylinder (Vogt et al., Proc. Natl Acad. Sci. USA, vol. 115, 2018, pp. 12674–12679). While in the $\varGamma = 2$ cylinder a single jump rope vortex occurs, the coexistence of four recirculating swirls is detected in this study. Their approach to the lid or the bottom of the convection box causes a temporal deceleration of both the horizontal velocity at the respective boundary and the vertical velocity in the bulk, which in turn is reflected in Nusselt number oscillations. The cellular flow regime shows remarkable similarities to properties commonly attributed to turbulent superstructures.

1. Introduction

Rayleigh–Bénard convection (RBC) is a classical problem in fluid dynamics that has been studied for more than 100 years and serves as a model system for thermally driven flows in nature and engineering (Bénard 1900; Chandrasekhar 1961; Ahlers, Grossmann & Lohse 2009; Chillà & Schumacher 2012). The RBC configuration considers a fluid exposed to a destabilizing vertical temperature gradient between two plane-parallel plates and is described by the Rayleigh number, the Prandtl number and the aspect ratio of the fluid volume:

(1.1a–c)\begin{equation} Ra = \frac{\alpha{g}\Delta TH^{3}}{\kappa\nu},\quad Pr = \frac{\nu}{\kappa},\quad\varGamma = \frac{L}{H}, \end{equation}

where $\alpha$, $\kappa$ and $\nu$ stand for the thermal expansion coefficient, thermal diffusivity and kinematic viscosity of the fluid. The symbols $g$ and $\Delta T$ denote the gravitational acceleration and the vertical temperature difference in a fluid layer of thickness $H$ and characteristic horizontal dimension $L$. Convection occurs when the fluid is sufficiently heated from below and cooled from above. Warm fluid rises towards the lid and cold fluid sinks to the bottom. This leads to large-scale flows with a spatial extent of the order of the height of the container $H$.

The flow fields show manifold patterns and can become quite complex. A large-scale circulation (LSC), which is also called the ‘wind of turbulence’ (Krishnamurti & Howard 1981; Niemela et al. 2001; Ahlers et al. 2009; Xi et al. 2009), can be observed in a wide range of parameters. There are opposing opinions in the literature as to whether the LSC at large $Ra$ develops directly from the steady flow patterns at small $Ra$ (Busse, Zaks & Brausch 2003; Hartlep, Tilgner & Busse 2003) or whether this is a turbulent mode independent of it (Krishnamurti & Howard 1981). The vast majority of studies to date have focused on the generic configuration of a cylinder with aspect ratio unity (Ahlers et al. 2009). Here, the LSC exists in the form of a single circulation roll. To a first approximation one can assume that the single-roll LSC has a vertical planar structure where the fluid elements mainly follow an approximately circular or elliptical path (Villermaux 1995; Funfschilling & Ahlers 2004; Brown & Ahlers 2009; Xi et al. 2009; Zhou et al. 2009). This structure usually exhibits distinct oscillations, which are attributed to torsional and sloshing modes (Funfschilling & Ahlers 2004; Sun, Xia & Tong 2005b; Brown & Ahlers 2009; Xi et al. 2009; Zhou et al. 2009; Stevens, Clercx & Lohse 2011). Torsion means that the transverse flow on both horizontal walls performs azimuthal oscillations, whereby a phase shift of $180^{\circ }$ can usually be observed between the upper half and the lower half (Brown & Ahlers 2009; Zürner et al. 2019). This is connected with the sloshing mode which describes a gradual horizontal displacement of the circulation roll (Brown & Ahlers 2009). The periodicity of both oscillation modes coincides with the turnover time $\tau _{to} = L_{LSC}/v_{LSC} \approx \pi H/v_{LSC}$, where $v_{LSC}$ and $L_{LSC}$ are the typical velocity magnitude and path length of the LSC, respectively (Zürner et al. 2019). Reorientations of the LSC occur on longer time scales, which can be caused by a rotation of the LSC plane or a cessation of the flow structure (Brown, Nikolaenko & Ahlers 2005). Azimuthal rotations of the LSC occur preferably in cylindrical geometries, because the LSC plane becomes very likely locked in a predominant orientation in the convection cell with rectangular or square cross-sections (Ahlers et al. 2009). In these geometries, the LSC preferentially aligns along one or the other diagonal. Orientation changes between the two diagonal planes are also observed, which result from a lateral rotation of the LSC (Bai, Ji & Brown 2016). During this, the LSC spends some finite time in a transition state in which the LSC plane is temporarily oriented parallel to a set of sidewalls (Foroozani et al. 2017).

Variations of the geometry of convection cell, in particular the aspect ratio, have a substantial impact on the large-scale structure of the flow. It turns out that the single-roll LSC represents only the special case of a large-scale flow for a cylinder with aspect ratio unity and the flow pattern becomes more complex the greater the distance from configurations with $\varGamma = 1$. Two or more rolls one above the other can exist in tall cylinders ($\varGamma < 1$) (Verzicco & Camussi 2003; Amati et al. 2005; Tsuji et al. 2005; Sun, Xi & Xia 2005a; Stringano & Verzicco 2006; Xi & Xia 2008), while in flat containers ($\varGamma > 1$) multiple rolls occur side by side or combinations of rolls and cell structures are observed (Hartlep et al. 2003; von Hardenberg et al. 2008; Bailon-Cuba, Emran & Schumacher 2010; Yanagisawa et al. 2010; Emran & Schumacher 2015; Pandey, Scheel & Schumacher 2018; Schneide et al. 2018; Stevens et al. 2018; Sakievich, Peet & Adrian 2016,  2020; Krug, Lohse & Stevens 2020). Recent studies in shallow fluid layers at very large aspect ratios up to $\varGamma = 128$ revealed the existence of coherent flow structures that are not affected by the lateral boundary conditions. These structures survive against the background of high-frequency turbulent fluctuations for time scales being considerably longer than the turnover time during which a fluid package covers a complete circulation within the structure (Emran & Schumacher 2015; Pandey et al. 2018; Schneide et al. 2018; Stevens et al. 2018; Krug et al. 2020). The horizontal extent of these turbulent superstructures exceeds the height of the fluid layer by several times. Pandey et al. (2018) suggested that these turbulent superstructures can be retraced up to flow patterns that are formed immediately after the onset of convection.

Compared with the large number of studies dealing with convection in the generic aspect ratio unity or for the case of $\varGamma \gg 1$, the region of moderate aspect ratios has been sparsely studied so far. Thus, it is still an open question as to how the transition from a single-roll LSC to turbulent superstructures proceeds with increasing $\varGamma$. In experiments performed in air ($Pr = 0.7$), du Puits, Resagk & Thess (2007) observed the breakdown of the single-roll LSC to an oscillatory two-roll structure when the aspect ratio exceeds a critical value of 1.68. Further increase of $\varGamma$ beyond 3.66 causes the emergence of an unstable multi-roll regime. Both thresholds are not necessarily universal, since a dependence on $Ra$ and $Pr$ must also be taken into account. The authors suggest the second value of $\varGamma$ as a lower limit above which the turbulent convection reaches a state which becomes unaffected by the influence of the lateral walls.

Turbulent RBC in moderate aspect ratios has also been addressed by a few studies using direct numerical simulations. Bailon-Cuba et al. (2010) investigated thermal turbulence in an air layer ($Pr = 0.7$) in a cylindrical geometry as a function of the aspect ratio in the range $0.5 \leqslant \varGamma \leqslant 12$ for $Ra$ between $10^{7}$ and $10^{9}$. The authors evaluated the heat transport in connection with transitions in flow patterns. The observed large-scale flow consists of multiple rolls or cellular structures with a pentagonal or hexagonal symmetry. These patterns appear to be similar to those observed at slightly supercritical conditions close to the onset of convection for sufficiently large aspect ratios. Sakievich et al. (2016) reported the formation of hub-and-spoke structures as long-living coherent structures in a cylindrical domain at $\varGamma = 6.3$ for $Pr = 6.7$ and $Ra = 9.6\times 10^{7}$. The flow patterns, which persist over about 600 free-fall times, appear to be similar to those obtained by Bailon-Cuba et al. (2010). The authors suggest an analogy to the ‘wind of turbulence’ occurring in RBC at low aspect ratios. A further study by Sakievich et al. (2020) provides a continuation of the characterization of this particular regime by Fourier modal decomposition and by an extension of the total run time of the simulations covering more than 3000 free-fall times. The authors showed that compared to the standard case at $\varGamma = 1$, the flow dynamics covers much longer time scales of the order of hundreds to thousands of free-fall times. This is explained by the fact that coherent structures with larger length scales can establish themselves in increasing-aspect-ratio domains.

A recent experiment conducted by Vogt et al. (2018a) in a cylinder with aspect ratio $\varGamma = 2$ revealed a completely new, and until then unexpected feature of the LSC. Their measurements in liquid gallium ($Pr = 0.027$) detected a strong fluctuating flow along a measuring line that was not within the plane of the LSC. This observation cannot be reconciled with the conventional image of a quasi-two-dimensional LSC plane that only performs torsional and sloshing oscillations. Accompanying numerical simulations in $\varGamma = 2$ and $\varGamma = \sqrt {2}$ revealed a complex three-dimensional flow structure resembling a twirling jump rope. Since the existence of this jump rope vortex (JRV) has not yet been demonstrated in any other geometry, it raises the question as to whether this new LSC mode is a general property for RBC flows in various geometries of different aspect ratios or whether it is a special feature for a small specific range of $\varGamma.$ It is known that at aspect ratios of $\varGamma \approx 2 \dots 3$ a transition range of the LSC from a single-roll regime to adjacent double rolls occurs (Bailon-Cuba et al. 2010). In the vicinity of such transitions, the formation of unconventional flow structures could possibly occur.

In this paper, we report a combined experimental and numerical work which considers the turbulent RBC in a cuboid container with square horizontal cross-section of aspect ratio 5 and continues a previous experimental study made by Akashi et al. (2019). We use the eutectic metal alloy GaInSn ($Pr = 0.03$) as working fluid. Low-$Pr$ convection is characterized by a high thermal diffusivity and an enhanced production rate of vorticity and shear which amplifies the small-scale intermittency and turbulence in the flow (Scheel & Schumacher 2017). The dominant influence of inertia in low-$Pr$ fluids qualifies them as a suitable object of investigation with respect to the formation and the dynamical behaviour of coherent flow structures in turbulent RBC. Several flow regimes were detected by Akashi et al. (2019) in the $Ra$ range $7.9 \times 10^{3} \leqslant Ra \leqslant 3.5 \times 10^{5}$. An increase of $Ra$ causes an increasing horizontal wavelength of the flow structure and a conversion from a four-roll pattern via transient four- and three-roll regimes to a cellular flow regime. This cellular structure is subject to pronounced regular oscillations, whose period corresponds to the turnover time, but proves to be quite stable over a long period of time. The evaluation of the temperature and velocity signals demonstrated that the flow in the investigated $Ra$ number range is turbulent. The previous study by Akashi et al. (2019) was not designed to uncover individual details of the cellular flow structure or to provide information about the changes that the structures undergo during the distinct oscillations. The motivation of the present study is to investigate the nature of these oscillations in more detail and to reveal the internal dynamics.

The paper is organized as follows. In § 2 we present a description of the experimental set-up and the numerical model. Section 3 reports on ultrasonic measurements of the flow structure and its temporal behaviour. Corresponding results obtained from numerical simulations are presented with respect to the three-dimensional flow structure, the dominating frequencies and the transport properties. A quantitative comparison of these data with the experimental findings confirms the reliability of our numerical approach. Section 4 is dedicated to the detailed analysis of the complex three-dimensional oscillations of the cellular pattern. Moreover, the consequences of the periodic changes of the flow field on heat transport are investigated.

2. Experimental methods and numerical scheme

2.1. Experimental set-up

The set-up employed for the experiments is almost identical to that used in previous studies (Tasaka et al. 2016; Vogt et al. 2018b; Akashi et al. 2019; Tasaka et al. 2021; Vogt et al. 2021; Yang, Vogt & Eckert 2021) to which the reader is referred for a detailed description. A noteworthy difference from our previous study (Akashi et al. 2019) is that additional measuring positions for the velocity sensors have been placed at the mid-height of the fluid container in order to better compare the flow structures found here with the results reported by Vogt et al. (2018a). Figure 1 shows schematics of the container and the arrangement of the ultrasonic measurement lines. The fluid layer has a square horizontal cross-section of $L \times L = 200 \times 200\ \textrm {mm}^{2}$ and a height $H$ of 40 mm, resulting in an aspect ratio of $\varGamma = 5$. In the vertical direction the fluid layer is bordered by two copper plates, whose temperatures, $T_{bot}$ and $T_{top}$, are kept constant by water circulation in branched channels inside the plates. The water temperatures are controlled by external thermostats. The vertical temperature difference $\Delta T=T_{bot}-T_{top}$ is varied between 3.0 and 15.5 K corresponding to the $Ra$ range from $6.7 \times 10^{4}$ to $3.5 \times 10^{5}$. The sidewalls are made of electrically non-conducting polyvinyl chloride. The entire fluid layer is encased in 30 mm closed-cell foam to minimize heat loss. The working fluid is the eutectic alloy $\textrm {Ga}^{67}\textrm {In}^{20.5}\textrm {Sn}^{12.5}$ that is liquid at room temperature. Temperature-dependent data of the thermo-physical properties can be found in Plevachuk et al. (2014). In particular, for a temperature of $25\,^{\circ } \textrm {C}$ the following values are reported: density $\rho = 6360\ \textrm {kg}\ \textrm {m}^{-3}$, thermal expansion coefficient $\alpha = 1.24 \times 10^{-4}\ \textrm {K}^{-1}$, kinematic viscosity $\nu = 3.4 \times 10^{-7}\ \textrm {m}^{2} \ \textrm {s}^{-1}$ and thermal diffusivity $\kappa = 1.03 \times 10^{-5} \textrm {m}^{2}\ \textrm {s}^{-1}$. This gives a Prandtl number of $Pr = 0.033$. The thermal diffusion time for the fluid layer, $t_\kappa = H^{2}/\kappa$, is $155$ s. Table 1 presents the $Ra$ investigated in the experiments together with the corresponding values for the free-fall velocity $u_{ff}= \sqrt {g\alpha \Delta TH}$ and the free-fall time $t_{ff} =\sqrt {H/(g\alpha \Delta T)}$ resulting from the respective temperature difference $\Delta T$.

Figure 1. Experimental set-up and arrangement of measurement lines: (a) top view and (b) side view, where light blue lines indicate the ultrasonic measurement lines.

Table 1. Parameters for the experiments including the Rayleigh number $Ra$, the Prandtl number $Pr$, the free-fall velocity $u_{ff}$ and the free-fall time $t_{ff}$.

Due to the limited possibilities for measuring fluid velocities in liquid metals, the majority of previous experimental studies of liquid metal convection were largely restricted to the description of heat transfer in the system and the measurement of temperatures and their fluctuations (Takeshita et al. 1996; Cioni, Cilberto & Sommeria 1997; Segawa, Naert & Sano 1998; Horanyi, Krebs & Müller 1999; Burr & Müller 2001). The development of the ultrasonic Doppler technique for flow measurements in liquid metals (Takeda 1987; Eckert & Gerbeth 2002) enables an experimental determination of flow structures in low-$Pr$ liquid metal convection (Cramer, Zhang & Eckert 2004; Mashiko et al. 2004; Yanagisawa et al. 2010,  2013; Tasaka et al. 2016; Zürner et al. 2019; Tasaka et al. 2021; Vogt et al. 2021; Yang et al. 2021). This powerful technique was successfully used in the previous study by Akashi et al. (2019) to identify the individual flow regimes and is also employed here. It provides instantaneous profiles of the velocity component projected onto each measurement line, $u_x(x,t)$ and $u_y(y,t)$, respectively. The spatial resolution of the velocity measurements is 1.4 mm in the direction of the ultrasonic beam line and about 5 to 8 mm in the direction perpendicular to it. The velocity resolution is about $0.5\ \textrm {mm} \ \textrm {s}^{-1}$ and a temporal resolution of 0.6 s is achieved. The measurement instrument (DOP 3010, Signal Processing SA) is applied in combination with ultrasonic transducers having a basic frequency of 8 MHz and a piezo-element with an active diameter of 5 mm. These sensors, which are mounted horizontally inside holes drilled into the sidewalls of the container, are in direct contact with the liquid metal. In general, velocity probes can be installed at sidewalls at different heights. In the present study, measurement results are obtained along three measurement lines (shown as light blue lines in figure 1), designated as $Middle \textit {1}$, $Middle \textit {2}$ and $Bottom$. $Middle \textit {1}$ and $Middle \textit {2}$ are at the mid-height of the container ($z = 0.5H$) and $Bottom$ is situated close to the bottom plate ($z = 0.25H$). The velocity profiles along each measurement line are acquired sequentially by multiplexing.

2.2. Numerical scheme

For the numerical simulation, a Boussinesq fluid in a cuboid container with a square horizontal cross-section and an aspect ratio of $\varGamma = 5$ is considered. Here, we used the non-magnetic version of a code originally developed for magnetohydrodynamic simulations in Yanagisawa, Hamano & Sakuraba (2015). Cartesian coordinates $(x, y, z)$ are used with the $z$ axis in the upward direction and the origin located at one of the bottom corners of the container. We solve the governing equations for the non-dimensional velocity ${\boldsymbol {u}}$, temperature $T$ and pressure $p$ as follows:

(2.1a,b)\begin{gather} \frac{\partial {\boldsymbol{u}}}{\partial t}={-}({\boldsymbol{u}}\boldsymbol{\cdot}\boldsymbol{\nabla}){\boldsymbol{u}}-\boldsymbol{\nabla} p+PrRa ({{T}}-\bar{T}(z)) {\boldsymbol{k}}_z+Pr\nabla^{2}{\boldsymbol{u}} , \quad \bar{T}(z)=1-z , \end{gather}

(2.2)\begin{gather}\boldsymbol{\nabla}\boldsymbol{\cdot}{\boldsymbol{u}}=0 , \end{gather}

(2.3)\begin{gather}\frac{\partial T}{\partial t}={-}({\boldsymbol{u}}\boldsymbol{\cdot}\boldsymbol{\nabla})T+\nabla^{2} T , \end{gather}

where ${\boldsymbol {k}}_z$ is the unit vector in the $z$ direction. Length and time are respectively non-dimensionalized by the layer thickness $H$ and the thermal diffusion time $t_\kappa =H^{2}/\kappa$. Accordingly, the velocity is normalized by $\kappa /H$. The temperature satisfies the isothermal boundary conditions at the top ($T_{top} = 0$) and bottom ($T_{bot} = 1$), while adiabatic conditions are assumed for the sidewalls. No-slip conditions are applied for the velocity at all boundaries. The second-order-accurate staggered-grid finite-difference method is used with a uniform grid interval in each direction. The time integration was carried out explicitly by means of the third-order Runge–Kutta method. This procedure modifies the equation of continuity to

(2.4)\begin{equation} \varepsilon^{2}{\partial p}/{\partial t}={-}\boldsymbol{\nabla}\boldsymbol{\cdot}{\boldsymbol{u}}, \end{equation}

where $\varepsilon$ is a small parameter from Chorin (1967). In our numerical simulations a value of $\varepsilon =0.002$ is used (Yanagisawa et al. 2015). The code is parallelized with MPI in combination with OpenMP.

To check the resolution of the simulation, we compared results from calculations made with different grid resolutions: $n_z = 80$, 128, 256 and $512$ (see table 2). For these pre-calculations, the parameters are set as $Pr = 0.025$, which is a representative value for liquid metals, and $Ra = 1.0\times 10^{5}$, $3.0\times 10^{5}$ and $6.0\times 10^{5}$. All simulations performed were able to reconstruct the same qualitative characteristics of the flow field: a cellular pattern with quasi-periodic oscillations. Quantitative comparisons are carried out using the values of the Nusselt number and the Reynolds number which were averaged both temporally and spatially. In our equations $\langle \cdot \rangle _{S, t}$ stands for the time–surface average on the top or bottom boundary and $\langle \cdot \rangle _{V, t}$ stands for the time–volume average. Nusselt numbers are calculated both at the top boundary $Nu_{top}$ and at the bottom boundary $Nu_{bot}$:

(2.5a,b)\begin{equation} Nu_{top}={-}\left \langle \left (\frac{\partial T}{\partial z} \right)_{z=1} \right \rangle_{S, t}, \quad Nu_{bot}={-}\left \langle \left(\frac{\partial T}{\partial z} \right)_{z=0} \right \rangle_{S, t}. \end{equation}

An alternative way to determine the Nusselt number is to integrate the upward heat transport in the entire volume (e.g. Pandey et al. 2018). For comparison we also use this definition for the Nusselt number $Nu_{vol}$:

(2.6)\begin{equation} Nu_{vol}=1 +\langle u_z T \rangle_{V, t}. \end{equation}

Table 2. Pre-tests of numerical simulations conducted with different numbers of grid points: the Rayleigh number $Ra$, the Prandtl number $Pr$, the number of grid points $n_x$, $n_y$ and $n_z$, the time-averaged Nusselt number $Nu$, the standard deviation of Nusselt number $Nu_{std}$, the number of grid points in the thermal boundary layer $\lambda _{\varTheta }$ estimated by the relation $Nu = 1/{(2\lambda _{\varTheta })}$, $n_{\lambda _{\varTheta }}$, the time-averaged Reynolds number $Re$, the standard deviation of Reynolds number $Re_{std}$ and the ratio of the Kolmogorov scale $\eta _{K}$ to the spatial resolution $l_z$.

When thermally balanced states are achieved after adequate time integrations, applying (2.5a,b) and (2.6) to our data yields a very good agreement. Therefore, in the further course of the paper only the values of $Nu_{top}$ will be presented as $Nu$.

The Reynolds number is calculated using the root-mean-square (r.m.s.) velocities $U_{rms}$ determined for the entire volume as characteristic velocity scale:

(2.7)\begin{equation} Re = \frac{U_{rms}}{Pr} \quad \mbox{with} \ U_{rms}=\sqrt{\langle {u_x}^{2}+{u_y}^{2}+{u_z}^{2} \rangle_{V, t}}. \end{equation}

The quantitative comparison of the pre-tests is shown in table 2. The agreement of the time-averaged Nusselt numbers $Nu$ is satisfactory within 0.4 %, and the Reynolds numbers $Re$ are almost the same within 0.8 %. As is shown below, quasi-periodicity is a dominant feature of the pattern; hence, fluctuations of the Nusselt number are also important. The fluctuations are evaluated in the form of the standard deviation and are listed in table 2 as $Nu_{std}$.

In addition, we checked the spatial resolution of our simulations by comparing them with the Kolmogorov scale $\eta _{K}$. Using the results of $Nu$ for given values of $Pr$ and $Ra$, $\eta _{K}$ can be estimated as (Scheel, Emran & Schumacher 2013)

(2.8)\begin{equation} \eta_{K} =\left(\frac{Pr^{2}}{(Nu-1)Ra}\right)^{{1}/{4}}. \end{equation}

Table 2 contains the ratio of $\eta _{K}$ to the spatial resolution $l_z= 1/n_z$. A ratio larger than one indicates a sufficient resolution in terms of direct numerical simulation; hence, these values are used as references. As shown in table 2, simulations with the highest grid resolution satisfy this criterion for the $Ra$ range considered here. Calculations of lower resolution are conducted in the case of long simulation times. Due to the high diffusivity of temperature, we assume that this is nevertheless sufficient to adequately evaluate $Nu$, $Re$, temperature distribution and flow pattern. Our independence grid study shows that the values of the quantities we focus on in this paper do not change significantly on a finer grid (see above), so we can assume their validity. Higher resolutions of the grid are necessary if details within the velocity boundary layers are to be examined, but this is outside the scope of this study.

Based on these comparisons and considering available computational resources, we performed simulations mainly with the resolution $n_z=128$. To check the long-term stability of the quasi-periodic behaviour and to perform long-time averaging over 100 thermal diffusion times (${\sim }300$ oscillations), we utilized the resolution of $n_z=80$. The parameters chosen for the numerical simulations in this study are listed in table 3. The numerical simulations cover the parameter range that is investigated in the experiments.

Table 3. Parameters for the numerical simulations of the experimental configuration including the Rayleigh number $Ra$, the number of grid points $n_x$, $n_y$ and $n_z$, the free-fall velocity $u_{ff}$ and the free-fall time $t_{ff}$.

3. Characterization of the large-scale flow structure

3.1. Three-dimensional cellular regime

Previous experiments carried out by Akashi et al. (2019) using the same set-up show a clear and reproducible dependence on $Ra$ with respect to the occurrence of different flow regimes. While the range of smaller $Ra$ is associated with various regimes of unsteady convection rolls, at $Ra$ above $6.5\times 10^{4}$ a new, stable flow pattern occurs, which is characterized by a high symmetry with respect to the square base area of the fluid container. This newly discovered flow pattern was designated as the cellular flow regime. We are aware that the term of a cellular flow structure in thermal convection is often associated with stationary cell patterns observed, for example, in the slightly supercritical state after the onset of convection. This is explicitly not the case here. Accompanying measurements of the temperature fluctuations indicate that the stage of developed thermal turbulence is reached in the $Ra$ range where the cellular regime occurs (Akashi et al. 2019). Our point of view is supported by results from Bailon-Cuba et al. (2010) who found extended convection rolls and pentagon-like cells in their direct numerical simulation data for a turbulent flow in a $\varGamma = 12$ cylinder at $Ra = 10^{7}$ and $Pr = 0.7$. After filtering out the small-scale turbulence, these resulting patterns of turbulent convection resemble the weakly nonlinear regime just above the onset of convection. However, the observation of the cellular regime by Akashi et al. (2019) is also surprising because the appearance of such cell-like structures is anticipated to be more prevalent in fluids with larger $Pr \gtrsim 1$, while the convection pattern at low $Pr$ is supposed to be dominated by roll-like structures (Breuer et al. 2004; Pandey et al. 2018). In order to analyse the cellular flow regime in more detail, the investigations in this study are focused on the $Ra$ range $6.7 \times 10^{4} \leqslant Ra \leqslant 3.5 \times 10^{5}$. The upper limit of $Ra$ is determined by the technical capabilities of the experimental equipment.

Figure 2(a) shows the spatio-temporal map of the velocity component $u_x$ at $Ra=1.2\times 10^{5}$ recorded by the ultrasonic sensor $Bottom$ at a height of $z = 0.25 H$. A negative velocity (blue) indicates a flow toward the transducers, while a positive velocity (red) is associated with the direction away from the sensor. In all results presented in the following, the velocity is given dimensionless with respect to the free-fall velocity $u_{ff}= \sqrt {g\alpha \Delta T H}$. In the convection cell considered here, the free-fall velocity reaches a value of $16.48\ \textrm {mm}\ \textrm {s}^{-1}$ at $Ra=1.2\times 10^{5}$. The dimensionless time is given in units of the free-fall time $t_{ff} =\sqrt { H/(g\alpha \Delta T)} = 2.43$ s.

Figure 2. (a) Spatio-temporal velocity map at $Ra=1.2\times 10^{5}$ recorded for $u_{x}$ along the $x$ direction (sensor $Bottom$ at $z = 0.25H$). (b) A schematic illustration of the cell structure with the sensor positions.

The velocity map obtained by the sensor $Bottom$ (figure 2a) reproduces a flow structure almost identical to those presented in the recent study by Akashi et al. (2019). The flow pattern divides the measured profile into two parts with respect to the centre of the container. The strength of the measured velocity is subject to strong quasi-periodic fluctuations, which at some points in time (e.g. $t/t_{ff} = 100, 150, 270$) even lead to short-term reversals of the flow direction. If the average behaviour of the flow is considered, then a fluid motion from the sidewall towards the centre is clearly dominant, which means that an ascending flow in the centre of the cell can be assumed. Based on the evaluation of the measured velocity profiles along six different measurement lines within a horizontal plane above the heated bottom, Akashi et al. (2019) obtained a rough reconstruction of the time-averaged velocity patterns in the convection cell. For the case of the cellular regime, they revealed a three-dimensional structure having upwelling flows at the centre and the four corners of the container as shown in figure 2(b) (Akashi et al. 2019). Once this flow pattern is established, it is stable for long measurement times of the order of a few thousand free-fall times and does not change spontaneously. Repeating the measurements several times, inverse flow patterns with a downward flow in the centre of the fluid container were also observed. These structures, which differ principally in the sign of the velocity, have similar probabilities of occurrence. Which variant is realized in the experiment apparently depends only on random asymmetries in the initial conditions.

Inspired by the work of Vogt et al. (2018a), here we additionally monitor the velocities in the horizontal centre plane of our container. If the convection pattern corresponds to a single circulation, which fills the space between the copper plates, one would expect that in the measuring plane at the mid-height the vertical component clearly dominates and only marginal parts of a horizontal flow can be found here. As shown in figure 3(b,c), it is surprising to observe pronounced flows also in the centre plane, their intensity only slightly below that detected along the measuring lines at $z = 0.25H$ and $z = 0.75H$. Moreover, these velocity plots are characterized by the same pronounced oscillations. In contrast to the measurement at $z = 0.25H$, however, an inversion of the flow direction occurs in every period where no dominance of one of the two directions can be determined. Figure 3 also contains a corresponding dataset from Vogt et al. (2018a) (see figure 3d,e). For better comparability, the time axis and the velocity are normalized with the free-fall time and the free-fall velocity, respectively. The comparison of both flow patterns shows striking similarities with respect to the manifestation of the dominant oscillations in terms of both the amplitude of the velocity and the time scale of the oscillations. This raises the question as to whether we are dealing with the same phenomenon in the $\varGamma =5$ box that is referred to as JRV and has already been observed in the $\varGamma =2$ cylinder. To answer this question, we take a closer look at the flow structure and its changes during the oscillations. Since the experiment reaches its limits here due to the small numbers of measuring sensors, we performed numerical simulations in parallel.

Figure 3. Spatio-temporal velocity distributions measured at the mid-height of the fluid container for $Ra = 1.2 \times 10^{5}$. (a) Sketch of the measuring line in the $\varGamma = 5$ box ($Pr = 0.033$). (b) Velocity pattern for the the full-length measuring line. (c) The same dataset as in (b) plotted for only one half of the measuring line. (d) Sketch of the measuring line in the $\varGamma = 2$ cylinder ($Pr = 0.027$). (e) Velocity pattern in the $\varGamma = 2$ cylinder. (Adapted from Vogt et al. 2018a.)

Figure 4(a,b) contains the numerical counterpart of the measurements as presented in figures 2 and 3. The qualitative agreement of the flow patterns proves the reliability of the numerical approach. The numerical velocity distributions are also dominated by strong oscillations, with periodic changes of direction in the centre plane of the container.

Figure 4. Results of the numerical simulations showing a reconstruction of the spatio-temporal velocity maps at $Ra=1.2\times 10^{5}$ in the cuboid container for (a) $u_x$ according to the measuring line of sensor $Bottom$ at $z = 0.25$ and (b) $u_y$ according to the measuring line of sensor $Middle \textit {1}$ at $z = 0.5$. The measurement positions are shown in figure 2(b).

Snapshots of temperature isosurfaces and velocity distributions in figure 5 provide a three-dimensional visualization of the convection pattern. Figures 5(a) and 5(b) display instantaneous isosurfaces of the temperature for the values $T = 0.2$ and $T = 0.8$. Figures 5(c) and 5(d) show the instantaneous distribution of the vertical velocity $u_z$ at different heights of $z = 0.75$ and $z = 0.25$. The colour indicates the direction of flow, meaning that blue (red) signifies a downwelling (upwelling) motion. The temperature surface appears to be smoother than the velocity fields, which is naturally due to the fact that low-$Pr$ fluids feature a high thermal diffusivity. Apart from that, temperature field and vertical velocity show an evident conformity. Previous studies of turbulent superstructures raised the question as to what extent the temperature distribution and the vertical velocity structure represent the same coherent structure (Pandey et al. 2018; Stevens et al. 2018; Krug et al. 2020). Sakievich et al. (2016) showed that the thermal updrafts and downdrafts match the three-dimensional field to a high degree in a $\varGamma =6.3$ cylinder. Further analysis by Sakievich et al. (2020) revealed a high correlation in the temporal dynamics of the low-order Fourier modes of temperature and vertical velocity. Based on this clear similarity and the findings of Krug et al. (2020), Sakievich et al. (2020) hypothesize similar principles of spatial organization of the coherent patterns in a $\varGamma =6.3$ cylinder and superstructures occurring in larger aspect ratios.

Figure 5. Results of the numerical simulations performed for $Ra=1.2\times 10^{5}$ representing snapshots of the temperature isosurface for (a) $T=0.2$ and (b) $T=0.8$, and the distribution of the vertical velocity $u_z$ at the height of (c) $z = 0.75$ and (d) $z = 0.25$. The colour scale of $u_z$ is the same as that in figure 4. See also supplementary movies 1–4 available at https://doi.org/10.1017/jfm.2021.996.

The temperature plots in figure 5 show that the cellular structure is characterized by an ascending flow in its centre and at the four corners of the container (see figure 5b), while the shape of the areas where the fluid descends resembles a rhombus (see figure 5a). It can also be seen that the zones of ascending warm fluids in the middle and the four corners are connected by diagonally running ridges (figure 5b). This rough characterization of the flow structure is consistent with the distributions of the vertical velocity component in figures 5(c) and 5(d). The appearance of finer structures in the velocity contours illustrates the character of inertia-dominated fluid turbulence and is in line with previous studies (Pandey et al. 2018; Vogt et al. 2018a; Krug et al. 2020; Sakievich et al. 2020) which demonstrated the occurrence of more small-scale velocity fluctuations compared with the temperature field. Thus, it can be stated that the results of the numerical simulations confirm the mean flow structure of the cellular regime suggested from the experiment based on a limited number of velocity sensors (Akashi et al. 2019). A striking property of this cellular convection is the occurrence of pronounced periodic oscillations, which consequently must be reflected in some form of periodic changes of the flow topology.

3.2. Dominant oscillation frequencies

The dominant frequency $f_{OS}$ of the periodic oscillations in the spatio-temporal velocity maps is determined by calculating the spatially averaged power spectral densities (PSDs) from velocity time series. The velocity time series are divided into multiple time intervals of 60 free-fall times. The spectra are calculated by fast Fourier transform for each time interval and for each measurement point. These interim results are finally averaged both spatially along the respective measuring line and temporally over the duration of the measurement.

Figures 6(a) and 6(b) present the PSDs obtained from the velocity maps along the measurement lines covered by the sensors $Bottom$ and $Middle \textit {1}$. The frequency is normalized by the thermal diffusion frequency $f_\kappa = \kappa /H^{2}$. The dominant oscillation frequency can be easily identified by the clear peaks in the PSDs. The same frequency is found at both measuring positions indicating a certain coherence of the flow structure. The spectra in both figures 6(a) and 6(b) show similar magnitudes over the whole frequency range and an identical slope in the inertial domain for $f/f_\kappa > 10^{1}$. There is a slight deviation between experiment and numerical simulation concerning the values of the normalized frequency $f/f_\kappa = 3.3$ and $f/f_\kappa = 3.1$, respectively.

Figure 6. Spatially averaged PSDs calculated from the spatio-temporal velocity maps for $Ra = 1.2\times 10^{5}$ (a) measured in the experiments by the sensors $Bottom$ and $Middle \textit {1}$ and (b) obtained from numerical simulations. (c) Oscillation frequency $f_{OS}$ as a function of $Ra$ from experimental and numerical data, where the solid line represents the least-squares fit of the experimental results, $f_{OS}H^{2}/\kappa = 0.016Ra^{0.46}$, and the dotted line shows the fitting of the results of numerical simulations, $f_{OS}H^{2}/\kappa = 0.011Ra^{0.48}$. (d) The period length of the oscillations normalized with the turnover time $\tau _{TO}$ versus $Ra$, where the the dotted line represents a value of one.

The oscillation frequency $f_{OS}$ increases with increasing $Ra$ as shown in figure 6(c). A power-law fit leads to a scaling of $f_{OS}/f_\kappa =0.016Ra^{0.46\pm 0.02}$ for the velocity measurements in the experiments and $f_{OS}/f_\kappa =0.011Ra^{0.48\pm 0.02}$ for the numerical data. The agreement is satisfactory, where the exponents agree within the error bars. The slight discrepancy may result from uncertainties in the material properties, which in turn could affect the exact determination of $Pr$ and $Ra$ in the experiment. In the paper by Plevachuk et al. (2014), from which the numerical values for the thermal diffusivity, viscosity and density of GaInSn were taken, the accuracy of these measured values is specified in the range from $\pm 1.5\,\%$ to $\pm 7\,\%$, which might be sufficient to explain the above mentioned deviations.

Previous results reported by Akashi et al. (2019) for experiments in the same container reported a scaling law of $f_{OS}/f_\kappa = 0.031Ra^{0.40\pm 0.02}$. The discrepancy can be explained by the fact that the range in which the power-law fit was calculated in Akashi et al. (2019) also includes smaller $Ra$, where unstable rolls occur instead of the cellular structure. This also indicates that the two separate flow patterns cause different exponents in the scaling law. Experimental (Tsuji et al. 2005; Vogt et al. 2018a; Zürner et al. 2019) and numerical (Scheel & Schumacher 2017) studies in cylindrical cells with $\varGamma = 1$ and $\varGamma =2$ come to smaller exponents. Zürner et al. (2019) found $f_{OS}/f_\kappa = (0.010\pm )\times Ra^{0.40\pm 0.02}$ at $Pr = 0.029$. The direct numerical simulation by Schumacher & Scheel (2016) at $Pr=0.021$ results in a scaling $f_{OS}/f_\kappa = (0.08\pm )\times Ra^{0.42\pm 0.02}$. Another scaling law, $f_{OS}/f_\kappa = (0.027\pm )\times Ra^{0.419\pm 0.006}$, is given by Vogt et al. (2018a).

Akashi et al. (2019) demonstrated that the ratio between the oscillation period $\tau _{OS} = 1/f_{OS}$ and the turnover time $\tau _{TO}$ is close to unity. The turnover time $\tau _{TO}$ describes the time required for a fluid parcel to complete a full circulation within the LSC structure. Except for the movement in the corners, the fluid elements follow in good approximation the path of an ellipse with the semi-axes $L/4$ and $H/2$. Using a simple approximation formula to determine the circumference of an ellipse, the turnover time can be calculated as follows:

(3.1)\begin{equation} \tau_{TO} = \frac{\pi \sqrt{2((L/4)^{2} + (H/2)^{2})}}{u_{rms}}. \end{equation}

The velocity $u_{rms}$ is calculated from the velocity data recorded by the sensor $Bottom$ as

(3.2)\begin{equation} u_{rms}=\overline{\sqrt{\frac{1}{L}\int^{L}_0u_x^{2}(x,t)\,{\textrm{d}}\kern0.06em x}}. \end{equation}

The ratio $\tau _{OS}$ to $\tau _{TO}$ is plotted versus $Ra$ in figure 6(d). All data points are in a close proximity to the value one.

The occurrence of characteristic frequency peaks is observed in many experimental studies (Castaing et al. 1989; Siggia 1994; Qiu & Tong 2001) and is apparently a characteristic feature of coherent large-scale structures occurring in thermal convection. Several studies construe these oscillations as an indication of the transport of individual plumes by the LSC (Cioni et al. 1997; Grossmann & Lohse 2000,  2002; Stevens et al. 2013). Another perspective is to associate these oscillations with the advection of deviations from an ideal two-dimensional LSC structure or other substructures by the LSC. This in a natural way explains an inherent relationship between the frequency of the oscillations and the turnover time. Thus, the torsion and sloshing phenomenon could also be understood as a deformation of the two-dimensional single-roll LSC that travels with a recurrence rate $1/\tau _{TO}$ throughout the fluid container (Zürner et al. 2019). A strong correlation between the turnover time and the periodicity of predominant oscillations has been demonstrated in cylinders with $\varGamma = 1$ (Brown & Ahlers 2009; Zürner et al. 2019) or $\varGamma = 2$ (Vogt et al. 2018a). In these geometries, where a single-roll LSC exists, the oscillatory behaviour is due to torsion and sloshing at $\varGamma = 1$ or to JRV dynamics at $\varGamma = 2$. In our case of the $\varGamma = 5$ box, a torsional movement of single-roll LSC can be excluded. Nevertheless, it is reasonable to assume that perturbations of the flow pattern or instabilities of the boundary layers are transported at the characteristic velocity of circulation within the coherent structures. This interpretation was also adopted by Cioni et al. (1997) who assumed that the dominating frequency of the oscillations reflected in their temperature signal is associated with the velocity of the global circulation and path length of the LSC. This assumption leads to the relation ${f_{OS}H^{2}}/{\kappa } \varpropto RePr$ (Cioni et al. 1997). In the next section it becomes evident that this applies also for the flow configuration studied here, since the Reynolds number and the oscillation frequency follow the same power law as a function of $Ra$.

3.3. Transport of momentum and heat

The Reynolds number $Re$ is the characteristic parameter to quantify the turbulent momentum transport in the convection cell. Our definition of $Re$ from the experiment relies on the r.m.s. velocity $u_{rms}$ determined from the linear velocity profile along the measuring line $Bottom$:

(3.3)\begin{equation} Re = \frac{u_{rms}H}{\nu}. \end{equation}

A simultaneous acquisition of all three velocity components in the entire container is difficult to achieve in the experiment. With the limitation to one velocity component on one single measuring line, $Re$ derived from this appears relatively arbitrary and the question arises as to what extent it is representative for the total flow and the momentum transfer. To check this, two $Re$ values are determined by means of the numerical simulations, one corresponding to the experiment with $u_{rms}$ determined on the $Bottom$ measurement line and the other with the r.m.s. velocity $U_{rms}$ averaged over the entire volume in accordance with (2.7). To compensate for the difference that all three components of the velocity were considered in the calculation of $U_{rms}$ instead of only one in $u_{rms}$, the corresponding values of $Re$ have been corrected by a factor of $\sqrt {1/3}$.

The scaling behaviour of the resulting $Re$ is plotted in figure 7. There is no significant difference between the two results of numerical simulation, which suggests that the time-averaged velocities determined on the $Bottom$ measurement line are representative for the global flow. A power-law fit reveals a scaling law of $Re = 3.43Ra^{0.44}$ for velocity measurements and $Re = 2.83Ra^{0.45}$ for the numerical simulations. The direct comparison between experiment and simulation reveals a slight difference in the scaling exponents, but the absolute values agree very well.

Figure 7. Scaling of the Reynolds number $Re$ versus the Rayleigh number $Ra$. Results of numerical simulations are calculated for two characteristic velocities $u_{rms}$ and $U_{rms}$ representing the r.m.s. velocity for the measuring line of sensor $Bottom$ and the entire volume, respectively. Least-squares approximation leads to $Re = 3.43Ra^{0.44}$ (grey solid line) for the experimental data, $Re = 2.83Ra^{0.45}$ for numerical simulations (grey dashed line) and $Re = 2.56Ra^{0.46}$ for the entire volume (grey dotted line).

A comparison with the previous study by Akashi et al. (2019) also reveals a slightly different scaling exponent, $u_{rms}(H/\kappa ) = 0.059Ra^{0.50\pm 0.02}$. In Akashi et al. (2019) also smaller $Ra$, where unsteady roll arrangements are observed but no cellular flow regime occurs, were included into the power-law fit. Thus, $Re$ in the cellular regime apparently shows a measurably more gradual increase with $Ra$ than in a roll regime. There are no differences between the current measurements and the results of Akashi et al. (2019), as long as only the range of the cellular regime is considered there.

Scaling exponents from $0.42$ to $0.53$ were reported by other experimental studies (Takeshita et al. 1996; Cioni et al. 1997; Yanagisawa et al. 2013,  2015; Vogt et al. 2018a; Zürner et al. 2019) and numerical simulations (Scheel & Schumacher 2017) for low-$Pr$ thermal turbulence. Yanagisawa et al. (2013) found a higher velocity magnitude and increased scaling exponent, $u_{rms}(H/\kappa ) = 0.08Ra^{0.53}$, in a square container with the same geometry filled with pure gallium corresponding to $Pr = 0.025$. Zürner et al. (2019) defined three different characteristic velocities: (1) the typical horizontal velocity magnitude near the plates $v_{LSC}$, (2) the typical vertical velocity magnitude of the LSC along the sidewall $v_{vert}$ and (3) the turbulent velocity fluctuations in the centre of the cell $v_{centre}$. This is reflected in slightly varying scaling exponents: (1) $v_{LSC} \propto Ra^{0.42\pm 0.03}$, (2) $v_{vert}\propto Ra^{0.42\pm 0.04}$ and (3) $v_{centre}\propto Ra^{0.46\pm 0.04}$ (Zürner et al. 2019). The smaller scaling exponent for the horizontal and vertical LSC is a result of probing different parts of the complex three-dimensional flow structure. It suggests that the $Ra$ dependence of the representative velocities depends on how and where it is measured in the convection cell. It seems plausible that the option of calculating the r.m.s. velocity over the entire cell can be considered as the most reliable method for evaluating the representative convection velocity. For instance, Scheel & Schumacher (2017) found a scaling $U_{rms}\propto Ra^{0.45\pm 0.01}$ in a cylinder with $\varGamma = 1$ at $Pr = 0.021$, which matches our results almost perfectly.

The heat transport of the three-dimensional cell structure is evaluated in numerical simulations by considering the Nusselt number dependence on the Rayleigh number, $Nu(Ra)$. Respective results are shown in figure 8. We find a scaling of $Nu= 0.14\times Ra^{0.26\pm 0.01}$ which is in agreement with the outcome of measurements and simulations in a cylindrical cell with $\varGamma = 1$ (Scheel & Schumacher 2017; Zürner et al. 2019). The laboratory experiments conducted by Zürner et al. (2019) for $Pr = 0.029$ provide $Nu\cong (0.12\pm 0.04)\times Ra^{(0.27\pm 0.04)}$ while the simulations by Scheel & Schumacher (2017) suggest $Nu\cong (0.13\pm 0.04)\times Ra^{(0.27\pm 0.01)}$ for $Pr = 0.021$. Recent experiments conducted in the same container by Vogt et al. (2021) show a scaling of $Nu = 0.166\times Ra^{0.25}$ and thus a satisfactory agreement with the result presented here.

Figure 8. Scaling of the Nusselt number $Nu$ versus $Ra$ obtained by numerical simulations, where the error bars represent the standard deviation of the temporal fluctuations and the grey line shows the fitted curve, $Nu = 0.14\times Ra^{0.26}$ (numerical simulation).

4. Oscillatory dynamics of the three-dimensional cell structure

4.1. Procedure of phase averaging

As shown by the flow measurements and the numerical simulations presented in figures 2, 3 and 4, the flow patterns are characterized by a variety of large- and small-scale fluctuations, which indicate the presence of both a periodically oscillating coherent flow filling the entire container and inertia-dominated turbulence on smaller scales. We are interested in the dynamics of the coherent flow structure and want to separate it from the background of the turbulent fluctuations by a suitable averaging method. Investigations of turbulent superstructures in very large aspect ratios are confronted with the same problem. Pandey et al. (2018) propose an averaging procedure where the averaging time to be used is in a range between the free-fall time $t_{ff}$ and the effective dissipation time $t_d = \max (H^{2}/\nu, H^{2}/\kappa )$. In the case of turbulent superstructures in very large aspect ratios, the patterns are slowly changing. The structure under consideration here undergoes significant oscillations with a constant periodicity. For this reason it is useful to apply a phase averaging to our numerical data. Such a procedure has already been suggested by Vogt et al. (2018a) to reveal the nature of the dynamic flow structures of the periodic oscillations of velocity and temperature in the $\varGamma = 2$ cylinder. The analysis covers one complete oscillation period $\tau _{OS} = 1/f_{OS}$, which in our case also corresponds to the turnover time for a respective fluid parcel (see figure 6). It is based on simulation results with a computation time comprising 16 oscillation periods with a total number of 5536 snapshots. The entire dataset is divided into 16 parts, each corresponding exactly to one oscillation period. This time period of one oscillation is again divided into 16 equidistant intervals, also called phases. Each interval or phase includes 21 snapshots of temperature, velocity and pressure. In the next step, these snapshots are averaged. Finally, for the distributions of temperature, velocity and pressure obtained for all oscillation periods are conditionally averaged with the respective distributions belonging to the same phase.

Vogt et al. (2018a) also had to adequately consider the orientation of the LSC in a convection cell for each time step, since the LSC symmetry plane in the cylindrical $\varGamma = 2$ cylinder wanders in the azimuthal direction. This did not have to be taken into account in our case, since the flow is approximately fixed in its orientation by the square geometry of the convection cell.

4.2. Variations of temperature and flow field during the oscillation period

The results presented in this section are restricted to $Ra = 1.2\times 10^{5}$ as an example. The flow at other $Ra$ shows a qualitatively similar behaviour. The numerical calculations of temperature and flow fields are performed for a duration of 16 oscillations which corresponds to about 300 free-fall times $t_{ff}$ and 5 thermal diffusion times $t_\kappa = H^{2}/\kappa$ with a grid resolution of $n_x \times n_y \times n_z = 400 \times 400 \times 80$. Figures 9(a) and 9(b) show the phase-averaged temperature distribution for the horizontal plane at the mid-height of $z = 0.5$ at two different phases. The temperature scale refers to the mean temperature $T_m = 0.5$ in the considered plane in such a way that the regions of the warm fluid marked by red coloration have a temperature above $T_m$, while colder fluid appears in blue coloration in the plots. We observe a distribution of temperature that resembles a checkerboard pattern, where the cold fluid occurs in the centre and the four corners of the container, while the warm fluid is located in centrally situated areas at the four sidewalls. This downwelling flow pattern is opposite to the upwelling one shown in the figure 5. The simulations reveal that the boundaries between the areas of warm and cold fluid which are identical to the isolines $T = T_m = 0.5$ are subject to significant deformation during the oscillations (see also supplementary movie 5). To quantify this behaviour in a simple way, we follow the displacements of some specially defined locations in the temperature distribution during one oscillation period. For this purpose we have selected the following positions: the intersection points of the $T_m$ isolines with the central lines of the centre plane along the $x$ axis (white triangles) and the $y$ axis (yellow triangles) as well as at the intersection points of the $T_m$ isolines with the sidewalls. The latter are marked by the grey and white squares in the $x$ direction and by the blue and white circles along the $y$ direction. Each of these marker points on the isolines can be assigned an $x$ coordinate $X_{T_m}$ and a $y$ coordinate $Y_{T_m}$. In the following, the movements of these points along the $x$ and $y$ centre lines and along the sidewalls are recorded. Our analysis takes into account that all marker points have only one degree of freedom and can change their position either exclusively along the $x$ or $y$ direction. The two snapshots in figures 9(a) and 9(b) show in a certain sense two extreme states in which the respective pairs of points have either a minimum or maximum distance to each other. The oscillation of the flow pattern is reflected accordingly in the coordinates of these markers.

Figure 9. Phase-averaged temperature distribution for $Ra = 1.2\times 10^{5}$ in the centre plane at (a) $t = 0$ and (b) $t = 0.5\tau _{OS}$; (c) time variations of the $x$ and $y$ positions of the mean temperature of $T_m = 0.5$ on the centre lines and the sidewalls of the container in the temperature distribution shown in (a,b) (numerical simulation).

This is demonstrated in figure 9(c) which displays the varying locations of the marker points during one oscillation period. A sinusoidal behaviour is detected, where the coordinates of each point accomplish a full oscillation in one cycle. It is noticeable that the selected pairs of markers always move in opposite directions; this concerns the white and yellow triangles on the centre lines as well as the points on the differently aligned side edges. If, for example, one considers the associated deformation of the fluid regions of higher temperature, this means that the approximately rectangular fields alternately retract to the area in front of the sidewall or extend into the central zone. For further consideration we have defined the phase shown in figure 9(a) as the starting point $t = 0$ of the oscillation period. Because of the symmetry of the configuration, both states in figures 9(a) and 9(b) are quasi-identical, since the marker points reach their extreme positions here. The distances between all marker points are nearly equal at intermediate times $t = 0.25\tau _{OS}$ and $t = 0.75\tau _{OS}$.

Figures 10(a)–10(d) show the cross-section of the phase-averaged temperature distribution, streamlines and velocity vectors in two vertical centre planes of the container at $y = 2.5$ and $x = 2.5$. The colour map represents the temperature of the fluid. The length of the velocity vectors indicates the magnitude of the velocity. The streamlines and the velocity vectors visualize two LSCs lined up side by side, where the fluid sinks in the centre of the cell and rises on the sidewalls. Size and location of the two dominating vortices are not constant but subject to permanent changes. The positions of the rotational axes of the vortices are not stationary, but move in the shown sectional planes. Naturally, this has significant consequences for the distribution of cold and warm fluid, which is reflected in the temperature colour plot (see also supplementary movie 6). The changes in the size of the respective zones of cold and warm fluids are clearly visible in the plots and confirm the dynamics already described in the context of figure 9. It should also be pointed out here that the temperature distribution does not only change in the centre of the convection box, but that significant alterations also occur near the bottom and the lid. We focus on this aspect in more detail later in § 4.4.

Figure 10. Phase-averaged temperature distribution with velocity vectors and streamlines for $Ra = 1.2\times 10^{5}$ in the cross-section of the centre lines of $y = 2.5$ (left column) and $x = 2.5$ (right column) at different phases: (a) $t=0$, (b) $t = 0.25\tau _{OS}$, (c) $t = 0.5\tau _{OS}$ and (d) $t = 0.75\tau _{OS}$. (e) Temporal change of representative velocity $u_{rms}$ calculated along the measurement lines of $z = 0.25$ and $z = 0.75$ in the cross-sections of $y = 2.5$ and $x = 2.5$ shown in (a–d) (numerical simulation).

For further evaluation, we consider the temporal behaviour of the r.m.s. velocities $u_{rms}$ recorded along two horizontal lines at $z = 0.25$ and $z = 0.75$ in both cross-sections. Time variations of $u_{rms}$ during one oscillation period are plotted as shown in figure 10(e). The velocities $u_{rms}$ monitored along the various horizontal lines reach their maximum (minimum) values at $t = 0.25\tau _{OS}$ or half a period later at $t= 0.75\tau _{OS}$. As already seen in the evaluation of the temperature distribution in figure 9(c), a sinusoidal behaviour of the velocities is observed. Furthermore, it is immediately noticeable that the velocities at the different heights in one plane and the velocities at the same height in different planes vary in opposite directions. This corresponds in an excellent way with the velocity measurements in the experiment (see Akashi et al. 2019). The periodic up and down of the velocity magnitude near the lid or bottom is related to an oscillatory shift of the centre of the dominating vortices in both vertical and horizontal directions. In figures 10(a)–10(d) it can be seen that smaller secondary vortices can form temporarily on the sidewalls if the axes of the larger vortices are located near the centre of the convection container. The vertical shift of the vortex positions lets the zone around the horizontal centre line alternately come into the effective range of the horizontal flow pointing towards the centre or towards the sidewalls, respectively. This explains the flow patterns with periodic reversals of flow direction found in the mid-height of the convection cell by both measurements and simulations (see figures 3b, 3c and 4b).

4.3. Delineation of the three-dimensional flow structure

In the previous section, it was shown that the location and properties of dominant vortices filling the entire height of the fluid layer change significantly during a period of oscillation.The position of a vortex centre can be determined by its axis of rotation, which is characterized by minimum pressure. The phase-averaged pressure distribution has been calculated for this purpose. Snapshots of the pressure fields in both vertical sectional planes are displayed in figures 11(a) and 11(b) for $t = 0.25\tau _{OS}$. Here, the positions of the minimum pressure are marked with dark blue coloured circles. The trajectories of these vortex axes during one oscillation period are presented in the figures 11(c) and 11(d). In addition to the information regarding the position of the pressure minimum in the respective sectional plane, the colouring of the circles in this diagram also indicates the time when the vortex centre was at which position. The trajectories of the axes of the LSCs form the central segment of a diagonal line running approximately between the lower corners of the convection cell and the centre of the upper plate, with the forward and return paths not exactly on the same line. It becomes also obvious that the two circulation centres always move in opposite directions in their respective sectional planes with an offset of half a period.

Figure 11. Phase-averaged pressure contour for $Ra = 1.2\times 10^{5}$ in the cross-section of (a) $y = 2.5$ and (b) $x = 2.5$ at the phase of $t = 0.25\tau _{OS}$. The minimum pressure points in the phase-averaged pressure contour during one oscillation period in the cross-sections of (c) $y = 2.5$ and (d) $x = 2.5$ (numerical simulation).

Figure 12 is intended to illustrate the movement of the circulation areas in three dimensions by means of the streamlines in the vicinity of the vortex centres. To visualize the third component in the vertical direction, the streamlines are coloured according to the vertical position of the structure. A first thing to notice is that four recirculating swirls exist in the convection cell. These four connected swirls are arranged as two pairs, with the axes of the respective pair aligned on average mainly along the $x$ and $y$ directions, respectively. Figures 12(a) and 12(b) show two snapshots at $t = 0.25\tau _{OS}$ and $t = 0.75\tau _{OS}$, respectively. These are the moments when the swirls are close to either the bottom or the top of the convection cell (see figures 10 and 11). A corresponding animation of these results (supplementary movies 7–9) during one oscillation period shows, in fact, that we are dealing here with the phenomenon of a three-dimensional flow structure that shows the same characteristics as JRV, which was observed and described first by Vogt et al. (2018a). The main difference from the geometry of a $\varGamma = 2$ cylinder, in which the dynamics of a JRV extends over almost the entire convection cell, is that four such structures exist in our $\varGamma = 5$ configuration, whose movements are restricted to each partial area of the fluid layer. Supplementary movie 10 shows that the phase-averaged temperature and velocity fields in the centre plane and the vertical cross-sections of the limited area are in very good agreement with the conditionally averaged temperature and velocity fields in the centre plane, the LSC symmetry plane and the perpendicular plane of the LSC symmetry plane in the cylinder of $\varGamma = 2$ (Vogt et al. 2018a). The direction of motion of the vortex centre is clockwise, whereas the direction of flow is counterclockwise, and vice versa. It is interesting to note that the revolution of these swirls in the convection cell occurs in the opposite direction to the main flow of the LSCs. The direction of motion of the vortex centre is clockwise, whereas the direction of flow is counterclockwise, and vice versa. Another peculiarity is that the trajectories of the circulation centres obviously do not follow a circle, as one would normally expect for a typical jump rope, but through a tilted ellipse-like trajectory as visualized in figure 11(c,d). One reason for this might be the fact that the limited height of the fluid layer imposes a geometric constraint on the dynamics of the JRV.

Figure 12. Three-dimensional visualization of phase-averaged streamlines at (a) $t = 0.25\tau _{OS}$ and (b) $t = 0.75\tau _{OS}$ for $Ra = 1.2\times 10^{5}$. The dashed lines mark the two centre planes $x = 2.5$ and $y = 2.5$, which also represent the symmetry axes of the flow pattern. The colour represents the vertical position of the structure in the container (numerical simulation).

4.4. Impact on heat transport

The temperature distributions presented in § 4.2 demonstrate that considerable redistribution of cold and warm fluid occurs in the course of the three-dimensional oscillations. In this section, we consider the extent to which this has measurable effects on the heat transport.

To verify the steadiness of the heat transport, Nusselt numbers are calculated from the phase-averaged temperature distributions. Figure 13 presents the phase-averaged $Nu$ numbers during one oscillation period. The results reveal that $Nu$ is not constant but rather shows an oscillatory behaviour with a verifiable difference between the maximum and minimum value. Moreover, at the first glance it is a very surprising finding that the period of the $Nu$ oscillation is only half as long as that of the fluctuations in the temperature field or the horizontal velocities.

Figure 13. Phase-averaged values of the Nusselt number $Nu$ during one oscillation period for for $Ra = 1.2\times 10^{5}$ (numerical simulation).

Figure 14 shows contour plots of phase-averaged local heat flux. These plots represent local heat flux calculated at the top boundary $F_{z, top}=- ({\partial T}/{\partial z} )_{z=1}$ and the bottom boundary $F_{z, bot}=- ({\partial T}/{\partial z} )_{z=0}$. Here, the average number of local heat flux for an entire plate means the Nusselt number of the convection. The red-coloured zones feature the local heat flux that is above the mean value of the entire area and period ($Nu \approx 3.0$), while local heat flux below the mean value is marked in blue colour. Here, the contrasting behaviour of the local heat flux at the bottom and the top is obvious which is caused by the flow that either impinges or detaches in the corresponding areas of the plates. The local heat flux distributions change during the oscillation period and agree to a certain degree with the temperature distributions presented in figure 9. According to figure 13, $Nu$ reaches a maximum value at time points $t = 0$ and $t = 0.5\tau _{OS}$, while minima occur at $t = 0.25\tau _{OS}$ and $t = 0.75\tau _{OS}$. The corresponding local heat flux distributions on the plates are geometrically similar. This is in agreement with the periodicity of the $Nu$ number oscillations.

Figure 14. Phase-averaged local heat flux at the top boundary $F_{z, top}$ (a–d) and at the bottom boundary $F_{z, bot}$ (e–h) for $Ra = 1.2\times 10^{5}$, at different phases: (a,e) $t=0$, (b,f) $t = 0.25\tau _{OS}$, (c,g) $t = 0.5\tau _{OS}$ and (d,h) $t = 0.75\tau _{OS}$. The contour interval is 0.5 (numerical simulation).

The heat transfer of the system is mainly determined by the conductivity and stability of the thermal boundary layers (Ahlers et al. 2009; Iyer et al. 2020). Our results presented in figure 10 suggest that the turbulent fluid motion might have a considerable impact on the extent of the thermal boundary layers. Indeed, our data also show that the instantaneous thermal boundary layer thickness at both the bottom and the top plates $\lambda _{\varTheta }$ oscillates with the same frequency as the JRV. Owing to the relation $Nu = 1/(2\lambda _{\varTheta })$, this is reflected accordingly in the changes of $Nu$ illustrated in figure 13.

Figure 15(a) compares time series of $Nu$ with time series of the horizontal velocity $u_{rms}$ obtained along the line $Bottom$ and the volume-averaged vertical velocity for the entire container $U_{z}$ obtained by numerical simulations. All data show quasi-periodic oscillations. The signals of $Nu$ and $U_{z}$ are synchronized with each other, i.e. they oscillate with the same frequency, while $u_{rms}$ exhibits twice the oscillation period. This is confirmed by the corresponding PSDs presented in figure 15(b). Both PSDs of $Nu$ and $U_{z}$ show dominant peaks at the same frequency, which is almost exactly twice the oscillation frequency $f_{OS}$. Here we see an effect that, on closer inspection, is already apparent in figure 10. Figure 10(e) reveals that the approach of the centre of the circulation to one of the horizontal interfaces suppresses the horizontal flow there, while the horizontal flow on the opposite side simultaneously reaches a maximum. The evaluation of the numerical data in the whole convection cell shows that the vertical flow component, on the other hand, becomes greatest when the centre of the circulation moves either from bottom to top or in the opposite direction. This vertical momentum transport thus additionally contributes to the heat transport.Accordingly, figure 15(a) demonstrates a simultaneous occurrence of maximal values for $Nu$ and $U_z$.

Figure 15. (a) Time variations of $Nu$, the volume average of the vertical velocity component $U_z$ and the horizontal velocity $u_{rms}$ obtained by the line $Bottom$ at $z = 0.25$ and (b) the PSDs calculated from the temporal fluctuations of $Nu$, $U_z$ and the PSD calculated from the horizontal velocity $u_{rms}$ obtained by the line $Bottom$ for $Ra = 1.2\times 10^{5}$ (numerical simulation).

5. Summary and conclusions

The present study combines flow measurements in a Rayleigh–Bénard set-up and corresponding numerical simulations. The flow in the cuboid container with an aspect ratio of $\varGamma = 5$ is fully turbulent at a Prandtl number of 0.03 and Rayleigh numbers in the range $Ra$ of $6.7\times 10^{4} \leqslant Ra \leqslant 3.5\times 10^{5}$. In this parameter range the thermally driven convection produces a coherent cellular flow structure that is extremely stable with respect to changes in topology and orientation. On the other hand, this cellular flow pattern is subject to distinct oscillations. The change in fluid velocity during the oscillations is quite dramatic: the local velocity value varies during one period between a maximum value of the measured velocity and a value close to zero; in some regions even a periodic flow reversal occurs.

Our analysis succeeded in clarifying the complex three-dimensional structure of the large-scale flow and identifying a multiple JRV structure that is the underlying nature of the oscillations. The JRV structure was originally discovered by Vogt et al. (2018a) in a $\varGamma = 2$ cylinder. As the name suggests, this is a vortex with a curved axis performing a cycling motion in three-dimensional space. Our results confirm the existence of such kind of flow structure also in a cuboid container with $\varGamma = 5$. This finding is most remarkable, since it shows that the structure of the JRV is not a rare phenomenon that occurs only in special geometrical arrangements. Instead, we showed that it is a more common property of turbulent convection.

The striking feature of the structure detected here is that we find four vortices at the same time, with the symmetry axes of the circling motion being aligned parallel to the sidewalls of the fluid container in each case. One full revolution of the JRV corresponds exactly to one periodic length of the oscillation of the cellular flow structure. The dynamics of the JRV is strongly correlated and characterized by the fact that the respective opposite vortices move towards the centre or the sidewall at the same time, while the other pair of vortices rotates with an offset of half a period. The JRVs interact significantly with the thermal boundary layers, the thickness decreasing each time the centres of the JRVs are in close proximity in the immediate vicinity of the respective copper plate. This also leaves detectable traces in the heat transfer in the form of an oscillation of the Nusselt number. Its periodicity is related to the dynamics of the JRVs. Furthermore, instead of following an almost circular path as observed in the $\varGamma = 2$ cylinder, the vortices move along an elongated and obliquely inclined ellipse. Thus, it is obvious that the movement of the jump rope is affected by the limited height of the fluid layer.

The previous study by Akashi et al. (2019) revealed an increase in the characteristic length scale of the flow patterns with increasing $Ra$. In the slightly supercritical region directly after the onset of convection, roll structures are formed whose horizontal wavelength corresponds approximately to the height of the fluid layer $H$. As $Ra$ grows, the number of rolls decreases, and thus the wavelength of the structures increases incrementally. Finally, there is a qualitative transition to the cellular regime considered here; the dimension of the cellular structure is ${\sim }5H$. An increase in the characteristic length scale of coherent large-scale velocity structures with $Ra$ has also been reported in earlier studies, in particular for turbulent superstructures (Fitzjarrald 1976; Hartlep et al. 2003; Pandey et al. 2018; Krug et al. 2020). Numerical studies identified the typical sizes of turbulent superstructures in large aspect ratios $\varGamma > 10$, which are of the order of 10$H$ for $Pr = 6.7$ (Busse 1994; Pandey et al. 2018) or $6H$ to $7H$ for $Pr = 0.7$ (Hartlep et al. 2003; Pandey et al. 2018; Stevens et al. 2018). Sakievich et al. (2016,  2020) refer to these ’natural’ dimensions of superstructures to explain the dominance of specific Fourier modes in a flat $\varGamma = 6.3$ cylinder at $Pr = 6.7$. From the data published by Pandey et al. (2018), it is evident that the horizontal wavelength of the superstructures in the range $Pr \leqq 10$ decreases when approaching smaller $Pr$. If we extrapolate the above estimates to the $Pr$ of 0.03 considered here, we achieve approximately a size of 5$H$, a value that is found for the cellular regime. However, due to the limitation imposed by the finite $\varGamma =5$ box, it is not possible to verify whether the velocity field cannot form larger structures.

In summary, the cellular flow regime we studied here exhibits some general properties that are also attributed to turbulent superstructures (Pandey et al. 2018; Stevens et al. 2018; Krug et al. 2020):

1. (i) The characteristic length scale of the flow structures increases with increasing $Ra$ as shown by Akashi et al. (2019).

2. (ii) The horizontal length scale is distinctly larger than the height of the flow domain.

3. (iii) The characteristic dynamic time scale of the evolving or altering large-scale flow structure is much longer than the free-fall time.

In this general respect, our study is comparable to the work of Sakievich et al. (2016,  2020), although the observed flow patterns differ from each other, which we believe is due to the differences in the shape of the convection vessels (cylinder versus box). In particular, it is not surprising that the flow in the box is restricted in its dynamics and does not form a rotational symmetry as shown by Sakievich et al. (2020). The cellular flow regime, instead, shows mirror symmetry with respect to the two vertical centre planes $x = 2.5$ and $y = 2.5$ of the convection box, respectively. This also reveals that the aspect ratios investigated in both studies (those of Sakievich et al. and ours) are still too small to produce flow structures that can develop unaffected by the sidewalls. New experiments with larger aspect ratios are needed to address this issue. Therefore, future work should consider further investigations in a wide parameter range of aspect ratios. In particular, it would be very interesting to find out to what extent the cellular structures studied in detail here also shape the flow in turbulent superstructures at large aspect ratios and whether such oscillatory dynamics with structures resembling JRVs can also be observed there.