A closure relationship between subgrid-scale (SGS) updraft–downdraft differences and resolvable-scale (RS) variables is proposed and tested for cloud-resolving models (CRMs), based on a data analysis of a large-eddy simulation (LES) of deep convection. The LES flow field is partitioned into CRM-RS and CRM-SGS using a cutoff scale that corresponds to a typical CRM grid resolution. This study first demonstrates the capability of an updraft–downdraft model framework in representing the SGS fluxes of heat, moisture, and momentum over the entire deep convection layer. It then formulates a closure scheme to relate SGS updraft–downdraft differences to horizontal gradients of RS variables. The closure is based on the idea that largest SGS and smallest RS motions are spectrally linked and hence their horizontal fluctuations must be strongly communicated. This relation leads to an SGS scheme that expresses vertical SGS fluxes in terms of horizontal gradients of RS variables, which differs from conventional downgradient eddy diffusivity models. The new scheme is shown to better represent the forward and backscatter energy transfer between CRM-RS and CRM-SGS components than conventional eddy-viscosity models.
With the increasing use of nonhydrostatic, kilometer-scale cloud-resolving models (CRMs) for weather and climate forecasts, the subgrid-scale (SGS) parameterization issue needs to be reexamined. For a CRM with a horizontal grid resolution of just a few kilometers, the resolvable-scale (RS) and the SGS are split right at or near the integral scale of vertical-velocity field, as pointed out by Moeng et al. (2010) based on a spectral analysis of a large-domain large-eddy simulation (LES). This argument is consistent with field observations that updraft–downdraft cores in deep cloud systems have an averaged diameter of 1–2 km (LeMone and Zipser 1980; Zipser and LeMone 1980; Khairoutdinov et al. 2009). The largest CRM-SGS and smallest CRM-RS motions are spectrally linked and strongly interact. They are also responsible for most of the moisture transport in deep convection systems (Moeng et al. 2010).
To examine the SGS transport process and its representation in CRMs, scientists at the National Science Foundation (NSF) Science and Technology Center for Multiscale Modeling of Atmospheric Processes (CMMAP) generated a benchmark LES, which used ~109 grid points (called Giga-LES) to simulate a tropical deep convection system (Khairoutdinov et al. 2009). The numerical domain of 204.8 km × 204.8 km × 27 km is large enough to resolve the mesoscale organization of a deep convection system, while the horizontal grid mesh of 100 m is fine enough to resolve individual clouds and energy-containing turbulent motions. The vertical grid stretches from 50 m near the surface to nearly 100 m at 18 km and to 300 m at the top of the domain. The simulated deep convection system reaches a quasi–steady state with its imposed large-scale environment after about 12 h of spinup process. This benchmark simulation has been used to examine the characteristics of the PBL under the influence of deep convection and to assess the eddy-diffusivity scheme commonly used in CRMs (Moeng et al. 2009), to investigate the assumed probability density function schemes for SGS models in CRMs (Bogenschutz et al. 2010), to study the effect of the largest SGS motions in the deep convection layer for applications to CRMs (Moeng et al. 2010), and to study the boundary layer moisture transport in CRMs (Moeng and Arakawa 2012). In those previous studies, the Giga-LES flow field was split into two components, CRM-RS and CRM-SGS, using either a low-pass smooth filter or by applying discretized area-averaging methods. The filter width or the averaging-area size mimics a CRM horizontal grid resolution.
There are two methods to retrieve SGS motions from a benchmark simulation. The traditional method adopted by the turbulence and engineering communities defines the SGS field as
(Wyngaard 2010), where c(x, y) is the total field of a variable c before filtering at a location x and y, is a filter-scale variable (i.e., a CRM-RS variable as implied in this paper), and G is a smooth-filter function (e.g., Gaussian or top hat). We call this smooth-filtering method. With smooth filtering, any SGS flux can be further decomposed into three components: Leonard, cross, and Reynolds terms [for their definitions, see Leonard (1974) and Germano (1986)]. The Leonard term depends only on the filter-scale field and hence can be analytically approximated in terms of filter-scale variables by Taylor expansion or a defilter process. This expression combined with the classical Smagorinsky model leads to the so-called mixed schemes (Zhou et al. 2001; Chow et al. 2005; Moeng et al. 2010). We will discuss this scheme later in the appendix.
On the other hand, from the point of view of meteorological modeling, SGS motions are often considered as deviations (or fluctuations) of the total field from grid-averaged (i.e., resolvable scale) value, which means
Here is a grid-averaged RS variable at the center of the averaging subdomain. We call this discretized-filtering method. With this method, the Leonard and the cross terms are exactly zero because a subdomain average of all subdomain fluctuations is zero. Thus, the two filtering procedures yield different SGS statistics.
Using the discretized-filtering method, Moeng and Arakawa (2012) found that a SGS updraft–downdraft representation (which will be briefly reviewed in section 2) can accurately represent the spatial distribution of the SGS moisture fluxes in the PBL for a typical CRM grid resolution. The spatial correlation between the SGS moisture fluxes diagnosed from their updraft–downdraft scheme and those retrieved directly from the benchmark SGS fluctuations is remarkably high, larger than 0.95, for the lower convection layer they examined.
The remarkable performance of the updraft–downdraft scheme prompts this study. The goal is to 1) expand the analysis of that updraft–downdraft scheme for the entire convection layer and also for heat and momentum transport, which will be shown in section 2; and 2) find a closure that can properly relate SGS updraft–downdraft properties to RS variables (in section 3). Section 4 discusses the energy transfer (or exchange) between CRM-RS and CRM-SGS and section 5 summarizes this work.
In this study, we continue to use the discretized filter because this definition of SGS is a better representative of SGS fluctuations in discretized atmospheric models like CRMs. We use the last recorded Giga-LES flow field (at hour 24 of the simulation) to examine the updraft–downdraft scheme and to come up with a closure scheme. We then choose a different record (in this case, the flow field at hour 20) to test the scheme for energy transfer between RS and SGS against the Giga-LES data. The results shown here are from a filter width of 4 km where the CRM-RS field is calculated by averaging over every 4 km × 4 km horizontal subdomain of the Giga-LES (over every 40 × 40 Giga-LES grid points) and the SGS field is retrieved using Eq. (2) for each subdomain centered at (x0, y0). Thus, the SGS motions defined here cover the range of scales from 4 km down to a few hundred meters.1
2. An updraft–downdraft scheme of SGS fluxes for the entire convection system
We first briefly review the updraft–downdraft scheme proposed by Moeng and Arakawa (2012). Their scheme is based on mass-flux representation and expresses SGS moisture fluxes for every CRM grid cell as
where the SGS mean updraft is computed as
and the SGS mean downdraft as
performed at every 4 km × 4 km subdomain; ws is the SGS vertical-velocity fluctuations; and Nu and Nd are the numbers of grid points of SGS updrafts and downdrafts, respectively, within each subdomain. The SGS mean moisture mixing ratio averaged over all updrafts and downdrafts, respectively, are
The SGS updrafts and downdrafts are defined regardless of cloudy or clear air, which is an important difference from the conventional mass-flux modeling where “updrafts” typically signify “cloud cores.” Our SGS motions do not distinguish between cloudy and environmental air.
Moeng and Arakawa (2012) found that the Eq. (3) can well represent the spatial distributions of τwq in the PBL (below 1 km) with A1 set to 0.4. The point-by-point spatial correlations of SGS fluxes computed from Eq. (3) with those retrieved directly from the Giga-LES are larger than 0.95 in the PBL. This remarkable success of Eq. (3) in representing local SGS fluxes in CRMs prompts us to further examine its use for the entire convection system and also for the liquid-water potential temperature flux and momentum fluxes τuw, τυw, and τww. We also need to check if this A1 value applies to the layer above the PBL and to the other SGS fluxes.
Before investigating the SGS scheme, we show the relative contribution of SGS fluxes to the total fluxes of q, θl, u, υ, and w in Fig. 1, for the filter width of 4 km. The solid curves show the total vertical fluxes, the dotted curves display the RS contribution, and dashed curves are the SGS contribution—all averaged over the entire convection domain. In the upper part of the deep convection layer, the RS contribution equals or larger than that of SGS. Below 5–6 km, the SGS fluxes dominate the transport process for q, θl, and w.
The thin dashed curves in Fig. 1 are computed from Eq. (3) with A1 set to 0.3 (the value will be discussed below) and with all SGS updraft–downdraft properties diagnosed from the Giga-LES via Eqs. (4), (5), (6), and (7). These updraft–downdraft fluxes agree well with the SGS fluxes (thick dashed curves) retrieved from the Giga-LES SGS field.
The horizontal-mean SGS fluxes shown in Fig. 1 are not the SGS statistics needed in the governing equations of CRMs. Local SGS fluxes, which differ by several orders of magnitudes at various locations in a convection system (e.g., Zipser 1977; Moeng et al. 2010), are the statistics needed in CRMs. Following Moeng and Arakawa (2012), we examine the performance of the updraft–downdraft scheme for local SGS fluxes of moisture, heat, and momentum by correlating the spatial distributions of the SGS fluxes between those retrieved from the benchmark and those diagnosed from the updraft–downdraft scheme. Figure 2a shows the correlation coefficients larger than 0.9 throughout the entire deep convection layer. One plausible reason for this remarkable outcome is the following: at a CRM grid cell where SGS ws and cs are well correlated (i.e., fluctuating together), its SGS flux τwc is expected to be large. At the same time, cu − cd, which is defined in Eqs. (6) and (7) according to the signs of ws, is expected to be maximized. Thus, larger SGS fluxes tend to occur at a CRM grid where the product of wu − wd and cu − cd is maximized.
Using a linear regression analysis, we estimate the value A1 as a function of height shown in Fig. 2b. The value A1 is nearly 0.3 throughout the deep convection layer, although slightly larger for the moisture flux near the surface (Moeng and Arakawa 2012). Setting A1 = 0.3 in Eq. (3) brings the ratios of the standard deviations close to 1 (Fig. 2c) for almost the entire convection layer.
The above analysis uses the updraft–downdraft properties diagnosed from the Giga-LES. For practical use in CRMs, these SGS updraft–downdraft properties need to be related to CRM-RS variables. We will address this closure problem next.
3. Relating SGS updraft–downdraft properties to CRM-RS variables
Most mass-flux types of models build their closure on the thermodynamic properties of just updrafts, by the use of an entraining-plume concept (e.g., Arakawa and Schubert 1974; Siebesma and Cuijpers 1995). For our purpose we find it easier to build a closure on the properties of updraft–downdraft differences. After all, it is the differences that are needed in Eq. (3). Also, dealing with their differences requires no modeling of the lateral mass exchange (entrainment/detrainment) between updrafts and downdrafts, which is a complicated process and difficult to model.
The SGS updraft–downdraft difference results primarily from the fluctuations of the largest SGS eddies (i.e., the most energy-containing SGS eddies). And because the largest SGS motions are spectrally linked to the smallest RS eddies (Moeng et al. 2010), it is reasonable to assume that horizontal fluctuations of the smallest RS motions and the largest SGS motions are strongly communicated. Based on this idea, we test the following closure relation:
where Δx and Δy are the differences of the RS variables w or c across the adjacent CRM grids in x and y; they represent the horizontal fluctuations of the smallest RS motions in this Giga-LES analysis. The right-hand side of Eq. (8) involves only CRM-RS variables.
Figures 3 and 4 display the horizontal distributions of the quantities computed from both sides of Eq. (8) for c = q, θl, u, υ, and w (with A2 set to 1.5, which will be discussed later) at z ~ 3 km. The spatial distributions are dominated by extreme events (where deep convection is active) of strong SGS fluctuations. SGS updraft–downdraft fluctuations are stronger where horizontal differentials of RS are larger.
For statistical measure, we compute their point-by-point spatial correlation coefficients, the slopes of linear regressions, and the ratio of standard deviations at all heights, shown in Fig. 5. The correlation coefficients are about 0.7 for the entire convection layer, except near the surface. The slopes of linear regressions (Fig. 5b) indicate that, in the bulk of the deep convection layer, A2 fluctuates about 1; it is larger for q, θ, w than for u and υ. Here for the purpose of a simple parameterization, we will assume A2 = 1.5 for all variables. With this constant, the ratios of standard deviations (Fig. 5c) lie between 0.5 and 1.5 below z ~ 7 km (which is about the height of the anvil-cloud base), though smaller above. (Near the surface, small-scale turbulent motions dominate the transport process; those small-scale eddies are either unresolved or poorly resolved in this Giga-LES. Thus, the analysis presented in this paper does not apply to the parameterization issue for small-scale turbulence near the surface.) Using this A2 value, Figs. 3 and 4 demonstrate the capability of the simple closure scheme in Eq. (8) in representing the spatial distributions of SGS updraft–downdraft differences, particularly for the extreme events. The performance in representing their magnitudes is less satisfactory, which may be improved by further tuning the A2 value.
This expression differs from the conventional eddy-diffusivity approach (i.e., downgradient K models) where vertical fluxes are expressed in term of vertical gradients of RS variables.
Expressing SGS vertical fluxes in terms of horizontal gradients of RS variables may seem unconventional but it is not a new concept. The mixed scheme (or gradient model) proposed in the turbulence community (Bardina et al. 1980; Piomelli et al. 1988; Zang et al. 1993; Chow et al. 2005; Lu and Porte Agel 2010) and that proposed by Moeng et al. (2010) for deep convection also includes an important term similar to Eq. (9). That term, often referred to as the Leonard term, is derivable only when a smooth-filter method in Eq. (1) is used to define SGS field. With the discretized filtering adopted in this study, the Leonard (or cross) term is exactly zero. Although the definitions of SGS fluctuations differ between the two filtering procedures, the net behavior of SGS transport should remain similar. In the appendix, we compare Eq. (9) to the mixed scheme described in Moeng et al. (2010) and offer a new interpretation to that previous study.
The SGS scheme given in Eq. (9) is not an eddy-viscosity type. Previous studies (Zang et al. 1993; Meneveau and Katz 2000) have suggested that non-eddy-viscosity type of SGS schemes cannot properly transfer energy to dissipation in LESs. Next we will examine the SGS scheme in representing such energy transfer in CRMs, using a different flow realization.
4. Energy transfer between CRM-resolvable and subgrid scales
Energy transfer between the RS and SGS components shows up in the resolvable-scale kinetic-energy budget. The SGS effect in the governing equation for CRM resolvable velocity is generally expressed as ∂τij/∂xj, where represents the SGS fluxes or stresses. Thus, in the CRM-RS kinetic energy budget, the SGS effect appears as , which can be further split into a diffusion term and a dissipation term: , where is the strain rate tensor defined as . The term τijSij thus represents energy transfer between CRM-RS and CRM-SGS components (Piomelli et al. 1991). Positive τijSij signifies energy source for the resolvable-scale kinetic energy with energy coming from SGS (i.e., upscale or backscattering transfer). Negative value indicates that energy is transferred downscale (forward) from RS into SGS, and thus representing SGS dissipation.
The top panel in Fig. 6 shows the horizontal distribution of τijSij at z ~ 3 km retrieved from the Giga-LES flow field at hour 20 of the simulation, which is a different flow realization from that used to determine the parameters A1 and A2. Warm colors show regions where energy flows upscale from SGS to RS, while cold colors indicate the other way around. At every CRM grid, chances of kinetic energy going up or downscale are nearly the same.
To see how the SGS scheme in Eq. (9) represents this energy transfer distribution, we expand its application to all component of wind stresses, as
Figure 7 shows the correlation coefficients between τij computed from Eq. (10) and those retrieved directly from Giga-LES. The coefficients are in the 0.6–0.9 range (except near the surface for τuw), indicating that the expanded application to τuu, τυυ, and τuυ is reasonable. Their scatterplots at z ~ 3 km are given in Fig. 8 where each symbol represents a 4 km × 4 km subdomain of the Giga-LES; the model represents the large SGS stress events better than the small stress events, although the slopes suggest that the A2 value may be tuned for various stress components. The scatterplot for τww looks the least satisfactorily, but Fig. 7 shows a correlation coefficient of about 0.8 for this component (short dashed curve) at z ~ 3 km. For the following demonstration, we keep A1 = 0.3 and A2 = 1.5.
Using the modeled τij to estimate the energy transfer rate yields a rather similar pattern compared to that retrieved from the Giga-LES, shown in the bottom panel in Fig. 6. This agreement is also demonstrated by the scatterplot given in Fig. 9. Most importantly, the sign of the transfer for the larger events is well captured. Note that the conventional Smagorinsky–Deardorff eddy-viscosity scheme can only produce negative τijSij; it is known as an “absolute dissipative” scheme (Piomelli et al. 1991; Sullivan et al. 2003). Thus, the SGS scheme given in Eq. (9) can represent the spatial distribution of energy transfer between CRM-RS and CRM-SGS much better than conventional eddy-viscosity models.
Unlike the horizontal velocity (where energy is concentrated at mesoscale), the vertical-velocity spectra peak at kilometer scales in the deep convection layer above the PBL (Moeng et al. 2010). A typical CRM thus splits its RS and SGS components at or near the w-spectral peak scale. A proper representation of the w component of energy transfer between RS and SGS in CRMs is thus crucial. The SGS term in the budget is . We compute these quantities using τuw, τυw, and τww retrieved from the Giga-LES and those computed from Eq. (10), and compare them in Fig. 10. The magnitude of the transfer rate is somewhat underestimated, but the sign is well captured.
The vertical profiles of correlation coefficients and the ratios of standard deviations between the retrieved and the modeled energy transfer rate at all heights are given in Fig. 11, for the total kinetic energy (solid curves) and the w-energy only (dotted curves). Below z = 9 km, the correlation coefficients are larger than 0.5. The ratio of standard deviations is close to 1 for total kinetic energy transfer rate, but is much larger for the w-energy transfer. Nevertheless, Fig. 10 shows that the sign of the transfer is still well captured, which is a much improvement over the conventional eddy-viscosity scheme.
5. Summary and discussion
Based on a large-domain large-eddy simulation (Giga-LES), we show that the updraft–downdraft model framework proposed by Moeng and Arakawa (2012) works extremely well in representing the SGS vertical fluxes of heat, moisture, and momentum in CRMs of deep convection. One plausible reason for this remarkable outcome is that well-correlated SGS fluctuations of ws and cs in a CRM grid cell not only lead to a larger SGS flux τwc but also maximize the difference between cu and cd. When ws and cs are poorly correlated, the SGS flux is negligibly small anyway. The point-by-point spatial correlation between the SGS fluxes and the updraft–downdraft differences (diagnosed from the Giga-LES) is larger than 0.9 throughout the entire deep convection layer. This prompts us to investigate a closure for the updraft–downdraft model framework for its potential use in CRMs.
We propose a closure scheme that relates SGS updraft–downdraft differences to horizontal gradients of CRM-RS variables, based on the following ideas: 1) SGS updraft–downdraft difference results primarily from the fluctuations of the largest (energy containing) SGS eddies, 2) the largest SGS eddies are spectrally linked to the smallest RS motions in CRMs, and thus 3) horizontal gradients of RS variables across a grid cell are likely to relate to the SGS updraft–downdraft differences in that grid cell. Point-by-point correlation coefficients between the quantities computed from the right- and left-hand sides of Eq. (8) are about 0.7 throughout the entire convection system, except near the surface (where the Giga-LES is too coarse to resolve dominant eddies).
This closure relation leads to a SGS scheme that represents CRM-SGS vertical fluxes by the form of horizontal gradients of CRM-RS variables, which is crucially different from the commonly used eddy-diffusivity models where vertical fluxes are related to vertical gradients of RS. Though unconventional, expressing SGS vertical fluxes in term of horizontal gradients of RS variables is not new; it has been advocated by the turbulence community to represent SGS stresses in LESs. Here we show that this unconventional scheme applies well to representing SGS fluxes in kilometer-scale CRMs of deep convection, as well as their kinetic energy transfer between RS and SGS components.
It is important to note that our analysis does not include motions smaller than a few hundred meters, which are underresolved in the Giga-LES. Also, this analysis does not consider numerical truncation that exists in actual CRM simulations. Like most SGS models, the proposed scheme relates SGS fluxes to local grid-scale resolvable variables, which strongly depend on numerical schemes adopted by CRMs. How the SGS scheme, numerical scheme, and near-grid-scale resolvable flow field respond to and affect each other in actual CRM simulations remains to be seen. It is particularly challenging for the vertical-velocity field where the energy peaks at or near CRM grid scales.
This material is based upon work supported by the National Science Foundation (NSF) Science and Technology Center for Multi-Scale Modeling of Atmospheric Processes (CMMAP) led by David Randall at the Colorado State University under Cooperative Agreement ATM-0425247. The benchmark LES was performed using the IBM BlueGene/L “New York Blue” supercomputer at the New York Center for Computational Sciences, which is a joint venture of Stony Brook University and Brookhaven National Laboratory. I would like to thank Marat Khairoutdinov (Stony Brook University) and Steve Krueger (University of Utah) for generating the benchmark simulation data, and Akio Arakawa (UCLA), David Randall (CSU), and Peggy LeMone and Peter Sullivan (NCAR) for lengthy discussions.
Comparing the Updraft–Downdraft and the Mixed Schemes
For discussion purposes, we briefly review the mixed scheme described in Moeng et al. (2010). The SGS vertical fluxes are expressed as
The first term appears to be a result of the use of a smooth filter in Eq. (1). With a smooth filter, SGS fluxes can be further decomposed into three terms: Leonard, cross, and Reynolds terms. The Leonard term involves only the resolvable-scale variables and with Taylor expansion can be approximately expressed in the form of the first term. The factor of 2 in front of the first term represents the cross term because Moeng et al. (2010) showed that the cross term has a spatial distribution and magnitudes similar to the Leonard term. This Leonard–cross term is similar to the SGS scheme shown in Eq. (9); they both express the SGS vertical fluxes in terms of the horizontal gradients of the RS variables [although continuous derivatives are used in Eq. (A1) because of its smooth RS and SGS flow fields].
The major difference is the coefficients of these horizontal gradient terms. Compared to A1 × A2 ~ 0.45 in Eq. (9) for τwq, the coefficient of the Leonard–cross term in Eq. (A1) is only ⅙, which is smaller by a factor of about 3. To look for this difference, we reexamine the spatial distribution and the magnitudes of the Leonard, cross, and Reynolds terms shown in Fig. 4 of Moeng et al. (2010). We find that in fact the Reynolds term also exhibits a spatial distribution very similar to the other two terms, and may as well be expressed in the same form as the Leonard–cross term. We also notice that the Reynolds term has a magnitude about 3 times larger than the Leonard or cross term. Adding the Reynolds term, the first term in Eq. (A1) then yields a coefficient that is close to 0.45.
The National Center for Atmospheric Research is sponsored by the National Science Foundation.
The numerical code SAM (the CSU's System for Atmospheric Model) that generates the Giga-LES uses a centered-finite differencing scheme for momentum equations, and hence tends to generate grid-scale noises. Thus, before analysis, we first perform 4 × 4 gridpoint averaging on the Giga-LES field to smooth out small-scale numerical noises. This presmoothing does affect the SGS updraft–downdraft mean properties, though not by a significant amount. Also, SAM adopts a staggered grid for various variables. For simplicity, we average all variables onto the w grid points first before performing the following analysis.