Parameterizing deep convection using the assumed probability density function method

Due to their coarse horizontal resolution, presentday climate models must parameterize deep convection. This paper presents single-column simulations of deep convection using a probability density function (PDF) parameterization. The PDF parameterization predicts the PDF of subgrid variability of turbulence, clouds, and hydrometeors. That variability is interfaced to a prognostic microphysics scheme using a Monte Carlo sampling method. The PDF parameterization is used to simulate tropical deep convection, the transition from shallow to deep convection over land, and midlatitude deep convection. These parameterized single-column simulations are compared with 3-D reference simulations. The agreement is satisfactory except when the convective forcing is weak. The same PDF parameterization is also used to simulate shallow cumulus and stratocumulus layers. The PDF method is sufficiently general to adequately simulate these five deep, shallow, and stratiform cloud cases with a single equation set. This raises hopes that it may be possible in the future, with further refinements at coarse time step and grid spacing, to parameterize all cloud types in a large-scale model in a unified way.


Introduction
Deep convective clouds are an integral part of the climate system, yet representing these clouds accurately in climate models is a challenge.Explicitly resolving these clouds in decades-long climate simulations is not yet computationally feasible; thus, it is beneficial to parameterize deep convection (e.g., Molinari and Dudek, 1992;Weisman et al., 1997;Warner and Hsu, 2000;Sato et al., 2009;Noda et al., 2010).Deep convective clouds are difficult to parameterize, however, partly because they produce copious ice and precipitation, which interact strongly with storm dynamics.For instance, convective rain may partially evaporate before reaching the ground, creating cold pools that spawn further convection (e.g., Liu et al., 1997;Tompkins, 2001;Khairoutdinov and Randall, 2006).
Deep convection has often been parameterized by the use of mass-flux schemes (e.g., Arakawa and Schubert, 1974;Kain and Fritsch, 1990;Donner, 1993;Zhang and Mc-Farlane, 1995;Arakawa, 2004;Kain, 2004;Gerard, 2007;Arakawa et al., 2011;Mapes and Neale, 2011;Bechtold et al., 2014).However, parameterization of deep convective clouds and turbulence may be aided by the use of the assumed probability density function (PDF) method.In the assumed PDF method, several lower-order subgrid moments are prognosed at each time step (e.g., Sommeria and Deardorff, 1977;Mellor, 1977;Smith, 1990;Lewellen and Yoh, 1993;Bony and Emanuel, 2001;Tompkins, 2002).The lower-order moments, along with an assumption about PDF shape, are used to estimate a joint PDF of subgrid-scale (SGS) moisture, temperature, vertical velocity, and hydrometeors.The subgrid PDF is then used to calculate various unclosed higher-order moments that appear in the equation set, thereby allowing the solution to be advanced another time step.
PDF parameterizations possess several advantages for the purpose of parameterizing deep convection.First, the use of a PDF that includes hydrometeor species and vertical velocity allows for a detailed representation of microphysics and dynamics, and of the coupling between them.For instance, a Published by Copernicus Publications on behalf of the European Geosciences Union.

R. L. Storer et al.: Parameterizing deep convection
PDF parameterization that prognoses lower-order moments is compatible with prognostic microphysics.The use of a prognostic microphysics scheme, in turn, allows for precipitation to be preserved between computational time steps, providing a better representation of hydrometeor collection (e.g., Posselt and Lohmann, 2009;Gettelman and Morrison, 2014).Representation of collection also benefits from the ability of PDF parameterizations to prescribe the correlations among hydrometeors, thereby estimating the degree of collocation among the hydrometeor types.Another benefit of a PDF parameterization is the detailed and rigorous treatment of the distribution of SGS vertical velocity, which in turn is important for parameterizing vertical turbulent transport.
Despite these advantages for parameterizing deep clouds, to date PDF parameterizations have been applied only to shallow clouds, except in higher-resolution cloud-resolving models (Cheng and Xu, 2006;Bogenschutz and Krueger, 2013) and except in Davies et al. (2013).Davies et al. (2013) includes a single-column model (SCM) simulation of a tropical deep convective case by the Cloud Layers Unified by Binormals (CLUBB) PDF parameterization.However, that study was an intercomparison of several SCMs and hence did not include details of the formulation or evaluation of CLUBB.CLUBB performed comparably to the other SCMs in Davies et al. (2013), but here we take the CLUBB parameterization and extend it in ways designed to better represent deep convection.In particular, the representation of the subgrid distribution of ice and precipitation is improved.
The modified version of CLUBB is then used to perform single-column simulations of three deep convective cases: a tropical case, a case of transition from shallow to deep convection, and a midlatitude case.To demonstrate that the modified version of CLUBB can still simulate shallow cloud cases, we also simulate a shallow cumulus layer and stratocumulus layer.To assess the quality of the single-column CLUBB simulations, we compare the deep convective cases to cloud-resolving simulations and the shallow cases to largeeddy simulations.
The remainder of this paper is organized as follows.In Sect.2, we will describe the single-column model CLUBB and the recent changes to it implemented in order to obtain improved simulations of deep convective clouds.In Sect.3, results from CLUBB for three deep convective cases will be compared to cloud-resolving model simulations.In Sect.4, we will show results of sensitivity tests to basic model setup, including reducing the vertical and time resolutions such as might be required to interface CLUBB with a climate model.In Sect.5, we will show results of two simulations of shallow clouds, demonstrating that this model framework can be used in simulations of both deep and shallow cloud types.In Sect.6, we will finish with discussion and conclusions.

Description of single-column model: CLUBB-SILHS
CLUBB is a single-column model of clouds and turbulence in the atmosphere (Golaz et al., 2002;Larson and Golaz, 2005;Larson and Griffin, 2013;Griffin and Larson, 2013).
CLUBB is designed to model a variety of cloud types and turbulence regimes in a unified way with a single equation set.CLUBB parameterizes both cumulus and stratocumulus clouds, and both dry convective and stably stratified boundary layers.The equation set consists of prognostic equations for various subgrid-scale higher-order moments.The moments include variances of vertical velocity (w 2 ), total water mixing ratio (r t 2 ), and liquid cloud water potential temperature (θ l 2 ).It also includes subgrid vertical turbulent fluxes of total water (w r t ) and liquid cloud water potential temperature, w θ l .(Here an overbar denotes a grid box average, and a prime denotes a subgrid deviation from a grid box average.)CLUBB's equation set is summarized in Appendix A. CLUBB has been coupled to the double-moment microphysics scheme of Morrison et al. (2009).The coupled model prognoses grid-mean number concentration and mixing ratios of rain, snow, graupel and cloud ice.The spatial variances of these hydrometeors are assumed to be proportional to the grid mean values squared.CLUBB's equation set contains various higher-order moments that are unclosed.CLUBB closes many of these moments using the assumed PDF method.That is, CLUBB assumes a functional form of the subgrid PDF and integrates over it in order to estimate the unclosed moments.CLUBB's subgrid PDF is a single, multivariate PDF with a double Gaussian (sum of two normals) distribution for vertical velocity w, total water mixing ratio r t , and liquid water potential temperature θ l (Larson et al., 2001) and a delta functionlognormal distribution for the hydrometeor variates (Griffin and Larson, 2014).The delta-lognormal PDF adds a delta function representing precipitation-free areas to the lognormal distribution of Larson and Griffin (2013).It is described further in Sect.2.2.1 below.
In order to drive microphysical processes, CLUBB has the option to call a microphysics scheme using a variant of traditional Monte Carlo sampling (McKay et al., 1979;Owen, 2003).This sampler is called the Subgrid Importance Latin Hypercube Sampler (SILHS; Larson et al., 2005;Larson and Schanen, 2013).SILHS draws one or more sample points per grid box and time step from CLUBB's joint subgrid PDF.SILHS collects a sample from each grid level into a subcolumn by assuming that correlation drops off exponentially with vertical distance and by using CLUBB's turbulent length scale to estimate the e-folding distance (Larson and Schanen, 2013).Each subcolumn is then fed into a microphysics scheme, one at a time.For each subcolumn, the microphysics scheme produces outputs, e.g., autoconversion rate.Then the resulting ensemble of autoconversion rates is averaged, yielding an estimate of the grid box mean autoconversion rate.
1. Construct the multivariate PDF of subgrid-scale variability (performed by CLUBB).CLUBB predicts a single, joint PDF for each grid box and time step that contains the following variates: heat content, moisture content, vertical velocity, and all relevant hydrometeor mixing ratios and number concentrations.

Draw subcolumns from the subgrid PDF (performed by SILHS).
A sample point is drawn from each grid level.These sample points, taken together, form a single vertically correlated profile, i.e., "subcolumn".The sampling is repeated until a set of subcolumns is constructed that collectively approximates the statistical properties of the grid column.
3. Feed subcolumns into physical parameterizations and compute process rate tendencies (performed by a microphysics scheme).Each subcolumn is fed into a microphysics scheme.For each subcolumn, the microphysics scheme computes a separate tendency.
4. Average microphysics tendencies from each subcolumn to form a grid box average profile.The set of tendencies is averaged to form a grid box mean tendency, which can then be fed back into the appropriate host model equations.The average may be a weighted average, with each subcolumn corresponding to different-sized area within the grid column.
The computational cost of CLUBB-SILHS is acceptable.CLUBB adds 20 % to the computational cost of version 5 of the Community Atmosphere Model (CAM5, Bogenschutz et al., 2013).SILHS's computational cost, when two sample points are used, is similar to that of CLUBB (Larson and Schanen, 2013).
Both CLUBB and SILHS are contained in a single SVN (Subversion) code repository that is freely available for noncommercial use at http://clubb.larson-group.com.After acquiring a free account the simulations presented here can be reproduced by checking out the branch located at http://carson.math.uwm.edu/repos/clubb_repos/branches/storer_et_al_paper_clubb/, compiling the source code, and running ./run_scm_all.bash.Complete details of all algorithms are available in the source code in the SVN repository, along with numerous code comments and a readme file that serves as a user manual.SAM (System for Atmospheric Modeling) output used for plotting and for prescribing radiative forcing can be found at http://www.larson-group.com/storer/.

Coupling subgrid variability to liquid and ice microphysics
In order to better simulate deep convection using CLUBB, we have generalized the treatment of subgrid variability.In particular, we have generalized CLUBB's multivariate PDF to include ice, incorporated the effects of ice and rain on subgrid buoyancy production, improved the vertical transport of hydrometeors, and allowed changes in the microphysics to affect scalar fluxes and variances.

Subgrid PDF of hydrometeors
CLUBB represents vertical velocity, liquid water potential temperature, and total water mixing ratio using a double Gaussian PDF.However, a double Gaussian is unbounded and hence is not well suited to representing hydrometeors.In this context, hydrometeors include the mixing ratios and number concentrations of cloud ice, snow, graupel, and rain.These quantities are nonnegative but often have small values.
In order to represent hydrometeors realistically and fully, a subgrid PDF should have two features.First, it should be multivariate, so that it can represent the collection of one hydrometeor species by another.Second, the subgrid PDF should include a hydrometeor-free region, so that hydrometeors do not suffer excessive evaporation or sublimation in otherwise clear regions.
The marginal subgrid PDF shape that we choose for the hydrometeor variates is a sum of a delta function and a multivariate lognormal.The delta function represents the portion of a grid box that is devoid of all hydrometeors (except liquid cloud droplets).The multivariate lognormal represents the hydrometeor-filled portion, i.e., the precipitation fraction.Details of the formulation of the lognormal are provided by the source code and in Larson and Griffin (2013) and in Griffin and Larson (2013).Although those papers used a multivariate lognormal to represent only the number concentration and mixing ratio of rain, here the multivariate lognormal is extended straightforwardly to represent the withinhydrometeor number concentrations and mixing ratios of all hydrometeors except liquid cloud water, which is diagnosed via saturation adjustment from total water.
At altitudes above the freezing level, the precipitation fraction is set equal to the fraction of a grid box that is supersaturated with respect to ice.Below the freezing level, the precipitation fraction is determined by a calculation similar to that of Morrison and Gettelman (2008).More details of the formulation are provided in Griffin and Larson (2014) and the source code.
The assumed PDF method requires that variances of each variate be provided.The within-hydrometeor standard deviation (σ i ) of the ith hydrometeor mixing ratio or number concentration, with mean µ i , is calculated as  where V i is the prescribed constant ratio between the variance and the mean squared of each variate.Values for V i used in the single-column simulations are provided in Table 1.
For multivariate PDFs, correlations among hydrometeor species must also be provided.In these simulations, the correlations are assumed to be constant within cloud throughout the course of the simulations.In order to assign values for these correlations, we examined the range of values found within cloud-resolving model simulations (which vary in space and time), and prescribed correlations that were representative.The correlations used for cloudy levels are shown in Table 2.As with the variances, the correlations are withinhydrometeor values.In future work, instead of prescribing correlations, one could attempt to diagnose correlations, following, e.g., Larson et al. (2011).

Incorporating the effects of latent heating of ice and rain on subgrid turbulence
Turbulence kinetic energy and turbulent fluxes of heat and moisture are influenced by the buoyancy of air parcels in local updrafts and downdrafts.The buoyancy variable in CLUBB is the virtual potential temperature θ v , where θ v includes the effects on buoyancy of water vapor and cloud droplet loading.In CLUBB, buoyancy affects turbulence through four buoyancy generation terms that appear on the right-hand side of prognostic equations for moments (see Eqs. A5, A9, A10, A11).If g/θ 0 denotes the acceleration due to gravity divided by a constant reference temperature, then (g/θ 0 )w θ v is a buoyancy term that generates (that is, contributes to the tendency of) w 2 , (g/θ 0 )w 2 θ v generates w 3 , (g/θ 0 )θ l θ v generates w θ l , and (g/θ 0 )r t θ v generates w r t .These buoyancy terms, in turn, are influenced by latent heating associated with phase change of hydrometeors and gravitational loading due to hydrometeors.
In prior versions of CLUBB, the only source of latent heating for the higher-order moments was condensation and evaporation of liquid cloud droplets, and the only source of loading was cloud liquid water mixing ratio.In other words, in the aforementioned four generation terms, the buoyancy perturbation θ v accounted only for variability due to the latent heating and water loading associated with cloud liquid water mixing ratio.(However, grid-mean values of θ l in prior versions of CLUBB were subject to the latent heating associated with rain and ice hydrometeors, where relevant, through microphysical calculations.) For the present simulations, we have incorporated a crude method to account for the latent heating and water loading associated with rain and ice hydrometeors in the four aforementioned buoyancy generation terms.Namely, we assume that for the purpose of calculating the four buoyancy generation terms, rain and all ice hydrometeors are assumed to be perfectly collocated in space with cloud liquid water mixing ratio.With this simple assumption, CLUBB's calculation of the four buoyancy generation terms remains unchanged, except that total condensate -including cloud liquid water, rain and ice mixing ratios -is input into the calculation in place of cloud liquid water.Conservation of heat and moisture for grid means still holds because this change only applies to the buoyancy generation terms.In the future, however, it would be desirable to examine the effects of this assumption.

Parameterizing transport of hydrometeors
Turbulent updrafts and downdrafts transport hydrometeors in the vertical, thereby acting to broaden the vertical extent of the hydrometeor profiles.In CLUBB, the tendency of turbulence to broaden the hydrometeor profiles is modeled by eddy diffusion.Although it is not obvious a priori that an eddy diffusion model is appropriate for cumulus layers, our large-eddy simulations indicate that eddy diffusion parameterizes the sign of hydrometeor transport satisfactorily in cumulus layers, but that the magnitude of the eddy diffusion is 1 or 2 orders of magnitude larger in cumulus than in stratocumulus layers (not shown).
CLUBB models the transport of a generic hydrometeor, r x , as where w denotes the vertical velocity, c = 0.75 a dimensionless constant, and z the altitude.The eddy diffusivity for hydrometeors, K, is given by where e is the subgrid turbulence kinetic energy and L is CLUBB's turbulent mixing length (Golaz et al., 2002;Larson et al., 2012).The two last factors in parentheses are designed to increase K by up to a factor of 100 in cumulus layers while leaving K small in stratocumulus layers.These two factors rely on the fact that spatial variability of rain and skewness of vertical velocity, Sk w , are larger in cumulus layers than in stratocumulus layers.

Allowing microphysical sources to influence subgrid distribution shape
Microphysical processes influence not only grid box means, but also the shape of the subgrid PDF.For instance, microphysical processes such as autoconversion and accretion occur preferentially in moister portions of a grid box.Such processes, in conjunction with sedimentation of precipitation, tend to remove moisture from the moistest part of a grid box, thereby diminishing the variance of total water in cloud layers (Khairoutdinov and Randall, 2002).In addition, evaporation of rain below cloud may at times occur in the coolest part of a grid box, thereby increasing variance of temperature below cloud.The increased subcloud temperature variance may, in principle, help initiate further convection.These effects are estimated using SILHS subcolumns.Each subcolumn is fed into the microphysics separately, and each produces a separate microphysical update to r t and θ l .Therefore, the sample variance among the set of subcolumns is different before and after the microphysics is computed.This change in variance can be added as a source term in CLUBB's prognostic equation for the scalar variance of r t or θ l .A similar calculation is made for the scalar fluxes of r t and θ l .
For a generic variable x, SILHS estimates the microphysical source as where Var(x) is the sample variance of x among subcolumns, which is calculated before and after the call to the microphysics, t denotes time, and t is the time between micro-physics calculations.Here, where n is the number of subcolumns, x i is the value of x in the ith subcolumn, x is the average over all subcolumns, and p i is the probability of choosing subcolumn i.The probability p i equals 1/n in the case of equally weighted subcolumns, and the probability weights must sum to 1.An analogous calculation is performed in order to update covariances between two variates.In our simulations, the following prognostic moments are updated by microphysical source terms: r t 2 , θ l 2 , r t θ l , w r t , and w θ l .The prognostic equations for moments that involve only w, such as w 2 and w 3 , do not contain microphysical source terms, even in theory.

Model configurations and case descriptions
The following analysis compares two models, a singlecolumn model (CLUBB-SILHS), and a 3-D cloudresolving/large-eddy model (SAM; Khairoutdinov and Randall, 2003) that provides reference simulations.
In order to assess the feasibility of unified parameterization, CLUBB-SILHS's configuration of all shallow and deep cloud layers is identical, other than case-specific options for droplet number concentration and radiative transfer that are described below.For all cloud cases, CLUBB-SILHS uses a 1 min computational time step and a 128-level stretched vertical grid, with a vertical grid spacing of approximately 100 m at an altitude of 1000 m.CLUBB-SILHS's microphysics is identical to the Morrison microphysics option in SAM, which is a two-moment scheme that predicts cloud wa-  ter, rain, cloud ice, snow, and graupel (Morrison et al., 2009).SILHS generates 16 sample subcolumns per time step.We compare each CLUBB-SILHS single-column simulation with a three-dimensional simulation performed by SAM.In order to help isolate model errors in CLUBB-SILHS, we configure SAM and CLUBB-SILHS identically in a number of aspects.Namely, SAM uses the same two-moment microphysics scheme, the same value of prescribed within-cloud cloud droplet number concentration, and the same treatment of radiative transfer.SAM is set up similarly in shallow and deep cases, except that SAM uses smaller grid spacing and time step for the shallow cases.More SAM options for all cases are described below in Table 3.

Deep convective cases
We choose to simulate three deep convective cases that have been studied in prior model intercomparisons, in part because doing so allows us to compare our 3-D simulations with those intercomparison results, building confidence in our 3-D simulations.We configure the deep convective simulations as per previous model intercomparisons of those cases.SAM simulations examined here are comparable to previous simulations and observations in general characteristics such as timing, precipitation, and liquid water path.We conclude, therefore, that they provide good reference simulations for evaluation of CLUBB-SILHS.For deep convective cases, both SAM and CLUBB-SILHS prescribe droplet number concentration to be 100 cm −3 , and both use the same 128-level vertical grid.SAM uses 1 km grid spacing in the horizontal.
The first deep convective case we present is from the Tropical Warm Pool International Cloud Experiment (TWP-ICE), which took place near Darwin, Australia, in early 2006 (May et al., 2008).The TWP-ICE case is a week-long simulation, beginning on 19 January and exhibiting multiple diurnal cycles of convective activity.Large-scale conditions are set up for the SAM and CLUBB simulations as in previous model intercomparison studies (Fridlind et al., 2012;Davies et al., 2013).Radiation in the SAM simulation is calculated by the Rapid Radiation Transfer Model (RRTM) (Iacono et al., 2008), and radiative heating values from SAM output are used to prescribe radiation in CLUBB-SILHS.The second case was taken from the Tropical Rainfall Measuring Mission Large-Scale Biosphere-Atmosphere (TRMM-LBA) experiment which took place in Brazil in early 1999 (Silva Dias et al., 2002).The LBA case occurred on 23 February, and consists of the transition from shallow to deep convection over land from early morning to early afternoon.Model configurations for SAM and CLUBB-SILHS follow a previous model intercomparison (Grabowski et al., 2006;Khairoutdinov and Randall, 2006).In this case, both SAM and CLUBB-SILHS use prescribed radiation, as specified in the intercomparison.
Third, we examine a midlatitude summertime deep convective case that occurred during the 1997 intensive observation period of the Atmospheric Radiation Measurement (ARM) program (ARM97).ARM97 is a 4-day simulation starting on 26 June that exhibits two weak precipitating convective events followed by a strong convective event.We configure the case as in a previous intercomparison (Xu et al., 2002;Khairoutdinov and Randall, 2003).As in TWP-ICE, SAM uses RRTM to compute radiation, and SAM output is used to prescribe radiation in the CLUBB-SILHS simulation.

Shallow cloud-layer cases
The shallow cases are configured as in previous intercomparison studies, including the appropriate options for radiation and cloud drop number concentration.The vertical grid spacing and time step for CLUBB-SILHS are identical to that used for the deep convective cases.SAM's configuration is summarized in Table 3.
The first shallow case we present is a 3-day simulation of drizzling trade-wind cumulus clouds observed during the Rain in shallow Cumulus over the Ocean (RICO) field campaign (Rauber et al., 2007).We configure RICO as in the intercomparison of vanZanten et al. (2011).As per the intercomparison specifications, we prescribe a cloud drop number concentration of 70 cm −3 in both CLUBB-SILHS and SAM, and radiation is turned off in both models.The second shallow case was taken from the Dynamics and Chemistry of Marine Stratocumulus (DYCOMS-II) field campaign (Stevens et al., 2003).We simulate the research flight 2 (RF02) case, which is a 6 h simulation of a nocturnal drizzling stratocumulus layer off the coast of California (Wyant et al., 2007).DYCOMS-II RF02 uses a cloud drop number concentration of 55 cm −3 in both CLUBB-SILHS and SAM simulations (Ackerman et al., 2009).Both models use a simplified analytic radiation model which ties radiative heating to the liquid water path (Stevens et al., 2005;Wyant et al., 2007;Larson et al., 2007).

Results: simulated time series and profiles
This section presents time series and profiles from our tropical deep convective case (TWP-ICE), our shallow-to-deep transitional case (LBA), and our midlatitude deep convective case (ARM97).For each of the three cases, we plot time series of the following vertically averaged quantities: liquid water path (LWP), rain water path (RWP), cloud ice water path (IWP), and snow water path (SWP).In addition, we plot profiles from selected time periods of liquid water potential temperature (θ l ), total water mixing ratio (vapor + cloud liq- uid water) (r t ), the fraction of a grid box that is occupied by liquid cloud water, and cloud liquid water mixing ratio (r c ).

Tropical deep convection: TWP-ICE
Time series of the four vertically integrated quantities from the TWP-ICE simulation are shown in Fig. 1.This was a convectively active period, and a strong diurnal signal can be seen in the time series.Because TWP-ICE is a tropical case, rainfall rate at the surface is constrained by prescribed largescale forcings.Instead, we show rain water path.Rain water path simulated by CLUBB-SILHS mimics that in SAM, as does snow water path.LWP in CLUBB-SILHS is too high compared to SAM, while cloud ice is too low, which may mean that the correlations between liquid and ice hydrometeors are not well tuned.Figure 2 depicts four profiles averaged over the heaviest precipitation event, which occurs during minutes 6000-8000.In this case, CLUBB-SILHS agrees well with SAM and exhibits only a few discrepancies.Liquid water potential temperature (θ l ) is slightly too cool in the mid and upper levels, and the liquid cloud fraction is slightly too large in the lower levels.Figure 3 shows mean profiles of the mixing ratios of hydrometeors, again averaged over the strongest convective episode.All of CLUBB-SILHS's profiles reproduce the corresponding SAM profiles results well.The CLUBB-SILHS simulation matches the CRM simulation well, aside from a lack of cloud in lower levels.

Transition from shallow to deep convection: LBA
LBA is a difficult case to simulate because it evolves substantially over a short (6 h) time period as it transitions from shallow to deep convection.Hence it is a challenge to simulate the timing of ice and precipitation formation.Time series from LBA are shown in Fig. 4. LWP in CLUBB-SILHS matches that in SAM well until the last hour of simulation time.In the last hour, LWP is depleted too much because rain forms more than an hour too late, and when it does form it produces excessive rain water path.The rain is delayed because precipitation efficiency is too low during the early period of weak forcing, possibly because of an inaccurate representation of the correlations among hydrometeors.The simulation overestimates cloud ice water path because the simulation produces too much cloud liquid aloft (Fig. 5).
Figure 5 shows profiles averaged over the last hour of the simulation, by which point the convection in SAM has reached a mature stage.CLUBB-SILHS's profiles of θ l and r t match those of SAM because none of the mean profiles has time to evolve during the 6 h simulation.However, CLUBB-SILHS produces too much cloud liquid in the upper levels.CLUBB-SILHS's hydrometeor profiles agree with those in SAM, with the exception of cloud ice mixing ratio (Fig. 6).Cloud ice is overpredicted because the excess liquid aloft leads to large rates of contact nucleation (not shown).The large numbers of nucleated ice crystals also lead to decreased snow production from autoconversion due to the smaller crystal sizes.

Midlatitude continental convection: ARM97
Time series plots for ARM97 are shown in Fig. 7.The timing and magnitude of LWP is reasonably captured in CLUBB-SILHS.In the first convective event, however, CLUBB-SILHS's cloud does not grow deep enough for the production of ice, and little rain is produced.The later convective episodes are better simulated, although CLUBB-SILHS's time series are too noisy.
Mean profiles for ARM97 are shown in Fig. 8, averaged over the third and strongest convective event (minutes 4320-5580).Similar to what was seen in TWP-ICE, θ l in CLUBB-SILHS is too cool aloft, and r t and liquid cloud fraction are too large through the lowest 6-7 km of the column.Hydrometeor profiles (Fig. 9) match fairly well to the SAM CRM (cloud resolving model) simulation, though CLUBB-SILHS does not produce enough ice aloft.The undesirably low precipitation efficiency seen in TWP-ICE is also evident in this case, as too much liquid is needed to produce the correct amount of rain.CLUBB-SILHS's simulation of the first event produces cloud liquid water but no precipitation.However, CLUBB-SILHS's simulation of the other two events match SAM well in timing and magnitude.

Examining the representation of dynamics
Deep convection can have important effects on large-scale dynamics through various feedbacks with microphysics (e.g., Liu et al., 1997;Tompkins, 2001;Khairoutdinov and Randall, 2006).As such, it is useful to examine the performance of CLUBB-SILHS with regard to dynamical fields such as vertical velocity.Figure 10 demonstrates some dynamic fields, the variance of vertical velocity and the turbulent flux of liquid water potential temperature and total water mixing ratio, in all three deep convective simulations.CLUBB-SILHS does a reasonable job of representing the variance of vertical velocity in all three cases.The turbulent flux of total water is quite similar in CLUBB-SILHS and SAM, though in ARM97 is a bit low.The turbulent flux of liquid water potential temperature is too low, particularly in TWP-ICE and ARM97.These cases are the most convectively active, so it is not surprising that the SCM would have the hardest time representing their dynamics.It is worth considering ways to improve the dynamics in future work, particularly by improving the coupling between microphysics and dynamics.Adding complexity to the formulation described in Sect.2.2.2 may be one way to accomplish this; however, currently, CLUBB-SILHS does a reasonable job in simulating the subgrid-scale dynamic variability.

Effects of our new methodology on deep convective simulations
Two of the changes to CLUBB described in Sect.2.2 have significant effects on the simulation of deep convection.First, improving the PDF for hydrometeors by including a nonunity precipitation fraction (Sect.2.2.1) -that is, allowing a region free of hydrometeors (cloud ice, snow, graupel, and rain) -improved LBA (Fig. 11).With the previous formulation, liquid cloud water increases to overly large values by the end of the simulation, and the excessive cloud water nucleates to form excessive cloud ice.However, no rain forms during the 6 h simulation.The addition of a hydrometeor-free region improves (increases) the precipitation efficiency of CLUBB-SILHS, allowing rain to form and remove excessive liquid cloud water.
Second, it turns out to be important for the simulations of deep convective cases to have strong turbulent transport of hydrometeors.The changes introduced to the eddy diffusivity calculation in Sect.2.2.3 allow for stronger values in convective cases (while not degrading shallow cases, as shown in Sect.6).Larger values of eddy diffusivity smooth the profiles of hydrometeors and transport hydrometeors farther upward, leading to closer agreement with SAM (see Fig. 12).The other two methodological changes -adding latent heating to subgrid moments (Sect.2.2.2) and adding microphysical source terms to prognosed second-order moments (Sect.2.2.4) -have more minor effects on the simulations.

Sensitivity to time step, number of samples, and vertical grid spacing
Aside from the errors explicitly mentioned above, the CLUBB-SILHS results presented here agree well with SAM.However, the configuration of CLUBB-SILHS used -with its fine vertical resolution, short time step, and numerous sample points -is computationally expensive.In this section, we present results that use configurations that are more computationally affordable.Simulations presented in prior sections use a 1 min time step.Figure 13 plots the effects on the deep convective simulations of increasing the time step to 5 min; the two lines plotted represent identical simulations except for the time step.A profile of liquid cloud fraction is shown alongside a time series of snow water path for each case in order to illustrate the changes that occur when changing the time step.The timing of convective events is little affected.However, the cases are too convectively active, leading to too much liquid cloud For each case, the variance of vertical velocity is shown on the left, the turbulent flux of liquid water potential temperature is in the middle, and the turbulent flux of total water mixing ratio is on the right (averaged over the same time periods as in previous figures).CLUBB-SILHS matches the CRM well in variance of vertical velocity, but underestimates the magnitude of some of the turbulent fluxes, particularly in the more convectively active cases.
fraction aloft.This may be partly related to an underprediction of subcloud evaporative cooling that occurs when prognostic microphysics schemes with sequential time splitting of sedimentation allow hydrometeors to sediment long distances between calls to the microphysical sources and sinks, such as evaporation (H.Morrison, personal communication, March 2014).
It is possible to use CLUBB-SILHS with any (even) number of subcolumns.Using more subcolumns leads to better sampling and hence more accurate estimates of averages of SGS variability (Larson and Schanen, 2013).However, feeding fewer subcolumns into the microphysics scheme reduces computational cost.We have presented results using 16 subcolumns, which is a compromise between accuracy and computational cost.Shown in Fig. 14 are the results of reducing the number of sample columns from 16 to 4. The results are similar, with the exception of LBA, which, with four sample points, does not produce surface precipitation, thereby permitting large amounts of cloud water to remain aloft.
Computational cost can also be reduced by coarsening the vertical resolution.To test the effects of this, Fig. 15 plots simulations of the deep convective cases, configured identically as before except that they use the 30-level stretched grid that is used in the CAM climate model (Neale et al., 2012).The simulations are adequate, but they are too convectively active and contain too much condensate.This is particularly apparent in LBA, which again has excessive cloud water aloft because the precipitation efficiency is too low.
That the results of the CLUBB-SILHS simulations degrade some when (vertical/temporal/sampling) resolution is decreased is not surprising.Such sensitivities are typical in all scales of models.For example, CAM has been shown to be sensitive to both horizontal and vertical resolution and the physics time step used (Bogenschutz et al., 2012;Reed et al., 2012;Williamson, 2013a, b;Wehner et al., 2014).The degradation of SCM simulations when brought to a resolution more realistic for including in a host model like CAM suggests that some more testing and tuning may be necessary; however, the simulations show initial promise.6 How does the generalized coupling to microphysics influence shallow cloud simulations?
Using the same configuration as used for the deep convective cases, we also simulate two shallow cloud cases.We find that these two cases are not degraded by the modifications introduced into CLUBB-SILHS in order to improve deep convective simulations.
Results from the RICO shallow cumulus case are shown in Fig. 16.CLUBB-SILHS simulates this trade-wind cumulus case reasonably well, although CLUBB-SILHS underes-Figure 13.CLUBB-SILHS simulations of TWP-ICE (above), LBA (middle), and ARM97 (below) showing sensitivity of simulations to the time step.For each case, liquid cloud fraction is shown on the left (averaged over the same time periods as in previous figures) and a time series of snow water path is shown on the right.The black line is the SAM CRM simulation, the red dashed line is the (default) CLUBB-SILHS simulation presented earlier in Sect.3, and the blue dotted line is CLUBB-SILHS with the same configuration as in the red line, but using a 5 min time step, rather than 1 min.With an increased time step, the CLUBB-SILHS simulations are too convectively active; however, the simulations are still reasonable when compared to the SAM CRM simulations.
timates cloud liquid water.However, the underestimate is not introduced by the modifications to CLUBB-SILHS.Rather, these errors are typical of some prior versions of CLUBB.
Results from the DYCOMS-II RF02 drizzling stratocumulus simulation are shown in Fig. 17.The mean fields are comparable to those in SAM, although, once again, cloud liquid and rain water path are underestimated in the CLUBB-SILHS simulation.The underestimate of liquid is mostly attributable not to the modifications made to CLUBB, but rather to the coarse resolution of CLUBB-SILHS versus SAM -in the 128-level grid, only 15 levels reside in Figure 14.CLUBB-SILHS simulations of TWP-ICE (above), LBA (middle), and ARM97 (below) showing sensitivity of simulations to the number of subcolumns utilized in SILHS.For each case, liquid cloud fraction is shown on the left (averaged over the same time periods as in previous figures) and a time series of snow water path is shown on the right.The black line is the SAM CRM simulation, the red dashed line is the CLUBB-SILHS simulation presented earlier in Sect.3, and the blue dotted line is CLUBB-SILHS with the same configuration as in the red line, but using only four subcolumns for SILHS, rather than 16.Results are degraded only slightly with a reduced number of subcolumns, except in the LBA case, where too much cloud forms aloft.
the lowest 1200 m of the domain.At higher resolutions, CLUBB-SILHS simulates greater liquid water content for the DYCOMS-II RF02 case (Griffin and Larson, 2013).
Although further study is necessary, the fact that reasonable results can be obtained for shallow convection with the same model configuration as for deep convection hints that a unified PDF parameterization may indeed be possible, with one equation set representing all turbulence and cloud types.(below) showing sensitivity of simulations to the number of vertical grid levels.For each case, liquid cloud fraction is shown on the left (averaged over the same time periods as in previous figures) and a time series of snow water path is shown on the right.The black line is the SAM CRM simulation, the red dashed line is the (default) CLUBB-SILHS simulation presented in Sect.3, and the blue dotted line is CLUBB-SILHS with the same configuration as in the red line, but using a 30-level vertical grid, rather than a 128-level grid.The CLUBB-SILHS simulations with decreased vertical resolution contain too much condensate, which is particularly apparent in the LBA case.

Conclusions
This paper presents single-column simulations of deep and shallow cloud layers.To simulate these clouds, our model contains one component (CLUBB) that estimates the subgrid PDF of clouds, turbulence, and hydrometeors, and a second component (SILHS) that draws samples from this PDF and feeds them into a microphysics scheme.This methodology provides a detailed representation of subgrid dynamics and hydrometeors, and also provides a mechanism to couple the two at the subgrid scale.important for deep convection, which has strong interactions between dynamics and microphysics.Most of the simulated fields are satisfactory, although there exist a few deficiencies.One is that rain forms 100 min too late in the LBA case of transition from shallow to deep convection (see Fig. 4).Another deficiency is that precipitation does not form in the first of the three convective events in ARM97 (see Fig. 7).Both these deficiencies suggest that CLUBB-SILHS produces insufficient precipitation when convection is weakly forced.Furthermore, when CLUBB-SILHS's time step is degraded to 5 min (see Fig. 13) or the number of vertical grid levels is degraded to 30 (see Fig. 15), the simulated updrafts are too strong, the surface precipitation is too weak, and too much cloud forms aloft.Again, these are symptoms of low precipitation efficiency.That is, the amount of precipitation produced in CLUBB-SILHS for an amount of condensate is too low in these cases -leaving too much cloud water remaining compared to the SAM simulations.
Nevertheless, CLUBB-SILHS produces quite acceptable overall results for both deep and shallow simulations.In particular, most profiles of clouds and hydrometeors match those of SAM acceptably, given the stochastic nature of convective precipitation.What ingredients in CLUBB-SILHS allow it to simulate deep convection?One is CLUBB-SILHS's detailed representation of and coupling between turbulence and microphysics and, in particular, ice microphysics.CLUBB-SILHS uses a delta-lognormal subgrid PDF of hydrometeors, which allows for the possibility of a hydrometeor-free region of a grid box and also allows the hydrometeors to fall preferentially through liquid cloud water.This, in turn, allows more precipitation to reach the ground, reducing cloud water aloft (see Fig. 11).Although the precipitation efficiency is still somewhat too low in the simulations, the deltalognormal PDF does improve precipitation efficiency.Another key ingredient is vertical turbulent transport of hydrometeors, which increases the altitude reached by snow and other hydrometeors.
These simulations suggest that -with improvements at coarse time steps and vertical grid spacings -it is feasible to develop a unified parameterization of deep convection, shallow convection, stratiform clouds, and turbulence.Such a parameterization would parameterize all cloud types in a largescale model with a single equation set.A unified parameterization would avoid the artificial categorization of clouds that is inherent in parameterization suites that use separate schemes for separate regimes.Use of a unified parameterization would also make it easier to ensure consistency of assumptions throughout a model.Greater consistency, in turn, would instill more confidence in a model's formulation.
Appendix A: CLUBB equations CLUBB includes predictive equations for the horizontal winds (u and v), total water r t (liquid + vapor), liquid water potential temperature θ l , and several higher-order moments.The following are the prognostic equations in the CLUBB model, provided for reference:  where R is the radiative heating rate, f the Coriolis parameter, u g and v g the geostrophic winds, and g is acceleration due to gravity.The equations are anelastic, where ρ s is the dry, static, base-state density and θ vs is the dry, base-state θ v , which both base-state variables change only with respect to altitude.In all of the equations, an overbar denotes a grid box average, and a prime denotes a subgrid deviation from that average.The subscript tol refers to a minimum threshold value and ls denotes large-scale forcing.Not shown are terms relevant for hole-filling, clipping, and damping.For more details on the various constants and closures used in this equation set, see Golaz et al. (2002) and Larson et al. (2012).

Figure 1 .
Figure 1.Time series of liquid water path (upper left), rain water path (upper right), ice water path (lower left), and snow water path (lower right) from the TWP-ICE deep convective simulation.Liquid water path is too high in the CLUBB-SILHS simulation, suggesting a low precipitation efficiency; however, the other hydrometeors follow SAM fairly closely.

Figure 2 .
Figure2.Profiles of liquid water potential temperature (θ l ) (upper left), total water mixing ratio (r t ) (vapor + liquid) (upper right), liquid cloud fraction (lower left), and cloud liquid water mixing ratio (lower right) from the TWP-ICE deep convective simulation.Profiles are averaged over minutes 6000-8000, capturing the strongest convective event.CLUBB-SILHS matches SAM well, with some small differences notable in θ l and liquid cloud fraction.

Figure 3 .
Figure3.Profiles of rain water mixing ratio (upper left), cloud ice mixing ratio (upper right), snow mixing ratio (lower left), and graupel mixing ratio (lower right) from the TWP-ICE deep convective simulation, averaged over the same time period as Fig.2(minutes 6000-8000).CLUBB-SILHS hydrometeor profiles match well with those simulated in SAM.

Figure 5 .
Figure 5. Profiles of liquid water potential temperature (θ l ) (upper left), total water mixing ratio (r t ) (vapor + liquid) (upper right), liquid cloud fraction (lower left), and cloud liquid water mixing ratio (lower right) from the LBA deep convective simulation.Profiles are averaged over the last hour of the simulation (minutes 301-360).The CLUBB-SILHS simulation matches the CRM simulation well, aside from a lack of cloud in lower levels.

Figure 6 .
Figure 6.Profiles of rain water mixing ratio (upper left), cloud ice mixing ratio (upper right), snow mixing ratio (lower left), and graupel mixing ratio (lower right) from the LBA deep convective simulation, averaged over the same time period as Fig. 5 (minutes 301-360).Profiles are consistent with what was seen in Fig. 4.

Figure 7 .
Figure 7. Time series of liquid water path (upper left), rain water path (upper right), ice water path (lower left), and snow water path (lower right) from the ARM97 deep convective simulation.CLUBB-SILHS's simulation of the first event produces cloud liquid water but no precipitation.However, CLUBB-SILHS's simulation of the other two events match SAM well in timing and magnitude.

Figure 9 .
Figure 9. Profiles of rain water mixing ratio (upper left), cloud ice mixing ratio (upper right), snow mixing ratio (lower left), and graupel mixing ratio (lower right) from the ARM97 deep convective simulation, average over the same time as Fig. 8 (minutes 4321-5580).The hydrometeor profiles match the CRM fairly well, though the ice hydrometeors are too few and do not reach high enough in the clouds.

Figure 10 .
Figure10.CLUBB-SILHS simulations of TWP-ICE (above), LBA (middle), and ARM97 (below) demonstrating the ability of CLUBB-SILHS to simulate dynamics for the deep convective cases.For each case, the variance of vertical velocity is shown on the left, the turbulent flux of liquid water potential temperature is in the middle, and the turbulent flux of total water mixing ratio is on the right (averaged over the same time periods as in previous figures).CLUBB-SILHS matches the CRM well in variance of vertical velocity, but underestimates the magnitude of some of the turbulent fluxes, particularly in the more convectively active cases.

Figure 11 .
Figure11.LBA simulations showing the sensitivity to the inclusion of the new precipitation fraction.The black line is the SAM CRM simulation, the red dashed line is the CLUBB-SILHS (default) simulation presented in Sect.3.2, which includes the precipitation fraction, and the blue dotted line is CLUBB-SILHS with the same configuration as in the red line, but with the option for the precipitation fraction set to false.Above are profiles of cloud liquid water mixing ratio (left) and cloud ice mixing ratio (right), averaged over the last hour of the simulation (minutes 301-360).Below are time series of liquid water path (left) and rain water path (right).Nonzero precipitation fraction is important for LBA, because it increases precipitation efficiency, allowing more rain to leave the atmosphere, thereby removing excess cloud liquid water aloft and reducing excessive ice formation and growth.

Figure 15 .
Figure15.CLUBB-SILHS simulations of TWP-ICE (above), LBA (middle), and ARM97 (below) showing sensitivity of simulations to the number of vertical grid levels.For each case, liquid cloud fraction is shown on the left (averaged over the same time periods as in previous figures) and a time series of snow water path is shown on the right.The black line is the SAM CRM simulation, the red dashed line is the (default) CLUBB-SILHS simulation presented in Sect.3, and the blue dotted line is CLUBB-SILHS with the same configuration as in the red line, but using a 30-level vertical grid, rather than a 128-level grid.The CLUBB-SILHS simulations with decreased vertical resolution contain too much condensate, which is particularly apparent in the LBA case.

Figure 17 .
Figure17.Results from the DYCOMS-II RF02 simulation of stratocumulus clouds.Above are time series of liquid water path (upper left) and rain water path (upper right).Below are mean profiles of liquid water potential temperature (middle left), total water mixing ratio (middle right), liquid cloud fraction (lower left), and cloud liquid water mixing ratio (lower right).Profiles are averaged over the last hour of the simulation (minutes 301-360).CLUBB-SILHS does a reasonable job of simulating these shallow stratocumulus clouds, though the amount of liquid is too low, which is related to the coarse vertical grid spacing.

Table 1 .
The prescribed ratio (V i ) between the within-precipitation variance and the mean squared of hydrometeors.Below-cloud values are the same values as used for cloudy levels.

Table 2 .
Hydrometeor correlations used for all five cases.The variables included in the correlation matrix are, from left to right or top to bottom: extended liquid water mixing ratio, orthogonal extended liquid water mixing ratio, vertical velocity, extended cloud drop number concentration, and the mixing ratios and number concentrations of rain, cloud ice, snow, and graupel.Shown are the within-cloud values.For below cloud another correlation table is used, and for our purposes it is identical except for the correlation between s and t.This correlation is set to 0.3 below cloud.The correlation matrix is symmetric, so the values below the diagonal are not filled in.

Table 3 .
Options for SAM simulations of the five cases.

www.geosci-model-dev.net/8/1/2015/ Geosci. Model Dev., 8, 1-19, 2015 Figure 12. TWP
-ICE (above), LBA (middle), and ARM97 (below) simulations showing the sensitivity to the boosted eddy diffusivity for convective cases.For each case, graupel mixing ratio is shown on the left and snow mixing ratio on the right.The black line is the SAM CRM simulation, the red dashed line is the CLUBB-SILHS simulation configured exactly as that presented earlier in Sect.3, and the blue dotted line is CLUBB-SILHS with the same configuration as in the red line, but without the new diffusivity calculation.The new formulation for eddy diffusivity smooths the hydrometeor profiles and transports hydrometeors farther aloft.