The inversion layer (IL) of a clear-sky, buoyancy-driven convective boundary layer is investigated using large-eddy simulations covering a wide range of convective Richardson numbers. A new model of the IL is suggested and tested. The model performs better than previous first-order models of the entrainment and provides physical insights into the main controls of the mixed-layer and IL growths. A consistent prognostic equation of the IL growth is derived, with explicit dependence on the position of the minimum buoyancy flux, convective Richardson number, and relative stratification across the inversion G. The IL model expresses the interrelationship between the position and magnitude of the minimum buoyancy flux and inversion-layer depth. These relationships emphasize why zero-order jump models of the convective boundary layer perform well under a strong inversion and show that these models miss the additional parameter G to fully characterize the entrainment process under a weak inversion. Additionally, the position of the minimum buoyancy flux within the new IL model is shown to be a key component of convective boundary layer entrainment. The new IL model is sufficiently simple to be used in numerical weather prediction or general circulation models as a way to resolve the IL in a low-vertical-resolution model.
Bulk models of the boundary layer, introduced by Ball (1960), Lilly (1968), and Betts (1973), have been useful tools to understand the dynamics of the convective boundary layer (CBL). The so-called zero-order models of the CBL assume a mixed layer of conserved variables capped by a jump corresponding to a sharp, infinitesimally thin inversion layer (IL). Even though zero-order models correctly predict the CBL growth in a linearly stratified fluid, those models cannot diagnose the relationship between the IL structure and CBL entrainment nor the CBL dynamics under conditions of a deep IL since, by construction, they assume an infinitesimally thin IL. The tight relationship between the structure of the IL and the entrainment rate becomes clear when one investigates the turbulence kinetic energy (TKE) budget (Zilitinkevich 1991; Fedorovich and Mironov 1995; Fedorovich et al. 2004). The depth of the IL defines the region of buoyancy consumption (Deardorff 1979; Hägeli et al. 2000) and is related to the CBL entrainment rate (Zilitinkevich 1991; Fedorovich and Mironov 1995; Fedorovich et al. 2004; Garcia and Mellado 2014).
In the IL, rapid vertical variations in the buoyancy fluxes can occur over a very short distance of a few tens of meters (Deardorff 1976; Randall 1980b, 1984). The IL is therefore unresolved in current generations of weather and climate models (Teixeira et al. 2008) and still remains a modeling challenge even for large-eddy simulations (LESs) (Stevens and Bretherton 1999; Bretherton et al. 1999; Stevens and Lenschow 2001; Sullivan and Patton 2011) so that direct numerical simulations (DNSs) of the IL may be needed (Mellado et al. 2010; De Lozar and Mellado 2013; Garcia and Mellado 2014). A better understanding of the dynamics of the IL is crucial, because the IL impacts (i) the convective inhibition and the transition from dry to shallow convection, as well as the magnitude of the cloud-base mass flux (Bretherton and Park 2009; Park and Bretherton 2009; Hohenegger and Bretherton 2011; Gentine et al. 2013a,b,c); (ii) the transition between stratocumulus and shallow cumulus (Randall 1980a); and (iii) the dynamics of the entrainment and the heat and moisture budgets of the boundary layer. The IL is therefore at the core of a correct representation of most boundary layer regimes and the transition between them.
To represent the finite depth of the IL, Betts (1974) introduced a first-order model of the dry CBL with a finite IL, where the profiles of conserved variables and their fluxes were linear. The entrainment at the mixed-layer top was derived assuming a constant potential temperature jump and constant inversion depth. Deardorff (1979) argued that first-order models of the boundary layer were not more useful than zero-order models, because in the former models, the mixed-layer top zm is assumed to match the minimum buoyancy flux height zi. This assumed correspondence between the mixed-layer top and minimum buoyancy flux height generates a singularity for the IL growth rate equation under strong inversions. VanZanten et al. (1999) reinvestigated Deardorff’s (1979) work but did not share his conclusions. Since Deardorff et al. (1980), the IL depth has generally been parameterized [e.g., defined as a function of the convective Richardson number (Ri)] to avoid the singularity pointed out by Deardorff (1979).
In existing first-order models based on Betts (1974), both the conserved variable and heat flux profiles are assumed linear in the inversion. It is clear that these profiles are incompatible: parabolic flux profiles should correspond to linear profiles of conserved variables. Deardorff (1979), followed by Fedorovich and Mironov (1995) and Fedorovich et al. (2004), proposed a more realistic representation of the IL, the so-called general structure model (GSM). In the GSM, consistent IL profiles and prognostic equations for the IL growth can be derived. Nonetheless, the structure of the IL in the GSM needs to be parameterized (Fedorovich and Mironov 1995; Fedorovich et al. 2004), which limits its applicability.
The objectives of the present paper are as follows: 1) to investigate the structure of the IL using LES and 2) to develop a minimal IL structure model, which does not require additional parameterization and can correctly predict the IL growth and structure. The focus is on the dry buoyancy-driven convective boundary layer; however, the results are expected to be extended to shear-driven and stratocumulus CBL. The paper is organized as follows: in section 2 we analyze previous bulk models; in section 3 we propose a structure of the IL based on LES results; in section 4 we develop the new inversion model that can resolve the problems of first-order CBL models, and we provide a derivation of the IL growth based on the TKE budget; in sections 5 and 6 we present the model results; and in section 7 we present conclusions of the study.
2. Discussion of current first-order models
The CBL is assumed to be dry (i.e., without moisture). Results are therefore presented in terms of potential temperature . The profiles of and of the vertical turbulent flux of , , in a first-order CBL model are depicted in Fig. 1 (Betts 1974; Deardorff 1979). In these models, is collocated with , and the flux is linearly decaying in the mixed layer. In the IL the potential temperature profile and the corresponding flux profile are linear, which is inconsistent, as discussed in the introduction.
The ensemble-averaging operator is denoted with an overbar. Neglecting horizontal advection and the radiation contribution, which are small in the dry CBL, and assuming a horizontally homogeneous CBL, the potential temperature conservation reads (neglecting overbars for averages for the mean variables) (e.g., Kim et al. 2006)
with being the turbulent vertical flux of potential temperature and being the horizontally averaged vertical velocity at height z. Taking the partial derivative in z of the potential temperature conservation law and inverting the time t and z derivative, we get the simple relationship (because the vertical gradient of vanishes in the mixed layer)
The profile of is therefore linear in the mixed layer. This is a well-known result of mixed-layer models (Lilly 1968).
In the IL, extending from zm to the top of the inversion h, if we take the second derivative of the conservation Eq. (1), we find that the third-order derivative of vanishes. A linear potential temperature profile therefore imposes a parabolic profile in the IL. With a parabolic flux profile, could be located above zm. We have evaluated the use of such a parabolic flux profile, and it generates an ill-defined solution when using the flux closure; therefore, it is not discussed further.
Soon after the introduction of the first-order model by Betts (1973, 1974), Deardorff (1979) pointed out that the representation of the inversion in such models was oversimplified for several reasons:
The observed maximum vertical gradient is generally much higher in observations than in first-order models. Indeed the profile exhibits a strong curvature in the IL with high gradient near the top of the inversion.
The minimum buoyancy level is located above the mixed-layer top.
The assumption that zm equals in first-order models generates an ill-defined solution: the prognostic equation of the IL leads to spurious results in the case of a sharp inversion.
Deardorff (1979) therefore concluded that first-order models were not any more useful than zero-order models. As a result, in subsequent uses of first-order models, the inversion depth has been parameterized (e.g., Deardorff et al. 1980; vanZanten et al. 1999; Pino et al. 2006; Kim et al. 2006).
To correct for the above-mentioned shortcomings of first-order models, Deardorff (1979) proposed a more realistic IL parameterization. Fedorovich and Mironov (1995) and Fedorovich et al. (2004) improved Deardorff’s GSM by introducing a self-similar representation of the buoyancy profile within the IL and evaluated the model using LES data. One of the issues with the GSM is that, in this model, the IL profile of buoyancy needs to be parameterized. How additional factors such as shear or radiation will affect such parameterization is unclear. We therefore try to derive an IL model that is sufficiently general so that additional components (i.e., shear, radiation, etc.) can be later added. We start with a simple dry shear-free CBL case; however, the model can easily be extended to include liquid water and radiation. To derive the model, we first investigate the typical IL structure in the shear-free CBL using LES data.
3. Large-eddy simulation analysis
a. Large-eddy simulation setups
In this analysis, we have used three shear-free CBL cases: (i) a weak-inversion case (Sullivan et al. 1998), (ii) a strong-inversion case (Sullivan et al. 1998), and (iii) a CBL growing against a linear stratification (Conzemius and Fedorovich 2006). All three simulations have been forced by a time-constant surface buoyancy flux and no-slip boundary conditions for the horizontal velocity and impermeability condition for the vertical component of the velocity. The three initial and final profiles are depicted in Fig. 2.
The cases have been run with MicroHH (C. van Heerwaarden et al. 2014, unpublished manuscript) (http://github.com/microhh). MicroHH is a three-dimensional computational fluid dynamics code that solves the filtered Navier–Stokes equations with the Boussinesq approximation applied on a staggered grid. The model uses a fully conservative second-order finite difference scheme in space (Morinishi et al. 1998) and a third-order low-storage Runge–Kutta scheme in time. The pressure is solved using fast Fourier transforms in the horizontal dimensions. The subfilter-scale diffusion has been modeled using a Smagorinsky scheme (Smagorinsky 1963) with a constant cs = 0.21 and a turbulent Prandtl number of ⅓. The first vertical level of the model is solved assuming the Monin–Obukhov similarity theory. At the top of the domain, a sponge layer is applied that dampens the propagating gravity waves in order to prevent reflection. The lateral boundary conditions are periodic. The time step of integration is found so that a Courant number of 1 is enforced in order to satisfy the Courant–Friedrichs–Lewy (CFL) condition.
The domain size in the case of the strong- and weak-inversion runs was 5120 m × 5120 m × 2048 m, with a grid spacing of 5.12 m × 5.12 m × 5.12 m; the domain size of the linearly stratified case was 9600 m × 9600 m × 3200 m, with a grid spacing of 12.5 m × 12.5 m × 6.25 m. All cases were started with zero values for the entire velocity vector and small random perturbations in the temperature field in the lowest 300 m that exponentially decay with height. The output time interval is 30 min in the strong-inversion case because the growth is slow and 180 s in the other cases. Average values in the LES are defined as horizontal mean values averaged over the output time step.
b. Mixed and inversion layer in large-eddy simulations
We here investigate the horizontally averaged LES profiles in order to gain insights on the structure of the IL. We especially focus on the following heights:
The minimum buoyancy flux height.
The intersection of the linear tangent at the maximum gradient of in the IL with the linear free-tropospheric profile [as in Garcia and Mellado (2014)]. This height is called h. This height is similar to—yet less noisy than—the height of the minimum second-order vertical derivative of potential temperature (Yamaguchi and Randall 2012), which corresponds to a rapid change in the vertical gradient of from the IL to the free-tropospheric value.
The height of vanishing buoyancy flux above zi, defining the top of the IL in the LES hLES. This height is notoriously difficult to properly retrieve from LES because of overdiffusion at the top of the inversion, insufficient spatial resolution, and limited domain size. We thus use a threshold for the magnitude of the flux (10−5 K m s−1) for all cases.
Figure 3 depicts the profiles of horizontally averaged potential temperature (top left), second vertical derivative of potential temperature (bottom left), sensible heat flux (top right) and vertical derivative of sensible heat flux (bottom right) for the weak-inversion case (Sullivan et al. 1998). The heights of interest are plotted on the temperature and flux profiles. It should be emphasized that the vertical resolution of the LES limits the exact retrieval of the heights.
A mixed layer is well defined below z0: the curvature is zero in the mixed layer and increases sharply above z0, as seen in Fig. 3. Using tank observations, Deardorff (1979) showed that the zero-crossover height of could lie slightly above the mixed layer, defined as the region of near-constant potential temperature. This is not the case in our LES. This difference could be due to either the measurement accuracy of Deardorff’s observations or to the limited LES resolution. Our LES experiments suggest that z0 is a good approximation of zm [as pointed out in Deardorff (1979); Fedorovich and Mironov (1995)], and we will thus assume that zm = z0 in the rest of the manuscript. Above z0, the potential temperature increases slightly up to zi. This increase is better seen from the curvature of the potential temperature profile. No spike is present in the curvature of the potential temperature at z0, emphasizing the continuity of the vertical gradient of potential temperature at the mixed-layer top z0.
In the weak-inversion case, is only very slowly decaying above h. Above this height the vertical gradient of is hardly distinguishable from the free-tropospheric gradient. This is due to the very small portion of the total profile present above h. Overall, h seems to be a good and reliable indicator of the IL top. We will use this height h as our definition of the LES-diagnosed IL top in the remainder of the manuscript. The advantage of this height is that it is easily diagnosed in the LES outputs, contrary to hLES. The noise in the LES data at the top of the inversions limits the exact retrieval of hLES (Mason 1989; Sullivan and Patton 2011).
Similar conclusions are reached in the strong-inversion case shown in Fig. 4 and in the case of the growth in a constant linear stratification shown in Fig. 5: the mixed layer is well defined below z0. Above z0, the potential temperature curvature increases slightly by up to 0.7 K in the strong-inversion case. Even in the strong-inversion case, h is close to hLES located about 30 m above. The curvature of potential temperature exhibits a drastic change above h, indicating a near discontinuity in the vertical gradient of at this height. In all cases, more than 99% of the TKE production by buoyancy flux is located below h (not shown), emphasizing that h is a physically meaningful definition of the top of the inversion based on the horizontally averaged profile.
4. New inversion-layer model
a. Model structure
The new inversion model structure is depicted in Fig. 6 and is based on the above LES analysis. The mixed layer in the model extends up to level z0, corresponding to the first zero crossing of . As highlighted earlier, in the mixed layer the profile is linear for consistency with the profile. The potential temperature and its vertical gradient are assumed continuous at z0 based on the LES profile analysis (section 3). For consistency, and its vertical derivative are continuous at z0: the vertical gradient of is thus assumed to be zero at z0. At the top of the IL (i.e., at h) the vertical derivatives of and are assumed to be discontinuous based on the LES diagnostics, which exhibited sharp gradient change above h. Consequently, the boundary conditions in h are as follows: is zero, and is continuously merging with the free-tropospheric profile. Those latter conditions are similar to the ones adopted at the upper IL interface in first-order models. The minimum value of occurs within the IL at zi. Like in the GSM parameterizations (Deardorff 1979; Fedorovich and Mironov 1995; Fedorovich et al. 2004), zi is a model diagnostic and is found as the height of the vanishing vertical derivative of (see discussion below).
To match the imposed IL boundary conditions, we use a minimal polynomial for the potential temperature profile. A second-order polynomial is thus used for the profile of in the IL. Given the continuity conditions, the profile of in the IL can be written as
with being the mixed-layer value, assumed to be vertically uniform; being the jump of across the inversion; and being the IL depth. Examples of the model profiles fitted on LES outputs are plotted in Fig. 7; the model is able to correctly characterize the structure of the IL, especially the strongest vertical gradient of near the top of the IL.
For consistency with the potential temperature profile in the IL, we use a third-order polynomial for . The profile can be expressed using the IL boundary conditions and the minimum flux in zi [see appendix A, Eq. (A4)]:
with (i.e., ) being the relative position of the minimum buoyancy flux within the inversion.
The new inversion model has several advantages over previous inversion models. Compared to first-order models, the buoyancy flux exhibits a minimum above the mixed-layer top; the maximum gradient of potential temperature is sharper than in first-order models, as observed (Deardorff 1979); and the profiles of and are mathematically consistent. Compared to general structure models (Deardorff 1979; Fedorovich and Mironov 1995; Fedorovich et al. 2004), the model does not require any added parameterization of the shape of the IL as a function of a stability parameter; only the specification (closure) of the minimum buoyancy flux is sufficient to close the system of equations similarly to first-order models (see section 4e below).
In addition to the imposed continuity conditions at the bottom and top of the IL, to close the system of equations of the mixed layer and IL, we need (i) an equation for the mixed-layer potential temperature budget, (ii) an equation for z0, (iii) an equation for h, and (iv) a closure for .
b. Mixed-layer budget
The buoyancy conservation [Eq. (1)] is integrated from the surface to z0 using Leibniz’s rule so that the mixed-layer conservation reads
c. Entrainment rate at the mixed-layer top
The total heat in the inversion is then simply
The heat budget in the inversion is (using Leibniz’s rule)
which leads to the following expression for the mixed-layer entrainment rate (see appendix A for a full derivation):
d. Entrainment at the inversion-layer top
At the CBL top (i.e., at level h) the vertical derivatives of potential temperature and are assumed to be discontinuous based on the LES analysis. The total derivative of can be decomposed into its partial derivatives, which can be evaluated either below (left derivative in h−) or above the inversion (right derivative in h+):
The continuity of the temperature at the inversion top h gives the entrainment rate (full derivation is provided in appendix A):
with being the potential temperature jump across the inversion, being the relative position of the minimum buoyancy flux, and being the inversion depth. This expression shares similarities with the entrainment derivation of Deardorff (1979) and vanZanten et al. (1999), except that the expression uses a reduced jump , similar to a jump computed half-way within the inversion. The expression also includes a flux divergence with the presence of the surface heating.
Using dimensionless scaling, the entrainment rate can be further simplified into
This expression resembles Sullivan et al.’s (1998): the first term in the brackets corresponds to the minimum buoyancy flux contribution, and the second term corresponds to the contribution of the heat storage in the inversion. We, however, point out the presence of a shape factor in the minimum buoyancy flux and storage term, which can strongly modify the effect of minimum buoyancy flux and heat storage on the CBL entrainment rate (see section 6). This term was not present in previously proposed first-order models (e.g., Betts 1974). This shape factor explains the decrease of the minimum buoyancy flux contribution growth, as well as the decrease of the storage term contribution at a large Richardson number noticed by Sullivan et al. (1998). Note the similarity between the mixed-layer and top of the inversion entrainment rates, with the presence of the minimum buoyancy flux, storage term, and shape of the IL profile . The presence of G in the denominator of the inversion top entrainment rate [Eq. (12)] shows that two-layer nonstratified fluid (G = 0) exhibits weaker entrainment at the PBL top compared to stratified fluid, in agreement with recent observations (Jonker and Jiménez 2014). The presence of G in the top of the IL entrainment also shows that the mixed-layer scaling parameters (height, convective velocity, and Richardson convective numbers) are not sufficient to fully characterize the IL growth. However, in the strong-inversion cases G ≪ 1, the IL growth does not depend on the relative stratification anymore, and mixed-layer convective scaling is sufficient to fully characterize the IL entrainment.
In the new IL model, the position of the minimum buoyancy flux in the IL (or similarly zi) is a diagnostic variable, defined as the height of vanishing vertical flux gradient: . Differentiating the flux with respect to zi using Eq. (4), we find
Expression (13) emphasizes the tight connection between the position of the minimum buoyancy flux, the magnitude of this flux, and the relative depth of the inversion. Expression (13) also shows the relationship between the relative inversion depth and the ratio .
with being the shape factor for the mixed-layer entrainment, and similarly,
with being the shape factor of the entrainment at the IL top.
A point to emphasize is that the minimum buoyancy flux’s position zi (or ) and the relative depth of the inversion are all connected, as emphasized in Eq. (13). The mixed-layer and top of the IL entrainment equations give two additional constraints on the relationship between the minimum buoyancy flux, , and the relative depth of the inversion. A single closure on the minimum buoyancy flux, as done in zero-order CBL models, should therefore be sufficient, as long as it is described by all parameters of the flow (convective mixed-layer parameters and relative stratification). Consequently, zero-order models include not only both the direct minimum buoyancy flux effect and the IL storage term without explicitly resolving the IL (Sullivan et al. 1998), but also the shape of the inversion (i.e., the position of the minimum buoyancy flux within the inversion ). Current first-order CBL models with a parameterized IL depth assume that the IL depth scales only with the mixed-layer convective parameters; our results show that this is not correct, and the additional parameter G is needed to fully describe the IL depth and entrainment.
The mixed-layer entrainment rate [Eq. (14)] does not exhibit any dependence on G (other than in front of the large-scale vertical velocity), contrary to the top of the IL entrainment rate [term 2 − G in the denominator]. This is consistent with recent findings regarding the IL structure in DNS (Garcia and Mellado 2014), which emphasize that the IL essentially consists of two layers: 1) a lower layer, the depth of which scales with the mixed-layer convective parameters (Ri, w*) and 2) an upper layer, the depth of which does not obey mixed-layer convective scaling and is instead obtained by the distance traveled by an overshooting parcel in the free troposphere. The upper-layer depth becomes very small compared to the lower layer in a strong-inversion case. Our theoretical results are in line with Garcia and Mellado (2014): the mixed-layer growth can be characterized by mixed-layer convective parameters [Eq. (14)] in the absence of subsidence, because is related to the minimum buoyancy flux and relative inversion depth [Eq. (13)]. Conversely, the top of the IL growth [Eq. (15)] needs the additional parameter G, the relative stratification across the inversion. In a strong-inversion case, so that both the mixed-layer and IL growths can be described by mixed-layer convective parameters only. It should therefore come as no surprise that zero-order bulk models of the CBL are able to accurately describe the entrainment process under a strong inversion, because the IL in this case should scale with mixed-layer convective parameters [in the absence of subsidence, as seen in the mixed-layer entrainment in Eq. (14)]. In the presence of a weak inversion and large G, zero-order models cannot accurately describe the entrainment process, because G is missing from the set of parameters used in the closure.
e. Minimum buoyancy with integrated TKE closure
The only unknown of the new IL model is . We turn to the TKE budget in order to define the closure. Previous closures have used a local TKE budget at zi (Tennekes 1973; Zeman and Tennekes 1977; Tennekes and Driedonks 1981; Driedonks and Tennekes 1984; Driedonks 1982; Pino et al. 2003) or the integral CBL TKE budget (Fedorovich 1995; Fedorovich and Mironov 1995; Fedorovich et al. 2004; Conzemius and Fedorovich 2006, 2007; Kim et al. 2006). Here we use the integral CBL TKE budget since it avoids the need for a parameterization of the transport and pressure terms: both terms are negligible when integrated over the entire CBL depth with neglected energy transport by waves above the CBL. Another advantage of the integral TKE budget is that the contribution of the storage term of TKE is much smaller than the other terms (Kim et al. 2006), whereas it can be a large contributor of the local TKE budget at zi (Zilitinkevich 1975).
In the assumed horizontally homogeneous CBL, the local TKE budget in the Boussinesq approximation at any height reads
with e being the TKE, B being the buoyant production of TKE, S being the shear production of TKE, T being the turbulent transport of TKE, P being the pressure fluctuation–induced TKE, being the TKE viscous dissipation rate, being the reference potential temperature, and being the reference hydrostatic air density. The normalized TKE budget is shown at different time steps of the LES in Fig. 8. The turbulent transport term vanishes when integrated over the entire CBL. The vertical integral of P is negligible in the absence of gravity wave transport (Fedorovich 1995), as well as the contribution of the storage term of TKE (Kim et al. 2006). After the initial transient stage, the integral TKE budget thus simplifies to
that is, the integral TKE production by buoyancy and shear is compensated by the integral viscous dissipation. In the case of a shear-free CBL, the shear term is zero so that Eq. (17) further simplifies to
In the mixed layer, the dissipation rate is typically nearly uniform vertically, only slightly decreasing with height in the mixed layer and then almost linearly decreasing to 0 at the PBL top (Fedorovich and Mironov 1995; Kim et al. 2006), as seen in Fig. 8. We therefore assume a uniform dissipation rate in the mixed layer and a linearly decaying dissipation rate in the IL so that the integrated dissipation rate over the entire CBL is
We do not take into account the larger TKE dissipation rate in the surface layer, because the surface layer is thin and is in quasi equilibrium (Zilitinkevich 1991) so that one may assume that the TKE is locally dissipated there (Charuchittipan and Wilson 2009). The integral of the TKE buoyancy production in the mixed layer is as follows (see appendix B for a complete derivation):
and in the IL, the production integral is
The minimum buoyancy flux is then related to the dissipation rate in the mixed layer, summing the buoyancy contribution in the mixed layer [Eq. (20)] and in the IL [Eq. (21)] and using the equality of the rate of change of TKE [Eq. (19)]:
For the closure, we use a CBL scaling of the dissipation rate with the cube of the mixed-layer convective velocity (Zilitinkevich 1991; Garratt 1994; Fedorovich and Mironov 1995). Figure 8 depicts this normalized mixed-layer dissipation rate. In this case, after quasi-steady state is reached, the mixed-layer dissipation term is estimated as
It should be emphasized that this closure [Eq. (24)] gives a relationship between (i) the minimum buoyancy flux, (ii) its , and (iii) . Equation (24), together with the diagnostic expression for the minimum flux height [Eq. (13)] gives two independent equations for , , and the minimum buoyancy flux. As described above, this interdependence of , , and the minimum buoyancy flux is implicitly imposed in zero-order models through the closure on the minimum flux only, since these models do not explicitly resolve the inversion layer. Here we have chosen a closure based on the minimum buoyancy flux, but any other closure based on or as a function of other CBL parameters could have been used in practice.
f. Inversion-layer growth rate
The first term in the largest bracket corresponds to the cooling induced by the negative buoyancy flux in the inversion. This term always increases the inversion depth. The second term corresponds to the change in heat storage in the inversion and is detrimental to the inversion growth. It is also important to stress the presence of a large-scale velocity term. Subsidence, < 0, decreases the IL depth as would be expected. Most importantly, G is needed to describe the growth of the IL and is present in neither zero-order models of the CBL, nor in parameterized first-order models (since the inversion is parameterized as a function of the mixed-layer convective parameters).
5. Model evaluation
We now evaluate our new IL model in three shear-free CBL cases described above: (i) a weak-inversion case (Sullivan et al. 1998), (ii) a strong-inversion case (Sullivan et al. 1998), and (iii) a growth against a linear stratification (Conzemius and Fedorovich 2006). The model is examined against vanZanten et al.’s (1999) first-order model, in which the IL depth is parameterized as a function of the mixed-layer convective parameters, and Deardorff’s (1979) first-order model, which has two independent prognostic equations for the growth of the mixed-layer and IL top. We again stress that Deardorff showed that this model led to spurious behavior in the presence of a strong inversion, since the inversion becomes unrealistically thin. The zero-order model results, with closure , are plotted as a baseline, even though this model cannot reproduce the IL and is therefore of limited use for our main purpose here: diagnosing the evolution of the IL.
In the weak-inversion case depicted in Fig. 9, Deardorff’s model correctly predicts the evolution of zi yet overestimates h. VanZanten et al.’s model strongly underestimates elevations of the mixed layer and IL. The zero-order model accurately predicts the time dependence of zi, even though it does not take G into account. The new inversion model accurately describes the mixed-layer and IL growths and favorably compares to LES outputs. All models perform well in terms of potential temperature (besides vanZanten et al.), as would be expected, since changes in the entrainment do not drastically impact the mixed-layer potential temperature (Gentine et al. 2013a).
In the strong-inversion case depicted in Fig. 10, vanZanten et al.’s model strongly underestimates the mixed-layer growth and overestimates the IL depth. The zero-order model performs very well and accurately reproduces zi. Deardorff’s model compares favorably to the LES in the first 2 h of the simulation, but the inversion depth is then strongly underestimated. In fact, at high Ri the inversion layer represented by this model collapses, leading to a singularity, as pointed out by Deardorff (1979). The new inversion model compares very favorably to the LES output. The mixed-layer growth is well captured, as well as the growth of the inversion depth.
The models are further evaluated in the case of the growth against a linear stratification in Fig. 11. The zero-order model again performs well. VanZanten et al.’s model underestimates the mixed-layer growth: the entrainment flux at the mixed-layer top is strongly underestimated in this model; the IL grows too quickly compared to the mixed-layer depth. Deardorff’s model underestimates the growth of the boundary layer. The new inversion model again performs well compared to the observations and accurately represents the IL growth. The new model offers important improvements compared to first-order models of the PBL, as it can accurately represent both the mixed-layer and IL growths in various conditions without any additional parameterization.
6. Inversion-layer structure with LES and inversion model
We now investigate the structure of the inversion layer and the relative position of the minimum buoyancy flux within the inversion . As seen in Fig. 12, the relative depth of the inversion decreases with the convective Richardson number. Deardorff’s (1983) parameterization of the inversion-layer depth as a function of the Richardson number, , overestimates the inversion-layer depth. Approaches based on a parcel theory approach (converting updraft buoyancy into kinetic energy consumption) (Stull 1973; Neggers et al. 2007) overestimate the dependence of the inversion depth on the Richardson number. The parameterization of Neggers et al. (2007) correctly performs at high Richardson numbers under a strong inversion. This emphasizes that, under a strong inversion, a parcel approach correctly represents the entrainment process and IL structure—this may reflect the fact that, under a strong inversion, the horizontally averaged IL represents an averaging of a locally varying interface with sharp edges separating a region of near-mixed-layer air from the free-tropospheric air (Lilly 2002a,b). Substantial departure from parcel theory is present at low Richardson numbers, showing that a parcel theory does not provide a good representation of the IL depth and entrainment in this regime. This may reflect the fact that under a weak-inversion, local stratification above the interfacial layer may substantially differ from the free-tropospheric value (Cohn and Angevine 2000). As discussed earlier, Garcia and Mellado (2014) showed, using DNS data, that the inversion layer may be considered as composed of two layers: a lower layer, the depth of which is related to the mixed-layer thickness, and an upper layer dominated by the penetrating thermals that are directly affected by the free-tropospheric stratification. The depth of this upper layer can be found based on a Lagrangian parcel overshooting in the free troposphere and is related to the Ozmidov length scale. Under a strong inversion, the upper layer becomes much thinner than the lower layer so that G (Ozmidov length scale) is not needed any more to characterize the IL. Garcia and Mellado (2014) also derived a parameterization of the relative IL depth: , with C = 0.3 and C1 = 0.76, plotted in Fig. 12. Our results are consistent with those recent DNS results and emphasize [Eqs. (14) and (25)] that the mixed-layer growth can be described by mixed-layer convective parameters, but the IL growth needs the additional parameter G, which is not accounted for in typical zero- and first-order CBL models.
The dependence of on the Richardson number is depicted in Fig. 12 (bottom panel). At a low Richardson number, the minimum buoyancy flux is located at 40% of the IL depth. At a low Richardson number, the position of the minimum buoyancy flux depends on the Richardson number. At a higher Richardson number, and in the strong-inversion case, asymptotically reaches ⅔. In addition, in the strong-inversion case, the IL depth is a near-constant fraction (¼) of the mixed-layer depth. The depth of the IL, which is defined as the layer between z0 and h, is nonnegligible, as assumed in mixed-layer models. If instead we used the region extending from zi to h for the IL definition, the IL depth would represent only 8% of z0 or 7% of zi. Under a weak stratification, the IL depth is a substantial fraction of the mixed-layer depth (35%–45%). The new inversion model can be used to understand how the entrainment modifies the structure of the IL and especially the position of the minimum buoyancy flux. First, the denominator of Eq. (13) imposes that < ⅔ (i.e., the minimum buoyancy flux must be located in the lower two-thirds portion of the IL). This theoretical asymptotical limit pertains under a strong inversion (i.e., small ). Our LES results in Fig. 12 indeed confirm that is always lower than this theoretical limit (⅔). Similarly, the equation of the entrainment at the IL top [Eq. (15)] imposes > ⅓. This result is confirmed by our LES data in Fig. 12.
Based on our LES data, the normalized minimum flux exhibits a high dependence for Ri < 15 (Fig. 13). For Ri > 15, the normalized flux is nearly constant with a mean value equal to 0.13, in good agreement with previous studies (e.g., Sullivan et al. 1998). The departure from a constant flux ratio and its dependence on the Richardson number has also been pointed out by Fedorovich et al. (2004). At a low Ri, the absolute LES minimum buoyancy flux is lower than the DNS value (Garcia and Mellado 2014) plotted in Fig. 13. The relative position of the minimum buoyancy flux height compared to the mixed-layer depth is also plotted on Fig. 13. Interestingly, this position is relatively constant across Ri values, on the order of 0.17, showing that the mixed-layer height provides a good scaling of the position of the minimum buoyancy flux, whereas the total IL depth cannot be characterized by mixed-layer scaling only.
The minimum buoyancy flux, the relative depth of the inversion, and the position of the minimum buoyancy are tightly coupled, as emphasized and diagnosed in the LES data and with the new inversion model in Fig. 14. The relative position of the minimum buoyancy flux in the inversion increases with , as shown in Fig. 14a, even though substantial spread of data is evident. Similarly, the relation between and is relatively loose, and no clear trend is seen (Fig. 14b). The relationship between , , and in Eq. (13) suggests that we should use instead of to evaluate the dependence between and . Indeed, using in lieu of substantially reduces the spread in the data (Fig. 14c). The model-inferred [Eq. (13)] correctly characterizes the LES-observed dependence of on , as seen in Fig. 14c. This again confirms the quality of the new inversion model predictions in a wide range of conditions. Some spread is present in the weak-inversion case but can be explained by the relatively small time step (180 s)—and hence noisy output—used in the analysis. The value exhibits a tight dependence on and decays almost linearly with (slope of −0.57). The presence of in could be expected to reduce the spread in the data, yet the near-linear decay emphasizes that and are coevolving. The relative position of the minimum buoyancy flux in the IL decreases sharply with the relative inversion depth (Fig. 14d). The modeled inversion structure is well characterized compared to LES data. It is generally assumed that the inversion layer is symmetric and that = 0.5, but in the weak-inversion case and constant stratification cases, is smaller than 0.5. In the strong-inversion case, the minimum buoyancy flux is located much higher within the IL ( 0.6).
The dependence of Y on the relative stratification (Deardorff 1979; Fedorovich and Mironov 1995) is investigated in the LES in Fig. 15. Contrary to Deardorff (1979) and Fedorovich and Mironov (1995), our LES results exhibit little dependence on G, aside from a slight decrease at high G. We point out that these differences with Deardorff (1979) and Fedorovich and Mironov (1995) are certainly explained by our different definition of the IL. Our definition of the base of the IL (first vanishing buoyancy flux) is similar to that of Deardorff (1979) and Fedorovich and Mironov (1995), but our IL top is always lower than theirs: indeed, in the LES diagnostic, we have used the intersection of the tangent at the maximum vertical gradient with the free-tropospheric profile (see section 3). We have used this alternate definition of the CBL top since the exact height of vanishing flux in the LES is highly variable and erratic in the LES outputs (Conzemius and Fedorovich 2006). As a consequence, our shape factor computation is different from previous studies. The mean shape factor diagnosed in the LES is 0.325, close to our theoretical value, Y = ⅓, assuming a parabolic potential temperature profile (see derivation in appendix A). There is a decrease of Y with G in the LES data in the weak-inversion and constant stratification cases, but the strong-inversion case does not seem to support this change of Y with G; Y is again close to ⅓ in this case. Equation (12) emphasizes that G must be less than 2 to avoid infinite entrainment. This limit is in good agreement with the maximum range of observed G values (0–1.8) based on a metadata analysis of LES and laboratory data (Fedorovich and Mironov 1995) and our range of observed values in Fig. 15.
We have proposed a new model for the clear-sky shear-free buoyancy-driven convective boundary layer that takes into account the curvature of the potential temperature profile in the inversion layer. To comprehend the structure of this layer, we have used large-eddy simulations in weak- and strong-inversion cases. The large-eddy simulations suggest that the following LES results can be applied to build a consistent model of the mixed layer and inversion layer:
The height of zero buoyancy flux and zero corresponding vertical gradient is not an accurate measure of the top of the inversion layer, since it depends on the LES resolution and numerics. Instead, we use the position of the intersection of the tangent at the steepest vertical gradient with the free-tropospheric profile or, similarly, the minimum curvature of the potential temperature profile h.
At h, the buoyancy flux and potential temperature vertical gradients can be considered discontinuous.
At the mixed-layer top, the buoyancy flux and potential temperature vertical gradients are continuous.
The proposed model resolves the inconsistencies of first-order models and does not require additional parameterization of the inversion structure, such as in the general structure model. The inversion-layer model is based on a second-order polynomial for the potential temperature profile and a third-order polynomial for the buoyancy flux profile. The main results are as follows:
The new inversion model can accurately prognosticate the rate of growth of the inversion and of the mixed layer in varying conditions (weak and strong inversions). The position of the minimum buoyancy flux can be deduced from the model and is located between the first and second third of the inversion depth, in agreement with the LES outputs. The minimum buoyancy flux height increases with the magnitude of the minimum buoyancy flux and asymptotically reaches two-thirds of the inversion-layer depth under a strong inversion. In turn, the inversion-layer depth decreases with increased magnitude of the minimum buoyancy flux.
Both the inversion depth and the position of the minimum buoyancy flux can be expressed in terms of the minimum buoyancy flux so that the CBL growth can be expressed as a single closure in terms of the minimum buoyancy flux.
The mixed-layer growth can be described by mixed-layer convective scaling (convective vertical velocity and Richardson number), whereas the inversion-layer growth needs the additional parameter G, the relative stratification across the inversion. This emphasizes that the set of parameters used in current zero- and first-order models to describe the entrainment is incomplete. The position of the minimum buoyancy flux seems, however, to scale well with the mixed-layer height, based on LES observations. When , the entrainment process and inversion-layer growth/depth can then be described by mixed-layer convective parameters only, as in zero-order models.
The relative stratification across the inversion is theoretically constrained to lie within the range 0–2, which is consistent with LES and laboratory observations. Two-layer nonstratified fluid (zero relative stratification) exhibits a smaller inversion-layer growth rate compared to stratified fluid.
This latter result emphasizes that zero-order models can correctly simulate the CBL growth even under weak inversions because they implicitly represent the flux profile curvature and depth of the inversion layer through their dependence on the minimum buoyancy flux. A correct closure on the minimum buoyancy flux is thus sufficient to represent (i) the magnitude of the entrainment flux, (ii) the storage term [as already emphasized by Sullivan et al. (1998)], and (iii) the shape of the inversion-layer profile. The shape of the inversion profile is a new component that provides additional insights in the entrainment process in the CBL.
The improved inversion model proposed here addresses and corrects the inconsistencies of first-order models pointed out by Deardorff (1979) and does not rely on a parameterization of the inversion-layer structure, as in general structure models. These results emphasize the impact of the inversion-layer vertical structure on the rates of growth of the mixed and inversion layers. The position of the minimum buoyancy flux is an important factor controlling the magnitude of the minimum buoyancy flux for the mixed-layer entrainment rate and for the inversion-layer growth. In climate models, we suggest that our new inversion-layer model can be used to refine the characterization of the inversion layer, avoiding nesting or the definition of an ambiguous layer [as in Park and Bretherton (2009)].
The authors thank Fabio D’Andrea and Alan K. Betts for stimulating discussions on the subject and three anonymous reviewers for their thorough comments. P. Gentine acknowledges support from Grant NSF-AGS 1035843, DOE-ASR ER65481, and NASA-NIP13-0095.
Derivation of Entrainment Velocities
We assume a quadratic profile for the potential temperature in the IL. To derive the profile, we use the following boundary conditions: (i) at the mixed-layer top (z = z0), in the inversion matches the mixed-layer value; (ii) the corresponding vertical gradient of is null; and (iii) is continuous at the top of the inversion and matches the free-tropospheric value. The profile in the inversion is thus
with being the mixed-layer potential temperature, being the jump of across the inversion, and being the IL depth.
The corresponding vertical turbulent flux of , , has to be a third-order polynomial, as discussed in section 4. The following conditions are assumed on the profile: (i) the flux vanishes at the mixed-layer top z0—this is, in fact, our definition of the mixed layer, as in Deardorff (1979) and Fedorovich (1995); (ii) the flux vanishes at the top of the inversion h; (iii) the flux vertical derivative is assumed to be continuous at the mixed-layer top and therefore equal to , assuming a linear flux profile in the mixed layer; and (iv) the flux is minimum at height zi within the IL, which is a model diagnostic (see below). Based on those assumptions, the profile of in the inversion can be derived:
with being the relative position of the minimum buoyancy in the inversion layer. It is obvious from the denominator of in Eq. (A2) that the upper bound of is ⅔ to avoid infinite flux in the inversion layer. Using Eq. (A2), we can express the minimum buoyancy flux as
Equation (A3) emphasizes the relationship between the minimum , the relative inversion depth, and the relative position of the minimum in the IL . In fact, only a single closure assumption is necessary, either on or on , to close the full set of equations, since is a prognostic variable. We here use a closure on , which is based on the vertically integrated TKE budget (see appendix B).
The total derivative of can be computed just above and below the top of the inversion, where the potential temperature is continuous, using the conservation equation [see Eq. (1)], neglecting radiation divergence and variations in the free-tropospheric profile. We use the definition of the total derivative as a function of partial derivatives just below (left derivative h−) or above the inversion (right derivative h+):
The local rate of change of , , is simply given by the conservation equation [see Eq. (1)]. Similarly, at the PBL top
in which is the entrainment velocity we(h) at the top of the boundary layer. Based on the assumed quadratic profile [Eq. (A1)], the potential temperature vertical gradient below the inversion top is
and the corresponding flux gradient is [using Eq. (A2)],
with . The entrainment rate of the boundary layer top is
Introducing the relative stratification (Deardorff 1979; Fedorovich and Mironov 1995), the convective Richardson number , Deardorff’s (1970) convective velocity , and the reference potential temperature under the Boussinesq approximation results in the following:
We now introduce the integral shape factor like Deardorff (1979):
With our parabolic profile [Eq. (A1)], Y = ⅓. The integral of the potential temperature in the inversion is thus
The heat budget in the inversion layer is (using Leibniz’s rule)
The first term on the left-hand side of the equation can be expanded using Eq. (A12):
The change in the potential temperature increment across the inversion is:
assuming a constant free-tropospheric profile. Assuming a constant Y and linear large-scale vertical velocity, Eq. (A13) then simplifies to
This can be rewritten using the convective velocity and Richardson number as
Using the entrainment velocity at the PBL top [Eq. (A10)], we can simplify this expression:
With our value of the shape factor Y = ⅓,
Integration of the TKE Budget
We assume that the storage term of the TKE equation is negligible compared to the other terms once the Ozmidov length scale is much smaller than the CBL depth (Dillon 1982). Zilitinkevich (1975) showed that the storage term could be an important term of the local TKE budget. In the inversion the storage is typically at least an order of magnitude larger than in the mixed layer (see Fig. 8). In our case, we are considering the vertically integrated TKE budget over the CBL depth, and the storage term is much smaller than other terms, as observed in the LES data. Only the dissipation term has to be parameterized to determine the closure of the TKE budget. Integrating the dissipation over the entire PBL gives, assuming uniform dissipation rate in the mixed layer and linear decay of the dissipation rate in the inversion layer,
The minimum buoyancy flux can now be expressed using the integrated PBL TKE budget:
with being the vertically integrated TKE buoyancy consumption (taken as negative) in the IL and being the TKE production contribution of waves, which is neglected here.
The integral of the TKE buoyancy production in the inversion is
and in the mixed layer
As the integrated TKE dissipation is , then the minimum buoyancy flux can be found as
or in dimensionless form,
Using the convective scaling for the dissipation rate, (Zilitinkevich 1991), gives
which simplifies to
This expression of the IL depth can be used in Eq. (13) to find the relationship or, conversely, .