2.1 Global angular impulse

The global angular momentum

serves as a measure of organised rotation (see for instance Molenaar et al. Reference Molenaar, Clercx and Van Heijst2004). Figure 1 shows a time series of the normalised angular momentum $L_{2D}/\overline{|L_{2D}|}$ . Two different regimes are observed. Blue areas correspond to periods of time where $L_{2D}$ changes sign spontaneously over time: positive (respectively negative) peaks in $L_{2D}$ alternate that are associated with a dominant counter-clockwise (respectively clockwise) central vortex. The blue areas consisting of a sequence of consecutive transitions is hereafter called the consecutive reversal (CR) regime. Outside this regime, the LSC is no longer well defined and one observes an extended cessation. Such complementary region is denoted here as the extended cessation (EC) regime.

Figure 1. Time evolution of $L_{2D}(t)/\overline{|L_{2D}|}$ shown for ( $Ra=5\times 10^{7}$ , $Pr=3$ ). Light blue areas correspond to a consecutive reversal (CR) regime, while blank areas correspond to an extended cessation (EC) regime. Some events (darker blue areas) may not be clearly assigned to the CR regime, see text. The two continuous lines correspond to the thresholds used by the filtering procedure. Value of the normalised standard deviation: $\unicode[STIX]{x1D70E}(|L_{2D}|)/\overline{|L_{2D}|}=0.499$ .

In order to differentiate in a precise manner the CR regime from the EC regime, a filtering algorithm has been devised which is modelled after (Lhuillier et al. Reference Lhuillier, Hulot and Gallet2013; Podvin & Sergent Reference Podvin and Sergent2015). We identify the set of consecutive times $r_{i}$ at which $L_{2D}$ changes sign. The time interval $[r_{i},r_{i+1}]$ is considered to be inside the CR regime if, during this interval, the value of $|L_{2D}|$ reaches at least once the threshold value $\overline{|L_{2D}|}+\unicode[STIX]{x1D70E}(|L_{2D}|)$ (light blue area in figure 1). A time interval where such threshold is not reached can be of two kinds corresponding to the darker blue areas or the white areas in figure 1. The first kind is sandwiched between two CR intervals and corresponds to a ‘rogue’ reversal, which likely belongs to the CR regime but has been filtered out by our criterion (the criteria for the selection of events in the CR regime is rather stringent as seen from the ‘rogue’ events displayed in figure 1). The second type, displayed in white, corresponds to an extended cessation.

For any interval $[r_{i},r_{i+1}]$ , its duration $\unicode[STIX]{x1D70F}_{1,i}\equiv r_{i+1}-r_{i}$ is also computed. When both intervals $[r_{i-1},r_{i}]$ and $[r_{i},r_{i+1}]$ are inside a CR regime, the duration $\unicode[STIX]{x1D70F}_{d,i}$ of the jump occurring around time $r_{i}$ between a clockwise and counter-clockwise central vortex or vice versa can be evaluated. It is computed by identifying the times located just before and just after time $r_{i}$ such that $|L_{2D}|$ reaches the threshold value $\overline{|L_{2D}|}-\unicode[STIX]{x1D70E}(|L_{2D}|)$ . $\unicode[STIX]{x1D70F}_{d,i}$ is simply the time lapse between these two events.

The change of the global angular momentum $L_{2D}$ may be better understood considering the relation directly obtained from the governing equation (2.2),

where $\boldsymbol{n}$ stands for the outwards unit normal vector to the domain boundary and $\text{d}l$ for a contour line differential element. Note that it is assumed that line integrals are performed in a counter-clockwise direction. The angular momentum thus evolves because of a bulk forcing term $M(t)$ known as the input torque (Molenaar et al. Reference Molenaar, Clercx and Van Heijst2004) and two boundary integral terms $I_{a}(t)$ and $I_{b}(t)$ . $I_{b}(t)$ is close, but not identical, to the integrated vorticity flux over the domain boundary. For a square cavity, the boundary term $I_{a}(t)$ simplifies to $I_{a}=(1/2)Pr\,Ra^{-0.5}\oint \unicode[STIX]{x1D714}\,\text{d}l$ . Vorticity on the boundary is related to the friction exerted by the fluid on the walls. This integral $I_{a}(t)$ is thus quantifying the friction along the boundary.

2.2 Mechanical energy balance

Other global quantities are useful to characterise at each time the instantaneous state of the system: the global kinetic energy $E_{kin}(t)\equiv (1/2)\int \boldsymbol{u}^{2}\,\text{d}x\,\text{d}y$ and the global potential energy $E_{pot}(t)\equiv -Pr\int y\unicode[STIX]{x1D703}\,\text{d}x\,\text{d}y$ . This latter quantifies the energy required to bring all fluid particles against gravity from their position at time $t$ to the reference level $y=0$ . However, one can introduce a more pertinent instantaneous quantity, namely the available potential energy (Sutherland Reference Sutherland2010), which is defined below. For a given time $t$ , the fluid is characterised by an instantaneous temperature field $\unicode[STIX]{x1D703}(x,y,t)$ with a lower (respectively upper) bound at $\unicode[STIX]{x1D703}_{min}$ (respectively $\unicode[STIX]{x1D703}_{max}$ ). Let us consider a one-to-one mapping $(x_{r}(x,y),y_{r}(x,y))$ from the square onto the square. This may be interpreted as a reordering of fluid particles inside the square cell. Rearranging modifies the temperature field leading to a new field $\unicode[STIX]{x1D703}_{r}(x,y,t)$ but this process is an adiabatic one i.e. $\unicode[STIX]{x1D703}_{r}(x_{r},y_{r},t)=\unicode[STIX]{x1D703}(x,y,t)$ . As a consequence, the probability distribution function (PDF) of temperature in the rearranged state is identical to the PDF in the instantaneous state. The potential energy of the rearranged state can be measured. Among the set of such mappings, there exists a subset which corresponds to the lowest potential energy. It is easy to understand that all mappings belonging to this subset have identical rearranged temperature field $\unicode[STIX]{x1D703}_{r}(y,t)$ , which does not depend on $x$ and monotonically increases with height (figure 2). This field characterises the background state. Although $y_{r}(x,y,t)$ is not a simple function of $(x,y)$ , the above remark implies that $y_{r}$ is a one-to one function of $\unicode[STIX]{x1D703}$ for a mapping in this subset. In practice (see Tseng & Ferziger Reference Tseng and Ferziger2001), the field $\unicode[STIX]{x1D703}_{r}(y,t)$ is computed by the following procedure. Firstly, the PDF of the instantaneous temperature at time $t$ , which is denoted by $P(\unicode[STIX]{x1D703})$ , is directly evaluated numerically within the interval $[\unicode[STIX]{x1D703}_{min},\unicode[STIX]{x1D703}_{max}]$ since $\unicode[STIX]{x1D703}(x,y,t)$ is known on the whole square box. Secondly, the conservation of temperature PDF with rearrangement imposes that $y_{r}(\unicode[STIX]{x1D703})$ be evaluated as a cumulative density function

This relation depends a priori on the domain geometry. For a square box of unit size, however, the proportionality is factor is reduced to unity. Lastly, function $\unicode[STIX]{x1D703}_{r}(y,t)$ is then obtained by a simple inversion of $y_{r}(\unicode[STIX]{x1D703},t)$ .

Figure 2. Temperature field $\unicode[STIX]{x1D703}(x,y,t)$ at a given instant $t$  (a), the corresponding background state $\unicode[STIX]{x1D703}_{r}(y,t)$  (b) and height $y_{r}(x,y,t)$  (c) for a square RB cell.

The background state is characterised by the lowest potential energy that can be reached by an adiabatic process starting from the instantaneous temperature field $\unicode[STIX]{x1D703}(x,y,t)$ . This quantity, called the background potential energy, is equal to $E_{bpot}\equiv -Pr\int y_{r}(x,y,t)\unicode[STIX]{x1D703}(x,y,t)\,\text{d}x\,\text{d}y$ . By a simple change of variable and using adiabaticity $\unicode[STIX]{x1D703}_{r}(x_{r},y_{r},t)=\unicode[STIX]{x1D703}(x,y,t)$ , one gets

The difference $E_{apot}(t)\equiv E_{pot}(t)-E_{bpot}(t)>0$ in potential energy between the instantaneous state and its background companion is called the available potential energy and represents the potential energy which could be effectively transformed from the instantaneous field $\unicode[STIX]{x1D703}(x,y,t)$ into motion (Lorenz Reference Lorenz1955; Winters et al. Reference Winters, Lombard, Riley and D’Asaro1995).

In analogy with $L_{2D}$ , the process may be better grasped by considering the evolution of the energies $E_{kin}$ , $E_{pot}$ , $E_{apot}$ , through some exact relations (see Winters et al. Reference Winters, Lombard, Riley and D’Asaro1995, Hughes et al. Reference Hughes, Gayen and Griffiths2013). For the kinetic energy $E_{kin}$ , the following relation holds ( $\unicode[STIX]{x1D626}_{ij}$  denotes the symmetric velocity gradient tensor)

The first bulk term $\unicode[STIX]{x1D6F7}_{y}(t)$ is a convective heat flux. More precisely, let us introduce the volume-averaged Nusselt number $Nu_{vol}\equiv Ra^{0.5}\langle v\unicode[STIX]{x1D703}\rangle _{vol}-\langle \unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D703}\rangle _{vol}$ . It is easily found that, for RB cells, $\unicode[STIX]{x1D6F7}_{y}(t)=Nu_{vol}(t)-1$ . The second bulk term $\unicode[STIX]{x1D716}(t)>0$ stands for the viscous dissipation rate. Finally one may write

The potential energy $E_{pot}$ verifies instead the relation

which contains the bulk term $Nu_{vol}$ and a boundary term

quantifying the conversion rate to $E_{pot}$ from external sources. More precisely, let us introduce the Nusselt number $Nu_{top}(t)\equiv -\int \unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D703}\,\text{d}x$ evaluated at the top $y=0.5$ as well as the Nusselt number $Nu_{bot}(t)\equiv -\int \unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D703}\,\text{d}x$ evaluated at the bottom plate $y=-0.5$ . For the present square cell, one easily verifies that $\unicode[STIX]{x1D6F7}_{b1}=(Nu_{top}+Nu_{bot})/2$ and consequently

Finally, the evolution equation for the available potential energy $E_{bpot}$ reads as

where the bulk term $\unicode[STIX]{x1D6F7}_{d}(t)$ quantifies the energy conversion rate due to diapycnal mixing. Since by definition $\unicode[STIX]{x2202}y_{r}/\unicode[STIX]{x2202}\unicode[STIX]{x1D703}>0$ , $\unicode[STIX]{x1D6F7}_{d}(t)$ is bound to be positive. The boundary term $\unicode[STIX]{x1D6F7}_{b2}$ provides the conversion rate from external sources. For the present RB cells, since $y_{r}(x,1/2,t)=-1/2$ and $y_{r}(x,-1/2,t)=1/2$ and because adiabaticity of side walls it is clear that $\unicode[STIX]{x1D6F7}_{b2}(t)=\unicode[STIX]{x1D6F7}_{b1}(t)$ . Finally by subtracting (2.9) by (2.12), one gets

5.1 Averaging procedure and generic reversal as function of $(Ra,Pr)$

Figure 5. From (a) to (i): normalised global angular impulse $L_{2D}/\overline{|L_{2D}|}$ , normalised kinetic energy $E_{kin}/\overline{|E_{kin}|}$ and normalised available potential energy $E_{apot}/\overline{|E_{kin}|}$ for $(Ra=5\times 10^{7},Pr=4.3)$ . (a,d,g) Each reversal cycle is centred and its time is rescaled by $\overline{\unicode[STIX]{x1D70F}}_{1}$ (only $10$ reversals are displayed and each colour is a different reversal); (b,e,h) each reversal cycle is centred and its time is rescaled by $\unicode[STIX]{x1D70F}_{1,i}$ ; (c,f,i) average value of rescaled curves obtained from the complete time series (thick lines) and curves corresponding to one standard deviation (dashed lines).

In order to identify similarities between different intervals in the CR regime, the following procedure is proposed to treat time series of global quantities such as $L_{2D}$ , $E_{kin}$ and $E_{apot}$ . Once the intervals $[r_{i},r_{i+1}]$ inside the CR regime are properly identified, all these intervals with $L_{2D}>0$ (respectively $L_{2D}<0$ ) are stacked together so that they have a common origin at $r_{i}$ (respectively $r_{i+1}$ ). If the time axis of each interval is rescaled by $\overline{\unicode[STIX]{x1D70F}_{1}}$ , one obtains figure 5(a,d,g). If the time axis of each interval is rescaled by its particular duration $\unicode[STIX]{x1D70F}_{1,i}$ , these curves display a consistent dynamical pattern (see figure 5 b,e,h). Note that, while figure 5 displays only 10 reversals to avoid cluttered graphs, these events displayed are considered as representative of the entire set. Obviously, all of the events inside the CR regime are taken into account in our procedure but ‘rogue reversal’ events are not. This procedure is similar to one used in the study of statistical properties of magnetic switches in the geodynamo problem (Valet et al. Reference Valet, Fournier, Courtillot and Herrero-Bervera2012; Lhuillier et al. Reference Lhuillier, Hulot and Gallet2013). For a sufficiently large number of recorded events, the average over these rescaled curves is expected to remove the noisy dynamics and to represent a generic reversal cycle (figure 5 c,f,i). This averaging technique, once applied to $E_{kin}$ and $E_{apot}$ , recovers the evolution of mechanical energies during the generic reversal cycle.

Figure 6. Curves corresponding to reversals for $Pr=3.0$ as obtained using the procedure described in figure 5. The average $L_{2D}/\overline{|L_{2D}|}$ is shown in thick lines and one standard deviation in dashed lines. From (a) to (c): ( $a$ )  $Ra=10^{7}$ , ( $b$ )  $Ra=5\times 10^{7}$ and ( $c$ )  $Ra=10^{8}$ .

Figure 7. Conditionally averaged fields $\boldsymbol{u}^{o}(x,y,t)$ and $y_{r}(\unicode[STIX]{x1D703}^{o}(x,y,t))$ at different instants during the generic reversal cycle for $(Ra=5\times 10^{7},Pr=3.0)$ . Fields are obtained as the ensemble average over 83 particular reversal cycles. Streamlines of the velocity field $\boldsymbol{u}^{o}$ are superposed over the colour map of field $y_{r}(\unicode[STIX]{x1D703}^{o})$ . Solid and dashed streamlines indicate the two senses of rotation.

Figure 6 shows the $L_{2D}/\overline{|L_{2D}|}$ curves for generic reversal cycle at $Pr=3.0$ : as we increase $Ra$ , the reversal cycle becomes more regular and the band representing the standard deviation narrows. The same averaging procedure can be applied to instantaneous temperature and velocity fields in order to obtain the evolution of a conditionally averaged temperature $\unicode[STIX]{x1D703}^{o}(x,y,t)$ and velocity $\boldsymbol{u}^{o}(x,y,t)$ fields during the reversal cycle (e.g. figure 7 for $(Ra=5\times 10^{7},Pr=3.0)$ ). Despite fluctuations between various realisations, dominant and persistent structures appear at specific times of the generic cycle.

5.2 Phases of reversal cycle

Figure 8. Generic reversal cycle for $(Ra=5\times 10^{7},Pr=3.0)$ . (a) System trajectory in phase space $(L_{2D}/\overline{|L_{2D}|},E_{kin}/\overline{|E_{kin}|},E_{apot}/\overline{|E_{kin}|})$ ; (b–d) average curves for $L_{2D}/\overline{|L_{2D}|}$ , $E_{kin}/\overline{|E_{kin}|}$ and $E_{apot}/\overline{|E_{kin}|}$ . Marks (i)–(x) and (a)–(e) indicate particular instants (see text for explanations). The accumulation phase is shaded and indicated by a blue colour bar. The release (respectively acceleration) phase is indicated by an orange (respectively green) colour bar. The same colour code is used in the phase space trajectory position.

In the phase space $(L_{2D}/\overline{|L_{2D}|},E_{kin}/\overline{|E_{kin}|},E_{apot}/\overline{|E_{kin}|})$ , let us consider the generic reversal cycle (figure 8). Consecutive instants (a)–(e) pinpoints particular dynamical times: $L_{2D}=0$ at instant (a); $E_{kin}$ reaches a local maximum at instant (b), $E_{kin}$ reaches a local minimum at instant (c); $E_{apot}$ reaches its minimum at instant (d); $|L_{2D}|$ reaches its maximum at instant (e). Points ( $\text{a}^{\prime }$ )–( $\text{e}^{\prime }$ ) are similar but correspond to an opposite rotation sign. Based on these instants, three successive phases are identified for the generic reversal cycle. They are called accumulation, release and acceleration.

Figure 9. (a,b) Generic reversal curves in the plane ( $L_{2D}/\overline{|L_{2D}|}$ , $E_{kin}/\overline{|E_{kin}|}$ ). Each curve corresponds to a different $Ra$ : $Ra=3\times 10^{7}$ (black line), $Ra=5\times 10^{7}$ (blue line) and $Ra=10^{8}$ (red line). (a) Corresponds to $Pr=3.0$ and (b) to $Pr=4.3$ . (c–e) Average curves during one cycle for $L_{2D}/\overline{|L_{2D}|}$ , $E_{kin}/\overline{|E_{kin}|}$ and $E_{apot}/\overline{|E_{kin}|}$ for $(Ra=5\times 10^{7},Pr=4.3)$ . Figures are displayed as in figure 8.

The accumulation phase is located between points ( $\text{e}^{\prime }$ ) and (a). It is characterised by a steady accumulation of $E_{apot}$ and a progressive decay of $|L_{2D}|$ and $E_{kin}$ . This phase ends when $E_{apot}$ reaches a maximum and $E_{kin}$ a minimum (figure 8). In terms of generic velocity field $\boldsymbol{u}^{o}(x,y,t)$ , a central vortex is present during this phase (see figure 7i–iii) until the global rotation switches signs at point (a) i.e. $L_{2D}=0$ (figure 7iv).

The release phase located between points (a) and (d), is defined by a sudden exchange from $E_{apot}$ to $E_{kin}$ . It can be split in three substeps. The first step from point (a) to point (b) contains a rapid increase of $E_{kin}$ to a maximum value and a rapid decrease of $E_{apot}$ . It corresponds to figure 7(v–vii). A second step follows from points (b) to (c) in which $E_{kin}$ suddenly decreases and $E_{apot}$ remains almost constant. This is associated with figure 7(viii,ix). After these two steps referred to as a rebound, a new increase of $E_{kin}$ is observed from points (c) to (d) concomitantly with a decrease of $E_{apot}$ until it reaches its minimum value.

Finally, the acceleration phase is located between points (d) and (e) and is characterised by an increase of $|L_{2D}|$ and $E_{kin}$ to peak values, whereas $E_{apot}$ remains almost constant. During this period, the flow reorganises gradually into a single dominant vortex (figure 7x).

For the $(Ra,Pr)$ considered inside the CR regime, we are able to recover a generic reversal cycle expressed in terms of the available mechanical energy (figure 9). Similarities between these curves for different $(Ra,Pr)$ suggest an equivalent underlying mechanism behind flow reversals, even if the intensity of the rebound decreases with $Pr$ . For $(Ra=5\times 10^{7},Pr=3.0)$ , the accumulation phase lasts longer (60 %), while the release and acceleration phases have shorter and similar durations (respectively 18 % and 22 % of the reversal cycle). For the range of $Ra$ considered these proportions are similar: for instance, the accumulation, release, and acceleration are observed to last 75, 13 and 12 % of the reversal cycle for $(Ra=10^{8},Pr=3.0)$ .

6.1 Time evolution of the available mechanical energy and flow structures

Figure 10. Instantaneous fields for a particular reversal cycle during the accumulation phase at regular time intervals at $(Ra=5\times 10^{7},Pr=3.0)$ . The first snapshot follows point ( $\text{e}_{p}^{\prime }$ ) while the last precedes point ( $\text{a}_{p}$ ). Streamlines are superposed either to vorticity $\unicode[STIX]{x1D714}$ (a), or to kinetic energy $1/2u_{i}u_{i}$ (b).

Figure 11. Instantaneous fields for a particular reversal cycle during the accumulation phase. They are identical to those of figure 10 but streamlines are superposed either to the field $Pr(y_{r}-y)\unicode[STIX]{x1D703}$ (a) or to the field $\unicode[STIX]{x1D735}y_{r}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}$ (b). Snapshots correspond to the same instants as in figure 10.

The spatial distribution of mechanical energy can be linked to flow structures observed during different phases of a reversal cycle. At the beginning of the accumulation phase (point ( $\text{e}_{p}^{\prime }$ )), a large diagonal vortex with small counter-rotating corner flows is observed (figure 10 a). On the one hand, it corresponds to a state of maximum kinetic energy condensed inside the central vortex (figure 10 b). On the other hand, the integrand $Pr\int (y_{r}-y)$ of $E_{apot}$ is mostly distributed inside the corner flows or along thin thermal boundary layers (figure 11 a). The field $\unicode[STIX]{x1D735}y_{r}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}$ is used to highlight the contour of thermal plumes (figure 11 b). This field is related to the spatial distribution of $\unicode[STIX]{x1D6F7}_{d}$  (see (2.12)). Small-scale thermal plumes are observed to be detached from the thermal boundary layers. They are then swept by the central vortex. Plumes are channelled into the corner flows, directly or after having been advected along the side walls. The progressive growth of the corner flows is illustrated on figures 10 and  11. This is consistent with previous observations by Sugiyama et al. (Reference Sugiyama, Ni, Stevens, Chan, Zhou, Xi, Sun, Grossmann, Xia and Lohse2010) that tied the time evolution of the corner-flows heights. The steady increase in $E_{apot}$ coincides with a deceleration of the central vortex and to a build-up of thermal energy inside the corner flows. Indeed, let us compute the contributions to $E_{apot}$ due to the region of thermal boundary layers (BL) and the bulk. The hot and cold BL are taken to have a constant thickness $\unicode[STIX]{x1D6FF}_{\unicode[STIX]{x1D703}}^{-1}=2\overline{Nu}$ . Contributions inside both BL (figure 12) amount to 30–40 % of $\overline{E_{apot}}$ and are fairly constant in time: standard deviation is less than 1 %. On the contrary, contributions to $E_{apot}$ from outside BL are directly influenced by the reversal cycle: a steady increase is observed during the accumulation phase until the release phase (figure 12). This suggests that the energy exchange observed during the release takes place only inside the bulk, while the boundary layers seem to be largely unaffected by the reversals.

Figure 12. Evolution of $E_{apot}$ contained within the thermal boundary layers (blue line) and outside (black line) at $(Ra=5\times 10^{7},Pr=3.0)$ .

Figure 13. Instantaneous fields for a particular reversal cycle during the early release phase (point ( $\text{a}_{p}$ )–( $\text{b}_{p}$ )) at regular time intervals $(Ra=5\times 10^{7},Pr=3.0)$ . Streamlines are superposed to the field $1/2u_{i}u_{i}$ . First snapshot follows point ( $\text{a}_{p}$ ), last precedes point ( $\text{b}_{p}$ ).

Figure 14. Instantaneous fields for a particular reversal cycle during the release phase at $(Ra=5\times 10^{7},Pr=3.0)$ . Snapshots correspond to points ( $\text{a}_{p}$ ) to ( $\text{d}_{p}$ ). Streamlines are superposed to either the field $y_{r}(x,z,t)$ (a) or to field $\unicode[STIX]{x1D735}y_{r}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}$ (b).

During the first part of the rebound (points ( $\text{a}_{p}$ ) to ( $\text{b}_{p}$ )), the opposing corner flows have become large and strong enough to deform and finally split the central vortex (figure 13): opposing corner flows then connect and form a single vortex (with opposite rotation with respect to the previous LSC). This allows the thermal energy stored inside corners to be rapidly released and transformed into kinetic energy. By the end of this exchange, the amount of thermal energy $E_{apot}$ stored outside the thermal boundary layers will be halved (figure 12). The complete rebound period (point ( $\text{a}_{p}$ )–( $\text{c}_{p}$ )) can be better illustrated by highlighting the role of the thermal plumes and boundary layers. To do so, one uses $y_{r}(x,y,t)$ which is a bijective function of temperature, and function $\unicode[STIX]{x1D735}y_{r}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}$ in figure 14. Once the vortex reconnection has taken place, the blobs of hot (respectively cold) fluid which are now inside the main central vortex, are allowed to travel upwards (respectively downwards) directly into the bulk. This results in the exchange of potential energy (contained inside small-scale structures) into kinetic energy (contained inside the central vortex). The newly formed circulation proceeds to rotate and is able to advect thermal blobs against the action of buoyancy forces (as seen between points ( $\text{b}_{p}$ ) and ( $\text{c}_{p}$ ) figure 14). Simultaneously as the thermals are released into the bulk, the surface separating hot and cold fluid increases which reinforces the mixing process (figure 14 b).

Figure 15. Instantaneous fields for a particular reversal cycle taken at regular intervals during the acceleration phase at $(Ra=5\times 10^{7},Pr=3.0)$ . The first snapshot follows point ( $\text{d}_{p}$ ); the last corresponds to point ( $\text{e}_{p}$ ). Streamlines are superposed either to vorticity $\unicode[STIX]{x1D714}$ (a), or to $\unicode[STIX]{x1D735}y_{r}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}$ (b).

The acceleration phase is illustrated on figure 15. Self-organisation of the LSC takes place, during which the kinetic energy and the angular impulse arrive to peak values. During that phase, the bulk contains less and less plumes: temperature becomes nearly homogeneous. By the end of this phase the flow settles to a large diagonal roll with small corner flows.

6.2 Mechanical energy transfer rates

Figure 16. Time evolution of the energy transfer rates. From (a) to (c) normalised $Nu_{vol}/\overline{Nu}$ , $(\unicode[STIX]{x1D716}+1)/\overline{Nu}$ , $\unicode[STIX]{x1D6F7}_{d}/\overline{Nu}$ , $\unicode[STIX]{x1D6F7}_{b1}/\overline{Nu}$ for $(Ra=5\times 10^{7},Pr=3.0)$ . Figure for $L_{2D}/\overline{|L_{2D}|}$ is given as a reference curve. (a) Corresponds the generic reversal: the curve shown by a black line corresponds to the average, while the hatched region represent one standard deviation from this average. (b) and (c) correspond to particular reversals displayed by an orange line alongside the standard deviation.

Figure 16 displays the evolution of energy transfer rates (given in (2.7) to (2.13)) for a generic reversal cycle as well as for two particular reversal cycles. More precisely, the bulk terms $Nu_{vol}=\unicode[STIX]{x1D6F7}_{y}+1$ , $\unicode[STIX]{x1D716}+1$ , $\unicode[STIX]{x1D6F7}_{d}$ and boundary term $\unicode[STIX]{x1D6F7}_{b1}=(Nu_{top}+Nu_{bot})/2$ normalised by $\overline{Nu}$ are presented. First, note that, while the time-averaged quantities $\overline{\unicode[STIX]{x1D716}}+1$ , $\overline{\unicode[STIX]{x1D6F7}_{d}}$ and $\overline{\unicode[STIX]{x1D6F7}_{b1}}$ converge to $\overline{Nu}$ , each corresponding term has a specific behaviour during the different phases of the reversal cycle. Let us describe each instantaneous transfer rates in turn.

Figure 17. Instantaneous fields for a particular reversal cycle taken at regular intervals during the early release phase at $(Ra=5\times 10^{7},Pr=3.0)$ . The first (respectively last) snapshot follows point ( $\text{a}_{p}^{\prime }$ ) (respectively precedes point ( $\text{b}_{p}^{\prime }$ )). Streamlines are superposed to field $2\unicode[STIX]{x1D626}_{ij}\unicode[STIX]{x1D626}_{ij}$ the value of which is given by the colour. Corresponding values of $\unicode[STIX]{x1D716}$ and $\unicode[STIX]{x1D6F7}_{d}$ below each snapshot.

The vertical heat flux $Nu_{vol}$ which measures the conversion from $E_{apot}$ to $E_{kin}$ (see (2.8) and (2.9)), is by far the term that is found the most fluctuating, notably during the release. During the first part of the rebound i.e. the interval between points ( $\text{a}_{p}$ ) and ( $\text{b}_{p}$ ) (respectively (a) and (b)) for the particular reversal (respectively for the generic curve), $Nu_{vol}$ reaches peak values which are several times larger than $\overline{Nu}$ . This is related to the release of thermal energy $E_{apot}$  (figure 8) and plumes (figure 14). Between points (b) and (c), the generic $Nu_{vol}$ abruptly decreases. In terms of a particular realisation, this is due to the rotation of the bulk acting against buoyancy forces (points ( $\text{b}_{p}$ ) and ( $\text{c}_{p}$ ) in figure 14) and may result in a negative heat transfer as seen in figure 16 for the particular reversal cycle. This is consistent with results from (Chandra & Verma Reference Chandra and Verma2013). During the acceleration and accumulation phases, $Nu_{vol}$ fluctuates less and slightly decreases.

In addition to $Nu_{vol}$ , the viscous dissipation rate $\unicode[STIX]{x1D716}$ governs the evolution of $E_{kin}$  (see (2.8)). Viscous dissipation increases briefly during the reorganisation periods: once during the release phase around point (b), and again during the rotation of the central vortex between points (c) and (d). A gradual increase in $\unicode[STIX]{x1D716}$ is observed during the acceleration phase, between points (d) and (e), followed by a progressive decay during the accumulation phase. For a particular reversal cycle, contributions to the viscous dissipation rate $\unicode[STIX]{x1D716}$ are located (see figure 17) primarily along the vertical side walls and along the horizontal plates where ascending (respectively descending) plumes collide.

In addition to $Nu_{vol}$ , $E_{apot}$ is governed by the mixing term $\unicode[STIX]{x1D6F7}_{d}$ (see (2.13)). As seen in fields $\unicode[STIX]{x1D735}y_{r}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D703}$ on figures 14 and 15, contributions to $\unicode[STIX]{x1D6F7}_{d}$ are distributed along thin filaments which trace the contour of thermal plumes. $\unicode[STIX]{x1D6F7}_{d}$ is thus affected by small scales and fluctuates during all the process.

However, fluctuations are slightly larger during the second part of the release phase (points ( $\text{b}_{p}$ ) to ( $\text{d}_{p}$ )). An increase on such mixing fronts during the rotation of the bulk leads to the corresponding increase of $\unicode[STIX]{x1D6F7}_{d}$ . Note that during the accumulation phase, the generic curve for $\unicode[STIX]{x1D6F7}_{d}$ slightly decreases. Finally, the forcing boundary term $\unicode[STIX]{x1D6F7}_{b1}$ does not fluctuate much on the whole cycle.

To summarise, the dominant mechanisms during the release phase are first the energy conversion $Nu_{vol}$ followed by the mixing $\unicode[STIX]{x1D6F7}_{d}$ during the rebound. The acceleration phase is characterised by the increase of the viscous dissipation $\unicode[STIX]{x1D716}$ . A constant decay of both dissipation $\unicode[STIX]{x1D716}$ and mixing $\unicode[STIX]{x1D6F7}_{d}$ processes is conversely observed during the accumulation phase.

6.3 Angular momentum transfer rates

Figure 18. Time evolution of the angular impulse transfer rates for $(Ra=5\times 10^{7},Pr=3.0)$ . From (a) to (c): evolution of $M/\overline{|M|}$ , $(M-I_{b})/\overline{|M|}$ , and $I_{a}/\overline{|M|}$ . Figure for $L_{2D}/\overline{|L_{2D}|}$ is given as a reference curve. (a) Corresponds to the generic reversal: the average curve is displayed by a black line and the hatched region represents one standard deviation from this average. (b) and (c) depict two particular reversal cycles.

Figure 19. Instantaneous fields for a particular reversal cycle taken at regular intervals during the accumulation phase at $(Ra=5\times 10^{7},Pr=3.0)$ . The first (respectively last) snapshot follows point ( $\text{e}_{p}$ ) (respectively precedes point ( $\text{a}_{p}^{\prime }$ )). Streamlines are superposed to field $0.5(\boldsymbol{x}\boldsymbol{\cdot }\boldsymbol{x})\unicode[STIX]{x2202}_{x}\unicode[STIX]{x1D703}$ (the value is provided by the colour code). Corresponding values of $M$ and $I_{a}$ are given below each snapshot.

Since the evolution of $L_{2D}$ is characteristic of flow reversals, it is of interest to examine the angular impulse transfer rates $M$ , $(M-I_{b})$ and $I_{a}$ of (2.4) (figure 18). In addition to points (a) to (e), we introduce points (f) and (g). These points are located inside the accumulation phase and coincide for point (f) with a change of sign for $M$ and for point (g) with a change of sign for $I_{a}$ .

Time evolution of the input torque $M$ has a maximum value during the rebound, followed by a local minimum near point (d), a slight increase during the acceleration phase and a monotonic decrease during the accumulation phase. The bulk term $M(t)$ and the boundary term $I_{b}(t)$ are well correlated and have similar orders of magnitude. This is why only the $M(t)-I_{b}(t)$ evolution is plotted.

Before point (e), the difference $M(t)-I_{b}(t)$ is larger than $I_{a}(t)$ , from (e) to (f) it is the opposite ( $\text{d}L_{2D}/\text{d}t$ changes sign). One may discriminate three time periods: from (a) to (f), $I_{a}$ opposes $L_{2D}$ (i.e. the central vortex), contrary to the torque $M$ or $(M-I_{b})$ ; from (f) to (g), $I_{a}$ opposes $L_{2D}$ similarly to $M$ or $(M-I_{b})$ ; and finally, from (g) to (a), $I_{a}$ contributes to $L_{2D}$ which is opposed by $M$ or $(M-I_{b})$ . During the accumulation phase, the integrand $0.5Pr[\boldsymbol{x}\boldsymbol{\cdot }\boldsymbol{x}]\unicode[STIX]{x2202}_{x}\unicode[STIX]{x1D703}$ of the input torque $M(t)$ is spatially distributed as follows (figure 19): it is negligible inside the bulk, and concentrated in two regions, i.e. along vertical side walls associated with the central vortex and inside corners flows. In these two regions they are overall of opposite sign. At point (f), $M(t)$ changes sign: the corner flows become dominant in the integrand of $M(t)$ and the overall torque is then opposing the central vortex (figure 18). Similarly, during the accumulation phase, the integrand of $I_{a}$ that is vorticity, is distributed in two sectors along the boundaries: one along the boundaries of the central vortex, and another of opposite sign along the top and bottom corners (figure 10). From point (g), $I_{a}$ has changed sign: this is related to the dominant contribution of the corner rolls.

As a consequence, it is noteworthy that the input torque $M(t)$ changes sign at point (f) long before the reversal time (point (a)). The interval from point (f) to point (g) can be seen as a transition period from a central vortex-dominated to a corner-rolls-dominated flow.

7.1 Linear stability approach

The significance of instant (g) where $I_{a}$ changes sign, is tentatively explained by two complementary approaches. In a first approach, one considers the following hypothesis: the change occurring in the accumulation phase near point (g), is due to a modification of the dynamics governing fluctuations around large-scale structures. The large-scale structures identified by fields $\unicode[STIX]{x1D703}^{o}(x,y,t)$ and $\boldsymbol{u}^{o}(x,y,t)$ are here obtained by an ensemble average over many realisations (see § 5.1). A linear stability analysis is thus defined. The base state at time $t_{o}$ is given by fields $\unicode[STIX]{x1D703}^{o}(x,y,t_{o})$ and $\boldsymbol{u}^{o}(x,y,t_{o})$ frozen at this particular time. The evolution of infinitesimal fluctuations $\unicode[STIX]{x1D703}^{\prime }(x,y,t)$ and $\boldsymbol{u}^{\prime }(x,y,t)$ is studied around this frozen base state. The linear stability analysis is performed by direct numerical simulations of the linearised Boussinesq equations,

Each linearised simulation starts with random disturbances of velocity and temperature fields and is computed for several hundred time units. The perturbation kinetic energy $\langle u_{i}^{\prime }u_{i}^{\prime }\rangle _{vol}$ or the square of the fluctuation temperature $\langle \unicode[STIX]{x1D703}^{\prime }\unicode[STIX]{x1D703}^{\prime }\rangle _{vol}$ are monitored in time. In figure 20, each curve is related to a different base state $t_{o}$ . Positive values of growth rate $\unicode[STIX]{x1D70E}$ are obtained in all cases.

For $t_{o}$ before transition point (g), the values of $\unicode[STIX]{x1D70E}$ are quite small ( ${\sim}10^{-3}-10^{-2}$ ) and the most amplified mode corresponds to a trace of the base state itself (see figure 20 c). This is possibly related to a slow variation in time of the large-scale flow. In contrast, once point (g) is reached, the growth rate $\unicode[STIX]{x1D70E}$ increases by a factor of 20, while a different most amplified mode appears. This amplified mode is reminiscent of the flow features observed in the bottom right corner of figure 10.

Figure 20. Linear stability analysis of the generic fields for ( $Ra=5\times 10^{7}$ , $Pr=3.0$ ). (a) Placement of the different values of $t_{o}$ with respect to the generic $L_{2D}$ and $E_{apot}$ curves. (b) Evolution of the energy of the normalised velocity and temperature perturbations for different base states (or equivalently $t_{o}$ ). Growth rate $\unicode[STIX]{x1D70E}$ measured for different values of $t_{o}$ : ▵ (green)  $\unicode[STIX]{x1D70E}=0.009$ , ▴ (green)  $\unicode[STIX]{x1D70E}=0.005$ , ▫ (cyan)  $\unicode[STIX]{x1D70E}=0.005$ , ▪ (magenta)  $\unicode[STIX]{x1D70E}=0.006$ , ○ (orange)  $\unicode[STIX]{x1D70E}=0.011$ , ● (orange)  $\unicode[STIX]{x1D70E}=0.250$ ; (c) streamlines corresponding to field $\boldsymbol{u}^{o}(x,y,t_{o})$ superposed to base state $\unicode[STIX]{x1D703}^{o}(x,y,t_{o})$ shown for two values of $t_{o}$ : ○ (orange) located between points (f) and (g) and ● (orange) located on point (g). (d) Disturbance field $\unicode[STIX]{x1D703}^{\prime }\unicode[STIX]{x1D703}^{\prime }(x,y,t)$ measured at the end of the curves corresponding to ○ (orange) and ● (orange).

7.2 Nonlinear approach with adiabatic boundary conditions

We have identified through a linear stability analysis a critical generic state around point (g). Let us now seek to relate this state to the ‘avalanche’ mechanism: a mechanism due to a localised accumulation of energy inside the fluid which is followed by a sudden transition. One could argue that, when localised accumulation is sufficient, a reversal event takes place even if the external thermal forcing is thereafter suppressed. A second approach, which is nonlinear, is based on this idea and confirms point (g) as a transition point associated with a threshold state.

First, let us consider the conditionally averaged fields ( $\unicode[STIX]{x1D703}^{o}(x,y,t_{o})$ and $\boldsymbol{u}^{o}(x,y,t_{o})$ ) obtained at time $t_{o}$ for the standard RB problem, see figure 7. Starting from such fields labelled by time $t_{o}$ , we perform a direct numerical simulation of the nonlinear Boussinesq equations changing boundary conditions from isothermal to adiabatic on the top and bottom plates. This effectively suppresses the external thermal forcing since $\unicode[STIX]{x1D6F7}_{b1}=\unicode[STIX]{x1D6F7}_{b2}=0$ and limits the amount of thermal energy contained inside the cavity. Note that different initial conditions i.e. different $t_{o}$ lead to different trajectories for the adiabatic problem.

Evolution of the normalised kinetic energy $E_{kin}/\overline{|E_{kin}|}$ and the normalised angular impulse $L_{2D}/\overline{|L_{2D}|}$ can be seen in figure 21, where each curve corresponds to different $t_{o}$  (shown as $\circ$ marks) inside the accumulation phase. Dashed (respectively solid) lines indicate negative (positive) values of $L_{2D}/\overline{|L_{2D}|}$ . For simulations preceding the transition point (g) a decay in both the angular impulse and the kinetic energy are observed. On the contrary, from point (g) on, a change from negative to positive values of $L_{2D}/\overline{|L_{2D}|}$ simultaneously as peak values in $E_{kin}/\overline{|E_{kin}|}$ are observed, before the ensuing decay. To illustrate this transition we follow the evolution of Fourier modes $\hat{v}_{mn}$ of the vertical velocity for two initial conditions before and after point (g) as seen on figure 21. For initial conditions preceding $I_{a}=0$ , an already weakened central vortex is unable to contain the corner flows in place but does not break (the mode $\hat{v}_{11}$ which is related to the monopole decreases but remains negative). On the contrary, for initial conditions following the transition point (g), Fourier coefficient $\hat{v}_{11}$ changes sign and thus indicates a LSC reversal, before the subsequent decay.

Figure 21. (a,c) Normalised angular impulse $|L_{2D}|/\overline{|L_{2D}|}$ (a,b) and normalised kinetic energy $E_{kin}/\overline{|E_{kin}|}$ (c,d) for $(Ra=5\times 10^{7},Pr=3.0)$ as a function of $t_{o}$ characterising the base state and time $(t-t_{o})$ . Each curve represents different initial conditions i.e. different $t_{o}$ , where solid lines indicate when a flow reversal takes place (see text). (b,d) Evolution of modal coefficients $\hat{v}_{mn}$ from two initial conditions: before and at point (g).