1. Introduction
Taylor's statistical theory of turbulence states that the turbulence is isotropic if the average value of any function of the velocity components, defined in relation to a given set of axes, is unaltered under axis rotation (Taylor Reference Taylor1935). However, the condition of isotropy is not satisfied for the energy-containing scales of turbulence, since no energy production can happen for isotropic turbulence due to its directional independence (Tennekes & Lumley Reference Tennekes and Lumley1972). Several metrics have been used to quantify turbulence anisotropy (see table1 for a brief review) and, out of those, one of the metrics to quantify the anisotropic signatures of the energy-containing motions at a point in the flow is the anisotropy Reynolds stress tensor (e.g. Krogstad & Torbergsen Reference Krogstad and Torbergsen2000). This tensor becomes zero in an isotropic turbulence and its anisotropy is quantified by using the invariants of
$b_{ij}$, an approach pioneered by Lumley & Newman (Reference Lumley and Newman1977), known as invariant analysis. The invariants of the anisotropy Reynolds stress tensor have been used extensively in the context of wall-bounded neutral flows to deduce the anisotropic characteristics of the energy-containing motions (e.g. Shafi & Antonia Reference Shafi and Antonia1995).
A brief summary of different approaches and metrics used to study anisotropy in a turbulent flow.

In convective turbulence, buoyant structures, such as thermal plumes, are the energy-containing motions that transport heat and drive the flow (Celani, Mazzino & Vergassola Reference Celani, Mazzino and Vergassola2001). These thermal plumes are well-organized structures of warm rising (warm-updraft) and cold descending (cold-downdraft) fluid, which generate ramp–cliff patterns in temperature time series when passing a thermal probe (Zhou & Xia Reference Zhou and Xia2002). Shang et al. (Reference Shang, Qiu, Tong and Xia2003) have shown that, in turbulent Rayleigh–Bénard convection, the time series of the instantaneous vertical heat flux associated with the thermal plumes displays intermittent characteristics. Intermittency is defined as a property of the turbulent signal which is quiescent for much of the time and occasionally bursts into life with unexpectedly high values more common than in a Gaussian signal (e.g. Davidson Reference Davidson2015). However, the effect of this intermittent heat transport on the anisotropic fluctuations in the velocity field of convective turbulence is not yet well understood, as acknowledged by Pouransari, Biferale & Johansson (Reference Pouransari, Biferale and Johansson2015). This problem is particularly relevant for the surface layer of a convectively driven atmospheric boundary layer, where the most prevalent coherent structures are the thermal plumes, and the heat transport characteristics associated with these plumes appear to be intermittent (e.g. Katul et al. Reference Katul, Albertson, Parlange, Chu and Stricker1994).
The previous works on the atmospheric surface layer (ASL) plumes have focused on the following: (a) deducing their detailed structures and dynamics (Wilczak Reference Wilczak1984); (b) identifying the coupling between the surface and air temperatures (Garai & Kleissl Reference Garai and Kleissl2011, Reference Garai and Kleissl2013); and (c) investigating the difference in the Monin–Obukhov similarity functions by conditioning on the updraft and downdraft motions (Li et al. Reference Li, Gentine, Mellado and McColl2018; Fodor, Mellado & Wilczek Reference Fodor, Mellado and Wilczek2019). However, some early investigators noted that in an unstable ASL there were certain intermittent bursts in the upward heat flux, persisting for approximately 10–20 s duration, which were associated with large downward momentum transport (Kaimal & Businger Reference Kaimal and Businger1970; Haugen, Kaimal & Bradley Reference Haugen, Kaimal and Bradley1971). They commented that the vertical velocity fluctuations associated with these heat-flux events could either transport momentum downwards in large bursts or transport it upwards. Businger (Reference Businger1973) termed these intermittent momentum bursts associated with the heat-flux events as ‘convection-induced stress’. Recently, Lotfy et al. (Reference Lotfy, Abbas, Zaki and Harun2019) also obtained the same result from a field experiment in an unstable ASL, where they observed that the persistent warm updrafts of 10–20 s duration were associated with a large amount of momentum flux in the downward direction. By investigating the large-eddy simulation results in convective conditions, Salesky & Anderson (Reference Salesky and Anderson2018) interpreted this phenomenon as a buoyancy-dominated scale modulation effect. They explained that, under highly convective conditions, the small-scale turbulence is excited in the updraft regions and suppressed in the downdraft regions, leading to intermittent periods of small-scale excitation in the momentum fluxes.
From the discussion above, it becomes apparent that, in an unstable ASL, the vertical velocity fluctuations associated with the coherent heat-flux events could transport large amount of momentum in intermittent bursts, in either the upward or downward direction. Since only the anisotropic part of the velocity fluctuations can carry momentum (Dey et al. Reference Dey, Ravi-Kishore, Castro-Orgaz and Ali2018; Könözsy Reference Könözsy2019), this indicates that the Reynolds stress anisotropy associated with these coherent heat-flux events must be different from the averaged whole flow. Therefore, studying the role of intermittent heat-flux events towards the anisotropy in the velocity fluctuations is of practical importance in the context of ASL turbulence. For a systematic investigation of this problem, invariant analysis of the anisotropy Reynolds stress tensor in an event-based framework is a well-suited approach.
The event-based approach in turbulence is based on the fact that coherent physical structures exist in a turbulent flow (e.g. Chapman & Tobak Reference Chapman and Tobak1985). Specifically, Narasimha et al. (Reference Narasimha, Kumar, Prabhu and Kailas2007) mentioned that, in this approach, the turbulent field can be expressed in terms of events, given that its types, magnitudes, arrival times, etc. are defined properly. The interest in the event-based description of turbulence started with the flow visualization studies of Kline et al. (Reference Kline, Reynolds, Schraub and Runstadler1967). They observed that the flow near the wall of a boundary layer was organized into streaks of high- and low-momentum fluid. Subsequently, the low-momentum streaks were seen to intermittently erupt away from the wall in a chaotic process named bursting. This accounted for much of the outward vertical transport of momentum and the production of turbulent kinetic energy in the boundary layer. A detailed review of different conditional sampling techniques to detect events in turbulence can be found in Antonia (Reference Antonia1981) and Wallace (Reference Wallace2016). The types of coherent structures whose signatures are associated with these events are reviewed in detail by Cantwell (Reference Cantwell1981), Robinson (Reference Robinson1991) and Jiménez (Reference Jiménez2018).
Sreenivasan, Antonia & Britz (Reference Sreenivasan, Antonia and Britz1979) first applied this event-based approach to investigate the effect on turbulence anisotropy associated with the coherent structures in a heated turbulent jet. Based on the premise that the fine structures were superposed on the large structures, Sreenivasan et al. (Reference Sreenivasan, Antonia and Britz1979) extracted the coherent ramp–cliff events in a heated turbulent jet, and then subtracted these patterns from the signal to get the superposed fluctuations. Their focus was to show that the skewness in the temperature gradient vanishes for the fine structures, thus confirming the local isotropy. Recently, following the work of Lozano-Durán, Flores & Jiménez (Reference Lozano-Durán, Flores and Jiménez2012), Dong et al. (Reference Dong, Lozano-Durán, Sekimoto and Jiménez2017) studied the connected regions of high-intensity momentum zones in three-dimensional simulations of homogeneous shear and channel flows and investigated the Reynolds stress anisotropy. They quantified anisotropy by the invariants of the Reynolds stress tensor within these high-intensity momentum zones along with their sizes; where the size was defined as the box diagonal of the parallelepiped that circumscribed these connected regions. Zhou & Xia (Reference Zhou and Xia2011) attempted to disentangle the role of thermal plumes on the velocity field in a Rayleigh–Bénard convection, by studying separately the anisotropy in the inertial subrange of the positive and negative vertical velocity increments. They showed that the negative increments at small separations deviated from the Kolmogorov scaling, which they attributed to the presence of the coherent structures such as thermal plumes.
The anisotropy directly associated with the intermittent occurrences of the coherent structures is regarded as a state-of-the-art theoretical and experimental problem (Pouransari et al. Reference Pouransari, Biferale and Johansson2015). To the best of our knowledge, very few studies have addressed this problem by adopting an event-based approach. This is particularly pertinent in the context of ASL turbulence, where there are no comprehensive studies to quantify anisotropy concomitant with the intermittent heat-flux events in convective conditions. The present study attempts to fill this gap, using a field-experimental dataset. Therefore, we define our objectives as follows:
(i) To investigate the detailed correspondence between the heat-flux events and turbulence anisotropy in an unstable ASL.
(ii) To formulate a structural description of the heat-flux events and to investigate whether they have any characteristic length scales associated with least anisotropic turbulence.
The present paper is organized in three different sections. In § 2 we describe the dataset and methodology to develop various statistical measures to quantify anisotropy associated with the heat-flux events. In § 3 we present and discuss the results, and in § 4 we conclude our findings and provide future directions for further research.
2. Data and methodology
We have used the dataset from the Surface Layer Turbulence and Environmental Science Test (SLTEST) experiment. The SLTEST experiment was conducted over a flat and homogeneous terrain at the Great Salt Lake desert in Utah, USA (40.14
$^\circ$N, 113.5
$^\circ$W), with the aerodynamic roughness length (
$z_0$) being
$z_{0}\approx 5~\textrm {mm}$ (Metzger, McKeon & Holmes Reference Metzger, McKeon and Holmes2007). The SLTEST site characteristics and the high quality of the dataset have been documented in detail in many previous studies (e.g. Hutchins & Marusic Reference Hutchins and Marusic2007). In this experiment, nine north-facing sonic anemometers (CSAT3, Campbell Scientific, Logan, USA) were installed on a 30 m tower approximately logarithmically at
$z=1.4$, 2.1, 3.0, 4.3, 6.1, 8.7, 12.5, 17.9, 25.7 m, levelled to within
${\pm }0.5^\circ$ from the true vertical. All CSAT3 sonic anemometers were synchronized in time and the sampling frequency was set at 20 Hz. The experiment ran continuously for nine days from 26 May 2005 to 3 June 2005.
2.1. Data processing
The data were divided into 30 min periods containing the 20 Hz measurements of the three wind components and the sonic temperature from all nine sonic anemometers. To select the 30 min periods for analysis, we followed the standard procedures listed below:
(i) The 30 min periods were selected from the fair weather conditions during the daytime periods with no rain.
(ii) The time series of all three components of velocity and sonic temperature were plotted and visually checked. No electronic spikes were found in the data (Vickers & Mahrt Reference Vickers and Mahrt1997).
(iii) The horizontal wind direction sector was limited to
$-30^\circ < \theta < 30^\circ$ (where
$\theta$ is the horizontal wind direction from North).(iv) The coordinate systems of all nine sonic anemometers were rotated in the streamwise direction by applying the double-rotation method of Kaimal & Finnigan (Reference Kaimal and Finnigan1994) for each 30 min period. The turbulent fluctuations in the wind components (
$u^{\prime }$,
$v^{\prime }$ and
$w^{\prime }$ in the streamwise, cross-stream and vertical directions, respectively), and in the sonic temperature (
$T^{\prime }$) were calculated after removing the 30 min linear trend from the associated variables (Donateo, Cava & Contini Reference Donateo, Cava and Contini2017).(v) Only those 30 min periods were chosen when the surface layer was unstable, i.e. the sensible heat flux was positive at all nine measurement heights, and the vertical variations in the 30 min averaged momentum and heat fluxes were less than 10 %.
Application of all these checks resulted in a total of 29 periods suitable for our analysis. For these periods,
$\sigma _{u}/\overline {u}$ was less than 0.2, so Taylor's hypothesis could be assumed to be valid (Willis & Deardorff Reference Willis and Deardorff1976). The Obukhov length (
$L$) was calculated for each of these 30 min periods as
where
$T_{0}$ is the surface air temperature, computed from the mean sonic temperature at
$z=1.4~\textrm {m}$,
$g$ is the acceleration due to gravity (
$9.8~\textrm {m}~\textrm {s}^{-2}$),
$H_{0}$ is the surface kinematic heat flux, computed as
$\overline {w^{\prime }T^{\prime }}$ at
$z=1.4~\textrm {m}$ (by the constant flux layer assumption),
$k$ is the von Kármán constant (0.4) and
$u_{*}$ is the friction velocity computed as
with
$\overline {u^{\prime }w^{\prime }}$ and
$\overline {v^{\prime }w^{\prime }}$ the streamwise and cross-stream momentum fluxes, respectively, computed at
$z=1.4~\textrm {m}$.
The range of
$-L$ values was between 2 and 20 m for these 29 periods suitable for our analysis. Since each 30 min period consisted of the nine level time-synchronized turbulence measurements from the CSAT3 sonic anemometers, a total of 261 combinations of the stability ratios (
$\zeta =z/L$) were possible for these selected periods. The entire range of
$-\zeta$ (
$12 \leq \zeta \leq 0.07$) was divided into six stability classes (Liu, Hu & Cheng Reference Liu, Hu and Cheng2011) and these were considered for the detailed analysis of the Reynolds stress anisotropy associated with the heat-flux events (table 2). We discuss the analysis methods in the following sections.
The six different stability classes formed from the ratio
$-\zeta =z/L$ in an unstable ASL flow, from highly convective (
$-\zeta >2$) to near-neutral (
$0<-\zeta <0.2$). The associated heights with each of the stability classes are also given.

2.2. Quadrant analysis
The quadrant analysis is a conditional sampling method of investigating the contributions to the turbulent transport of scalars and momentum in terms of the organized eddy motions present in the flow (Wallace Reference Wallace2016). The four different quadrants of the
$u^{\prime }\text {--}w^{\prime }$ and
$T^{\prime }\text {--}w^{\prime }$ planes are defined in table 3. In the
$T^{\prime }\text {--}w^{\prime }$ (
$u^{\prime }\text {--}w^{\prime }$) quadrant plane, the warm updrafts (I) (ejections (II)) and cold downdrafts (III) (sweeps (IV)) are the down-gradient motions. On the other hand, the remaining two quadrants represent the counter-gradient motions generated due to the turbulent swirls in the flow (Gasteuil et al. Reference Gasteuil, Shew, Gibert, Chilla, Castaing and Pinton2007).
The four quadrants of
$u^{\prime }\text {--}w^{\prime }$ and
$T^{\prime }\text {--}w^{\prime }$ in an unstable ASL.

In the quadrant analysis method applied to the ASL, the momentum- or heat-flux fractions and time fractions from each quadrant of
$u^{\prime }\text {--}w^{\prime }$ or
$T^{\prime }\text {--}w^{\prime }$ are reported over smooth and rough surfaces (McBean Reference McBean1974; Antonia Reference Antonia1977; Narasimha et al. Reference Narasimha, Kumar, Prabhu and Kailas2007; Zou, Zhou & Sun Reference Zou, Zhou and Sun2017). The flux fractions (
$F_{f}$) and time fractions (
$T_{f}$) for each quadrant (
$X$) are evaluated as
\begin{equation} \left. \begin{array}{ll} {}(F_{f})_{X} = \displaystyle\dfrac{\displaystyle\sum[(w^{\prime}x^{\prime})I_{X}]} {\displaystyle\sum{w^{\prime}x^{\prime}}} & \quad (x=u,T) , \\ {}(T_{f})_{X} = \displaystyle\dfrac{\displaystyle\sum I_{X}}{N} & \quad ({X}=\textrm{I},\textrm{II},\textrm{III},\textrm{IV}), \end{array}\right\} \end{equation}where
\[ I_{X} = \begin{cases} 1 & \text{if } \{w^{\prime},x^{\prime}\} \in {{X}}{,}\\ 0 & \text{otherwise}{,} \end{cases} \]
and
$N$ is the total number of points in a run.
However, following Chowdhuri & Burman (Reference Chowdhuri and Burman2020), we extend the quadrant analysis method to study the anisotropy Reynolds stress tensor in relation to the heat-flux events occurring in the
$T^{\prime }\text {--}w^{\prime }$ quadrant plane. We normalize
$w^{\prime }$ and
$T^{\prime }$ by their respective standard deviations, and use the symbol
$\hat {x}$ to denote the turbulent fluctuations in
$x$ normalized by its standard deviation (
$\hat {x}=x^{\prime }/\sigma _{x}$, where
$x$ can be
$u$,
$w$ or
$T$). Before describing the methodology, we give a short description of the anisotropy Reynolds stress tensor.
2.2.1. Anisotropy Reynolds stress tensor
The anisotropy Reynolds stress tensor is widely used to express the anisotropy in the energy-containing motions (Pope Reference Pope2000), and is defined in the Cartesian tensor notation as
where
$i= 1$, 2 and 3 denote the streamwise, cross-stream and vertical directions,
$q$ is the turbulent kinetic energy, and
$\delta _{ij}$ is the Kronecker delta. Note that
$b_{ij}$ is a symmetric and traceless tensor, bounded between
$-1/3 \leq b_{ij} \leq 2/3$, and equal to zero for isotropic turbulence. From the Cayley–Hamilton theorem, the two invariants
$\xi$ and
$\eta$ of
$b_{ij}$ are defined as
where
$\xi$ represents the topology of the anisotropy Reynolds stress tensor and
$\eta$ represents the degree of isotropy.
The different realizable anisotropic states of turbulence are defined based on the values of
$\xi$ and
$\eta$ and are represented on the
$\xi \text {--}\eta$ plane, known as the anisotropy-invariant map (Choi & Lumley Reference Choi and Lumley2001). The anisotropy Reynolds stress tensor (
$b_{ij}$) has three limiting anisotropic states based on the shape of the energy distribution in the three principal axes associated with the three eigenvalues and eigenvectors of
$b_{ij}$, also known as the componentiality of turbulence (Kassinos, Reynolds & Rogers Reference Kassinos, Reynolds and Rogers2001). These three limiting states of
$b_{ij}$ are 1-component anisotropy (rod-like energy distribution,
$b_{1 { c}}$), 2-component anisotropy (disk-like energy distribution,
$b_{2 { c}}$), and 3-component isotropy (spherical energy distribution,
$b_{3 { c}}$), represented in the principal axes coordinate system as
\begin{align} b_{1 { c}} = \begin{bmatrix} 2/3 & 0 & 0 \\ 0 & -1/3 & 0 \\ 0 & 0 & -1/3 \end{bmatrix}, \quad b_{2 c} = \begin{bmatrix} 1/6 & 0 & 0 \\ 0 & 1/6 & 0 \\ 0 & 0 & -1/3 \end{bmatrix}, \quad b_{3 { c}} = \begin{bmatrix} 0 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} . \end{align} An alternative to the anisotropy-invariant maps is the barycentric map introduced by Banerjee et al. (Reference Banerjee, Krahl, Durst and Zenger2007), where each realizable anisotropic state of
$b_{ij}$ is written as a linear combination of the three limiting states
$b_{1 { c}}$,
$b_{2 { c}}$ and
$b_{3 { c}}$ as
where the coefficients
$C_{1 { c}}$,
$C_{2 { c}}$ and
$C_{3 { c}}$ are the three corresponding weights associated with the three limiting states, defined as
\begin{equation} \left.\begin{array}{c} C_{1 { c}} = e_{1}-e_{2},\\ C_{2 { c}} = 2(e_{2}-e_{3}),\\ C_{3 { c}} = 3e_{3}+1, \end{array}\right\} \end{equation}with
where
$e_{1}$,
$e_{2}$ and
$e_{3}$ are the three eigenvalues of
$b_{ij}$ in the order
$e_{1} > e_{2} > e_{3}$. Note that these three coefficients
$C_{1 { c}}$,
$C_{2 { c}}$ and
$C_{3 { c}}$ are bounded between 0 and 1. In the extreme case, one of the coefficients taking the value 0 signifies that the particular limiting state associated with that coefficient does not exist. Similarly, 1 signifies only that that particular limiting state exists while the other two states are non-existent (Banerjee et al. Reference Banerjee, Krahl, Durst and Zenger2007). Given the linearity in the construction of the barycentric map, it provides a non-distorted visualization of anisotropy (Radenković, Burazer & Novković Reference Radenković, Burazer and Novković2014). Banerjee et al. (Reference Banerjee, Krahl, Durst and Zenger2007) defined the coefficient
$C_{3 { c}}$ as the degree of isotropy, such that, the higher the value of
$C_{3 { c}}$, the more the anisotropic state of
$b_{ij}$ is dominated by the 3-component isotropy. The anisotropic states of
$b_{ij}$ can be represented by the
$RGB$ colour map of Emory & Iaccarino (Reference Emory and Iaccarino2014) as
\begin{equation} \begin{bmatrix} R \\ G \\ B \end{bmatrix} = C_{1 { c}} \begin{bmatrix} 1 \\ 0 \\ 0 \\ \end{bmatrix} + C_{2 { c}} \begin{bmatrix} 0 \\ 1 \\ 0 \\ \end{bmatrix} + C_{3 { c}} \begin{bmatrix} 0 \\ 0 \\ 1 \\ \end{bmatrix} , \end{equation}such that the 1-component anisotropy is red, the 2-component anisotropy is green and the 3-component isotropy is blue. All other states within the barycentric map are linear combinations of these three colours.
Since in this study we will be using the barycentric map to visualize the anisotropic states of
$b_{ij}$, some details about its construction are appropriate here. The barycentric map is spanned by a Euclidean domain where the three limiting states of
$b_{ij}$ are placed at the three vertices of an equilateral triangle having the coordinates
$(0, 0)$ for the 2-component anisotropy,
$(1, 0)$ for the 1-component anisotropy, and (
$1/2, \sqrt {3}/2$) for the 3-component isotropy (Stiperski & Calaf Reference Stiperski and Calaf2018). This is graphically illustrated in figure 1. The coordinate system (
$x, y$) of the barycentric map is defined as
\begin{align} x &= C_{1 { c}}x_{1 { c}}+C_{2 { c}}x_{2 { c}}+C_{3 { c}}x_{3 { c}}\nonumber\\ &= C_{1 { c}}+\frac{C_{3 { c}}}{2} \end{align}and
\begin{align} y &= C_{1 { c}}y_{1 { c}}+C_{2 { c}}y_{2 { c}}+C_{3 { c}}y_{3 { c}}\nonumber\\ &= \frac{\sqrt{3}}{2}C_{3 { c}}, \end{align}
such that the distance from the base of the equilateral triangle is directly proportional to the degree of isotropy (
$C_{3 { c}}$) of
$b_{ij}$. At the centroid of the barycentric map (
$1/2, \sqrt {3}/6$), from (2.11) and (2.12) it can be shown that
This barycentric map can also be divided into three equal regions, R1 (right-third portion), R2 (left-third portion) and R3 (top-third portion), by drawing three perpendicular bisectors from the centroid of the map (figure 1). From (2.11) and (2.12), along with the constraint defined in (2.9), these perpendicular bisectors can be represented mathematically as
Subsequently from symmetry it follows that: in the region R1,
$C_{1 { c}} > \{C_{2 { c}}, C_{3 { c}}\}$; in the region R2,
$C_{2 { c}} > \{C_{1 { c}}, C_{3 {c}}\}$; and in the region R3,
$C_{3 { c}} > \{C_{1 { c}}, C_{2 { c}}\}$. Therefore, in each of these three regions (R1, R2 or R3), the anisotropic state of
$b_{ij}$ is dominated by a particular limiting state associated with its coefficient (
$C_{1 { c}}$,
$C_{2 { c}}$ or
$C_{3 { c}}$).
An example of a barycentric map spanned by an equilateral triangle is shown to graphically illustrate the anisotropic states of
$b_{ij}$, using the
$RGB$ colour map of Emory & Iaccarino (Reference Emory and Iaccarino2014). The three vertices of the equilateral triangle represent the three limiting states with coefficients
$C_{1 { c}}$,
$C_{2 { c}}$ or
$C_{3 { c}}$ being equal to 1. At the sides opposite to the vertices, any one of these coefficients is 0, which indicates the absence of that particular anisotropic state associated with it. The black circle is the centroid of the equilateral triangle where
$C_{1 { c}}=C_{2 { c}}=C_{3 { c}}=1/3$. The three black lines are the three perpendicular bisectors that divide the equilateral triangle into three equal regions: R1 (right-third portion), R2 (left-third portion) and R3 (top-third portion). In each of these three regions, the anisotropic state of
$b_{ij}$ is dominated by a particular limiting state associated with its coefficient
$C_{1 { c}}$,
$C_{2 { c}}$ or
$C_{3 { c}}$.

2.2.2. Representation of anisotropy on
$\hat {T}\text {--}\hat {w}$ quadrant plane
To study the detailed correspondence between the anisotropic states of
$b_{ij}$ and the heat-flux events occurring in the
$\hat {T}\text {--}\hat {w}$ quadrant plane, we first linearly bin
$\hat {T}$ and
$\hat {w}$ into a uniform
$50\times 50$ grid for each run belonging to a particular stability class. The width of each grid is defined as
We choose the maximum (
$\hat {T}_{max}$,
$\hat {w}_{max}$) and minimum (
$\hat {T}_{min}$,
$\hat {w}_{min}$) values over all runs from a particular stability class to ensure the same grid for individual runs. Subsequently, we find the points lying between {
$\hat {T}_{bin}(m)<\hat {T}<\hat {T}_{bin}(m)+\textrm {d}\hat {T}$,
$\hat {w}_{bin}(n)<\hat {w}<\hat {w}_{bin}(n)+\textrm {d}\hat {w}$}, where
$1\leq m \leq 50$,
$1\leq n \leq 50$, and
$\hat {T}_{bin}(m)$ and
$\hat {w}_{bin}(n)$ are the edges of a particular
$(m,n)$ grid. For these points, we construct the anisotropy Reynolds stress tensor at
$(m,n)$ grid as
\begin{align} &\langle b_{ij} \mid\{\hat{T}_{bin}(m)<\hat{T}<\hat{T}_{bin}(m)+\textrm{d}\hat{T}, \hat{w}_{bin}(n)<\hat{w}<\hat{w}_{bin}(n)+\textrm{d}\hat{w}\} \rangle \nonumber\\ &\quad =\frac{{\left(\displaystyle\sum{{u^\prime}_{i} {u^\prime}_{j}}\right)_{m,n}}}{{\left(\displaystyle\sum{{u^\prime}_{i} {u^\prime}_{i}}\right)_{m,n}}} - \frac{1}{3}\delta_{ij}, \end{align}conditioned on the heat-flux events occurring between
and assign it to the value
$\{\hat {T}_{bin}(m),\hat {w}_{bin}(n)\}$.
In (2.16), the terms
${(\sum {{u^\prime }_{i} {u^\prime }_{j})}_{m,n}}$ are the contributions to the Reynolds stress tensor from each
$(m,n)$ grid. Accordingly, the trace of
$b_{ij}$ from (2.16) is zero due to the kinetic energy term
${(\sum {{u^\prime }_{i} {u^\prime }_{i})}_{m,n}}$ appearing in the denominator. Note that this kinetic energy is the energy contained in each
$(m,n)$ grid, rather than the total kinetic energy over the whole 30 min period. This formulation is similar to the scale decomposition of
$b_{ij}$, where at each scale the anisotropic Reynolds stress tensor is normalized by the kinetic energy contained in that scale to make it trace-free (Yeung & Brasseur Reference Yeung and Brasseur1991; Liu et al. Reference Liu, Yuan, Mei, Sun, Liu and Wang2017; Brugger et al. Reference Brugger, Katul, De Roo, Kröniger, Rotenberg, Rohatyn and Mauder2018).
To assess the frequency of occurrences of these heat-flux events, we compute the joint probability density function (j.p.d.f.) between
$\hat {T}$ and
$\hat {w}$ as
where
$N_{m,n}$ is the number of points lying in
$(m,n)$ grid and
$N$ is the total number of points in a 30 min run (36 000 for SLTEST data). Following Nakagawa & Nezu (Reference Nakagawa and Nezu1977), we also calculate the bivariate Gaussian j.p.d.f. for each grid as
\begin{align} &G(\hat{T}_{bin}(m),\hat{w}_{bin}(n))\nonumber\\ &\quad =\frac{1}{2 {\rm \pi} \sqrt{1-{R^{2}_{wT}}}} \exp{\left[-\left(\frac{{\hat{T}^2_{bin}(m)}-2R_{wT}\hat{T}_{bin}(m)\hat{w}_{bin}(n) +{\hat{w}^2_{bin}(n)}}{2(1-{R^{2}_{wT}})}\right)\right]}, \end{align}
where
$R_{wT}$ is the correlation coefficient between
$w$ and
$T$ (
$\overline {w^{\prime }T^{\prime }}/\sigma _{w}\sigma _{T}$).
If the three eigenvalues of
$b_{ij}$ (as defined in (2.16)) are
$e_{1b}$,
$e_{2b}$ and
$e_{3b}$, respectively, with
$e_{1b}>e_{2b}>e_{3b}$, we can calculate the degree of isotropy for
$(m,n)$ grid as
and the
$RGB$ colour map of its anisotropic states as
\begin{align} \begin{bmatrix} R \\ G \\ B \\ \end{bmatrix} &=\langle C_{1 { c}} \mid\{\hat{T}_{bin}(m),\hat{w}_{bin}(n)\} \rangle \begin{bmatrix} 1 \\ 0 \\ 0 \\ \end{bmatrix} + \langle C_{2 { c}} \mid\{\hat{T}_{bin}(m),\hat{w}_{bin}(n)\} \rangle \begin{bmatrix} 0 \\ 1 \\ 0 \\ \end{bmatrix} \nonumber\\ &\quad +\langle C_{3 { c}} \mid\{\hat{T}_{bin}(m),\hat{w}_{bin}(n)\} \rangle \begin{bmatrix} 0 \\ 0 \\ 1 \\ \end{bmatrix}, \end{align}with
Since we construct the same linear grid values of
$\hat {T}$ and
$\hat {w}$ for all runs belonging to a particular stability class, we take the average of the j.p.d.f.,
$C_{3 { c}}$, and the
$RGB$ colour matrices over all the individual periods. This averaging is necessary since it reduces the variability that exists from one run to another, due to the chaotic nature of turbulence. For a particular stability range, we can thus plot the averaged two-dimensional matrices of
$P(\hat {T}_{bin}(m),\hat {w}_{bin}(n))$,
$\langle C_{3 { c}} \mid \{\hat {T}_{bin}(m),\hat {w}_{bin}(n)\} \rangle$, and the
$RGB$ colour maps for each
$(m,n)$ grid of the
$\hat {T}\text {--}\hat {w}$ quadrant plane. While presenting the results in § 3.2, these averaged metrics are referred to as being associated with
$(\hat {T},\hat {w})$, without explicitly mentioning that these are the binned values. It is worth noting that the results obtained from this method are almost insensitive to the choice of the grid size. We verified this by changing the grid sizes of
$\hat {T}$ and
$\hat {w}$ by a factor of 2 (
$25\times 25$ and
$100\times 100$) and repeating the calculations, with no appreciable change being noticed in the results (not shown).
By performing the binning exercise in
$\hat {w}$ and
$\hat {T}$ as discussed above, we mask any time dependence, and hence no information can be obtained about the time scales of the associated heat-transporting events. It is thus interesting to formulate a description of the distribution of Reynolds stress anisotropy associated with different time scales of these heat-transporting events. To extract that information in addition to the quadrant analysis, we turn our attention to persistence analysis.
2.3. Persistence analysis
In non-equilibrium systems, persistence is defined as the probability that the local value of a fluctuating field does not change sign up to a certain time (e.g. Majumdar Reference Majumdar1999). The concept of persistence has earlier been used by Chamecki (Reference Chamecki2013) to study the non-Gaussian turbulence in canopy flows. He showed that an asymmetric velocity distribution inside the canopy can have very different persistent time scales for ejection and sweep events. Chamecki (Reference Chamecki2013) also noted that the persistent time is equivalent to the inter-pulse periods between the subsequent zero crossings of the turbulent signal (Sreenivasan, Prabhu & Narasimha Reference Sreenivasan, Prabhu and Narasimha1983; Kailasnath & Sreenivasan Reference Kailasnath and Sreenivasan1993; Bershadskii et al. Reference Bershadskii, Niemela, Praskovsky and Sreenivasan2004). We can apply this definition of persistence to the joint fluctuations in vertical velocity and temperature, to characterize the distribution of the time scales of the heat-flux events from four different quadrants of
$T^{\prime }\text {--}w^{\prime }$.
In order to implement our method, we choose the time series of
$w^{\prime }$ and
$T^{\prime }$ from any 30 min period belonging to a specific stability class (table 2), and conditionally sample the events occurring in the four different quadrants of the
$T^{\prime }\text {--}w^{\prime }$ plane. The events conditionally sampled from each quadrant of
$T^{\prime }\text {--}w^{\prime }$ (I, II, III or IV) can persist either as a single pulse or as a block of many consecutive pulses with a certain duration
$T_{B}$, before switching to another quadrant. The duration
$T_{B}$ is computed as the number of points residing within a single block, multiplied by the sampling interval of 0.05 s. In figure 2 we provide a graphical illustration of this method by showing a segment of a time series belonging to a particular stability range (
$-\zeta =9$), sampled from the warm-updraft (I) and cold-downdraft (III) quadrants. The shaded blocks in figure 2(a,b) represent warm updrafts (red) and cold downdrafts (blue) respectively, which persist for around 10–20 s duration. Associated with these blocks, we also show the horizontal velocity fluctuations (
$u^{\prime }$ and
$v^{\prime }$) in figure 2(c,d).
A 60 s long section of a time series of
$w^{\prime }$,
$T^{\prime }$,
$u^{\prime }$ and
$v^{\prime }$ for
$-\zeta = 9$. The red and blue shaded regions show two particular blocks of heat-flux events corresponding to warm-updraft and cold-downdraft quadrants, respectively, which persist for a time
$T_{B}$ of around 10–20 s.

We convert the block duration (
$T_{B}$) to a streamwise length by using Taylor's hypothesis, that is, multiplying
$T_{B}$ with the mean wind speed (
$\overline {u}$) computed over the 30 min period. We then scale
$T_{B}\overline {u}$ with a relevant length scale. The possible candidates as the relevant length scales in an unstable ASL are the measurement height
$z$ and the boundary layer depth
$z_{i}$. However,
$z_{i}$ was not measured directly at SLTEST and hence an alternative large-eddy length scale
$\lambda$ was used by Chowdhuri, McNaughton & Prabha (Reference Chowdhuri, McNaughton and Prabha2019), where
$\lambda$ was computed as the peak wavelength of the horizontal velocity spectrum at
$z=25.7$m. This was based on the observation that the large-scale structures contribute directly to the horizontal velocity spectrum in the ASL (Kaimal et al. Reference Kaimal, Wyngaard, Haugen, Coté, Izumi, Caughey and Readings1976; Panofsky et al. Reference Panofsky, Tennekes, Lenschow and Wyngaard1977; Banerjee et al. Reference Banerjee, Katul, Salesky and Chamecki2015). As discussed by Chowdhuri et al. (Reference Chowdhuri, McNaughton and Prabha2019), a model spectrum of the form
is fitted to the streamwise velocity (
$u$) spectrum, where
$\kappa$ is the streamwise wavenumber and
$a$ and
$b$ are the best-fit constants. By maximizing (2.22) with respect to
$\kappa$,
$\lambda$ is evaluated as
$4{\rm \pi} b/3$. The other details and the rationale behind the computation of
$\lambda$ can be found in Chowdhuri et al. (Reference Chowdhuri, McNaughton and Prabha2019).
The spectrum or the scalewise distribution of the normalized streamwise lengths of the blocks (
$T_{B}\overline {u}/\ell$, where
$\ell$ can be either
$z$,
$\lambda$ or the combination of the two) can be at least a few decades wide, given the large variation in
$T_{B}$, ranging from a minimum of 0.05 s (sampling interval) to a few seconds. For the blocks associated with each
$T^{\prime }\text {--}w^{\prime }$ quadrant, we thus logarithmically bin their scaled streamwise lengths (
$T_{B}\overline {u}/\ell$) into 60 bins, where the minimum and maximum are chosen over all the 30 min periods that fall within a particular stability class. Below we discuss the method to compute the Reynolds stress anisotropy associated with these blocks of different normalized streamwise lengths.
2.3.1. The distribution of the Reynolds stress anisotropy
For any particular
$T^{\prime }\text {--}w^{\prime }$ quadrant, we collect all the blocks of the heat-flux events having their normalized streamwise lengths between
where
$(T_{B}\overline {u}/\ell )_{bin}\{m\}$ is the logarithmically binned value,
$\textrm {d}\log (T_{B}\overline {u}/\ell )$ is the bin width and
$m$ is the index of the bin (
$1\leq m \leq 60$). The bin width is defined as
We construct the anisotropy Reynolds stress tensor associated with these blocks as
\begin{align} & \langle b_{ij}\mid[(T_{B}\overline{u}/\ell)_{bin}\{m\}<(T_{B}\overline{u}/\ell)< (T_{B}\overline{u}/\ell)_{bin}\{m\}+\textrm{d} \log(T_{B}\overline{u}/\ell)]\rangle \nonumber\\ &\quad =\frac{{\displaystyle\sum{{u^\prime}_{i} {u^\prime}_{j}}}}{\displaystyle\sum{{u^\prime}_{i} {u^\prime}_{i}}} - \frac{1}{3}\delta_{ij}, \end{align}
and assign it to a streamwise size of
$(T_{B}\overline {u}/\ell )_{bin}\{m\}$. In (2.24), the terms
$\sum {{u^\prime }_{i} {u^\prime }_{j}}$ are the contributions to the Reynolds stress tensor from all the blocks having their sizes between
$(T_{B}\overline {u}/\ell )_{bin}\{m\}$ and
$(T_{B}\overline {u}/\ell )_{bin}\{m\}+\textrm {d}\log (T_{B}\overline {u}/\ell )$.
Similar to § 2.2.2, we calculate the three coefficients associated with 1-component anisotropy, 2-component anisotropy and 3-component isotropy (
$C_{1 { c}}$,
$C_{2 { c}}$ and
$C_{3 { c}}$) of
$\langle b_{ij}\mid (T_{B}\overline {u}/\ell )_{bin}\{m\} \rangle$ as
\begin{equation} \left.\begin{array}{l} \langle C_{1 { c}} \mid(T_{B}\overline{u}/\ell)_{bin}\{m\} \rangle = {\tilde{e}}_{1b}-{\tilde{e}}_{2b} {,}\\ \langle C_{2 { c}} \mid(T_{B}\overline{u}/\ell)_{bin}\{m\} \rangle = 2({\tilde{e}}_{2b}-{\tilde{e}}_{3b}) {,}\\ \langle C_{3 { c}} \mid(T_{B}\overline{u}/\ell)_{bin}\{m\} \rangle = 3{\tilde{e}}_{3b}+1, \end{array}\right\} \end{equation}
where
${\tilde {e}}_{1b}$,
${\tilde {e}}_{2b}$ and
${\tilde {e}}_{3b}$ are the three eigenvalues of
$\langle b_{ij}\mid (T_{B}\overline {u}/\ell )_{bin}\{m\} \rangle$ with
${\tilde {e}}_{1b}>{\tilde {e}}_{2b}>{\tilde {e}}_{3b}$.
Since we construct the same logarithmic grids of
$(T_{B}\overline {u})/\ell$ for all the runs belonging to a particular stability class, we take the average of these three coefficients over all these periods to reduce the run-to-run variability.
2.3.2. Probability and flux distributions
The probability density function (p.d.f.) of the normalized streamwise lengths of the blocks belonging to any particular
$T^{\prime }\text {--}w^{\prime }$ quadrant is calculated as
where
$N_{b}$ is the number of blocks lying between
and
$N_{tot}$ is the total number of blocks detected over a 30 min period (from the same quadrant). The heat and momentum fluxes within these blocks are defined as
\begin{align} & \langle w^{\prime}x^{\prime} \mid [(T_{B}\overline{u}/\ell)_{bin}\{m\}< (T_{B}\overline{u}/\ell)< (T_{B}\overline{u}/\ell)_{bin}\{m\}+\textrm{d} \log(T_{B}\overline{u}/\ell)] \rangle \nonumber\\ &\quad = \frac{\sum{w^{\prime}x^{\prime}}}{N \times \textrm{d} \log(T_{B}\overline{u}/\ell)} \quad (x=u,T){,} \end{align}
where
$N$ is the number of samples in a 30 min run for the SLTEST data. These heat- and momentum-flux distributions are scaled by the product of the standard deviations
$\sigma _{w}\sigma _{T}$ and
$\sigma _{u}\sigma _{w}$, respectively. When these scaled flux distributions from (2.27) are integrated over the whole spectrum of
$(T_{B}\overline {u})/\ell$, the results show the strength of the coupling between
$w^{\prime }$ and
$T^{\prime }$ (
$u^{\prime }$) from each quadrant of
$T^{\prime }\text {--}w^{\prime }$.
Similar to § 2.3.1, we take the average of the p.d.f.s and the heat- and momentum-flux distributions over all the 30 min periods belonging to a particular stability class. While presenting the results in § 3.3, these averaged distributions of the degree of isotropy, probability and fluxes are referred to as being associated with
$(T_{B}\overline {u})/\ell$ only, without explicitly mentioning that these are the binned values. Apart from that, the amount of spread between the individual 30 min runs for a particular stability class is computed as one standard deviation from the ensemble average and shown as error bars. Similar to quadrant analysis, the results obtained from this method have also been verified for sensitivity to the choice of the number of bins.
3. Results and discussion
We begin by discussing the general characteristics of the Reynolds stress anisotropy with the change in the stability ratio
$-\zeta$. We also highlight the correspondence between the intermittent nature of turbulent heat transport and the Reynolds stress anisotropy. By presenting the relevant results, this correspondence is further investigated in detail, complemented with the quadrant and persistence analyses of the heat-flux events. The possible physical interpretations of these results are also discussed.
3.1. The characteristics of Reynolds stress anisotropy with stability
Figure 3(a) shows the anisotropic states of the Reynolds stress tensor for the 30 min averaged flow plotted on the barycentric map (see (2.10), (2.11) and (2.12)) with the stability ratio
$-\zeta$. The variations of the three associated coefficients (
$C_{1 { c}}$,
$C_{2 { c}}$ and
$C_{3 { c}}$) with
$-\zeta$ are provided in figure S1a in the supplementary figures (available at https://doi.org/10.1017/jfm.2020.471). As shown in figure 1, the barycentric map is spanned by an equilateral triangle, which can be divided into three regions R1, R2 and R3, where the anisotropic states of
$b_{ij}$ are dominated by 1-component anisotropy, 2-component anisotropy and 3-component isotropy, respectively. As evident from figure 3(a), the anisotropic states of
$b_{ij}$ move towards the region R3 from the region R2 as
$-\zeta$ approaches the local free convection limit (
$-\zeta >1$). This implies that the anisotropic state of
$b_{ij}$ is more dominated by the 3-component isotropy as the surface layer becomes highly convective. The reason for this is that, in a highly convective surface layer, the turbulent kinetic energy (TKE) is mainly generated in the vertical direction through a buoyancy production term, while in the horizontal direction the production of TKE due to shear is almost negligible. However, the pressure–strain correlation in highly convective conditions efficiently redistributes the TKE generated in the vertical to the horizontal direction, thus driving the turbulence to be dominated by the 3-component isotropy (McBean & Elliott Reference McBean and Elliott1975; Zhuang Reference Zhuang1995; Bou-Zeid et al. Reference Bou-Zeid, Gao, Ansorge and Katul2018). The effectiveness of pressure–strain correlation in redistributing the TKE in highly convective conditions is related to the covariance between the pressure and vertical velocity fluctuations, as detailed in the physical model of McBean & Elliott (Reference McBean and Elliott1975). On the other hand, for small values of
$-\zeta$ (
$0<-\zeta <0.2$), the anisotropic state of
$b_{ij}$ is dominated by the 2-component anisotropy, as the blue shaded points in figure 3(a) remain concentrated within the region R2. From table 2, it is clear that the near-neutral stability class (
$0<-\zeta <0.2$) corresponds to the lowest three levels of the SLTEST experiment (
$z=1.4$, 2.1, 3.0 m), where due to the blocking of the ground the vertical velocity fluctuations are suppressed. Therefore, the turbulence very close to the ground is in a 2-component anisotropic state dominated by the horizontal velocity components. This is in agreement with the studies by Krogstad & Torbergsen (Reference Krogstad and Torbergsen2000) and Ali et al. (Reference Ali, Hamilton, Cortina, Calaf and Cal2018).
Plots of the (
$a$) anisotropic states of
$b_{ij}$ on the barycentric map (see (2.11) and (2.12)), (
$b$) the degree of isotropy
$C_{3 { c}}$ (see (2.8)) versus the wind shear (
$\partial {\overline {u}}/\partial z$), and (
$c$) scaled vertical velocity and temperature standard deviations (
$\sigma _{w}/u_{*}$ and
$\sigma _{T}/T_{*}$),
$C_{3 { c}}$ and the correlation coefficient between the vertical velocity and temperature (
$R_{wT}$) versus the stability ratio (
$-\zeta$). In panel (
$c$), the left
$y$ axis is logarithmic, the right
$y$ axis is linear and the
$x$ axis is reversed such that the
$-\zeta$ values proceed from large to small. The thick blue and red lines denote the local free convection scalings for
$\sigma _{w}/u_{*}$ and
$\sigma _{T}/T_{*}$. The colour bar at the bottom corresponds to both panels (a) and (b), showing the stability ratios as
$\log _{10}(-\zeta )$.

Apart from the anisotropic states of
$b_{ij}$, we can also evaluate its degree of isotropy
$C_{3 { c}}$ to quantify how close the turbulence is towards the 3-component isotropy. From figure 3(b) we note that for
$0<-\zeta <0.2$, strong anisotropic turbulence (
$C_{3 { c}} \approx 0.1$) is associated with large wind shear (
$\partial {\overline {u}}/\partial {z}$). The wind shear is approximated using the finite-difference scheme (Stull Reference Stull1988; Arya Reference Arya2001) as
where
$z_{2}$ and
$z_{1}$ are the two adjacent levels from the SLTEST data with
$z_{2}>z_{1}$ and
$z_{m}=(z_{2}+z_{1})/2$. However, in the limit of local free convection (
$-\zeta > 1$), the effect of wind shear is weak and the turbulence is more dominated by 3-component isotropy (
$C_{3 {c}}>1/3$). This result is consistent with the observations of Stiperski & Calaf (Reference Stiperski and Calaf2018), where they found that, in an unstable surface layer as we approach
$z \to 0$ (associated with small values of
$-\zeta$), the anisotropic characteristics of turbulence are dominated by strong wind shear. On the contrary, as the local free convection is approached (
$-\zeta > 1$), the effect of wind shear weakens and the turbulence becomes less anisotropic. This is in agreement with Jin, So & Gatski (Reference Jin, So and Gatski2003), where they showed both analytically and from numerical simulations that, in a buoyant shear flow, the effect of increase (decrease) in buoyancy (shear) was to drive the turbulence towards isotropy.
To investigate this further, figure 3(c) shows the scatter plot of the scaled vertical and temperature standard deviations (
$\sigma _{w}/u_{*}$ and
$\sigma _{T}/T_{*}$) along with the correlation coefficient between
$w$ and
$T$ (
$R_{wT}$) and the degree of isotropy (
$C_{3 { c}}$), against the stability ratio
$-\zeta$. Note that the temperature scale (
$T_{*}$) is defined here as
$H_{0}/u_{*}$ with the omission of the negative sign, to keep the quantity
$\sigma _{T}/T_{*}$ positive. The local free convection scalings for
$\sigma _{w}/u_{*}$ and
$\sigma _{T}/T_{*}$ are given as
where the coefficients are fitted from the data and match well with the values reported by Wyngaard, Coté & Izumi (Reference Wyngaard, Coté and Izumi1971). It is interesting to note that, after
$-\zeta < 0.5$, the local free convection scaling does not hold for
$\sigma _{w}/u_{*}$, but it extends for
$\sigma _{T}/T_{*}$. Khanna & Brasseur (Reference Khanna and Brasseur1997) explained this as: the buoyancy-induced motions contribute more to the temperature fluctuations than the shear-induced motions.
In an unstable surface layer, the horizontal velocity variances depend on the global stability ratio
$-z_{i}/L$, rather than on
$-\zeta$ (Monin & Yaglom Reference Monin and Yaglom1971; Panofsky Reference Panofsky1974; Panofsky et al. Reference Panofsky, Tennekes, Lenschow and Wyngaard1977; Wyngaard Reference Wyngaard2010). Therefore, the variation in degree of isotropy (
$C_{3 { c}}$) with
$-\zeta$ is mainly determined by the strength of the vertical velocity fluctuations (
$\sigma _{w}$), decreasing from
$C_{3 { c}} \approx 0.6$ to
$C_{3 { c}} \approx 0.1$ as
$\sigma _{w}/u_{*}$ decreases with
$-\zeta$ (figure 3c). From figure 3(c) and supplementary figure S1b, we note that
$w^{\prime }$ is more strongly coupled to
$T^{\prime }$ than to
$u^{\prime }$ in the local free convection (
$R_{wT}\approx$ 0.65 and
$R_{uw}\approx -0.05$). However, with decrease in
$-\zeta$, the correlation coefficient between
$w^{\prime }$ and
$u^{\prime }$ increases (
$R_{uw}\approx -0.25$) whereas it decreases between
$w^{\prime }$ and
$T^{\prime }$ (
$R_{wT}\approx 0.4$). This is also reflected in the transport efficiencies of heat (
$\eta _{wT}$) and momentum (
$\eta _{uw}$), defined as
\begin{equation} \eta_{wx}=\frac{\left({\displaystyle\sum{w^{\prime}x^{\prime}}}\right)_{{\mathit{down}\text{-}\mathit{gradient}}}+ \left({\displaystyle\sum{w^{\prime}x^{\prime}}}\right)_{{\mathit{counter}\text{-}\mathit{gradient}}}} {({\sum{w^{\prime}x^{\prime}}})_{{\mathit{down}\text{-}\mathit{gradient}}}}, \end{equation}
where
$x$ can be either
$u$ or
$T$ (Li & Bou-Zeid Reference Li and Bou-Zeid2011). From supplementary figure S1b, it is evident that in local free convection
$\eta _{uw}\to 0$ whereas
$\eta _{wT}$ almost approaches a constant value of 0.9. However, with decrease in
$-\zeta$,
$\eta _{uw}$ increases to
${\approx } 0.6$ and
$\eta _{wT}$ decreases to
${\approx } 0.75$. We next investigate the p.d.f.s of
$T^{\prime }$,
$w^{\prime }$ and
$u^{\prime }$ to establish a correspondence between the turbulence anisotropy and its transport characteristics.
Figure 4(a) and (b) show the skewness and kurtosis of the p.d.f.s of
$T^{\prime }$,
$w^{\prime }$ and
$u^{\prime }$ (
$\overline {{x^{\prime }}^3}/\sigma _{x}^3$ and
$\overline {{x^{\prime }}^4}/\sigma _{x}^4$, where
$x=u,w,T$) along with the degree of isotropy (
$C_{3 { c}}$). The associated p.d.f.s are shown in supplementary figure S2. For a perfect Gaussian distribution, the skewness and kurtosis have values of 0 and 3, respectively. Physically, the skewness is associated with the asymmetry in the p.d.f.s whereas the kurtosis is related to intermittency (Davidson Reference Davidson2015).
The scatter plot of the (
$a$) skewness and (
$b$) kurtosis of the temperature, vertical velocity and streamwise velocity fluctuations (
$T^{\prime }$,
$w^{\prime }$ and
$u^{\prime }$) are shown against
$-\zeta$. The red, blue and pink open circles denote
$T^{\prime }$,
$w^{\prime }$ and
$u^{\prime }$, respectively, with their skewness and kurtosis being plotted on the left-hand
$y$ axis. The black stars show the degree of isotropy (
$C_{3{c}}$, see (2.8)) with its values being plotted on the right-hand
$y$ axis. The thick horizontal black lines denote the values of 0 and 3, which are the skewness and kurtosis for the Gaussian distribution.

From figure 4(a) and (b), it is clear that the skewness and kurtosis of the temperature fluctuations are strongly non-Gaussian (
${\approx } 1.5$ and 5, respectively) in the local free convection limit (
$-\zeta > 1$). The strong non-Gaussian nature of temperature fluctuations in highly unstable conditions is remarkably consistent with the previous studies in the ASL (Chu et al. Reference Chu, Parlange, Katul and Albertson1996; Garai & Kleissl Reference Garai and Kleissl2013; Lyu et al. Reference Lyu, Hu, Liu, Xu and Cheng2018). Similar behaviour has also been observed in turbulent Rayleigh–Bénard convection experiments of Adrian, Ferreira & Boberg (Reference Adrian, Ferreira and Boberg1986). The strong non-Gaussianity in
$T^{\prime }$ in highly unstable conditions is caused due to the intermittent bursts associated with warm updrafts, interspersed with relatively more frequent quiescent cold downdrafts bringing well-mixed air from aloft (Adrian et al. Reference Adrian, Ferreira and Boberg1986; Chu et al. Reference Chu, Parlange, Katul and Albertson1996). However, the skewness and kurtosis of
$T^{\prime }$ become closer to Gaussian (0.5 and 3, respectively) for the near-neutral stability (
$0 <-\zeta < 0.2$). The close-to-Gaussian characteristics of the
$T^{\prime }$ p.d.f.s in a near-neutral ASL are in agreement with Chu et al. (Reference Chu, Parlange, Katul and Albertson1996) and with the pipe flow experiment of Nagano & Tagawa (Reference Nagano and Tagawa1988), where temperature behaved more like a passive scalar.
On the other hand, the p.d.f.s of
$u^{\prime }$ remain near-Gaussian for all values of
$-\zeta$, with its skewness and kurtosis approaching 0 and 3, respectively. For
$w^{\prime }$, the skewness stays almost constant at 0.4 to 0.5 for all values of
$-\zeta$, implying the consistent upward transport of vertical kinetic energy (Chiba Reference Chiba1978; Hunt, Kaimal & Gaynor Reference Hunt, Kaimal and Gaynor1988). However, the kurtosis for
$w^{\prime }$ increases from 3 to 4 as
$-\zeta$ decreases. This observation is consistent with Chu etal. (Reference Chu, Parlange, Katul and Albertson1996), where they found that the kurtosis in
$w^{\prime }$ increased from 3.12 in highly unstable conditions to 3.77 in near-neutral conditions. Chiba (Reference Chiba1984) postulated that this increase in the kurtosis of
$w^{\prime }$ at small
$-\zeta$ values is related to the increasing importance of the small-scale eddies near the ground. However, Hong et al. (Reference Hong, Choi, Ishikawa and Kim2004) hypothesized it to be related to the low-speed streaks, initiating inactive and active turbulence interactions with increasing intermittency.
We note that the degree of isotropy (
$C_{3 { c}}$) also decreases in a similar way as the skewness and kurtosis of the temperature fluctuations approach a near-Gaussian distribution with decrease in
$-\zeta$ (figure 4). Katul et al. (Reference Katul, Hsieh, Kuhn, Ellsworth and Nie1997) demonstrated that the temperature skewness was directly related to the difference in the time fractions (
${\rm \Delta} T_{f}$) of the warm-updraft and cold-downdraft events (asymmetry) as
where
$Q_{3}=\overline {{T^{\prime }}^3}/\sigma _{T}^3$, by assuming that the time fractions spent in the counter-gradient quadrants of the
$T^{\prime }\text {--}w^{\prime }$ plane could be ignored. This implies that the asymmetry in the distributions of the warm-updraft and cold-downdraft events associated with the skewness of the temperature fluctuations, has a strong correspondence with the anisotropy in the Reynolds stress tensor. In supplementary figure S3, we show the heat-flux fractions (
$F_{f}$) and the time fractions (
$T_{f}$) associated with each quadrant of the
$T^{\prime }\text {--}w^{\prime }$ plane. It indicates that, in highly unstable conditions (
$-\zeta > 1$), the warm updrafts carry more heat flux even though they occur for less time than the cold downdrafts (see supplementaryfigure S3c).
The same observation can also be made from figure 5(a), where the strong non-Gaussianity in the temperature fluctuations in highly unstable conditions (
$-\zeta >2$) introduces a large asymmetry in the p.d.f.s of the scaled heat flux (
$P(\hat {w}\hat {T})$). The intermittent bursts associated with warm updrafts characterized by large kurtosis carry more heat flux than predicted by the distribution if
$\hat {w}$ and
$\hat {T}$ were both standard Gaussian random variables. According to Krogstad (Reference Krogstad2013), the p.d.f. of the product of two standard Gaussian random variables
$\hat {x}$ and
$\hat {y}$ can be expressed as
where
$\mathrm {K}_{0}(|\hat {x}\hat {y}|)$ is the modified Bessel function of the second kind. However, this strong non-Gaussianity is not felt in the p.d.f.s of the scaled momentum flux (
$P(\hat {u}\hat {w})$) because the probability distributions of the
$u$ and
$w$ fluctuations are closer to Gaussian compared to those for temperature (figures 5(b) and S2).
The p.d.f.s of the scaled (
$a$) heat flux (
$P(\hat {w}\hat {T})$) and (
$b$) momentum flux (
$P(\hat {u}\hat {w})$) are shown for the six different classes of
$-\zeta$ as indicated in the legend on the right. The thick black curves denote the modified Bessel function of the second kind, which corresponds to the p.d.f. of
$\hat {w}\hat {x}$ (
$x=u,T$), if
$\hat {u}$,
$\hat {w}$ and
$\hat {T}$ were all standard Gaussian random variables (see (3.5)). The grey shaded portions show the hyperbolic hole, defined as
$|\hat {x}\hat {w}|=1$ (
$x=u,T$).

In a nutshell, from figures 4 and 5 one can infer that the characteristics of the turbulent heat transport in an unstable surface layer are strongly (weakly) non-Gaussian for highly (feebly) convective conditions, associated with less (more) anisotropic turbulence. However, till now we have presented the anisotropic characteristics of the averaged flow, which comprises the heat-flux events from all four quadrants of
$T^{\prime }\text {--}w^{\prime }$. Given the asymmetric and intermittent nature of turbulent heat transport, it is thus imperative to employ event-based analysis to investigate ‘whether the strongest (weakest) heat-flux events are associated with less (more) anisotropic turbulence’. Therefore, we turn our attention towards the quadrant analysis to deduce the anisotropic characteristics of the Reynolds stress tensor, associated with the heat-flux events from the four different quadrants.
3.2. Quadrant analysis of Reynolds stress anisotropy
From quadrant analysis, we study the detailed correspondence between the Reynolds stress anisotropy and the heat-flux events of varying intensities with their frequency of occurrences. Figure 6 shows the
$RGB$ colour map computed by (2.20) with the superposed contours of degree of isotropy (see (2.19)) on the
$\hat {T}\text {--}\hat {w}$ quadrant plane. From the
$RGB$ colour map, the anisotropic states of the Reynolds stress tensor in the red, green and blue shaded regions of the
$\hat {T}\text {--}\hat {w}$ quadrant plane are dominated by 1-component anisotropy, 2-component anisotropy and 3-component isotropy, respectively. We also include the hyperbolic hole, defined as
$|\hat {T}\hat {w}|=1$, to identify the strong heat-flux-producing events that lie in the region outside of it (Smedman et al. Reference Smedman, Högström, Hunt and Sahlée2007). The six different panels in figure 6(a)–(f) correspond to the six different stability classes as mentioned in table 2.
The quadrant maps of the degree of isotropy (
$\langle C_{3{c}} \mid (\hat {T},\hat {w}) \rangle$, see (2.19)) plotted on the
$\hat {T}\text {--}\hat {w}$ quadrant plane are shown for six different classes of the stability ratios: (
$a$)
$-\zeta >2$, (
$b$)
$1<-\zeta <2$, (
$c$)
$0.6 < -\zeta <1$, (
$d$)
$0.4<-\zeta <0.6$, (
$e$)
$0.2< -\zeta <0.4$ and (
$f$)
$0<-\zeta < 0.2$. The anisotropic states of
$\langle b_{ij}\mid (\hat {T},\hat {w}) \rangle$ are represented by the
$RGB$ colour map such that the red, green and blue shaded regions of the quadrant plane are dominated by 1-component anisotropy, 2-component anisotropy and 3-component isotropy, respectively (see (2.20)). The thick pink lines denote the hyperbolic hole
$|\hat {T}\hat {w}|=1$. Quadrants I and III represent the warm updrafts and cold downdrafts, whereas quadrants II and IV represent the cold updrafts and warm downdrafts.

From figure 6(a), we notice that, in highly convective conditions (
$-\zeta >2$), the anisotropic states of the Reynolds stress tensor for strong heat-flux events (
$|\hat {T}\hat {w}|>1$) are mostly dominated by either 3-component isotropy or 1-component anisotropy (indicated by blue and red, respectively). However, for weak heat-flux events (
$|\hat {T}\hat {w}|<1$) the anisotropic states of the Reynolds stress tensor are dominated by 2-component anisotropy (indicated by green). This implies that the influence of the three limiting states of the Reynolds stress tensor are associated with specific heat-flux events, residing within the red (1-component anisotropy), blue (3-component isotropy) and green (2-component anisotropy) regions of the
$\hat {T}\text {--}\hat {w}$ quadrant plane. We also notice from figure 6(a) that the zones of 3-component isotropic states (blue regions) reside mainly within the warm-updraft (
$C_{3{c}}\approx 0.5$) and cold-downdraft (
$C_{3{c}}\approx 0.4$) quadrants. On the other hand, from figure 7(a) we note that the j.p.d.f. contours between
$\hat {T}$ and
$\hat {w}$ depart significantly from the bivariate Gaussian distribution (see (2.18)) in highly convective conditions. By comparing the features in figure 6(a) with the j.p.d.f. contours in figure 7(a), we observe that the 1-component anisotropy zones (red regions) are associated with extremely low-probability events of very high heat fluxes, located well beyond the hyperbolic hole (
$|\hat {T}\hat {w}|\gg 1$).
The contour maps of the j.p.d.f.s between
$\hat {T}$ and
$\hat {w}$ (
$P(\hat {T},\hat {w})$ as thick black lines, see (2.17)) and the bivariate Gaussian distribution (
$G(\hat {T},\hat {w})$ as dotted red lines, see (2.18)) are shown for six different classes of the stability ratios: (
$a$)
$-\zeta >2$, (
$b$)
$1<-\zeta <2$, (
$c$)
$0.6<- \zeta <1$, (
$d$)
$0.4<-\zeta <0.6$, (
$e$)
$0.2< -\zeta <0.4$ and (
$f$)
$0< -\zeta <0.2$. The same
$RGB$ colour map and hyperbolic hole from figure 6 are shown here too.

However, as
$-\zeta$ becomes smaller, the j.p.d.f. contours become progressively close to bivariate Gaussian distribution (figure 7a–f), with the green regions (2-component anisotropy) being systematically more prominent (figure 6a–f). On the other hand, the blue regions (3-component isotropy) become systematically less visible (figure 6a–f). This is consistent with figure 3(a), where the anisotropic states of the Reynolds stress tensor become progressively more dominated by 2-component anisotropy as the near-neutral stability is approached. Furthermore, this is also in agreement with figure 4, where highly anisotropic turbulence is associated with an almost symmetrical distribution of the warm updrafts and cold downdrafts in near-neutral stability, due to the small values of skewness in
$T^{\prime }$ (see (3.4)).
It is interesting to note that the 1-component anisotropy indicated by the red regions in the
$\hat {T}\text {--}\hat {w}$ quadrant plane does not appear to have a signature in the 30 min averaged Reynolds stress anisotropy (figure 3a). This is because this anisotropic state is associated with highly intermittent low-probability events of very high heat fluxes. In addition to that, the Reynolds stress anisotropy is dominated by the 3-component isotropic state specifically for those heat-flux events which reside within the blue regions of the
$\hat {T}\text {--}\hat {w}$ quadrant plane (figure 6). This observation is non-trivial and this outcome would not be possible without an event-based description. Since the approaches based on time-averaged statistics would predict that higher convective conditions (high heat fluxes) are associated with less anisotropic turbulence, this analysis shows that the connection between the intensity of the heat flux and turbulence anisotropy is more intricate than that. Therefore, one can ask ‘whether there are any characteristic sizes of the heat-flux events associated with least anisotropic turbulence’. However, the quadrant analysis does not give information about the time scale or size of the heat-flux events. We thus focus our attention on persistence analysis to investigate the anisotropic states of the Reynolds stress tensor associated with the streamwise sizes of the heat-fluxevents.
3.3. Persistence analysis of Reynolds stress anisotropy
We employ persistence analysis to characterize the streamwise sizes of the heat-flux events from each quadrant of
$T^{\prime }\text {--}w^{\prime }$. This is achieved by converting the persistent time
$T_{B}$ to a streamwise length
$T_{B}\overline {u}$ from Taylor's hypothesis. We begin by discussing the persistence p.d.f.s to highlight the physical characteristics of these heat-flux events and the aspect of non-Gaussianity. Along with that, we also investigate the anisotropic states of the Reynolds stress tensor associated with these heat-flux events of different sizes. The spread in the averaged plots is shown as the error bars, computed as one standard deviation from the ensemble mean for a particular stability class.
3.3.1. Persistence p.d.f.s of heat-flux events
Figure 8(a)–(f) show the p.d.f.s of the normalized streamwise sizes (
$(T_{B}\overline {u})/z$) for the heat-flux events occurring in each quadrant of
$T^{\prime }\text {--}w^{\prime }$, corresponding to the six different stability classes (table 2). We choose to normalize the streamwise sizes by
$z$, under the assumption that these heat-flux events are associated with the thermal plumes which grow linearly with height (Tennekes & Lumley Reference Tennekes and Lumley1972). The associated histograms of
$(T_{B}\overline {u})/z$ for the heat-flux events from each quadrant are also shown in figure 17(a)–(f) (see the Appendix). Typically, for the warm-updraft and cold-downdraft quadrants we encounter 100–200 heat-flux events corresponding to the large sizes
$(T_{B}\overline {u})/z>4$.
The persistence p.d.f.s of the heat-flux events are shown for the six different stability classes: (
$a$)
$-\zeta >2$, (
$b$)
$1<-\zeta <2$, (
$c$)
$0.6<-\zeta <1$, (
$d$)
$0.4<-\zeta <0.6$, (
$e$)
$0.2<-\zeta <0.4$ and (
$\,f$)
$0<-\zeta <0.2$. The markers for different quadrants are explained in the legend in panel (a). A distinct power-law of exponent
$-0.4$ is shown as a thick grey line in all the panels. The thick black lines correspond to the log-normal distribution. The error bars in all the panels show the existing spread between individual 30 min runs for each of the stability classes, computed as one standard deviation from the ensemble mean.

The most distinct feature we notice from the highly convective (
$-\zeta > 2$) stability class (figure 8a) is that the persistence p.d.f.s of the warm-updraft and cold-downdraft events collapse with a power law of an exponent
$-0.4$,
which approximately extends up to
$(T_{B}\overline {u})/z \approx 1$. A similar power law was reported by Chamecki (Reference Chamecki2013) for the persistent p.d.f.s of
$u$ and
$w$ fluctuations smaller than the integral time scale in a plant canopy. Apart from that, Yee et al. (Reference Yee, Kosteniuk, Chandler, Biltoft and Bowers1993) and Katul et al. (Reference Katul, Albertson, Parlange, Chu and Stricker1994) also documented a power-law behaviour in the p.d.f.s of the small sizes of the heat-flux bursts from an unstable ASL. Additionally, we also note that the p.d.f.s of
$(T_{B}\overline {u})/z$ for the counter-gradient events are in close agreement with the p.d.f.s of the down-gradient events for small values of
$(T_{B}\overline {u})/z< 0.2$. Beyond those sizes, the p.d.f.sof the counter-gradient events drop faster than the down-gradient events. This implies that these counter-gradient events have a statistical tendency to occur in smaller sizes and do not persist for a long time. This is in agreement with the simulations of Dong et al. (Reference Dong, Lozano-Durán, Sekimoto and Jiménez2017), where they found that the p.d.f.s of counter-gradient and down-gradient momentum events of small sizes agreed with each other and diverged for larger sizes.
However, this power-law segment systematically disappears as we approach the near-neutral stability (
$0 <-\zeta < 0.2$) and gets replaced by a log-normal distribution (figure 8f). This is broadly in agreement with Sreenivasan & Bershadskii (Reference Sreenivasan and Bershadskii2006), where they commented that, for an active scalar such as temperature in highly convective turbulence, the p.d.f.s of the inter-pulse periods followed a power law. Conversely, in a shear-driven turbulence when the temperature behaved more like a passive scalar, the p.d.f.s followed a log-normal distribution.
For
$(T_{B}\overline {u})/z> 1$, the p.d.f.s of the warm-updraft and cold-downdraft events significantly differ from each other in the highly convective case (
$-\zeta > 2$, figure 8a). However, they systematically agree with each other as the near-neutral stability (
$0<-\zeta < 0.2$) is approached (figure 8a–f). As we will show later, this is related to the asymmetry in the distributions of the warm-updraft and cold-downdraft events due to strong non-Gaussianity in temperature fluctuations in a highly convective surface layer. As discussed by Chamecki (Reference Chamecki2013), these large values of
$(T_{B}\overline {u})/z$ are exponentially distributed according to a Poisson-type process, which could be studied by considering the cumulative distribution functions (CDFs) of
$(T_{B}\overline {u})/z$. In general, the CDFs are comparatively smoother than the p.d.f.s, thus yielding a more robust fit for the exponential distribution. The CDF (
$F[(T_{B}\overline {u})/z]$) is defined as
\begin{equation} F [(T_{B}\overline{u})/z ]= \int_{{[(T_{B}\overline{u})/z]}_{max}}^{[(T_{B}\overline{u})/z]} P [(T_{B}\overline{u})/z ]\,\textrm{d}\log[(T_{B}\overline{u})/z]. \end{equation} Figure 9 shows the CDFs of the heat-flux events from the four quadrants of
$T^{\prime }\text {--}w^{\prime }$. We note that the power-law region is not seen clearly in the CDFs, as they converge to 1 for the small sizes
$(T_{B}\overline {u})/z<$ 1. For the large sizes (
$1 \leq (T_{B}\overline {u})/z \leq 10$), we plot the CDFs in a log-linear coordinate system (see the insets in figure 9) such that the exponential decay,
in such plots would appear as a straight line with a slope of
$-k$. From the insets in figure9, we notice that, for larger values of
$(T_{B}\overline {u})/z$,
$F[(T_{B}\overline {u})/z]$ indeed decays exponentially according to (3.8). We also find that the slopes corresponding to the warm-updraft and cold-downdraft events are significantly different from each other (
$k=1.3$ and
$k=0.7$, respectively) for highly convective stability (figure 9a). However, these two slopes become systematically close to each other as the near-neutral stability is approached (
$k=0.65$, figure 9f). On the other hand, with stability, no appreciable change in the slope is observed for the counter-gradient events. This difference in the slopes for warm-updraft and cold-downdraft events is linked to the non-Gaussianity in temperature fluctuations in a highly convective surface layer. Sreenivasan et al. (Reference Sreenivasan, Prabhu and Narasimha1983) mentioned that the long intervals (large
$(T_{B}\overline {u})/z$) are a consequence of large-scale structures passing the sensor and the short intervals (small
$(T_{B}\overline {u})/z$) are a consequence of the nibbling small-scale motions superposed on the large-scale structures. From that perspective, we expect that the non-Gaussian characteristics of the warm-updraft and cold-downdraft events might be related to the large-scale structures.
Same as figure 8, but the CDFs are shown instead of the p.d.f.s.: (
$a$)
$-\zeta >2$, (
$b$)
$1<-\zeta <2$, (
$c$)
$0.6<-\zeta <1$, (
$d$)
$0.4<-\zeta <0.6$, (
$e$)
$0.2<-\zeta <0.4$ and (
$f$)
$0<-\zeta <0.2$. In each panel, the inset shows the enlarged area between
$1\leq (T_{B}\overline {u})/z \leq 10$ (grey shaded region), where the CDFs are plotted on log-linear axes to indicate the exponential decay as a straight line. The markers are explained in the legend in panel (f). The equations related to the exponential decay are shown in the legend of each panel.

To summarize, from figures 8 and 9 we have observed that the warm-updraft and cold-downdraft events having sizes
$(T_{B}\overline {u})/z<$ 1 are scale-invariant owing to a power-law dependence in the highly convective stability. This scale-invariant property disappears systematically as the near-neutral stability is approached. Apart from that, the effect of non-Gaussianity (Gaussianity) appears mostly at the sizes
$(T_{B}\overline {u})/z> 1$ in a highly (weakly) convective surface layer, possibly associated with the large-scale structures (Sreenivasan et al. Reference Sreenivasan, Prabhu and Narasimha1983). We will revisit this while investigating the linkage between the persistence p.d.f.s and the degree of isotropy of the Reynolds stress tensor in § 3.3.4. Next we discuss the anisotropy characteristics of the Reynolds stress tensor associated with these heat-flux events of different sizes.
3.3.2. The degree of isotropy, heat- and momentum-flux distributions
We begin by discussing the amount of heat flux associated with the normalized streamwise sizes
$(T_{B}\overline {u})/z$ (see (2.27)), corresponding to the six different stability classes. From figure 10(a) and (b), we note that the
$z$-scaling of the streamwise sizes of the heat-flux events collapses the scaled heat-flux peak positions at
$(T_{B}\overline {u})/z \approx 2$ and
$(T_{B}\overline {u})/z \approx 2.5$ for the warm-updraft and cold-downdraft quadrants, respectively. We also observe that the heat-flux events from the counter-gradient quadrants contribute insignificantly to the total heat flux (figure 10c and d). This result is in agreement with the heat-flux fractions shown in supplementary figure S3. To infer whether the least anisotropic turbulence is associated with the peak positions of the heat flux, we investigate the distributions of the degree of isotropy (
$C_{3 { c}}$, see (2.25)) associated with these events.
The heat-flux distributions (see (2.27)) are plotted against the normalized sizes
$(T_{B}\overline {u})/z$ of the heat-flux events corresponding to (a) warm-updraft, (b) cold-downdraft, (c) cold-updraft and (d) warm-downdraft quadrants. In the upper panels, the black arrows indicate the collapsed position of the peaks of the heat-flux distribution associated with the warm updrafts and cold downdrafts. The different colours represent the six different stability classes as indicated in the legend on the right.

From figure 11, we note that the down-gradient heat-flux events corresponding to warm-updraft and cold-downdraft quadrants are associated with larger values of
$C_{3 { c}}$, compared to the counter-gradient heat-flux events from cold-updraft and warm-downdraft quadrants. Therefore, we may infer that the counter-gradient events which carry significantly less heat flux are associated with more anisotropic turbulence than the warm-updraft and cold-downdraft events. Apart from that, we observe that there is a critical size of warm-updraft and cold-downdraft events associated with the maximum value of
$C_{3 { c}}$ and this critical size is larger for the cold downdrafts compared to the warm updrafts. Also, the maximum value of
$C_{3 { c}}$ decreases systematically as the near-neutral stability is approached. This is in agreement with our previous observations for the averaged flow, where the degree of isotropy systematically decreased from highly convective to near-neutral stability (figures 3 and 4). Moreover, the heat-flux events corresponding to warm updrafts are associated with relatively less anisotropic turbulence than the cold downdrafts, as the values of
$C_{3 { c}}$ are larger in general. However, the peak positions of
$C_{3 { c}}$ associated with warm-updraft and cold-downdraft events do not match with the peak positions of the heat-flux distributions (figure 11a and b). This mismatch is more apparent for the warm-updraft events than the cold downdrafts. The deviation in the peak positions of
$C_{3 { c}}$ and heat-flux distributions complements the results from the quadrant analysis (figures 6 and 7), where we found that the large-heat-flux events do not necessarily relate to the least anisotropic turbulence.
The distributions of the degree of isotropy (
$C_{3 { c}}$) are plotted against the normalized sizes
$(T_{B}\overline {u})/z$ of the heat-flux events corresponding to (a) warm-updraft, (b) cold-downdraft, (c) cold-updraft and (d) warm-downdraft quadrants. In the upper panels, the black arrows indicate the collapsed position of the peaks of the heat-flux distribution associated with the warm updrafts and cold downdrafts (figure 10).

From the definition of isotropy (Könözsy Reference Könözsy2019), there are two possible reasons contributing to the Reynolds stress anisotropy associated with the sizes of the warm-updraft and cold-downdraft events:
(i) The amplitudes of the horizontal velocity fluctuations exceed the vertical velocity fluctuations.
(ii) The vertical velocity fluctuations contribute substantially to the upward or downward transport of streamwise momentum.
To investigate the first of the two aforementioned reasons, in figure 12 we show the ratios of the vertical and horizontal velocity variances associated with the normalized sizes
$(T_{B}\overline {u})/z$ of the warm-updraft and cold-downdraft events. The velocity variances associated with
$(T_{B}\overline {u})/z$ from each quadrant of
$T^{\prime }\text {--}w^{\prime }$ are defined similarly as in § 2.3.2, such that
\begin{align} & \langle {x^{\prime}}^2 \mid [(T_{B}\overline{u}/z)_{bin}\{m\}< (T_{B} \overline{u}/z)< (T_{B}\overline{u}/z)_{bin}\{m\}+\textrm{d} \log(T_{B}\overline{u}/z)] \rangle \nonumber\\ &\quad =\frac{\sum{{x^{\prime}}^2}}{N \times \textrm{d} \log(T_{B}\overline{u}/z)} \quad (x=u,v,w) {,} \end{align}
where the symbols carry their usual meaning. Since in isotropic turbulence the three velocity variances in
$x$,
$y$ and
$z$ directions are equal to each other, it follows that
\begin{equation} \frac{\langle {w^{\prime}}^2 \mid [(T_{B}\overline{u})/z] \rangle}{\langle {u^{\prime}}^2 \mid [(T_{B}\overline{u})/z] \rangle + \langle {v^{\prime}}^2 \mid [(T_{B}\overline{u})/z] \rangle}=0.5. \end{equation}The distributions of the degree of isotropy are plotted against the normalized sizes
$(T_{B}\overline {u})/z$ of the heat-flux events from (
$a$) warm-updraft and (
$b$) cold-downdraft quadrants. The ratio between the vertical and horizontal velocity variances are plotted against the normalized sizes
$(T_{B}\overline {u})/z$ of the heat-flux events from (
$c$) warm-updraft and (
$d$) cold-downdraft quadrants. The horizontal black arrows in panels (
$c$) and (
$d$) indicate the ratio 0.5.

From figure 12(c) and (d), we note that, for the sizes of the warm-updraft and cold-downdraft events
$(T_{B}\overline {u})/z < 2$ and
$(T_{B}\overline {u})/z < 2.5$, respectively, the horizontal velocity variances dominate, since the variance ratio is smaller than 0.5. At the peak position of the heat flux, the variance ratio indeed becomes closer to 0.5 for the highly convective stability and then systematically decrease as the near-neutral stability is approached (figure 12c and d). However, for the highly convective stability the maximum in the degree of isotropy associated with these events occurs at relatively smaller sizes than the peak positions of the heat-flux distribution (figure 12a and b). Therefore, the disagreement in the peak positions of heat flux and degree of isotropy might be related to the second reason associated with momentum transport.
Figure 13 shows the distributions of the heat and momentum fluxes associated with the warm-updraft and cold-downdraft events of different sizes,
$(T_{B}\overline {u})/z$ (see (2.27)). It is clear that the heat-flux peak position associated with warm-updraft events (
$(T_{B}\overline {u})/z \approx 2$) corresponds to a significant amount of down-gradient momentum in a highly convective surface layer (figure 13a and c). However, for the peak position
$(T_{B}\overline {u})/z \approx 2.5$ associated with the cold-downdraft events, the momentum transport is rather erratic in nature (figure 13b and d). The association of highly erratic momentum transport with the cold downdrafts has been observed in the numerical simulations of Li et al. (Reference Li, Gentine, Mellado and McColl2018) and in the observations of Chowdhuri & Prabha (Reference Chowdhuri and Prabha2019). Salesky & Anderson (Reference Salesky and Anderson2018) interpreted this as: under highly convective conditions, the small-scale turbulence is excited in the updraft regions and suppressed in the downdraft regions, leading to intermittent periods of small-scale excitation in the momentum fluxes.
The scaled heat- and momentum-flux distribution plotted against the normalized streamwise sizes
$(T_{B}\overline {u})/z$ of the heat-flux events corresponding to the warm-updraft (a and c) and cold-downdraft (b and d) quadrants.

Summarizing these observations, we note that there is a characteristic size of the warm-updraft and cold-downdraft events associated with least anisotropic turbulence, which does not scale with
$z$. On the other hand, the sizes of the warm-updraft and cold-downdraft events which carry the maximum heat are found to scale with
$z$. The mismatch in the peak positions of heat flux and degree of isotropy is related to the fact that the warm-updraft events which carry the maximum amount of heat are also associated with significant down-gradient momentum transport (figure 13a and c). However, for the cold-downdraft events these two peak positions almost coincide (figures 11b and 12b). This might be related to inefficient momentum transport associated with the cold-downdraft events, unlike the warm updrafts (figure 13c and d).
So far, we have focused on the degree of isotropy (
$C_{3 { c}}$) of the Reynolds stress tensor associated with the heat-flux events of different sizes, to quantify how close the turbulence is to 3-component isotropy. However, apart from 3-component isotropy, there are 1-component and 2-component anisotropic states whose dominance is described by the other two coefficients
$C_{1 { c}}$ and
$C_{2 { c}}$, respectively. Therefore, we can ascertain the entire anisotropic states of the heat-flux events of different sizes, by investigating the distributions of all three coefficients such as:
$C_{1 { c}}$,
$C_{2 { c}}$ and
$C_{3 { c}}$. We focus our attention on the down-gradient events, since for these events there is an intricate relation between the heat-flux intensity and the degree of isotropy.
3.3.3. The anisotropic states of the Reynolds stress tensor
Figure 14 shows the three coefficients associated with the three limiting states of the anisotropy Reynolds stress tensor to describe the anisotropic states of the warm-updraft and cold-downdraft events (see (2.25)). These three coefficients describe the corresponding weights associated with each of the three limiting states of the anisotropy Reynolds stress tensor, such as:
$C_{1 { c}}$ is related to 1-component anisotropy (red lines in figure14),
$C_{2 { c}}$ is related to 2-component anisotropy (green lines in figure 14) and
$C_{3 { c}}$ is related to 3-component isotropy (blue lines in figure 14). Note that from (2.25) the sum of the three anisotropy coefficients should be 1 for each
$(T_{B}\overline {u})/z$ value. Nevertheless, the graphs shown in figure 14 are ensemble-averaged over all the runs from a particular stability class. Owing to such ensemble averaging, the sums of the three coefficients may differ from 1 for some
$(T_{B}\overline {u})/z$ values. This occurs because some of the bins of
$(T_{B}\overline {u})/z$ values may remain empty for some specific runs which construct the ensemble.
The three coefficients (
$C_{1 { c}}$,
$C_{2 { c}}$ and
$C_{3 { c}}$) associated with the three limiting states of the anisotropy Reynolds stress tensor are plotted against the normalized streamwise sizes
$(T_{B}\overline {u})/z$ of the heat-flux events corresponding to warm-updraft and cold-downdraft quadrants: (
$a$)
$-\zeta >2$, (
$b$)
$1<-\zeta <2$, (
$c$)
$0.6<-\zeta <1$, (
$d$)
$0.4<-\zeta <0.6$, (
$e$)
$0.2<-\zeta <0.4$ and (
$f$)
$0<-\zeta <0.2$. The markers associated with these three coefficients are explained above panel (a).

For the highly convective stability (figure 14a), we note that for the warm-updraft events the maximum in
$C_{3 { c}}$ is located at
$(T_{B}\overline {u})/z \approx 0.5$. Moreover, we also observe that, for the warm-updraft events smaller than this size (
$(T_{B}\overline {u})/z < 0.5$), the values of the coefficient
$C_{2 { c}}$ exceed the other two coefficients. On the other hand, for the sizes of warm-updraft events larger than
$(T_{B}\overline {u})/z > 0.5$, the coefficient
$C_{1 { c}}$ is the largest amongst the three and its peak position coincides with the heat-flux peak position. This implies that the anisotropic states of the Reynolds stress tensor associated with the warm-updraft events smaller (larger) than the critical size
$(T_{B}\overline {u})/z \approx 0.5$ are dominated by 2-component (1-component) anisotropy.
However, for the cold-downdraft events from highly convective stability (figure 14a), the coefficients
$C_{1 { c}}$ and
$C_{2 {c}}$ dominate over
$C_{3 { c}}$ for all sizes, albeit
$C_{2 { c}}$ being the largest for sizes
$(T_{B}\overline {u})/z < 1$. For sizes
$(T_{B}\overline {u})/z > 1$, the coefficient
$C_{1 { c}}$ is the largest and its peak position almost coincides with the maximum heat flux associated with cold-downdraft events. Interestingly enough, we also find that, as the near-neutral stability is approached (figure 14a–f), the coefficient
$C_{2 { c}}$ systematically becomes the largest amongst the three for most of the sizes of warm-updraft and cold-downdraft events. Nonetheless, at larger sizes (in the order of the heat-flux peak positions) there remains a tendency for the coefficient
$C_{1 { c}}$ to dominate the anisotropic state of the Reynolds stress tensor.
The results from figure 14 are in accordance with the results from the quadrant analysis in figures 6 and 7. For the highly convective stability, we note that the blue regions (dominated by 3-component isotropy) in the anisotropy contour maps (figure 6a) broadly correspond to the critical sizes of warm updrafts and cold downdrafts (
$(T_{B}\overline {u})/z \approx 0.5$ and
$(T_{B}\overline {u})/z \approx 1$) where
$C_{3 { c}}$ values are maximum. For the sizes smaller (larger) than this, the anisotropic states of the Reynolds stress tensor are dominated by 2-component (1-component) anisotropy, shown as green (red) regions in figure 6(a). However, for the near-neutral stability (figures 6f and 14f) the dominance of 2-component anisotropy is associated with almost all sizes of the warm-updraft and cold-downdraft events, except at the larger sizes where there is a signature of 1-component anisotropy.
Hitherto, from analyses presented in figures 11–14 we have found that the least anisotropic turbulence is associated with particular sizes of warm-updraft and cold-downdraft events. These sizes do not exactly correspond to the peak positions of the heat-flux distribution and also do not scale with
$z$. With stability (highly convective to near-neutral) this particular size changes from
$(T_{B}\overline {u})/z \approx 0.5$ to
$(T_{B}\overline {u})/z \approx 2$ for the warm-updraft events and from
$(T_{B}\overline {u})/z \approx 1$ to
$(T_{B}\overline {u})/z \approx 3$ for the cold-downdraft events. From the persistence p.d.f.s and CDFs of these events presented in figures 8 and 9, we have noted that there is a power-law behaviour associated with sizes
$(T_{B}\overline {u})/z < 1$, followed by an exponential decay (Poisson-type process) for sizes
$(T_{B}\overline {u})/z > 1$. In the following section, we present results to investigate whether there is any correspondence between these p.d.f.s and the critical sizes of warm-updraft and cold-downdraft events associated with least anisotropic turbulence.
3.3.4. The linkage between degree of isotropy and persistence p.d.f.s
Before discussing anisotropy, to highlight non-Gaussianity we convert the p.d.f.s in figure 8 to a distribution about the time fractions (
$T_{f}$) spent in each quadrant of
$T^{\prime }\text {--}w^{\prime }$, by presenting the same in a premultiplied form. If, from a particular quadrant of
$T^{\prime }\text {--}w^{\prime }$, a number
$N_{tot}$ of blocks are being detected, with each
$N_{i}$th block containing
$n_{i}$ points, then we can write
\begin{equation} \sum_{i=1}^{N_{tot}}{N_{i}n_{i}} \propto T_{f}, \end{equation}
where
$T_{f}$ is the time fraction spent in that particular quadrant. Since the probability of finding a block containing
$n_{i}$ points is
$N_{i}/N_{tot}$, from (2.26) we can write
Therefore (3.11) can be expressed as
\begin{equation} \int_{({T_{B}\overline{u}}/{z})_{min}}^{({T_{B}\overline{u}}/{z})_{max}} \left(\frac{T_{B}\overline{u}}{z}\right) P\left[\left(\frac{T_{B}\overline{u}}{z}\right)\right] \textrm{d}\log\left(\frac{T_{B}\overline{u}}{z}\right) \propto T_{f}. \end{equation}From (3.12) we can also write
\begin{equation} \int_{({T_{B}\overline{u}}/{z})_{min}}^{({T_{B}\overline{u}}/{z})_{max}} \left\{ \left[\frac{T_{B}\overline{u}}{z} P\left(\frac{T_{B}\overline{u}}{z}\right)\right]_\textrm{III}- \left[\frac{T_{B}\overline{u}}{z} P\left(\frac{T_{B}\overline{u}}{z}\right)\right]_\textrm{I} \right\} \textrm{d}\log\left(\frac{T_{B}\overline{u}}{z}\right) \propto {\rm \Delta} T_{f}, \end{equation}
where the subscripts I and III refer to the warm-updraft and cold-downdraft quadrants and
${\rm \Delta} T_{f}$ is the difference in the time fractions spent in those quadrants. From (3.4) we know that
${\rm \Delta} T_{f} \approx \overline {{T^{\prime }}^3}/\sigma _{T}^3$, given the assumption that the time fractions spent in the counter-gradient quadrants could be neglected. Since from figure 8 we have noticed that the persistence p.d.f.s of the counter-gradient events decrease faster than the down-gradient events for the large sizes, we may rewrite (3.13) as
\begin{equation} \int_{({T_{B}\overline{u}}/{z})_{min}}^{({T_{B}\overline{u}}/{z})_{max}} \left\{ \left[\frac{T_{B}\overline{u}}{z} P\left(\frac{T_{B}\overline{u}}{z}\right)\right]_\textrm{III}- \left[\frac{T_{B}\overline{u}}{z} P\left(\frac{T_{B}\overline{u}}{z}\right)\right]_\textrm{I} \right\} \textrm{d}\log\left(\frac{T_{B}\overline{u}}{z}\right) \propto \frac{\overline{{T^{\prime}}^3}}{\sigma_{T}^3}. \end{equation} Figure 15(a)–(f) show the premultiplied p.d.f.s of
$(T_{B}\overline {u})/z$ corresponding to the warm updrafts and cold downdrafts for the same six different stability classes, along with the degree of isotropy. Upon close inspection, we note that these premultiplied p.d.f.s can be divided into two regions, which approximately intersect at
$(T_{B}\overline {u})/z \approx 1$. The first region extends up to
$(T_{B}\overline {u})/z \approx 1$, where the premultiplied p.d.f.s of the warm-updraft and cold-downdraft events collapse with a power law in highly convective stability (
$-\zeta > 2$). This power-law region progressively diminishes as the near-neutral stability is approached (figure 15a–f). The second region extends beyond
$(T_{B}\overline {u})/z \approx 1$, where these premultiplied p.d.f.s are widely separated in highly convective stability, while agreeing with each other in near-neutral stability (figure 15a–f). In the premultiplied form we can relate the difference in the values between the warm updrafts and cold downdrafts to the non-Gaussianity through (3.14). Therefore, we claim that the effect of non-Gaussianity (Gaussianity) in a highly (weakly) convective surface layer is only felt through those warm-updraft and cold-downdraft events having sizes
$(T_{B}\overline {u})/z>1$. This also explains why at sizes
$(T_{B}\overline {u})/z> 1$, the p.d.f.s and CDFs of the warm-updraft and cold-downdraft events differ most for the highly convective stability (figures 8a and 9a).
The log–log plots of the premultiplied p.d.f.s of
$(T_{B}\overline {u})/z$ (see (3.12)) corresponding to the heat-flux events from the warm-updraft and cold-downdraft quadrants are shown for the six different stability classes: (
$a$)
$-\zeta >2$, (
$b$)
$1<-\zeta <2$, (
$c$)
$0.6<-\zeta <1$, (
$d$)
$0.4<-\zeta <0.6$, (
$e$)
$0.2<-\zeta <0.4$ and (
$f$)
$0<-\zeta <0.2$. In all the panels, the right-hand
$y$ axis is linear and used to represent the distribution of the degree of isotropy (
$\langle C_{3 { c}}\mid [(T_{B}\overline {u})/z] \rangle$) associated with the warm-updraft and cold-downdraft events. The thick grey line shows the same power law as in figure 8, but owing to premultiplication the exponent changed to
$+0.6$. The grey shaded region represents
$(T_{B}\overline {u})/z< 1$, and the black line denotes the value of 1. The markers are explained in the legend in panel (f).

From figure 15 we can also compare the distribution of the degree of isotropy between the warm-updraft and cold-downdraft events. We find that, with the change in stability, the peak positions of the degree of isotropy shift systematically from the region
$(T_{B}\overline {u})/z<$ 1 to the region
$(T_{B}\overline {u})/z> 1$. From figure 9, we noted that,for sizes
$(T_{B}\overline {u})/z> 1$, the characteristics of the warm-updraft and cold-downdraft events might be related to the passing of the large-scale structures over the measurement points.
To get a preliminary insight into this systematic shift, we empirically investigated the distributions of the degree of isotropy associated with the warm-updraft and cold-downdraft events by normalizing their streamwise sizes with a mixed length scale. This mixed length scale is a geometric mean of two length scales such as the large-eddy length scale
$\lambda$ (see (2.22)) and
$z$, represented as
$z^{0.5}\lambda ^{0.5}$. This has been discovered in the context of event-based analysis, where Rao, Narasimha & Badri Narayanan (Reference Rao, Narasimha and Badri Narayanan1971) showed that the frequency of the burst events in a turbulent boundary layer scaled with a mixed time scale, involving both inner and outer variables. Similarly, Alfredsson & Johansson (Reference Alfredsson and Johansson1984) found that the governing time scale of the near-wall region of a channel flow was a mixture of outer and inner scales. They interpreted this as a sign of the interaction of outer and near-wall flows. This mixed scale has been reviewed in detail by Buschmann, Indinger & Gad-el Hak (Reference Buschmann, Indinger and Gad-el Hak2009) and Gad-el Hak & Buschmann (Reference Gad-el Hak and Buschmann2011). Recently, McKeon (Reference McKeon2017) noted that this mixed length scale can be derived from first principles through matched asymptotic expansions, a theory proposed by Afzal (Reference Afzal1984). Afzal (Reference Afzal1984) showed that, by matching the inner and outer expansions of the Reynolds shear stress, an intermediate layer could be formulated for wall-bounded turbulent flows where the appropriate length scale was the geometric mean of the inner and outer length scales.
Figure 16 shows that, by normalizing the streamwise sizes of the warm-updraft and cold-downdraft events by the mixed length scale, one could reasonably collapse the peak positions of the degree of isotropy at
$(T_{B}\overline {u})/(z^{0.5}\lambda ^{0.5}) \approx 0.08$ and
$(T_{B}\overline {u})/(z^{0.5}\lambda ^{0.5}) \approx 0.15$, respectively. From figures 8 and 15, we have found that the warm-updraft and cold-downdraft events having sizes
$(T_{B}\overline {u})/z< 1$ are scale-invariant owing to a power-law dependence in the highly convective stability. This scale-invariant property disappears systematically as the near-neutral stability is approached. Apart from that, the effect of non-Gaussianity (Gaussianity) appears mostly at the sizes
$(T_{B}\overline {u})/z> 1$ in a highly (weakly) convective surface layer. Therefore, this mixed length scaling to collapse the peak positions of the degree of isotropy may suggest that the least anisotropic turbulence might be associated with an interaction between two different physical processes. One of these processes might be related to scale invariance, while the other with non-Gaussianity, associated with the warm-updraft and cold-downdraft events. This is at present a conjecture, which needs to be verified from theoretical arguments. Recently Tong & Ding (Reference Tong and Ding2020) have proposed a matched asymptotic expansion for the convective surface layer, to derive the scaling of the mean velocity profile. By following their footsteps, along with the line of reasoning developed by Afzal (Reference Afzal1984), it might be possible to derive this mixed length scale from first principles for convective surface layer turbulence. However, this is beyond the scope of the present article. We present our conclusions in the next section.
The distributions of the degree of isotropy (
$C_{3 { c}}$) for the
$z$- and mixed length-scaled sizes of the heat-flux events from warm-updraft and cold-downdraft quadrants (
$(T_{B}\overline {u})/z$ and
$(T_{B}\overline {u})/(z^{0.5}\lambda ^{0.5})$) are shown in the top and bottom panels, respectively. The scale
$\lambda$ is the large-eddy length scale obtained from (2.22). In the top panels, the grey shaded region represents
$(T_{B}\overline {u})/z< 1$ and the black line denotes the value of 1. In the bottom panels, the black arrows indicate the peak positions of
$C_{3 { c}}$ corresponding to the mixed length-scaled sizes of the heat-flux events from warm-updraft and cold-downdraft quadrants.

Same as in figure 8, but the histograms are shown for the heat-flux events from each quadrant.
$(\textit{a}) \ - \zeta > 2, \ (\textit{b}) \ 1 < -\zeta < 2, \ (\textit{c}) \ 0.6 < -\zeta < 1, \ (\textit{d}) \ 0.4 < -\zeta < 0.6, \ (\textit{e}) \ 0.2 < $
$-\zeta < 0.4 \hbox{ and } (\,\textit{f})\ 0 < -\zeta < 0.2.$

4. Conclusions
We report novel comprehensive results of Reynolds stress anisotropy associated with intermittent heat transport in an unstable ASL, from the SLTEST experimental dataset. We adopt an event-based description of the heat-transporting events occurring intermittently and persisting over a wide range of time scales. The Reynolds stress anisotropy is quantified by using a metric called degree of isotropy, computed from the smallest eigenvalue of the anisotropy Reynolds stress tensor. The important results from this study can be broadly summarized as follows:
(i) The anisotropic state of the Reynolds stress tensor evolves from being dominated by 2-component anisotropy to being dominated by 3-component isotropy as the stability changes from weakly to highly convective. The degree of isotropy of the Reynolds stress tensor is governed by the strength of the vertical velocity fluctuations, which preferentially couple with the temperature fluctuations. These temperature fluctuations exhibit strong (weak) non-Gaussian characteristics in a highly (weakly) convective surface layer.
(ii) The Reynolds stress anisotropy in an unstable surface layer is strongly related to the asymmetric and intermittent nature of heat transport, associated with non-Gaussianity in the temperature fluctuations.
(iii) By adopting an event-based approach, it is found that not all the heat-flux events are associated with the same anisotropic state of turbulence. The anisotropic states associated with highly intermittent large-heat-flux events are dominated by 1-component anisotropy; whereas, the anisotropic states associated with more frequent but weak heat-flux events are dominated by 2-component anisotropy. On the other hand, the anisotropic states associated with moderate-heat-flux events which lie between these two extremes are dominated by 3-component isotropy.
(iv) There is a critical size associated with the organized heat-flux events (warm updrafts and cold downdrafts) which corresponds to the maximum value of the degree of isotropy (i.e. least anisotropic turbulence). By investigating the anisotropic states of the Reynolds stress tensor, it is found that, in a highly convective surface layer, the warm-updraft and cold-downdraft events smaller (larger) than this critical size are associated with anisotropic states dominated by 2-component (1-component) anisotropy. However, in a near-neutral surface layer, the anisotropic states are mostly dominated by 2-component anisotropy, regardless of the sizes of the warm-updraft and cold-downdraft events.
(v) This critical size associated with least anisotropic turbulence does not scale with
$z$. However, the
$z$ scaling is successful in collapsing the peak positions of the heat-flux distribution associated with the sizes of the warm-updraft and cold-downdraft events. This disagreement occurs because the sizes of the warm-updraft events corresponding to maximum heat flux are also associated with a significant amount of streamwise momentum. This causes a drop in the degree of isotropy associated with their sizes.
Note that the findings from this study should be verified from field experiments in an unstable ASL flow conducted over rough surfaces and in complex terrains. Our preliminary investigation shows that this critical size probably scales with a mixed length scale
$z^{0.5}\lambda ^{0.5}$, where
$\lambda$ is the large-eddy length scale. We propose a conjecture that this mixed length scaling may reflect an interaction between two different physical processes, one of which may be associated with scale invariance and the other with the non-Gaussianity in turbulence. The verification of this conjecture is beyond the scope of the present study. An inevitable limitation of this study is the unavailability of three-dimensional velocity and temperature information. Owing to this constraint, the intermittent heat-flux events and the associated Reynolds stress anisotropy cannot be connected to the three-dimensional topology of the coherent structures in convective turbulence. In the future, we would address this problem through large-eddy or direct numerical simulations. This study also raises a few important questions, which deserve future attention:
(i) Is there a theoretical framework to explain the mixed length scale in convective turbulence?
(ii) What is the physical connection between the event-based (related to flow structures) and scale-based (related to harmonic analysis) description of turbulence anisotropy?
Acknowledgments
The Indian Institute of Tropical Meteorology (IITM) is an autonomous institute fully funded by the Ministry of Earth Sciences, Government of India. S.C. and S.K. gratefully acknowledge the Director, IITM, Dr Thara Prabhakaran, Dr P. Mukhopadhyay, Dr Shivsai Ajit Dixit and Mr Pramit Kumar Deb Burman for the constant encouragement and stimulating discussions during the research. T.B. acknowledges the new faculty start-up funding from the Department of Civil and Environmental Engineering, the Henry Samueli School of Engineering, University of California, Irvine. The authors would like to sincerely thank the three anonymous reviewers and the editor whose insightful comments helped to improve the quality of the article. The authors are also immensely grateful to DrKeith G. McNaughton for kindly letting them use the SLTEST dataset for this research, pointing out some very important references on mixed length scaling, and providing deep insightful feedback on the article. The computer codes used in this study are available upon request from the corresponding author.
Declaration of interests
The authors report no conflict of interest.
Supplementary figures
Supplementary figures are available at https://doi.org/10.1017/jfm.2020.471.
Appendix. Histograms of the heat-flux events
In figure 17(a)–(f), we show the histograms of the heat-flux events from each quadrant of
$T^{\prime }\text {--}w^{\prime }$, corresponding to the six different stability classes. The number of events (
$n$) shown in figure 17(a)–(f) is computed after considering all 30 min runs from a particular stability class (e.g. 55 number of 30 min runs for
$-\zeta > 2$, amounting to 27.5 hours of observation). For each stability class, the total number of heat-flux events counted over all sizes
$(T_{B}\overline {u})/z$ from each quadrant are given in table 4. Typically, for the warm-updraft and cold-downdraft quadrants, we encounter more than 100–200 of heat-flux events corresponding to the sizes
$(T_{B}\overline {u})/z>4$. For the counter-gradient quadrants, the total number of heat-flux events corresponding to large sizes
$(T_{B}\overline {u})/z>1$ is also more than 100, although the histograms decrease faster than the down-gradient quadrants. This implies that these counter-gradient events have a statistical tendency to occur in smaller sizes and do not persist for a long time. Therefore, the mean statistics shown in figures 8–16 for the heat-flux events from all four quadrants have been averaged over more than 100–200 events for the streamwise sizes
$(T_{B}\overline {u})/z>1$. We note that many statistics textbooks (e.g. Ross Reference Ross2014) as well as the seminal paper by Student (Reference Student1908) consider that a sample size of more than 30 is enough for ensuring the statistical convergence of the mean to the actual population mean, from the weak law of large numbers. Nevertheless, we also performed the Student's
$t$-test to ensure the statistical significance of the ensemble mean. For the large values of
$(T_{B}\overline {u})/z$, based on the sample size of around 100–200 events, the margin of error in the ensemble mean computed over these samples is approximately 7 %–10 % with a confidence level of 95 %.
The total number of heat-flux events from each quadrant of
$T^{\prime }\text {--}w^{\prime }$ are tabulated for all sizes
$(T_{B}\overline {u})/z$, corresponding to each stability class as shown in table 2. The number in parentheses denotes the number of 30 min runs in each stability category.









































































































































































































