A Stochastic Closure for TwoMoment Bulk Microphysics of Warm Clouds: Part II, Validation
Submitted on 2 May 2016 to
Meteorology and Atmospheric Physics
Abstract
The representation of clouds and associated processes of rain and snow formation remains one of the major uncertainties in climate and weather prediction models. In a companion paper (Part I), we systematically derived a two moment bulk cloud microphysics model for collision and coalescence in warm rain based on the kinetic coalescence equation (KCE) and used stochastic approximations to close the higher order moment terms, and do so independently of the collision kernel. Conservation of mass and consistency of droplet number concentration of the evolving cloud properties were combined with numerical simulations to reduce the parametrization problem to three key parameters.
Here, we constrain these three parameters based on the physics of collision and coalescence resulting in a “region of validity.” Furthermore, we theoretically validate the new bulk model by deriving a subset of the “region of validity” that contains stochastic parameters that skillfully reproduces an existing model based on an a priori droplet size distribution by Seifert and Beheng (2001). The stochastic bulk model is empirically validated against this model, and parameter values that faithfully reproduce detailed KCE results are identified.
Furthermore, sensitivity tests indicate that the stochastically derived model can be used with a time step as large as 30 seconds without significantly compromising accuracy, which makes it very attractive to use in medium to long range weather prediction models. These results can be explored in the future to select parameters in the “region of validity” that are conditional on environmental conditions and the age of the cloud.
MSC 35L02 and 65C20 and 65C35 and 65M08 and 65M30 and 65Z99 and 86A10
O O1.0em,0.5ex m m m
1 Introduction
Climate and numerical weather prediction models are large computer codes based on a discretization of the governing equations for fluid dynamics and thermodynamics for the coupled atmosphereoceanearth system. Due to limitations in computing resources only a small part of the spectrum of scales involved in this interactive complex system are explicitly modelled and many small scale processes are represented via various recipes known as parameterizations. The parameterization of cloud and precipitation processes constitute a major uncertainty in climate and weather prediction models. For example, reliable determinations of the onset of rainfall and the radar reflectivity in clouds are crucial for the accuracy of models that rely on the parameterizations of cloud microphysical processes [14]. A large size gap exists between the length scales on which the cloud microphysical processes occur and the grid cell size of global and region climate models [3]. Bulk parameterizations are simple equations that represent these processes in a computationally efficient manner [4].
In warm clouds, bulk parameterizations can include terms representing nucleation, collision and coalescence, and precipitation [7]. The simplest models evolve an initial droplet size distribution (DSD) using only collision and coalescence terms: selfcollection (cloud and rain), autoconversion, and accretion. Auto conversion involves droplets of a size that are more strongly influenced by turbulence then the other collision and coalescence processes [18], and it is the only process that is present in all of the equations that evolve both rain and cloud aggregates. In a bulk parameterization of warm clouds, for example, a common technique was to derive or ‘fit’ the auto conversion term for rain mixing ratio to data, and then make the remaining auto conversion terms to be functions of the rain mixing ratio term [8, 16, 4]. In Part I, we broke with that tradition by independently deriving all auto conversion terms while preserving conversation of mass and accuracy of number [2]. The independent derivations of all auto conversion terms permitted further control of the auto conversion parameter. This new technique is based on a systematic partitioning of the droplet distribution spectrum into cloud and rain aggregate means plus sets of pointwise random fluctuations. This stochastic representation, which is in essence a new way of taking into account uncertainty in the DSD, is used as a basis to close high order terms to derive a two moment parameterization for bulk cloud microphysics directly from the kinetic coalescence equation (KCE).
Kessler is widely regarded as the first to parameterize a collision and coalescence process in a cloud microphysical context. In a 1969 paper he used a Heaviside function to terminate autoconversion when a critical threshold of cloud liquid water content was reached [11]. Since then there have been many parameterizations; some continuing to only model auto conversion [11, 12] and others modelling a full suite of one and two moment bulk parameterizations [8, 16, 4, 7]. The development of an auto conversion only model included one that softened the termination of auto conversion [19]. Others introduced a dependence on number concentration [11] and [12]. A few bulk parameterizations have used fitted models, one to DNS [4] and another to detailed schemes [8]. Similarly to the new stochastic model, Seifert and Beheng (2001, hereafter SB01) used the kinetic collection equation (KCE) to derive a twomoment bulk cloud parameterization of collision and coalescence in the form of an ODE system which evolves the mixing ratios and number concentrations of both cloud and rain aggregates [16]. The KCE is an integrodifferential equation which evolves the droplet spectrum by summing over all possible collisions leading to two integrals modelling respectively the gain and loss of droplets of every possible size:
(1)  
where is the drop mass, is the number concentration, a density function, and is the collisioncoalescence kernel so that is the ratedensity of collisioncoalescence between two droplets of mass and , respectively.
Khairoutdinov and Kogan (2000) [8], Seifert and Beheng (2001) [16], and Franklin (2008) [4], derived expressions for the autoconversion process, and for the accretion process, as they affected rain mixing ratio. They used conservation of mass to deliver expressions for loss of cloud mixing ratio due to these processes, and used the equation for mean droplet mass,
(2) 
to deliver expressions for other terms in their parameterizations. The autoconversion parameterizations for the remaining evolved quantities (, , ) are simply functions of the rain mixing ratio terms.
Deriving a model independent of the DSD is important because any assumed distribution may not be universally appropriate for differing cloud types [21, 20]. Thus a climate model utilizing only the new stochastic bulk model as a microphysical parameterization scheme, could represent different cloud types, ages, and turbulent phenomena or intensities by simply selecting appropriate stochastic parameter values to be used in each climate model cell.
As already mentioned, in lieu of applying any particular droplet size distribution (DSD), the stochastic bulk parameterization represented number concentration density and mixing ratio density as the sum of the state space mean and a pointwise fluctuation, and used a 2D Taylor expansion to represent the collision kernel [2]. Further simplifications and a significant reduction in the number of parameters was achieved in the light of detailed numerical simulations based on the KCE by invoking two independent order of magnitude arguments: (i) the mean and standard deviation of the stochastic processes, and (ii) the temporal fluctuations of the evolving quantities [2]. Moreover, the resulting twomoment bulk cloud parameterization is independent of the underlying collision kernel and yet is consistent with existing models [8, 16, 4], and in particular the SB01 model can be recovered exactly with an appropriate choice of parameters. The new stochastic model also independently derived autoconversion and accretion terms for each of the four evolved quantities while preserving conservation of mass and consistency of number [2]. This new stochastic model is, to the authors knowledge, unique among existing stochastic models in that it is a stochastic differential equationbased method (‘SDEbased’ stochastic method) rather than being based on sampling an assumed probability distribution (‘samplebased’ stochastic method). See Part I for a brief list of ‘samplebased’ stochastic methods used to model cloud microphysical processes [2].
The paper is organized as follows. In Section 2, we recall the stochastic bulk parameterization model equations and use physical constraints associated with the processes of collision and coalescence to constrain the parameters and define a “region of validity.” In Section 3, we apply the piecewise polynomial kernel used in [16] to the stochastic bulk parameterization and conduct validation tests based on detailed simulations of KCE. In particular, we do a term by term matching of the SB01 and the stochastically derived model equations and define the region in parameter space delimited by the applicability of the SB01 model. The “region of validity,” in parameter space, based on collision and coalescence physics, and the region defined by the SB01 model are not identical but they are not disjoint either. A thorough exploration of these regions in parameter space is conducted for two typical initial cloud distributions, one representing a polluted/continental cloud and one representing a clean/oceanic cloud and the sets of parameters that reproduce the KCE results are identified for each distribution. As we will see in Section 4, while the KCE results are reproduced by a set of parameters within the intersection of the two regions for the case of a clean cloud, in the case of a polluted cloud the matching parameters are found outside the SB01 model limits. A concluding discussion is given in Section 5.
2 Controlling the Autoconversion Parameter
In Part I, we systematically decomposed the number concentration density and mixing ratio density into statespace means and sets of pointwise fluctuations within cloud and rain droplet aggregates, summed the fluctuations over defined intervals and considered each instantaneous sum to be a point in a stochastic process. Closure of the kinetic collection equation was achieved by representing higher moments as functions of the pointwise fluctuations. The stochastic processes were defined as a mean and a standard deviation, and order of magnitude arguments eliminated the standard deviation terms. The means of the stochastic processes are parameters. Consistency of number and conservation of mass were used to reduce the degrees of freedom of the stochastic parameterization. Number concentration and mixing ratio are related by Equation 2. Writing the four stochastic bulk rate equations using three independent parameters yields
(3)  
The physical meanings and effects of these parameters are
controls the strength of cloud selfcollection relative to autoconversion  
controls the strength of autoconversion  
controls the strength of rain selfcollection 
The derivations and details which resulted in the stochastic parameterization presented in Equation 3 are given in the first paper. The values of each of the three stochastic parameters could be acquired from the results of a binpair interaction detailed method.
The ‘SDEbased’ stochastic model is general enough for any kernel, but here we validate it by applying the piecewise polynomial kernel given in Equation 8, and graphically acquire values of the stochastic parameters after constraining the stochastic parameters using conservation of mass, consistency of number, and the definition of auto conversion. We compare it to sievert and Beheng’s model using Bott’s detailed method as a benchmark. The values for results driven by a hydrodynamic kernel and by two turbulent kernels are presented in a third paper.
2.1 Physical Limits: Consistency of Number and Conservation of Mass
Since the kernel () and number concentration () are nonnegative, they are omitted from the following expressions. Cloud selfcollection and autoconversion each decrease cloud number concentration, thus: . Autoconversion will cause rain number concentration to increase, thus: . Both of these conditions result in . Because cloud selfcollection does not affect rain number concentration, the cloud loss will be at least twice the rain increase: . This results in . Auto conversion will decrease cloud mixing ratio and cause rain mixing ratio to increase, thus: . Thus we have
(4) 
The bounds on the two parameters related to cloud selfcollection and auto conversion due to accuracy of number and conservation of mass can be represented on a single graph shown in the picture on the left in Figure 1. The lower bound on , or the strength of autoconversion, is dependent on the distance of the mean cloud radius to the threshold radius: the closer they are together the stronger the auto conversion rate and a smaller mean cloud radius demands a weaker autoconversion rate. This specific relationship can be quantified by examining the coefficients of the terms in Equation 3.
Additionally, rain selfcollection decreases rain number concentration and is constrained by the following inequality:
2.2 Autoconversion and the Threshold Droplet Size
There is a need to constrain the auto conversion parameter [20]. Franklin (2007) and Khairoutdinov and Kogan (2000) have constrained it in an adhoc manner using best fit methods to DNS and to a detailed, bin based method, respectively. Neither their constraints nor their parameters have physically based meanings [5, 8]. In the stochastic bulk parameterization, the auto conversion parameter is further constrained using the relationship between the mean cloud radius and the threshold radius, and the definition of autoconversion; namely, that precollision, and postcollision, droplets must be smaller than and larger than the threshold radius, respectively. All of the mass in an auto conversion collision moves from cloud to rain.
The effect that autoconversion has on the rain spectrum is described by the first terms on the rhs side of the second and fourth equations in Equation 3. These are related to number and mixing ratio, respectively. Therefore, given , the equation requiring that postcollision droplets have a mass greater than the mass of the threshold droplet is
Letting and solving for gives the equation for the upper boundary of the triangular region in the graph on the right in Figure 1:
(5) 
The lower boundary is constructed from the requirement that precollision droplets are each smaller than the threshold droplet. This is accomplished algebraically by considering that the number of cloud droplets lost to autoconversion is twice the number rain droplets produced by the same process. The number of cloud droplets lost is given by doubling the first term on the rhs if the second differential equation in Equation 3 and proceeding similarly as for the upper boundary.
which yields the lower bound on the triangular region in the graph on the right in Figure 1:
(6) 
Using as is depicted in the right hand graph in Figure 1, we get for the lower and upper bounds, respectively:
When a cloud is young or heavily polluted, we can expect to be large and the values for the stochastic parameters to be near the upper right corner of the region in the plane. This demands strong cloud selfcollection and weak autoconversion, which is to be expected when the cloud droplets are primarily small. By contrast in a mature cloud or one with primarily large droplets, will be smaller and the region giving values for will be different.
The strength of cloud selfcollection is determined from the relationship between the first terms in the first two differential equations in Equation 3. When the value of is chosen such that the first term in the equation for is zero (i.e. ), then none of the droplets subtracted in the equation for creates rain droplets and thus all of them create cloud droplets and cloud self collection is maximized. Consequently, autoconversion in nonexistent and the region of valid values for vanishes at the intersection of the lines for max cloud selfcollection and no autoconversion (i.e. and ).
Alternately, when the value of is chosen such that all of the droplets subtracted in the evolution of in Equation 3 produce rain droplets in the second equation, (i.e. ), then cloud selfcollection is nonexistent and auto conversion is maximized. In other words, the coefficient of the first term in is twice the coefficient of the first term in because two cloud droplets produce one rain droplet in the process of autoconversion.
2.3 Constraint of the Autoconversion Parameter
The constraining of the auto conversion parameter is possible because the four auto conversion terms are derived independently. The independent derivations allow the ratios of mixing ratio and number to be compared to the threshold mass. This is not possible when one auto conversion term is a function of another auto conversion term. The equations for the two sloped sides of the triangle in Figure 1 appear because of the independent derivations of the number concentration and mixing ratio equations. The third side is resultant from the two cloud droplet terms for number being derived independently.
Conservation of mass and consistency of number create the same bounds as the upper right vertex of the triangle in Figure 1 which is the location on the graph representative of strong cloud selfcollection and weak (to nonexistent) auto conversion. Stochastic parameter values in this area of the triangle model a polluted cloud with many small droplets. Consistency of number constrains the maximal and minimum value of the cloud selfcollection parameter. Conservation of mass only constrains the maximal value (weakest strength) of the auto conversion parameter.
The maximum strength of auto conversion is limited by the lower vertex of the triangle. Auto conversion is maximal when the mean cloud radius is close to the threshold radius which gives minimal values to in the two vertices on the left side of the triangle which collaborates some of Franklin’s results [4]. Alternatively, when cloud selfcollection is strong, the lower bounds on the value of the auto conversion parameter is restricted further than the apparent bounds given by the perimeter of the triangle: a strong cloud selfcollection parameter () disallows taking values from the lower area of the triangle. Thus is constrained to be ‘near’ to 0.5, i.e. modelling a weaker auto conversion rate.
3 The Stochastic Model with a PieceWise Kernel
The stochastically derived model (3) is general and can be adapted for any kernel, but here we validate it for the case of the piecewise polynomial kernel in (8). This particular kernel was used by SB01 to derive a similar set of twomoment equations based on the averaging the KCE with a specified droplet distribution to close high order terms. For the sake of completeness, we also derive a region of validity for the stochastic model in (3) by matching it to the SB01 model based on their assumptions of parameters, namely, the shape parameter of the assumed distribution, , and the time scale for autoconversion and accretion processes. We call this region the SB01 region of validity in contrast to the stochastic region of validity derived above based on the physical constraints motivated by the properties of the collision and coalescence process of cloud droplets. We then use detailed simulations of KCE as benchmarks to find parameter values that are physical based on two measures (the time at which 50% of cloud mass is converted into rain and the size of the mean rain droplet at that time) within the regions of validity of the stochastically derived model, with this specific piecewise polynomial kernel, visavis the performance of the SB01 model. The cases involving a hydrodynamic kernel and two turbulent kernels will be considered in a future publication.
3.1 Seifert and Beheng’s Model (2001)
SB01 developed a parameterization that models the effects of collision and coalescence on the evolution of four bulk cloud microphysical properties. Their parameterization was subsequently used in 2010 in a 1D rain shaft model that included droplet breakup and condensation and evaporation [17].
In [16] the centres of cloud and rain aggregates: and are systematically used, and a gamma distribution and are assumed as initial conditions to close the twomoment problem. They derived autoconversion and accretion expressions as functions of their tuning parameters: the gamma distribution shape parameter and a nondimensional time parameter :
(7)  
where is a threshold mass corresponding to a radius of m. The functions and represent internal time scales for the corresponding processes that are derived empirically using simulated data with [16].
SB01[16] used a piecewise combination of Golovin’s sum of mass kernel and Long’s quadratic polynomial kernel, referred to here simply as the piecewise polynomial kernel.
(8) 
where cm gs, and cm gs. Using abbreviated notation, their bulk model is
(9)  
The parameters and are used in seven of eight terms in their bulk parameterization.
The terms in Equation 9 are labeled by the collision processes that they model. Auto conversion occurs in each of the four evolving quantities and accretion occurs in three equations. The coefficients of seven of the eight terms are dependent on an adhoc shape parameter , and all are dependent on one or more of six other parameters (, , , , , ).
Although SB01 omitted from the accretion terms in their final parameterizations, this factor is included in in Equation 7 for completeness, and it appears in precisely the same terms in the stochastic model once the same piecewise polynomial kernel is applied. These eight coefficients, excluding the collision kernels, are equated with the coefficients of the corresponding terms in the stochastic parameterization.
3.2 Reproduce Existing Parameterization
We reproduce a bulk model driven by a piecewise kernel. Applying a different kernel, say a hydrodynamic or a turbulent one, would yield a different parameterization. The piecewise polynomial kernel shown in Equation 8 is applied to the stochastic bulk parameterization Equation 3. This analytic kernel introduces the centroids of the cloud and rain aggregates into the stochastic parameterization which appear in the accretion terms. The stochastic bulk parameterization of cloud microphysical collision processes driven by the LongGolovin kernel is:
(10)  
A termbyterm comparison of Equations (9 and 10) reveals an exact match regarding the dependencies on aggregate number concentrations ( and ) and mixing ratios ( and ). The coefficients of corresponding terms are equated to determine the range of the stochastic parameters as functions of the range of values that Seifert and Beheng gave to their parameters, and .
3.3 Select Stochastic Parameter Values
We make a systematic comparison by matching corresponding terms. In both Seifert and Beheng’s parameterization and the stochastic parameterization, the auto conversion terms for cloud and rain mixing ratio are equal with the opposite sign. The same is true for the cloud and rain accretion terms. Consequently, with regard to mixing ratio, only the auto conversion and accretion terms that evolve rain mixing ratio are compared because the cloud terms have been restricted by conservation of mass. The stochastic term for rain selfcollection is compared with the corresponding term in the Seifert and Beheng parameterization. The following calculations use a threshold radius of m and a mean cloud radius of m. The auto conversion term for rain number concentration contains a stochastic parameter () that is a different parameter than the one for mixing ratio (), and so the rain number coefficients are compared. The first terms in both cloud number concentration equations combine the processes of auto conversion and cloud selfcollection. These are compared; however, they contain the same stochastic parameter () as in rain number concentration. The more meaningful range of values for is derived from cloud number concentration. The physically permissible values of the stochastic parameters that reproduce Seifert and Beheng’s bulk parameterization are collected tabularly and displayed graphically in Section 3.4.
Auto conversion, mixing ratio
The coefficients for auto conversion of rain (and cloud) mixing ratio are set equal
and solved for the stochastic parameter:
The shape parameter fraction is monotonic over the interval , and the time parameter function has a maximum at . The bounds that Seifert and Beheng set for and and the maximum of are used to determine the range of the stochastic parameter governing the strength of auto conversion:
(11) 
This range is displayed graphically (blue patterned region) in Figure 2.
Auto conversion, rain number concentration
Equating the corresponding autoconversion terms for rain number concentration gives
Solving for the stochastic parameter gives:
(12) 
The stochastic parameter governing the strength of cloud selfcollection relative to auto conversion is constrained by Equation 12 to values less than but close to the maximum permissible value: The lower bound is consistent with the value of when given the constraint that a rain droplet produced during an auto conversion process must be larger than the threshold mass .
Auto conversion, cloud number concentration
Equating the corresponding autoconversion terms for cloud number concentration gives:
For any value of specified by Seifert and Beheng (), the bounds on the stochastic parameter is readily calculated: . The derivation of their auto conversion term for rain umber concentration uses the assumptions that and . These assumptions may be true during the early stages of a cloud’s lifetime, but are less likely to be true in a mature cloud. Therefore, the lower bound ( derived in this section will be used to constrain this stochastic parameter. Since the upper bound based on cloud number concentration terminates auto conversion, the upper bound () based on rain number concentration will be used.
(13) 
This range is displayed graphically (blue patterned region) in Figure 2.
Accretion
The none of the three accretion terms in the stochastic bulk parameterization contain any stochastic parameters. Except for the sign, the coefficients for the three accretion terms are identical:
The nondimensional time function is an adhoc parameter that was included after completing their derivations, and the shape parameter is a byproduct of the assumed droplet gamma distribution. Seifert and Beheng removed the shape parameter function from the accretion terms in their final parmeterization scheme. There are no arbitrary or adhoc parameters in the stochastic bulk rate parameterization. Consequently, the only difference between the stochastic parameterization and their ‘adjusted’ parameterization is the nondimensional time function which was added after derivations were complete. This function serves to slow accretion during the early and late stages of a cloud’s lifetime and accelerate accretion during the middle stages, maximally when .
A first order truncation of the collision kernel in the stochastic parameterization may give stochastic parameters in the accretion terms which may serve the same purpose as their nondimensional time function. However, all stochastic parameters are associated with concrete physical meanings. Trading the higher order Taylor expansion terms of the collision kernel for the simplicity of the ‘SDEbased’ stochastically derived set of coupled ODE’s and omitting the adhoc nondimensional parameter inclusions, all accretion terms are identical.
Rain selfcollection, number concentration
Equating the corresponding rain selfcollection terms for number concentration gives:
(14) 
3.4 Regions of Validity in Parameter Space
The stochastic parameterization reproduces Seifert and Beheng’s parameterization when , , and are restricted to values identified in Equations (11, 13, and 14). These restrictions are formed by making the stochastic parameters to be functions of Seifert and Beheng’s parameters: and . Equations (11 and 13) bound the values of and and create a 2D space of parameter values which is displayed as the patterned (blue) in Figure 2. Values within this region control the strength of cloud selfcollection and auto conversion. This rectangular region intersects the triangular (solid, red) region formed using the physical constraints from Sections 2.1 and 2.2.
This intersection is a pentagonal region with vertices (counterclockwise from lower left in the graph on the right of Figure 2 at (0.49066, 0.20094), (0.49533, 0.20094), (0.49997, 0.49808), (0.49997, 0.49829), and (0.49995, 0.49829). This pentagon is most clearly depicted in Figure 4 when discussing sensitivity in Section 4.3. The bounds of values of and in this pentagonal region are listed in Table 1. Three of the sides of the pentagon are apparent in the graph on the right in Figure 2. The graphs in Figure 2 and the values in Table 1 are derived using a threshold radius of m and an initial mean cloud radius of m which is representative of a terrestrial or ‘polluted’ cloud. The equations for the two sloped sides of the solid, red triangle in Figure 2 are given as upper and lower bounds of in Table 1.
Parameter  Lower  Condition  Upper 

0.49066  —  0.49997  

0.20094  3215.5  

6431.5  

0.49829  

=0.25  

[black, ultra thick]Top MarkCenter MarkBottom Mark
For the ‘clean’ cloud, an initial mean cloud radius of m is used with the same threshold radius. The bounds of the stochastic parameters for a ‘clean’ or marine cloud when the evolution is driven by the LongGolovin kernel is shown graphically in Figure 3 and numerically in Table 2. In a ‘clean’ cloud, the ratio forces a sufficiently smaller upper bound on (weaker cloud selfcollection and stronger auto conversion) such that the intersecting region is a quadrilateral. The vertices (counterclockwise from lower left in the graph on the right of Figure 3 of this foursided figure are (0.34320, 0.72496), (0.4216, 0.72496), (0.49300, 0.39063), and (0.49300, 0.44531). Three of the sides of the quadrilateral are apparent in Figure 3. This quadrilateral is most clearly depicted in Figure 5 when discussing sensitivity in Section 4.3. The equations for the two sloped sides of the solid, red triangle are given as upper and lower bounds of in Table 2.
Parameter  Lower  Condition  Upper 

0.34320  —  0.49300  

0.7250  7.81253.40625  

15.6257.3125  7.81253.40625  

=0.25  

4 Numerical Simulations
The intersection of the rectangular (blue), ‘SB subset region,’ and triangular (red), ‘region of validity’ of parameter values form polygons. Parameter values within these polygons both represent physically real conditions and are a function of Seifert and Beheng’s parameters. We examine the evolution of aggregate number concentration and mass by the stochastic bulk model by (i) varying the cloudself collection and auto conversion parameters within the intersection of the ‘region of validity’ and the ‘SB subset region,’ (ii) observing the time required to convert 50% of cloud droplet mass to rain droplets within this intersection, and (iii) comparing the “50% time” to results from the Linear Flux Method [1]. We show the effects of macroscopic cloud properties within the redblue intersection of by calculating the initial rain droplet radius and the radius at the time of 50% conversion. Finally, sensitivity tests are performed w.r.t the time step and the rain selfcollection parameter in Section 4.3.
4.1 Model Setup
We perform simulations for two cloud types: polluted and clean. The polluted clouds have an initial mean cloud droplet radius of 10 m which is 239 droplets per cm given the typical liquid water content of 1 g cm. Clean clouds have an initial mean cloud droplet radius of 16 m and a corresponding droplet concentration of 58.3 cm. The third stochastic parameter , regulating rain selfcollection is fixed at 0.25. All simulations use 3 order AdamsBashforth with Improved Euler and AB2 for initialization.
The stochastic model is validated against results from a detailed evolution of the cloud droplet spectrum. We use the Linear Flux Method by Bott (1998) with 70 bins that double every second bin. The smallest droplet radius on the discretized spectrum is 1 m and the largest is 3.251 mm.
4.2 Results
Three types of results are shown for each cloud type. We present contour graphs of the 50% conversion time that span relevant portions of the intersection of the the region of validity and the ‘SB subset region.’ The relevant portions are where changes in the 50% conversion time as a function of the cloud selfcollection and auto conversion parameters are noticeable. This occurs where cloud selfcollection is stronger and auto conversion is weaker, where and are near by slightly less than 0.50. Within this same region in the () space we present contour graphs showing how the mean rain radius, at the 50% conversion time, changes as the stochastic parameters are varied. Additionally we present another type of contour graph that shows how the mean rain radius at the first time step varies as a function of (). The third result is a set of curves of massovertime graphs comparing the conversion of cloud to rain mass by the stochastic model, Seifert and Beheng’s model, and the detailed Linear Flux Method [1].
4.2.1 50% Mass Conversion: Contour Graphs
Contour graphs depicting the time needed to convert 50% of cloud mass to rain mass are shown for simulations with an initial mean cloud radius () of 10 m (polluted) and 16m (clean). The polluted simulation requires a cloud selfcollection phase to grow the mean cloud droplet size large enough for auto conversion to have a significant effect. Therefore, the cloud selfcollection parameter can be expected to be strong and the auto conversion parameter to be weaker as is shown by the location of the thinlined black box in the upper corner of the lefthand graph in Figure 4.
Both graphs in Figure 4 are blowups of the graphs in Figure 2. In the graph on the left, the solid (red) region contains values of the cloud selfcollection and auto conversion parameters that adhere to consistency of number and conservation of mass. The patterned (blue) region contains values which reproduce Seifert and Beheng’s bulk parameterization. This graph is not to scale, but identifies where the contour map is located relative to the graphs in Figure 2. The two black line segments in the contour map correspond to the edge of the patterned (blue) area within the black box in the graph on the left. Stochastic parameter values below and to the left of the black lines will reproduce Seifert and Beheng’s bulk parameterization given the the bounds on their parameters and . Stochastic parameter values which are above and to the right of the black lines do not reproduce their bulk parameterization.
The contour map (on the right) in Figure 4 shows how strongly the evolution is affected by small changes in the auto conversion parameter in this small corner in the region of permissible values. The black triangle identifies the stochastic parameter values (0.49997, 0.49808) and the red diamond is (0.49997, 0.49904). These parameters convert 50% of cloud droplet mass to rain droplet mass in 1824 and 1944 seconds, respectively. The black diamond and the red triangle have the same %50 mass conversion time as the corresponding marks of the same colour. This is depicted in the contour graph by noting that contour lines are horizontal, i.e. mass conversion by the polynomial kernel is independent of the cloud selfcollection parameter. The evolution curves corresponding to these parameters are in Figure 8.
The red marks are located outside of the region of stochastic parameter values that will reproduce Seifert and Beheng’s bulk parameterization. The stochastic parameters are constant throughout the evolution. Seifert and Beheng’s nondimensional time parameter changes as the cloud matures. The auto conversion parameter can be considered as an average of possible auto conversion parameter values during the maturation of a cloud. The function is concave. Jensen’s inequality shows that the concave function of a mean gives an over estimate [10], and thus unduly increasing the strength of the auto conversion parameter. When the incidental strengthening of the auto conversion parameter is mitigated by relaxing the restrictions imposed by the bounds of and on the values of the stochastic parameters, the stochastic model produces results consistent with detailed results using stochastic parameters identified by the (red) diamond.
The contour graphs for simulations with an initial mean cloud radius of 16m is shown in Figure 5 with the thin black boxes and lines, and the solid (red) and patterned (blue) regions serving the same purpose as in Figure 4. These graphs detail the regions presented in Figure 3. The black triangle, and square, identify the stochastic parameter pairs (0.493013, 0.390833) and (0.486060, 0.391093), respectively. These each convert 50% of mixing ratio to rain mixing ratio in 880 seconds. The evolution curves of both pairs are shown in Figure 9. The red diamond identifies the stochastic parameter values (0.482682, 0.229417) and converts 50% of mixing ratio to rain mixing ratio in 726 seconds.. A cloud with a mean cloud radius of 16m has more liquid water content contained larger cloud droplets. A stronger (smaller valued) auto conversion parameter is expected, and is seen when comparing the results from the simulations with an initial droplet radius of 10m.
4.2.2 Mean Rain Radius: Contour Graphs
The mean rain radius is an important complement to the 50% conversion time metric in a zerodegree parcel cloud model. The 50% conversion time quantifies how quickly cloud water is redistributed from being contained in cloud droplets to forming rain droplets. However, in a parcel model simulation the newly formed rain droplets are essentially “still in the cloud,” and have not started precipitating. This is a significant difference between a parcel model and a 1D rain shaft model. By analyzing the conversion time in conjunction with mean rain droplet radius at 50% conversion time, we can get an idea of the onset of precipitation for particular stochastic parameter values relative to other stochastic parameter values and relative to other bulk and detailed models. A 1D model would be appropriate when comparing to observational data. However, the inclusion of observational data is beyond the scope of this work.
A larger mean rain droplet radius at the time of 50% conversion will lead to a quicker onset of precipitation, and a smaller mean rain droplet radius at that time indicates a slower onset of precipitation. The mean rain radius at the 50% time is shown in Figure 6 for both the polluted (10 m) case and the clean (16 m) case. The boundary of the ‘SB subset region’ and the same triangular and square markers, identifying the location of particular parameter values, as shown in Figures 4 and 5 are replicated in the rain radius contour plot. These two contour plots also show that there is a greater dependence, of the mean rain radius at the 50% time of conversion, on the cloud selfcollection parameter than on the auto conversion parameter.
The initial (at the first time step) mean rain radius is shown in the two contour plots in Figure 7. In both cases, the upper boundary has an initial mean rain radius of 40 m which is equal to the threshold radius. This boundary separates the precollision droplets that form a rain droplet during the auto conversion process from ones that form a cloud droplet during this particular collision process. Since the latter contradicts the definition of auto conversion, this upper boundary separates physically allowable stochastic parameter values from ones that are not allowable. At the lower boundary, the initial mean rain radius is 50.4 m. This boundary separates the precollision cloud droplets that form a rain droplet during the auto conversion process from ones that formed by pairs of rain droplets during this particular collision process. Since the latter contradicts the definition of auto conversion, this lower boundary also separates physically allowable stochastic parameter values from ones that are not allowable. The range of initial rain radii is important when using particular turbulent kernels in the stochastic bulk model.
For instance, a tabulated turbulent kernel by [13] has the greatest enhancement when smaller cloud droplets are coalescing with larger cloud droplets. A turbulent kernel such as this would be best served by stochastic parameters along the upper boundary (m) since the strongest enhancement of collisions occurs with the smallest of droplets, the predominate size of newly produced rain droplets be approximately equal to the size of the threshold droplet: m.
In contrast, a parameterized kernel by [4] has stronger enhancement factors when droplets of nearly the same size coalesce and the greatest enhancements are when the largest of the cloud droplets coalesce with each other. This type of turbulent kernel would be best represented by the stochastic bulk model when using parameters along the lower boundary (m). Two cloud droplets with radii slightly less than that of the threshold droplet would produce a rain droplet with a radii of .
4.2.3 Evolution of Mass Conversion
Cloud and rain mixing ratios are evolved using three methods: (i) the new ‘SDEbased’ stochastic bulk parameterization, (ii) the bulk parameterization by Seifert and Beheng (2001), and (iii) a detailed scheme using Bott’s Linear Flux Method (LFM) from 1998. The results from the detailed method are taken as the benchmark. Mixing ratios in polluted clouds are evolved for 45 minutes, and clean clouds (larger initial mean droplets) for 30 minutes. The stochastic parameter pairs chosen for these simulations are taken from the marker locations on the contour plots in Figures 4 and 5.
(, )  second  

seconds  m)  m)  
(0.49994008, 0.49808267)  1824  628.2  40.0 
(0.49997003, 0.49808213) 
1824  791.5  50.4 
(0.49997003, 0.49904107) 
1944  773.9  40.0 
(0.49998501, 0.49904053) 
1944  975.2  50.4 
Bulk: Seifert & Beheng (2001) 
2253  151.3  39.6 
Detailed: Linear Flux Method 
1985  

For the polluted cloud case, the stochastic parameter pairs used to generate the two stochastic curves in Figure 8 are in the second and third rows of Table 3. The first and fourth rows show simulations that generate the same 50% conversion time, but have different mean rain radii at the 50% time. The stochastic parameters that produce a 50% conversion time of 1944 seconds were specifically chosen to match the detailed results of 1985 seconds. Parameters that exactly reproduced the time of 1985 seconds we’re available. However, the evolution curves were indistinguishable. So, we chose values that highlighted how closely the curves matched. Seifert and Beheng’s conversion rate was considerably slower. However, the retarded growth by Seifert and Beheng’s parameterization is consistent with Franklin’s (2008) analysis that their model of auto conversion underestimates the auto conversion rate [4].
Referring to Figure 4, we see that the intersection of the ‘SBsubset region’ and the ‘region of validity’ is a truncation of the triangle due to the limits imposed by the values of and that Seifert and Beheng chose. When recreating their model with the stochastic parameterization, the space of stochastic parameter values was created using the maximums and minimums of and . However, changed during the evolution while the stochastic parameters remained constant. This arbitrary limitation calls for the relaxation of the truncation on the triangle due to the bounds on and imposed by and because the function is concave and by Jensen’s inequality overestimates the auto conversation terms [10]. To compensate for this acceleration, values of and which inhibit the evolution are allowed. These values are stochastic parameter pairs of (, such that and are allowed.
Choosing stochastic parameters that adhere to the threshold droplet restriction, i.e. contained within the (solid, red) triangle in Figure 2, but not necessarily adhering to the derived restrictions on the auto conversion parameter, produces a bulk evolution curve more reflective of the detailed benchmark as shown in Figure 8. Figure 4 plots the locations of the coordinates ( ) used to generate the results of the mean stochastic parameterization in Figure 8 relative to the linear boundaries given by the relationship of auto conversion and the droplet threshold.
(, )  second  

seconds  m)  m)  
(0.482682, 0.229417)  726  132.4  50.4 
(0.493013, 0.390833) 
880  167.4  50.4 
(0.486060, 0.391093) 
880  132.9  40.0 
Bulk: Seifert & Beheng (2001) 
678  72.3  39.6 
Detailed: Linear Flux Method 
721  

For the clean cloud case, the stochastic parameter pairs used to generate the two stochastic curves in Figure 9 are in the first two rows of Table 4. The evolution curve, ‘Stoch A,’ is detailed in the first row of Table 4 and closely matches the 50% conversion time of the detailed simulation. For the initial 57 minutes the conversion of cloud water to rain water is excessive. At approximately 400 seconds the rate of conversion, as indicated by the slope of the respective graphs, is less than the rate for the detailed simulation. The slower conversion rate of the stochastic parameterization persists for the duration of the evolution.
When m cloud mass is closer to the cloudrain threshold. To accurately model this initial condition a more aggressive auto conversion parameter is needed. Consequently, values of the auto conversion parameter need to be smaller and are situated more towards the centre of the quadrilateral region. The preferred stochastic value relative the region on which one can take values is shown in Figure 5. This is in contrast to a more polluted initial condition such as when m. In this polluted situation the cloud mass is further from the cloudrain threshold and is more accurately modelled by a weaker auto conversion parameter, and weaker to the extant that the concavity bias of Seifert and Beheng’s nondimensional time parameter places an excessive restriction on the cloud selfcollection and auto conversion stochastic parameters as is shown by the diamond (red) marker in Figure 4 being outside the patterned (blue) region.
4.3 Sensitivity Tests
Two types of sensitivity tests were performed. The effect that the time step had on the time required to convert 50% of cloud droplets to rain droplets and the rain droplet radii at that time. The third sensitivity test examines the effect that the rain number parameter has on these metrics. Both tests were performed using an initial mean cloud radius of 10 m, and cloud selfcollection and auto conversion parameters of =(0.49997,0.49808).
As expected, the rain number stochastic parameter only effects the evolution of rain number. When is varied from 0.10 to 0.40, the mean rain radius varies by 12.7%.
seconds  (m)  
0.40  1824  742.7 
0.30 
1824  775.0 
0.25 
1824  791.6 
0.20 
1824  808.5 
0.10 
1824  843.3 

timestep  

seconds  seconds  (m) 
120  1920  797.6 
60 
1860  818.0 
30 
1830  795.9 
10 
1830  800.1 
1 
1824  791.6 
0.1 
1823.6  791.0 

Varying the time step from 0.1 seconds to 30 seconds has a negligible effect on the time required to convert half of the initial cloud mass to rain droplets (column two in Table 6). The mean rain radius at the time that half of the liquid water content is converted to rain varies by 4.9m (0.62%) between time steps of 0.1 and 30 seconds. No advantage is gained by reducing the time step smaller than the 1 second used in the simulations in this paper. It may be possible to relax the time step to 30 seconds without significantly compromising accuracy. It should be noted that at 50 minutes of evolution the mean rain radius was greater than half a centimetre indicating that the omission of rain removal and droplet breakup were effecting the results.
5 Conclusion
Assumptions and simplifications are necessary to produce computationally affordable parameterizations that represent cloud microphysical processes. Parameterizations of cloud microphysical processes over the past fortyfive years have made assumptions regarding the droplet size distribution [16, 11, 12]. Many parameterizations depend on adhoc parameters [8, 16, 11, 12, 4]. The stochastic bulk parameterization of collision and coalescence studied here does not rely on any particular distribution of the droplet size spectrum, but rather assumes that a distribution exists and it has a mean. All of the parameters in the stochastic parameterization have physical meaning, and their values can be recovered from data. The parameters represent the first moments of time series of aggregates of fluctuations and of product fluctuations. These are fluctuations from the state space mean of mixing ratio and number over defined portions of the droplet size spectrum. Recovering values of that data can be done graphically because the stochastic bulk parameterization sufficiently constraints the auto conversion parameters. Values can also be recovered from detailed simulations, which can be deterministic or stochastic, such as Bott’s computationally efficient binbased method [1] or the more computationally demanding GillespieSeeelberg algorithm [6, 15] as will be shown in a subsequent paper.
The stochastic bulk parameterization of collision and coalescence is independent of any assumed droplet size distribution and yet retains the flexibility of applying any collision kernel to the resultant set of coupled ordinary differential equations. The flexibility of using any collision kernel in the stochastic parameterization comes at the expense of retaining only the zeroth order term in the 2D Taylor expansion of the collision kernel centred at either the cloud, or rain, mean mass. This flexibility requires that the value of the collision kernel be computed, or retrieved from a lookup table, while the stochastic parameterization is being utilized by a climate model.
Stochastic parameters can be chosen that reflect particular cloud characteristics like type (maritime vs. continental and clean vs. polluted) and age. This is particularly important because now a single bulk model of collision and coalescence can be used without the restriction of needing the same droplet size distribution for every cell in a climate model. The new stochastic bulk model gives the flexibility of choosing different stochastic parameters for different cloud types. Also, different cells in a climate model can use stochastic parameters associated with various collision kernels depending on the intensity and type of turbulence found in the clouds in that climate model cell.
We used the relationship of auto conversion to threshold droplet mass to define ‘regions of validity’ for the auto conversion () and the cloud selfselection () parameters such that the selection of values for these parameters would cause the auto conversion terms in the stochastic bulk model to produce a rain droplet when two cloud droplets coalesced. This restriction of the auto conversion parameter was only possible because each of the four occurrences of an auto conversion term in the stochastic bulk model was derived independently of every other auto conversion term. To the best of the authors’ knowledge, no other bulk cloud microphysics model of collision and coalescence was able to accomplish this physicsbased restriction on the auto conversion term because every other bulk model only derived one auto conversion term and then made the other three auto conversion terms to be functions of the derived term.
Other parameterizations have been restricted to a specific kernel [16, 11, 12, 4]. The kernel in Franklin’s parameterization was dependent on the turbulent strength, but limited in the range of turbulent kinetic energy used in that kernel. The stochastic bulk parameterization can be used with a variety of kernels without any further derivations. The value of the kernel at the mean cloud radius and mean rain radius is used in the stochastic bulk parameterizations. Seifert and Beheng’s [16] parameterization used a kernel that contained both mean radii, and could be directly applied to the resultant stochastic derivation. This application changed the powers of the evolving quantities in the stochastic equations, and the resultant powers of the evolved quantities matched the ones in Seifert and Beheng’s differential equations.
Moreover, physical constraints based on the physics of collision and coalescence processes were used to define a “region of validity” for the stochastic parameters. We matched both models term for term and solved for the stochastic parameters and as functions of their adhoc parameter. We defined an ‘SB subset region’ in the 2D () space by applying the bounds that they gave for their parameters. Their bulk model is restricted to a piecewise polynomial collision kernel, while ours retains the flexibility to use a variety of collision kernels. Both models had precisely the same dependence on cloud and rain mixing ratio and number concentration.
For given values of the three stochastic parameters, the new stochastic parameterization produced results that closely matched results from detailed simulations that used Bott’s LFM when the same initial conditions were used. This was done for each of polluted and clean cloud types. A simple metric of the time required to convert half of the cloud mixing ratio to rain mixing ratio was used to validate the results of the stochastic bulk model. It is demonstrated here that while adequate parameter values can be found for the stochastic model, they can depend heavily of the initial mean cloud radius, i.e, the environmental conditions or the cloud’s age. Moreover, for the two experiments considered here, despite the utilization of an ad hoc adjustment function to correct of excessive accretion and autoconversion rates, Seifert and Beheng’s model either overestimates or underestimates the time at which 50% of the cloud mass is converted into rain. Thus, the stochastic parameters that are found to faithfully reproduce the detailed KCE results can sometimes lies outside the ‘SB subset region.’ This is essence demonstrates the usefulness of the stochastic bulk parameterization. As alternative to their adhoc auto conversion and accretion adjustment functions, Markov jump processes can be used to sample the “region of validity,” with the jump processes conditional on environmental conditions (e.g. clean v.s. polluted, marine v.s. continental) and the age of the cloud [9].
The sensitivity of the ”50% conversion time” was examined as a function of the 2D  space. This conversion time metric was sensitive to small changes in the auto conversion parameter near the values that would reproduce the evolution curve produced by the detailed method, but showed very little sensitivity elsewhere in the domain of physically permissible values. The conversion time metric showed little sensitivity to the cloud selfcollection parameter anywhere on this domain for a given value of the rain selfcollection parameter. The mean rain radius at the time of 50% conversion was sensitive to the value of the rain selfcollection parameter. Both the conversion time and the mean rain radius at 50% conversion time showed very little sensitivity to changes in the time step in the range 0.130 seconds. Subsequent research in the development of the stochastic bulk parameterization will apply hydrodynamic and turbulent kernels for the purpose of developing a methodology that models cloud microscopic processes using only cloud macroscopic information.
The intersection of the “region of validity” and the “SB subset region” defines a set of values for the stochastic parameters of cloud selfcollection and auto conversion that both adheres to physical constraints and reproduces the bulk model by Seifert and Beheng (2001). We validated the stochastic bulk model by finding stochastic parameter values that very closely reproduced the 50% mass conversion time and the evolution curves. This was done for the cases of a polluted cloud (239 cm) and a clean cloud (58.3 cm). The model was further validated for both cloud types by noting that the mean rain radius at the first time step had physically meaningful bounds related to the process of auto conversion, i.e. the bounds on the rain radius at the first time step was dependent on the fact that the process of auto conversion produces one rain droplet from two cloud droplets.
Acknowledgements This research is part of D. Collins’s Ph.D. thesis. The research of B. Khouider is partly supported by a grant from the Natural Sciences and Engineering Research Council of Canada. D. Collins’s fellowship is partly funded through this grant.
References
 [1] Bott, A.: A flux method for the numerical solution of the stochastic collection equation. J. Atmos. Sci. 55, 2284–2293 (1998)
 [2] Collins, D., Khouider, B.: A stochastic closure for twomoment bulk microphysics of warm clouds: Part I, derivations. in preparation (2016)
 [3] Devenish, B., Bartello, P.: Review article: Droplet growth in warm turbulent clouds. J. Q. R. Meteorol. Soc. 138, 1401–1421 (2012)
 [4] Franklin, C.: A warm rain microphysics parameterization that includes the effects of turbulence. J. Atmos. Sci. 65, 1795–1816 (2008)
 [5] Franklin, C., Vaillancourt, P., Yau, M.K.: Statistics and parameterizations of the effect of turbulence on the geometric collision kernel of cloud droplets. J. Atmos. Sci. 64, 938–954 (2007)
 [6] Gillespie, D.T.: An exact method for numerically simulating the stochastic coalescence process in a cloud. J. Atmos. Sci. 32, 1977–1989 (1975)
 [7] Grabowski, W.W., Wang, L.: Growth of cloud droplets in a turbulent environment. Annu. Rev. Fluid Mech. 45, 293–324 (2013)
 [8] Khairoutdinov, M., Kogan, Y.: A new cloud physics parameterization in a large eddy simulation model of marine stratocumulus. Mon. Wea. Rev. 128, 229–243 (2000)
 [9] Khouider, B., Biello, J., Madja, A.J.: A stochastic multicloud model for tropical convection. Comm. Math. Sci. 8(1), 187–216 (2010)
 [10] Larson, V.E., Wood, R., R., F.P.: Systematic biases in the microphysics and thermodynamics of numerical models that ignore subgrid scale variability. J. Atmos. Sci. 58, 1117–1128 (2001)
 [11] Liu, Y., Daum, P.: Parameterization of the autoconversation process. Part I: Analytical formulation of the kesslertype parameterizations. J. Atmos. Sci. 61, 1539–1548 (2004)
 [12] Liu, Y., Daum, P., McGraw, R., Wood, R.: Parameterization of the autoconversation process. Part II: Generalization of the sundqvisttype parameterizations. J. Atmos. Sci. 63, 1103–1109 (2006)
 [13] Pinsky, M., Khain, A., Krugliak, H.: Collisions of cloud droplets in a turbulent flow. part v: Applications of detailed tables of turbulent collision rate enhancement to simulation of droplet spectra evolution. J. Atmos. Sci. 65, 357–374 (2008)
 [14] Pruppacher, H., Klett, J.: Microphysics of Clouds and Precipitation. Kluwer Academic Publishers, Boston (1997)
 [15] Seesselberg, M., Trautmann, T., McGraw, R., Thorn, M.: Stochastic simulations as a benchmark for mathematical methods solving the coalescence equation. Atmos. Res. 40, 33–48 (1996)
 [16] Seifert, A., Beheng, K.D.: A double moment parameterization for simulating autoconversion, accretion, and selfcollection. J. Atmos. Sci. 5960, 265–281 (2001)
 [17] Seifert, A., Stevens, B.: Microphysical scaling relations in a kinematic model of isolated shallow cumulus clouds. J. Atmos. Sci. 67, 1575–1590 (2010)
 [18] Shaw, R.: Particle turbulent interactions in atmospheric clouds. Annu. Rev. Fluid Mech. 35, 183–227 (2003)
 [19] Sundqvist, H.: A parameterization scheme for nonconvective condensation including prediction ofcloud water content. Quart. J. Met. Soc 104, 667–690 (1978)
 [20] Wood, R.: Drizzle in stratiform boundary layer clouds. Part II: Microphysical aspects. J. Atmos. Sci. 62, 3034–3050 (2005)
 [21] Zawadzki, I., Fabry, F.: The development of drop size distribution is light rain. J. Atmos. Sci. 51(11001114) (1994)