2.1. LES of the ABL

The ABL LES were designed to systematically control the global stability parameter $-z_i/L$ from neutral ($-z_i/L = 0$) to moderately convective ($-z_i/L \sim \text {--}10$) while maintaining quasi-stationary periods of fully developed ABL turbulence over which to collect statistics. To this end, we systematically increase surface heat flux while keeping the imposed geostrophic wind vector fixed. Quasi-stationarity is maintained through the capping inversion design. The base computational domain was designed to capture 5–10 integral scales in the horizontal, which is common in the LES literature (Schmidt & Schumann Reference Schmidt and Schumann1989; Moeng & Sullivan Reference Moeng and Sullivan1994; Khanna & Brasseur Reference Khanna and Brasseur1998; Sullivan & Patton Reference Sullivan and Patton2011) when modelling neutral and unstable ABLs. As we will see in subsequent sections, some of the ABL stability states considered in this work show surprisingly large coherence lengths. Therefore, we assess the influence of extended domains (up to $\sim 22 z_i$) on the statistics underlying the results and conclusions from this study. This additional analysis is described in appendix A.3 and discussed further in § 7. To minimize the influence of the upper boundary, the vertical domain is maintained at least twice the capping inversion depth.

The filtered momentum, energy (potential temperature) and continuity equations upon which the LES of the ABL is developed are

(2.1)\begin{gather} \frac{ \partial \boldsymbol{u}^{r}}{ \partial t}+ \boldsymbol{\nabla} \boldsymbol{\cdot} (\boldsymbol{u}^{r}\boldsymbol{u}^{r})^{r} ={-}\frac{1}{\rho_0} \boldsymbol{\nabla} p^{*}- \boldsymbol{\nabla} \boldsymbol{\cdot} {\boldsymbol{\tau}}_a^{sfs}+\frac{\boldsymbol{g}}{\theta_0} (\theta^{r}-\theta_0)+\boldsymbol{f}\times (\boldsymbol{V}_g - \boldsymbol{u}^{r}), \end{gather}

(2.2)\begin{gather} \frac{ \partial \theta^{r}}{ \partial t}+ \boldsymbol{\nabla} \boldsymbol{\cdot} (\boldsymbol{u}^{r} \theta^{r})^{r} ={-} \boldsymbol{\nabla} \boldsymbol{\cdot} {\boldsymbol{\tau}}_{\theta}^{sfs}, \end{gather}

(2.3)\begin{gather} \boldsymbol{\nabla} \boldsymbol{\cdot}\boldsymbol{u}^{r} = 0, \end{gather}

where LES resolved-scale variables are indicated with superscript $r$, $\boldsymbol {u}^{r}$ is the resolved velocity vector, $p^{*}$ is the effective resolved fluctuating pressure (see below) and $\theta ^{r}$ is the resolved potential temperature. The sub-filter-scale (SFS) ‘stress’ tensor ${\boldsymbol {\tau }}^{sfs}=( \boldsymbol {u}^{r}\boldsymbol {u}^{s}+ \boldsymbol {u}^{s}\boldsymbol {u}^{r}+\boldsymbol {u}^{s}\boldsymbol {u}^{s})^{r}$ and temperature flux vector ${\boldsymbol \tau }^{sfs}_{\theta }=( \theta ^{r}\boldsymbol {u}^{s}+ \theta ^{s}\boldsymbol {u}^{r}+\theta ^{s}\boldsymbol {u}^{s})^{r}$ must be modelled, where superscript $s$ indicates SFS variables. Typically, ${\boldsymbol {\tau }}^{sfs}$ is separated into its isotropic and deviatoric parts, with its isotropic part included in $p^{*}$ (see below) and its deviatoric part $\boldsymbol {\tau }^{sfs}_a$ modelled using the one-equation eddy viscosity closure described in Moeng (Reference Moeng1984) (see also Deardorff Reference Deardorff1980; Sullivan, McWilliams & Moeng Reference Sullivan, McWilliams and Moeng1996; Khanna & Brasseur Reference Khanna and Brasseur1997):

(2.4a)\begin{equation} \boldsymbol{\tau}_a^{sfs}={-}2\nu_T\boldsymbol{\mathsf{s}}^{r}, \quad \nu_T=C_kl\sqrt{e}, \quad \varDelta=(\varDelta_x,\varDelta_y,\varDelta_z)^{1/3}, \quad C_k=0.08, \end{equation}

where

(2.4b)\begin{equation} l = \left\{ \begin{array}{@{}ll} \varDelta, & \partial\theta^{r}/ \partial z \leq 0 \text{ (unstable)}, \\ l_s=0.76e^{1/2}\left(\dfrac{g}{\theta_0}\dfrac{ \partial\theta^{r}}{ \partial z}\right)^{{-}1/2}, & \partial\theta^{r}/ \partial z > 0 \text{ and } l_s<\varDelta \text{ (stable)}. \end{array} \right. \end{equation}

Here $\boldsymbol{\mathsf{s}}^{r}$ is the resolved strain-rate tensor, $\varDelta _{x,y,z}$ are the cell dimensions for uniform hexahedral grid cells, $l$ and $\varDelta$ represent the SFS and grid length scales, respectively, and $e$ is the ‘SFS turbulent kinetic energy’, solved using a prognostic equation as in Moeng (Reference Moeng1984). In all simulations $C_k=0.08$. The SFS temperature flux vector $\boldsymbol {\tau }_{\theta }^{sfs}$ is also modelled with an eddy diffusivity closure (Deardorff Reference Deardorff1980; Moeng Reference Moeng1984):

(2.5)\begin{equation} {\boldsymbol \tau}_{\theta}^{sfs}=\nu_{\theta} \boldsymbol{\nabla}\theta^{r}, \quad \nu_{\theta}=\nu_T\left[1 + \frac{2l}{\varDelta}\right]. \end{equation}

The momentum equation (2.1) adopts the Boussinesq approximation for gravitational force with the background state variables $\rho _0$ and $\theta _0$ taken as the initial surface values. The local resolved pressure $p^{r}=p^{r\,\prime }+P$ is separated into ensemble mean ($P$) and fluctuating parts, with the fluctuating part included in effective pressure, $p^{*}=p^{r\,\prime }-({\rho _0}/{3})\, \textrm {tr} \{ \boldsymbol {\tau }^{sfs} \}$. Mean pressure gradient $\boldsymbol {\nabla } P$ is written in terms of a specified ‘geostrophic wind’ vector $\boldsymbol {V}_g$, defined through the Ekman balance between pressure and Coriolis force:

(2.6)\begin{equation} \boldsymbol{\nabla}_h P\equiv \rho_0 \,\boldsymbol{f}\times \boldsymbol{V}_g, \end{equation}

where $\boldsymbol {f}=f\hat {\boldsymbol {e}}_z$ is twice the projected Earth rotation-rate vector. In the canonical ABL, $\boldsymbol {V}_g$ and ${ \boldsymbol {\nabla }} P$ are horizontal and $\boldsymbol {V}_g$ describes the mesoscale wind, which is set to a nominal value of $10\ \textrm {m}\,\textrm {s}^{-1}$ for this study. In our simulations the Coriolis parameter $f\equiv |\,\boldsymbol {f}|$ is set to $10^{-4}\ \textrm {s}^{-1}$ as representative of the mid-latitudes. The specification of $\boldsymbol {V}_g$ in concert with surface heat flux $Q_0$ establishes the stability state $-z_i/L$ of the boundary layer after reaching quasi-equilibrium.

2.1.1. The quasi-equilibrium state

The ABL is capped by a stabilizing inversion layer in potential temperature, established in part through specification of the initial condition in the potential temperature field. The initial potential temperature profile is constant up to the initial capping inversion height $z_{i_o}$ and increases from $z_{i_o}$ to the top of the computational domain with constant slope. The rate at which the capping inversion depth grows with time ($\textrm {d}z_i/\textrm {d}t$) increases with surface heat flux $Q_0$, leading to concerns with potential influence on stationarity, and therefore equilibrium. In addition, if $z_i$ were to increase to ${\sim }2/3$ the domain height, it would become necessary to increase the extent of the vertical domain and modify the grid, which might impact subtle transitions in ABL structure. To minimize these concerns, the initial capping inversion strength, quantified by the vertical gradient in mean potential temperature within the capping inversion layer, $[\textrm {d}\theta ^{r}/\textrm {d}z]_{CI}$, was carefully chosen to minimize the rate of growth in $z_i$ at the highest $Q_0$ (and $-z_i/L$) while remaining consistent with observations of capping inversion strengths in the daytime atmosphere. In the current study, the magnitude of $[ \textrm {d}\theta ^{r}/\textrm {d}z ]_{CI}$ averaged over the capping inversion layer was ${\simeq }87\ \textrm {K}\,\textrm {km}^{-1}$ with a standard deviation of $37\ \textrm {K}\,\textrm {km}^{-1}$ over 14 simulations (table 1). For comparison, the LES study of Sullivan & Patton (Reference Sullivan and Patton2011) adopted a capping inversion temperature gradient of $\simeq 80\ \textrm {K}\,\textrm {km}^{-1}$. To estimate realistic ranges of capping inversion strengths, we analysed data from the open literature (Sorbjan Reference Sorbjan1996; Sullivan et al. Reference Sullivan, Moeng, Stevens, Lenschow and Mayor1998; Hennemuth & Lammert Reference Hennemuth and Lammert2006) as well as from a database collected from aircraft flights in Indiana (Davis Reference Davis2017). Using these, we estimated the capping inversion strengths in the daytime ABL to lie between $9.8\ \textrm {K}\,\textrm {km}^{-1}$ and $120\ \textrm {K}\,\textrm {km}^{-1}$ with a mean of $56\ \textrm {K}\,\textrm {km}^{-1}$ and standard deviation of $40\ \textrm {K}\,\textrm {km}^{-1}$, indicating that the capping inversion strengths in our LES experiments are within observational ranges. Additionally, the depths of capping inversions in our LES were in the range of those analysed from available field data (see Grabon et al. Reference Grabon, Davis, Kiemle and Ehret2010).

Table 1. Simulation parameters over the quasi-stationary analysis periods.

In figure 1 we plot the ratio of capping-inversion-growth velocity to buoyancy-driven mixed layer velocity scale $w_*$ for the capping inversion strength $[ \textrm {d}\theta ^{r}/\textrm {d}z ]_{CI}$ adopted in our computations at our highest $Q_0$ ($-z_i/L=3.84$), demonstrating that the worst case $\textrm {d}z_i/\textrm {d}t$ is well over two orders of magnitude smaller than the dominant characteristic velocity scale (consistent with Sullivan & Patton (Reference Sullivan and Patton2011)), thus maintaining quasi-stationarity. The fully developed quasi-equilibrium period for analysis was chosen from the simulation window over which the friction velocity ($u_*$), buoyancy velocity scale ($w_*$) as well as the net mass flow rate through the domain become nearly constant with time. This quasi-equilibrium is realized after ${\sim }200$ eddy turnover times, $\tau _u\equiv z_i/u_*$. The time periods over which statistics were calculated in this analysis varied between $60\tau _u$ and $80\tau _u$ for the different $-z_i/L$.

Figure 1. Time evolution of the ratio of capping-inversion-growth velocity to the buoyancy-driven mixed layer velocity scale at the highest $-z_i/L\simeq 3.84$ for the chosen initial capping inversion strength.

The dynamical system was discretized within a rectangular computational domain $5\ \text {km} \times 5\ \text {km} \times 2\ \text {km}$ using a uniformly spaced geometric grid of dimension $324 \times 324 \times 96$. The pseudo-spectral algorithm was applied in the horizontal with periodic boundary conditions; second-order finite differencing was applied in the vertical on a staggered grid. Dealiasing by truncation in the horizontal using the 2/3 rule creates an effective grid of $216 \times 216 \times 96$ points with cell aspect ratio of ${\simeq }1$ on which the resolved scales evolve. Note that effective grid resolution (as opposed to geometric grid resolution) is quoted in all figures. A series of $14$ simulations were carried out to equilibrium with systematically varying surface heat flux for fixed geostrophic wind velocity magnitude of $10\ \textrm {m}\,\textrm {s}^{-1}$ and direction (aligned with $x$). Over the periods of data collection, the vertical extent of the domain is $2z_i$–$3z_i$, so that the boundary layer was resolved in the vertical with 31–51 points and the surface layer with 8–12 points. Time advance was performed using third-order Runge–Kutta integration.

For each simulation, the initial potential temperature distribution was constant in $z$ up to the initial capping inversion layer of thickness ${\simeq }200\ \text {m}$ centred at $z\simeq 640\ \text {m}$ with temperature difference of 30 K ($[\textrm {d}\theta ^{r}/\textrm {d}z]_{CI}=0.144\ \textrm {K}\,\textrm {m}^{-1}$). Above this layer, the initial vertical temperature gradient was set to $0.003\ \textrm {K}\,\textrm {m}^{-1}$. On the upper boundary, the mean temperature gradient is fixed at the initial value and the vertical gradients of the horizontal velocity components, the SFS viscosity and the vertical velocity were all fixed at zero. On the lower surface, the temperature flux $Q_0=\langle w\theta \rangle _0$ was prescribed at the $14$ values given in table 1. A lower boundary condition on total surface fluctuating shear stress vector was applied similarly to the nonlinear model described in Moeng (Reference Moeng1984), which was also applied in the studies by Khanna & Brasseur (Reference Khanna and Brasseur1997, Reference Khanna and Brasseur1998), Brasseur & Wei (Reference Brasseur and Wei2010) and Sullivan & Patton (Reference Sullivan and Patton2011). The surface stress model introduces a roughness parameter $z_0$ in a modelled relationship for the ratio of mean horizontal velocity magnitude at the first grid level $U_1$ to the friction velocity $u_*=|\langle \boldsymbol {u}_hw \rangle |^{1/2}_0$. This modelled relationship for $U_1/u_*$ is based on a continuous function of $z/z_0$ given by the empirical form of Monin–Obukhov law-of-the-wall (LOTW) scaling developed by Paulson (Reference Paulson1970) for the surface layer. The roughness parameter $z_0$ is set to $0.16\ \text {m}$ in our study.

2.1.2. LES grid design and overshoot

Grid design for LES is complicated by the interaction between resolution, grid cell aspect ratio and SFS parametrization. Specifically, an issue inherent to LES of the rough-surface ABL is the problem of the near-surface ‘overshoot’, an overprediction of mean velocity gradient at the first few LES grid cells adjacent to the surface. The overshoot is characterized as a peak in the vertical mean velocity gradient after normalization by inertial LOTW scaling in $z$ (see Piomelli (Reference Piomelli2008) and Brasseur & Wei (Reference Brasseur and Wei2010, BW10), and references therein, for detailed discussions of this issue). This long-standing problem with LES prediction has been diagnosed by BW10 as arising from the inconsistency between the inertia dominance requirement in the LOTW scaling argument for the inertial surface layer and non-physical frictional content within the SFS stress model. The strength of the spurious frictional content depends on the design of the grid (specifically, vertical resolution and grid cell aspect ratio). The base effective grid (i.e. after dealiasing in the horizontal) of dimension $216 \times 216 \times 96$ used in this study was designed to eliminate this spurious overshoot after following through a systematic procedure involving numerous LES of different grids (nearly $80$ simulations across the various $-z_i/L$). In appendix A.1 we describe the application of the BW10 framework to systematically remove the overshoot from all simulations. This analysis is extended in appendix A.2 to the design and application of grids to study impacts of varying resolution.

In total, 14 simulations (figure 2 and table 1) were carried out to quasi-equilibrium and analysed with instability state parameter $-z_i/L$ from $0$ to $3.84$ corresponding to $Q_0$ from $0$ to $0.050\ \textrm {K}\,\textrm {m}\,\textrm {s}^{-1}$ with fixed geostrophic wind speed of $10\ \textrm {m}\,\textrm {s}^{-1}$. The stability parameter $-z_i/L$ varies approximately linearly with $Q_0$ (figure 2a), largely due to surface friction velocity $u_*$ and capping inversion height $z_i$ varying relatively little with $-z_i/L$ (figure 2b). In contrast, the Obukhov length $L$ and buoyancy velocity scale $w_*\equiv ( gQ_0z_i/\theta _0 )^{1/3}$ vary approximately as $1/Q_0$ and $Q_0^{1/3}$, respectively.

Figure 2. (a) Variation of instability parameter with surface heat flux. (b) Variation of simulation parameters with instability state, $-z_i/L$.

2.2. Analytical methods applied to the LES data

To study transitions in the structure of the most energy-dominant turbulence motions within the ABL, we quantify two-point correlations and coherence lengths. These converged statistics are estimated by leveraging horizontal homogeneity and quasi-stationarity to average over horizontal planes and the ${\sim }70\tau _u$ analysis window over which the capping inversion depth $z_i$ grew by ${\sim }10\,\%$ for the most convective instability state ($-z_i/L = 3.84$), less for less unstable states. Some second- and third-order moments are also quantified to characterize variances, covariances and terms in the budget equations for turbulent kinetic energy. Unless otherwise specified, all statistics are calculated with $x$ and $y$ rotated to align with the mean velocity vector in that horizontal plane such that the rotated coordinate system $(x^{r}_1, x^{r}_2, x_3) \equiv (x_{rot}(z), y_{rot}(z), z)$ varies with $z$.

The two-point velocity–velocity correlation on horizontal planes in rotated coordinates is given by

(2.7)\begin{equation} {R}_{ij,k}(r_k,z)=\langle u^{r}_i(\boldsymbol{x}^{r}_h,z)u^{r}_j(\boldsymbol{x}^{r}_h+r_k\boldsymbol{e}^{r}_k,z) \rangle, \end{equation}

where $u^{r}_i$ are the fluctuating velocity components in the rotated coordinates and there is no sum on $k$. Homogeneity in the horizontal implies that the correlation is independent of starting point $\boldsymbol {x}^{r}_h$ on the horizontal plane and that the ensemble averages can be replaced by averages over horizontal planes at fixed z. In (2.7) $r_k$ is restricted to the horizontal plane $(k = 1,2)$, all in the rotated coordinates. Thus, $R_{ij,k}$ depends only on $z$, in addition to $r_k$. The two-point velocity–velocity correlation in the $z$ direction is given by

(2.8)\begin{equation} {R}_{ij,3}(r_3,z_0)=\langle u_i(\boldsymbol{x}_h,z_0)u_j(\boldsymbol{x}_h,z_0+r_3) \rangle, \end{equation}

where $u_i$ and $\boldsymbol {x}_h \equiv (x_1, x_2)$ are the fluctuating velocity components and horizontal location in unrotated coordinates. In the inhomogeneous $z$ direction, the two-point correlation depends on the starting point, $z_0$ (chosen as $0.1z_i$) and $r_3=z-z_0$ is the distance in the vertical direction. In our analysis, the angle brackets in (2.7) and (2.8) imply averaging both in space over the horizontal $z$ planes and in time over at least $60$ integral time scales $\tau _u$ during the quasi-equilibrium period.

The corresponding integral scales are given by

(2.9)\begin{equation} {L}_{ij,k}(z) =\int_0^{\infty} \frac{R_{ij,k}(r_k,z)}{R_{ij,k}(0,z)}\,\textrm{d}r_k \end{equation}

on a horizontal plane and by

(2.10)\begin{equation} {L}_{ij,3}(z_0) =\int_0^{\infty} \frac{R_{ij,3}(r_3,z_0)}{R_{ij,3}(0,z_0)}\,\textrm{d}r_3 \end{equation}

in the vertical direction. The integral scale in the vertical is estimated using non-rotated coordinates. We have confirmed that the vertical integral scales are mostly insensitive to the rotation of the horizontal coordinates. In practice, the integral in (2.9) is calculated using Fourier decomposition in the horizontal, which effectively ends the integration of the two-point correlation at half the horizontal domain width. To calculate (2.10), the integral extends to the capping inversion, $z_i$. Further discussion on the calculation of integrals is given in § 7 and appendix A.3.

3.1. Trends in coherence lengths

The production of streamwise velocity fluctuations originates statistically from the interaction between mean shear and Reynolds shear stress. Furthermore, stability theory (Lilly Reference Lilly1966; Brown Reference Brown1980) and studies of homogeneous turbulent shear flow (Lee et al. Reference Lee, Kim and Moin1990; Brasseur & Wang Reference Brasseur and Wang1992; Brasseur & Lin Reference Brasseur and Lin2005) show that streamwise coherence is created by the dynamical interactions between mean shear rate $S$ and the underlying turbulence fluctuations, an effect that strengthens with normalized shear rate, $S\tau _u$. Therefore, one anticipates that the streamwise coherence length of streamwise velocity fluctuations ($L_{11,1}$) in the surface layer will be highest in the purely shear-dominated neutral state ($-z_i/L = 0$) and that any addition of surface heating (i.e. $-z_i/L > 0$) will interfere with this mechanism and reduce streamwise coherence. If this is the case, then the streamwise coherence length of streamwise velocity fluctuations, $L_{11,1}$, should decrease monotonically with increasing $-z_i/L$ from the neutral state.

Surprisingly, this turns out not to be the case, as shown in figure 3(a), where $L_{11,1}$ is plotted against $-z_i/L$ at different distances from the surface. Our simulations predict that streamwise coherence lengths of streamwise fluctuations $u^{\prime }$ actually increase with increasing $-z_i/L$ as surface heat flux $Q_0$ is added to a previously unheated neutral canonical ABL. This response is strong and occurs with small increases in $Q_0$ and instability state, and is observed through the entire boundary layer from near the surface ($z/z_i=0.1$) to near the capping inversion ($z/z_i=0.7$). The anticipated reduction in streamwise coherence length with increase in $-z_i/L$ does occur, but not until $L_{11,1}$ has first increased by over a factor of 2 everywhere. The doubling of horizontal coherence length $L_{11,1}$ initiates at $-z_i/L \simeq 0.30$ and peaks at $-z_i/L = 0.433$ after the addition of very small levels of surface heating (from figure 2(a) and table 1, $Q_0 = 0.0055\ \textrm {K}\,\textrm {m}\,\textrm {s}^{-1}$ at $-z_i/L = 0.433$). The sudden increase in $L_{11,1}$ from $-z_i/L \simeq 0.30$ and $\simeq 0.43$ is followed by the expected monotonic decrease.

Figure 3. Variation of streamwise coherence lengths of (a) streamwise velocity fluctuations ($L_{11,1}$) and (b) vertical velocity fluctuations ($L_{33,1}$) with the global instability parameter $-z_i/L$. Each curve corresponds to fixed height within the boundary layer as indicated. The three vertical lines in each panel correspond to the following important transition states: $-z_i/L=0.433$ (dashed), $-z_i/L=1.080$ (dot-dashed) and $-z_i/L=2.86$ (dotted).

A similar trend is also observed in the streamwise coherence length of vertical velocity fluctuations ${w^{\prime }}$ ($L_{33,1}$), but with significant differences (figure 3b). In particular, a sudden increase in coherence length occurs at $-z_i/L=0.433$ for $L_{33,1}$ as it does in $L_{11,1}$, albeit with less variation when $-z_i/L\leq 0.304$. Also, the variation is much less pronounced at $z_i/z=0.1$, where vertical velocity fluctuations are damped at all $-z_i/L$. More importantly, however, whereas $L_{11,1}$ decreases monotonically when $-z_i/L$ increases above the critical value of $0.433$, $L_{33,1}$ continues to increase to a peak at $-z_i/L=1.08$, after which it decreases monotonically. It is interesting that the normalized horizontal coherence lengths at the highest values of $-z_i/L$ simulated are all roughly in the range $0.6z_{i}\text {--}0.7z_i$ for both horizontal and vertical turbulence fluctuations at all levels in the ABL, with the exception of $L_{33,1}$ at the highly damped location $z/z_i = 0.1$. Also note that, whereas both $L_{11,1}$ and $L_{33,1}$ individually increase to a peak coherence length, the total relative increase in streamwise coherence length is, overall, much larger for vertical fluctuating velocity due to its smaller values in the neutral regime. These observations will be discussed in more detail in § 4, where we analyse the existence of a ‘critical’ instability at $-z_i/L\simeq 0.43$. This intriguing and non-monotonic transition in the coherence lengths of the energy-containing turbulent velocity fluctuations ($L_{11,1}$ and $L_{33,1}$) also manifests in the organization and overall structure of the energy-dominant turbulent eddies as $-z_i/L$ increases, thereby suggesting a complex interplay between energy-containing buoyancy and shear-driven turbulent motions.

In appendices A.2 and A.3 we describe an analysis of the impacts of grid resolution and horizontal domain resolution, respectively, on the results presented in figure 3. Overall, we find that, although resolution and domain size impacts are not negligible, the key features identified here are accentuated with increase in resolution or horizontal domain size (on appropriately designed grids). These issues are further discussed in § 7.

3.2. Trends in single-point statistics

Figure 4 shows the changes in select terms in the Reynolds stress tensor and their budget equations with systematic variation in stability state. In general, the variances and covariances show less dramatic variations with $-z_i/L$ (figure 4) as compared to the relative changes in streamwise coherence length with $-z_i/L$ (figure 3). However, all the different single-point statistical metrics show perceptible changes occurring at $-z_i/L=0.433$ (highlighted by the dashed vertical line), indicating a critical transition at this instability state. For example, the streamwise fluctuating velocity variance (figure 4a) increases to a peak at $-z_i/L=0.433$ over the entire boundary layer before gradually decreasing to smaller values at higher $-z_i/L$, not dissimilar to the variation of $L_{11,1}$ with $-z_i/L$ (figure 3a). The transition in streamwise variance production near the surface (figure 4e) also shows a peak at $-z_i/L=0.433$ followed by a nearly linear decrease.

Figure 4. Transition of turbulence statistics with the global instability parameter $-z_i/L$ at fixed vertical locations in the ABL, as specified by the line types and symbols defined in panel (a): (a) streamwise velocity variance, $\langle {u^{\prime }}^{2} \rangle$; (b) vertical velocity variance, $\langle {w^{\prime }}^{2} \rangle$; (c) rate of production of vertical velocity variance by buoyancy force (proportional to vertical potential temperature flux), $2\langle w^{\prime }\theta ^{\prime }\rangle g/\theta _0$; (d) vertical turbulent momentum flux, $\langle u^{\prime }w^{\prime } \rangle$; (e) shear production rate of streamwise velocity variance, $-2\langle u^{\prime }w^{\prime }\rangle ({\textrm {d}\langle u\rangle }/{\textrm {d}z})$; and (f) time rate of change of horizontal velocity variance due to vertical turbulent transport, $-({\textrm {d}\langle w^{\prime }u^{\prime }u^{\prime }\rangle }/{\textrm {d}z})$. Here ($u^{\prime },v^{\prime },{w^{\prime }}$) are the fluctuating parts of resolved velocity in a local frame of reference with $x$ aligned with the mean velocity vector at that height, $z$. The two (three) vertical lines correspond to the key transition states: $-z_i/L=0.433$ (dashed), $-z_i/L=1.080$ (dot-dashed) and $-z_i/L=2.86$ (dotted).

Away from the surface, the magnitude of shear-induced production rate (figure 4e) decreases while its transition with $-z_i/L$ shows a sharp decrease near $-z_i/L=0.433$ before asymptoting at the more unstable states. The lack of a peak in shear-induced streamwise variance production away from the surface does not support the observed peaks in the streamwise variance at these heights (figure 4a). To explain this, consider vertical turbulent transport in the streamwise variance budget (figure 4f). This figure shows that streamwise velocity fluctuations generated closer to the surface are transported to the mixed layer, on average, by vertical velocity fluctuations. The extent of this transport grows sharply around $-z_i/L \simeq 0.43$ and peaks at $-z_i/L \simeq 1.08$. These observations suggest that the streamwise variance structure across the entire ABL is affected by a combination of production from shear closer to the surface and transport of turbulence by buoyancy-driven vertical velocity fluctuations away from the surface, in a manner dependent on the stability state parameter $-z_i/L$. Furthermore, the non-monotonic variation of shear-driven streamwise variance with $-z_i/L$, in spite of constant geostrophic wind forcing, is partially a consequence of buoyancy and shear-rate contributions to pressure–strain-rate correlation.

Both vertical velocity variance (figure 4b) and its production from buoyancy (figure 4c) grow nearly linearly with $-z_i/L$ throughout the ABL, a trend consistent with the near-linear increase in surface heat flux, $Q_0=\langle w^{\prime }\theta ^{\prime } \rangle _0$, which is specified at the surface to control the changes in global stability parameter in our computations. However, we observe slight decreases in $\langle w'^{2} \rangle$ near $-z_i/L=0.433$ in the mixed layer region that can be attributed to mechanisms other than production. The transition in Reynolds stress $\langle u'w' \rangle$ with $-z_i/L$ (figure 4d), whose evolution is tied to the vertical variance through the mean gradient production term $\langle w'^{2} \rangle \partial \langle u \rangle /\partial z$, also grows monotonically with $-z_i/L$ similar to $\langle w'^{2} \rangle$, but with sharper changes around $-z_i/L = 0.433$.

In summary, almost all of the statistical measures shown in figure 4 display intriguing transition around $-z_i/L = 0.433$ although the specific details differ among the variables. The mechanics of this transition depends not only on turbulence production, but also on turbulence transport processes driven by buoyancy forces, especially around $-z_i/L = 1.08$. These observations suggest that the interesting responses to increasing $-z_i/L$ in the two-point statistical measures are dynamically correlated with the responses in single-point statistics (figure 4). This is particularly true of the streamwise coherence lengths of streamwise velocity fluctuations (figure 3a) and vertical velocity fluctuations (figure 3b), which peak at $-z_i/L = 0.433$ and $-z_i/L=1.08$, respectively.

4.1. The existence of a critical transition in statistical coherence

We argue that the extreme sensitivity shown in figure 3 between very small increases in heating rate and ABL turbulence structure suggests critical transition dynamics. The streamwise coherence lengths of both streamwise and vertical turbulence velocity fluctuations ($L_{11,1}$ and $L_{33,1}$) more than double across the ABL with extremely small increases in instability state parameter at $-z_i/L$ values considered to be ‘near-neutral’. In the mixed layer, in particular, $L_{11,1}$ increases from $0.7z_i$ to $2.1z_i$, a factor of 3 increase, while $L_{33,1}$ increases from roughly $0.16z_i$ to $0.70z_i$, over a four-fold growth from neutral to critical states. These dramatic increases in coherence lengths at critical transition occur with much smaller relative increases in horizontal variances (figures 4a and 4b).

Figure 5 shows that, at the critical stability state, the entire ABL structure changes dramatically. Here we plot the horizontal integral scales of horizontal (figure 5a) and vertical velocity (figure 5b) fluctuations at fixed stability states across the boundary layer. In the current discussion, we focus on the transition between the blue and red curves. The blue curves indicate relatively minor changes in horizontal integral-scale structure in the near-neutral precritical (subcritical) states, culminating with the line with symbols at $-z_i/L=0.304$ just before the sudden transition in structure represented by the red curves. Whereas $L_{11,1}$ changes relatively little within the boundary layer in the precritical (subcritical) period, $L_{33,1}$ remains virtually unchanged. However, with very small changes in stability state parameter from $-z_i/L =0.304$ (blue symbolled curve) to the critical state $-z_i/L=0.433$ (red symbolled curve), the streamwise integral scales of both horizontal and vertical velocity fluctuations increase sharply over the entire boundary layer. The value of $L_{11,1}/z_i$ changes from roughly 1 to over 2 in the lower mixed layer ($z/z_i \simeq 0.34$), while $L_{33,1}/z_i$ increases from a small fraction to order $1$ in the mid mixed layer ($z/z_i \simeq 0.44$). As in figure 3, figure 5 shows that $L_{11,1}$ reaches a sharp peak at the critical state ($-z_i/L\simeq 0.433$) while $L_{33,1}$ continues to grow more gradually beyond an initial large jump at the critical state to a peak at $-z_i/L\simeq 1.08$, a post-critical process that will be discussed in § 5 in context with the formation of large-scale rolls. Interestingly, the sudden transition in $L_{11,1}$ and $L_{33,1}$ at the critical state is much stronger in the mixed layer than either near the ground or nearer to the capping inversion, where blockage in the vertical causes reductions in correlation lengths.

Figure 5. Transition in streamwise coherence lengths over the vertical extent of the ABL at fixed instability states, $-z_i/L$. (a) Streamwise coherence length of fluctuating horizontal velocity ($L_{11,1}$). (b) Streamwise coherence length of vertical fluctuating velocity ($L_{33,1}$). The line colours are chosen to be consistent with discussions in the text on the changes in ABL coherence structure surrounding special stability states: blue indicates subcritical stability states near neutral, while red indicates the rapid transition from $-z_i/L=0.304$ (blue line with symbols) to the critical state $-z_i/L=0.433$ (red line with symbols). Yellow lines indicate the transition from the critical state to the ‘maximum coherence state’ $-z_i/L=1.08$ dominated by large-scale rolls (yellow line with symbols), followed by green lines indicating further transition to the ‘moderately convective’ states when $-z_i/L\gtrsim 2.86$ (green line with symbols).

Figure 5 shows a mild subcritical response of horizontal coherence length to slight increases in surface heating, but only in the streamwise, not the vertical, velocity fluctuations. However, in the transition to the critical state (red curves), the entire streamwise coherence of the boundary layer changes dramatically, with maximum streamwise coherence lengths peaking in the lower mixed layer for horizontal velocity fluctuations ($L_{11,1}$) and in the mixed layer for vertical fluctuations ($L_{33,1}$). This sudden change at the critical state is primarily in the coherence structure of the ABL with relatively minor impact on fluctuation levels (variances in figures 4(a) and 4(b) show at most $10\,\%$ change).

Such abrupt restructuring of the entire ABL structure suggests the initiation of a global instability from the interaction between shear-driven turbulent motions and low-level buoyancy-driven motions. Additionally, in the transition from subcritical ($-z_i/L \lesssim 0.304$) to supercritical ($-z_i/L \gtrsim 0.433$), streamwise coherence lengths change from being roughly uniform over the ABL to peaking in the mixed layer away from the surface (figure 5) – suggesting a dramatically altered dynamical state through sudden linking of the lower and mid boundary layer.

To quantify this linkage, we plot in figure 6 the coherence lengths in the vertical direction from two-point correlations in $z$ referenced to $z/z_i=0.1$ (integral length scales are therefore normalized on $0.9z_i$). Interestingly, whereas in the neutral state the horizontal and vertical velocity fluctuations have the same vertical coherence length (roughly $0.2z_i$), these vertical scales diverge with the addition of low levels of surface heat flux in the subcritical state. In particular, the vertical coherence length of horizontal fluctuations ($L_{11,3}$) increases by $56\,\%$ at the critical state while the average vertical coherence length of vertical velocity concentrations changes little. As the streaks suddenly extend in the horizontal by a factor of 3 in the lower mixed layer, they fatten in the vertical by roughly $50\,\%$. At the same time, the concentrations of vertical velocity do not change in scale in the vertical, but extend in the horizontal by over a factor of 4 (!) with slight increase in surface heat flux. We will show in the next section that, as $-z_i/L$ continues to increase from critical, the vertical scale of horizontal fluctuations decreases post-critical while that of vertical velocity initiates a monotonic increase until they match.

Figure 6. Transitions of vertical coherence lengths with instability parameter $-z_i/L$ calculated vertically from level $z=0.1z_i$. The solid, dashed and dotted curves correspond to vertical coherence length of horizontal, spanwise and vertical velocity fluctuations, respectively. The three vertical lines correspond to transition states: $-z_i/L=0.433$ (dashed), $-z_i/L=1.080$ (dot-dashed) and $-z_i/L=2.86$ (dotted).

We conclude that the initiation of the dramatic transition in ABL structure occurs primarily in the streamwise coherence of both horizontal and vertical velocity fluctuations in the lower (horizontal) and mid (vertical) surface layer, with a fattening of the horizontal fluctuations (streaks) in the vertical. The implication is that the critical state corresponds to an increase in communication between the lower and mid ABL through the stimulation of an overlap in the lower-to-mid mixed layer between concentrations of horizontal and vertical fluctuations. This takes place as the horizontal coherence lengths of both increase dramatically – without the stimulation of thermals in the vertical. However, whereas transition to the critical state does not stimulate a change in vertical coherence of vertical fluctuations, $L_{33,3}$, it does initiate a monotonic growth (of $L_{33,3}$ discussed in §§ 5 and 6). It appears that the dramatic restructuring of the entire boundary layer around the critical transition state is not due to a sudden creation of non-local communication between the lower and upper ABL through buoyancy-driven vertical motions. Rather, the critical state initiates a supercritical ABL that is now susceptible to the generation of buoyancy-driven thermals with continued increases in surface heating.

4.2. Role of near-surface streaks in the ‘critical’ transition in turbulence structure

A bifurcation-like critical transition that manifests as an abrupt change in the streamwise coherence lengths of both streamwise and vertical velocity fluctuations, as well as in vertical coherence lengths of streamwise velocity fluctuations, also initiates a supercritical ABL where shear-driven low-speed streaks below align with concentrations of high-temperature and vertical velocity fluctuations above – a process culminating in the formation of large-scale rolls (§ 5).

In figure 7 we superpose an isosurface of vertical velocity at $+2\sigma _w$ ($\sigma _w$ is the standard deviation of vertical velocity fluctuations) over a near-surface plane ($z/z_i=0.1$) of isocontours of fluctuating velocity in the unrotated $x$ direction. These plots show the ABL from above over the $5\ \textrm {km}\times 5\ \textrm {km}$ horizontal domain at different $-z_i/L$. Whereas the red isosurface shows localized coherent thermal updraughts in the mixed layer, the plane of isocontours shows the coherent shear-generated turbulence structures in the surface layer – ‘low-speed streaks’ in blue and the less coherent ‘high-speed streaks’ in red.

Figure 7. Simultaneous visualization of isocontours of $x$-component velocity on the $z/z_i=0.1$ plane in unrotated coordinates (the underlying blue–green–red isocontours) and a single isosurface of vertical velocity (dark red) at a ${w^{\prime }}$ level of $2\sigma _w$ ($\sigma =$ standard deviation). The isocontours of $u^{\prime }$ vary from $-2\sigma _u$ (dark blue) to $+2\sigma _u$ (bright red); the solid green colour is close to zero. The dark blue and blue isocontours identify coherent ‘low-speed streaks’ generated by mean shear rate, which is dominant in the surface layer. The bright red/orange isocontours identify ‘high-speed streaks’. The dark red isosurface identifies concentrations of strong buoyancy-driven updraughts in the mixed layer. From these images, one can identify the spatial relationships between shear-driven low-speed streaks below and updraughts above, as a function of stability state. The transitional states discussed in the text are indicated with bold in the following: (a) $-z_i/L=0.000$; (b) $-z_i/L=0.304$; (c) $-z_i/L=\boldsymbol {0.433}$ (critical); (d) $-z_i/L=0.744$; (e) $-z_i/L=0.916$; (f) $-z_i/L=\boldsymbol {1.08}$ (maximally coherent roll state); (g) $-z_i/L=1.42$; (h) $-z_i/L=\boldsymbol {2.86}$ (transition to the canonical moderately convective ABL); and (i) $-z_i/L={3.84}$.

Figures 7(a) and 7(b) correspond to subcritical instability states. At these states, the regions of concentrated vertical velocity visually display much less coherence than do the low-speed streaks below, with relatively little apparent alignment between the updraughts and streaks. In contrast, the turbulence at critical and supercritical instability states (figures 7c and 7d) indicate higher levels of strength and organization in vertical velocity in the mixed layer. These developing vertical updraughts not only are more elongated in the horizontal, consistent with the dramatic increase in streamwise integral scales of horizontal and vertical velocity fluctuations (figure 3), but also are much more obviously aligned with the underlying low-speed streaks in the surface layer. Indeed, comparison of the supercritical state (figure 7d) with the subcritical state (figure 7b) clearly shows a strong change in both the structure of updraughts – they become more concentrated, well defined and coherent along the streaks – and more obvious alignment with the stronger low-speed streak structures below. Khanna & Brasseur (Reference Khanna and Brasseur1998) argued that an increased alignment between mixed layer updraughts above and surface layer streaks below results from the near-surface concentration of higher-temperature fluid within shear-generated, highly coherent, low-speed streaks, concentrating buoyancy force there and driving flow vertically within sheet-like regions of concentrated vertical velocity fluctuations. In this way, the low-speed streaks operate as ‘seeds’ to updraughts when buoyancy force drives vertical motions in the mixed layer post-critical.

In summary, we see from both observation (figure 7) and quantification (figure 3) that the streamwise coherence of surface layer elongated streaks suddenly increases with increasing surface heating (and $-z_i/L$) at the critical state ($-z_i/L = 0.433$). In the post-critical states, observations up to $-z_i/L =1.08$ (figures 7d–7f) indicate that the low-speed streaks become increasingly well defined and coherent visually even as the streamwise coherence length $L_{11,1}$ decreases with increasing $-z_i/L$ (figure 3a).

5.1. The transition from critical to ‘maximum coherence’ stability states

We have learned that the transition between the neutral and moderately convective states in the equilibrium ABL is non-monotonic, with a critical stability state occurring at 1 %–2 % of peak midday flux and sudden changes in streak and updraught structure. Specifically, the streaks elongate in the horizontal and thicken in the vertical with increase in $-z_i/L$ such that the aspect ratio $L_{11,1}/L_{11,3}$ in the upper surface layer changes from nearly $5$ in the neutral ABL to ${\sim }9$ in the critical ABL state (figures 3a, 5a and 6). The updraughts elongate, but do not thicken over this critical transition, with the aspect ratio $L_{33,1}/L_{33,3}$ increasing from ${\sim }1$ to ${\sim }5$ (figures 3b, 5b and 6) in the mid mixed layer, suggesting amalgamation into ribbon-like structures. The centre of these elongated updraughts lies just over that of the thickened streaks below (figure 5). Interestingly, the fact that the vertical extent of the updraughts has not changed significantly at the critical state suggests that the updraughts are not yet sensitive, in the vertical, to direct forcing from buoyancy. However, the concentrations of horizontal velocity fluctuations have thickened in the vertical to now partially overlap with updraughts (figure 5). We argue that this restructuring at the critical state with sudden spatial overlap between the upper margins of the streaks and lower margins of the updraughts yields a dynamically altered ABL in how it responds to further increases in surface heat flux post-critical as discussed below.

Continued increases in surface heating (and buoyancy force) post-critical appear to strengthen the coherence of the shear-driven horizontal fluctuations within low-speed streaks (figures 7d–7f) in spite of the decrease in $L_{11,1}$ (figure 3a). The decrease in $L_{11,1}$ post-critical appears to be driven by increasingly incoherent concentrations of the high-speed horizontal velocity fluctuations in between increasingly coherent shear-driven low-speed streaks. The nonlinear interplay between surface-heating-induced buoyancy force and shear-driven streak coherence manifests differently in the dynamically altered supercritical ABL as compared to the subcritical ABL. Whereas horizontal coherence lengths of both horizontal and vertical fluctuations respond mildly to slight increases in surface heat flux subcritical, visualizations (figure 7) and quantifications (figure 5) show that further increases in $-z_i/L$ post-critical, lead to vertical extension of the ribbon-like concentrations of updraughts. This vertical extension is in response to increases in buoyancy force (figure 6), even as the updraughts continue to extend in the horizontal (figure 3b). This is particularly the case away from the surface where vertical velocity is not damped, suggesting that the critical state initiates a new susceptibility of updraughts to the influence of buoyancy.

The visualizations in figures 7(c)–7(f) indicate that the critical state initiates a strong spatial relationship between the horizontally extending and vertically growing updraughts above, and the fattened low-speed streaks below. The supercritical regime is characterized by spatial correlation between streaks and updraughts – turbulence structures created by fundamentally different physical mechanisms that strengthen with increasing surface heat flux and $-z_i/L$. This leaves the supercritical states receptive to buoyancy and buoyancy-induced vertical motions that increasingly couple the lower and upper atmosphere with increasing $-z_i/L$. We argue that this post-critical regime with spatial correlation between streaks and updraughts (figure 7) is the initiation of the Khanna & Brasseur (Reference Khanna and Brasseur1998) mechanism, whereby low-speed streaks act as seeds for updraughts.

Figure 3 shows the existence of two well-defined peaks in horizontal coherence length: $L_{11,1}$, which peaks at the critical state $-z_i/L=0.433$ and declines monotonically with $-z_i/L$ post-critical; and $L_{33,1}$, which peaks at $-z_i/L = 1.08$, a ‘maximum coherence state’ for updraughts. The horizontal scale of the updraughts increases by approximately $70\,\%$ in the mixed layer from the critical to maximum coherence states (figures 3b and 5b), an increase of approximately 8 times the scale at the initial neutral ABL state. The change in the relative streamwise coherence scale of updraughts to streaks, $L_{33,1}/L_{11,1}$, is shown in figure 8(a). Over the range of supercritical states $-z_i/L\sim 0.433\text {--}1.08$, where the horizontal scale of the updraughts peaks, the horizontal scale of the streaks ($L_{11,1}$) has decreased to become only ${\sim }30\,\%$ larger than that of the updraughts at the maximum coherence state. Figure 8(a) shows that the streamwise coherence length of updraughts ($L_{33,1}$) approaches that of the streaks ($L_{11,1}$) with increasing $-z_i/L$ post-critical, with a major reduction in the rate of increase in $L_{33,1}/L_{11,1}$ with $-z_i/L$ taking place as $-z_i/L$ increases from below to above the maximum coherence state, $-z_i/L=1.08$.

Figure 8. Transition in the ratio of coherence lengths of vertical to horizontal velocity fluctuations with stability state $-z_i/L$ at fixed heights: (a) streamwise coherence lengths and (b) vertical coherence lengths referenced to $0.1z_i$. The three vertical lines in grey correspond to the key transition states: $-z_i/L=0.433$ (dashed), $-z_i/L=1.08$ (dot-dashed) and $-z_i/L=2.86$ (dotted). The horizontal dotted line at 1.0 identifies the transition states at which the scales of horizontal and vertical velocity fluctuations match.

The maximum coherence state is characterized not only by a peak in horizontal scale of the updraughts (figure 3b) and by a sudden change in the relative rates of increase in horizontal scale of updraughts versus streaks with $-z_i/L$ (figure 8a), but also by the changes in the vertical scales of updraughts versus streaks. As shown in figure 6, increased buoyancy force extends the updraughts vertically ($L_{33,3}$) while the vertical thickness of the streaks ($L_{11,3}$) reduces. Other than the neutral state where vertical velocity fluctuations are driven indirectly by mean shear, there is one stability state where the vertical correlation length of buoyancy-induced updraughts matches the vertical scale of the shear-driven streaks: the maximum coherence state $-z_i/L=1.08$ (figure 6), with strong vertical overlap between them. However, as surface heat flux and $-z_i/L$ increase beyond the maximum coherence state, the thickness of streaks continues to decrease and that of the updraughts continues to increase, as shown in figure 8(b), where the ratio of vertical length scales, $L_{33,3}/L_{11,3}$, is plotted against $-z_i/L$. Similar to the ratio of horizontal scales shown in figure 8(a), the rate of increase in $L_{33,3}/L_{11,3}$ with $-z_i/L$ suddenly decreases at the maximum coherence state where the ratio passes through 1, and eventually saturates at $-z_i/L\sim 3$ (§ 6.2).

5.2. Large-scale rolls at the maximum coherence state

A key characteristic of the post-critical equilibrium ABL state is the progressive increase in the spatial overlap between streaks and updraughts (figure 7) with increases in instability. In this new supercritical ABL state the vertical motions not only are more receptive to buoyancy, but also thicken in the vertical to increasingly couple the lower and upper ABL with increases in surface heat flux. A manifestation of this coupling is ‘large-scale rolls’, a well-known structural feature of the ABL consisting of $z_i$-scale rows of highly elongated updraughts between broader downdraught regions in the presence of mean streamwise flow, creating overall streamwise helical patterns that move warmer air to the capping inversion and cooler air to the ground. Large-scale rolls are often described in the literature in relationship to parallel lines of clouds near the capping inversion that are created by the condensation of moisture carried aloft within the coherent sheet-like updraughts (LeMone Reference LeMone1973; Etling & Brown Reference Etling and Brown1993; Weckwerth et al. Reference Weckwerth, Wilson, Wakimoto and Crook1997, Reference Weckwerth, Horst and Wilson1999).

One can study the changes in ABL streamwise coherence structure across the boundary layer (figure 5) as it transitions from the critical to the maximum coherence states to derive insights into the formation of rolls. The updraughts (figure 5b) transition from subcritical to critical (blue to red symbolled curves) by the sudden development of a strong well-defined peak in streamwise coherence length in the mid mixed layer ($z/z_i\approx 0.44$). In the transition to maximum coherence state (from red to yellow symbolled curves) the peak broadens over the mixed layer as the updraughts amalgamate in the horizontal and thicken in the vertical. By contrast, after the horizontal scales of the horizontal fluctuations ($L_{11,1}$) increase at the critical state and peak in the lower mixed layer ($z/z_i\approx 0.34$), as shown by the red symbolled curve in figure 5(a), the transition to the maximum coherence state $-z_i/L=1.08$ involves a complete restructuring of streamwise coherence structure across the ABL. The transition is to one that, we argue, reflects the formation of maximally coherent large-scale rolls. In particular, the peak in $L_{11,1}$ in the critical state collapses to a local minimum while two new peaks in $L_{11,1}$ appear near the upper and lower ABL, where vertical motions are suppressed. A smaller peak ($L_{11,1}/z_i\approx 1.56$) appears in the upper surface layer ($z/z_i\approx 0.22$) while a larger peak ($L_{11,1}/z_i \approx 1.70$) appears just below the capping inversion ($z/z_i\approx 0.88$). The upper peak, we argue, arises from strong vertical velocities, induced in updraughts highly receptive to buoyancy, that are forced horizontal by blockage at the capping inversion. This reflects the maximum coherence state, where the streamwise coherence length of updraughts ($L_{33,1}$) is maximum in the mixed layer. In this way the large-scale roll state is distinguished by its cylindrical topology – exceptionally strong horizontal coherence of vertical velocity in the mixed layer combines with exceptionally strong horizontal coherence of horizontal velocity fluctuations ($L_{11,1}$) in the upper boundary layer. This change in morphology of ABL coherent structure initiates the existence of roll eddies and corresponds to the matching of vertical coherence length of updraughts with that of the streaks (figures 6 and 8b).

Figure 7 visually illustrates the transition to large-scale rolls in the approach to the maximum coherence state. At stability states below and above $-z_i/L=1.08$, many small-scale concentrations of low-speed horizontal velocity fluctuation are observed between the dominant updraughts that are aligned with the dominant low-speed streaks. The maximum coherence state (figure 7f) is characterized by exceptionally coherent updraughts that spatially overlap exceptionally coherent and extended low-speed streaks, with few low-speed horizontal fluctuations in between. Fluid particles that originate at the lower margins of the exceptionally strong updraughts (red isosurfaces) move vertically, on average, as they move longitudinally with the mean flow. At the capping inversion the fluid particles are forced horizontally towards the centreline between neighbouring updraughts, generating a local maximum in $L_{11,1}$. Between the updraughts, fluid particles are forced downward within less coherent downdraughts, where they are again forced horizontal by surface blockage, creating another local peak in $L_{11,1}$ and completing the helical structure of the large-scale rolls.

Because the lateral separation distance between coherent thermal updraughts scale on $z_i$, the helical streamline patterns are expected to extend over multiple $z_i$ in the streamwise direction. Interestingly, while the streamwise coherence lengths of the vertical and horizontal velocity fluctuations peak at $1z_{i}\text {--}2z_i$ in the maximum coherence state, the visual streamwise extent of the coherent large-scale roll structure, in both figure 7 and in meteorological images of cloud-topped large-scale rolls, extend well outside the simulation domain of $7z_{i}\text {--}8z_i$. Cloud visualizations (Weckwerth et al. Reference Weckwerth, Wilson, Wakimoto and Crook1997, Reference Weckwerth, Horst and Wilson1999) suggest that the coherent large-scale rolls often extend over tens of $z_i$. This is consistent with simulations performed over an extended domain as reported in appendix A.3, where the rolls appear to extend beyond $O(10z_i)$ at the maximum coherence state.

In summary, as heat flux and $-z_i/L$ continue to increase post-critical ($-z_i/L\gtrsim 0.433$), the dynamically altered supercritical ABL, with updraughts becoming more receptive to buoyancy force, undergoes a dramatic, albeit monotonic, change in its coherence structure towards a ‘maximum coherence state’ at $-z_i/L=1.08$. This special state is characterized by maximally coherent large-scale rolls realized as $L_{33,1}$ approaches $L_{11,1}$ in the mixed layer and when $L_{33,3}$ approximately matches $L_{11,3}$ in response to a strong spatial overlap between the lower margins of the updraughts and upper margins of the streaks, creating exceptionally strong coupling across the depth of the ABL. Thus, the formation of rolls is a manifestation of the Khanna & Brasseur (Reference Khanna and Brasseur1998) mechanism at a state with maximal global response to buoyancy force in the mixed and upper layers of the ABL, localized to exceptionally coherent low-speed streaks in the upper surface layer.

6.1. The subcritical stability states ($0\leq -z_i/L< 0.371$)

The subcritical stability states may be traditionally described as ‘near neutral’ – an equilibrium ABL with finite but sufficiently small levels of surface heat flux such that the ABL has properties ‘close’ to the neutral state. In principle, near-neutral states are characterized by negligible role of buoyancy force, with temperature effectively a passive scalar. The transition to temperature as an active scalar initiates near the critical stability state in the range $-z_i/L\sim 0.371\text {--}0.433$ (figure 5), with very low levels of surface heat flux that increasingly destroys the near-neutral character of the ABL.

From neutral through $-z_i/L\lesssim 0.304$, however, our simulations produce an ABL that appears to be roughly neutral. Interestingly, however, the streamwise velocity fluctuations – but not the vertical velocity fluctuations – respond mildly to the extremely low subcritical levels of heat addition to the neutral ABL. Figure 4(a) indicates small changes in the variance of streamwise velocity. Figures 3(a), 5(a) and 6 show small changes in both the streamwise and vertical coherence structure of the streamwise velocity fluctuations across the boundary layer, while figures 3(b), 5(b) and 6 indicate negligible subcritical response in the vertical velocity fluctuations. These results suggest the existence of mild, but measurable, impacts of extremely low-level surface heating on the streamwise velocity fluctuations – just below a critical level of surface heat flux that causes complete restructuring of ABL turbulence and a fundamental change from passive to active temperature response. It is interesting that this low-level response is entirely in the streamwise velocity fluctuations; vertical fluctuations respond only at and post-critical.

6.2. Transition to the ‘moderately convective’ stability state

The ‘maximum coherence state’ associated with highly coherent large-scale rolls (see § 5) is quantitatively characterized by a peak in streamwise coherence length of the amalgamated updraughts at $-z_i/L = 1.08$ and preceded by a sudden change in ABL receptivity to surface heating at the critical state. Both of these transitions in ABL character are visually identifiable by the structure of horizontal versus vertical fluctuating velocity concentrations and the spatial relationships between them.

However, there is a third transition in ABL turbulence structure around $-z_i/L = 2.86$ that is characterized quantitatively by asymptotes and peaks in scale relationships (figures 6 and 8) and qualitatively by changes in updraughts versus streak structure (figure 7). This occurs as increases in $-z_i/L$ move the ABL into the ‘moderately convective’ states described in previous studies (Moeng & Sullivan Reference Moeng and Sullivan1994; Khanna & Brasseur Reference Khanna and Brasseur1998). Figure 6 shows that the vertical coherence lengths of both vertical and horizontal fluctuating velocities asymptote at $-z_i/L > 2.86$, suggesting a regime whereby additional increases in surface heat flux do not cause much change in their vertical scales. Post-critical, the vertical scale of updraughts monotonically increases while that of the streaks decreases monotonically (figure 6). The ratio $L_{33,3}/L_{11,3}$ passes through 1 at the maximum coherence state and continues to grow to a plateau at $-z_i/L > 2.86$ (figure 8b), where the vertical scale of updraughts is approximately $56\,\%$ greater than that of streaks.

The horizontal coherence lengths of the updraughts versus streaks, by contrast, approach one another with the ratio of horizontal scales $L_{33,1}/L_{11,1}$ (figure 8a) hovering near 1 in the mixed layer at $-z_i/L\gtrsim 2.86$. Near the ground, where vertical motions are suppressed ($z/z_i=0.1$), the horizontal scales of both vertical and horizontal fluctuations reduce slowly when $-z_i/L$ exceeds 2.86 (figure 3) so that the ratio $L_{33,1}/L_{11,1}$ plateaus with streaks nearly twice as long in the horizontal as updraughts, in contrast to the mixed layer.

Visualizations in figures 7(f)–7(i) clearly show that the horizontal coherence of streaks and updraughts breaks down at higher $-z_i/L$ beyond the maximum coherence state. In addition, the spatial correlation between streaks and updraughts, which starts at the critical state $-z_i/L=0.433$ and strengthens up to maximal spatial correlation at $-z_i/L=1.08$, reduces in the transition from the maximum coherence state to the moderately convective limit. Presumably, this process will continue into the fully convective state as the impacts of buoyancy on vertical motions fully dominate impacts of shear and the rolls transition to turbulent Rayleigh–Bénard convection cells. What distinguishes the moderately convective ABL is the continued importance of shear rate and, correspondingly, a roll-like structure, even as convection increasingly dominates and roll coherence diminishes.

Analysis of the ABL streamwise coherence structure (figure 5) shows that the horizontal scales of both vertical and horizontal velocity fluctuations reduce monotonically at all heights as the stability state transitions from $-z_i/L=1.08$ (yellow curves) to $-z_i/L=3.84$ (green curves). However, in comparing the moderately convective states (green curves) with the neutral state (blue curves), the streamwise coherence lengths of horizontal velocity fluctuations in the moderately convective state are below that of the streaks in the neutral state. By contrast, the horizontal scale of the vertical fluctuations remains well above that in the neutral state, by a factor of 3 in the mixed layer. This difference between the neutral and moderately convective states reflects the two different sources of vertical velocity fluctuation, namely from buoyancy (post-critical) and indirectly from shear (in the neutral). Correspondingly, in the transition to moderately convective ABL, the buoyancy-driven mixed layer velocity scale $w^{*}$ progressively dominates the shear-driven surface layer velocity scale $u^{*}$ (figure 2b).