Relationship between drizzle rate , liquid water path and droplet concentration at the scale of a stratocumulus cloud system

Relationship between drizzle rate, liquid water path and droplet concentration at the scale of a stratocumulus cloud system O. Geoffroy, J.-L. Brenguier, and I. Sandu GAME/CNRM (Météo-France, CNRS) 42, av. Gaspard Coriolis 31057 Toulouse Received: 12 December 2007 – Accepted: 8 January 2008 – Published: 26 February 2008 Correspondence to: O. Geoffroy (geoffroy@knmi.nl) Published by Copernicus Publications on behalf of the European Geosciences Union.


Introduction
The formation of liquid precipitation in convective clouds involves a suite of microphysical processes driven by the turbulent cloud circulation: production of supersaturation, activation of cloud condensation nuclei (CCN) where supersaturation occurs, droplets growth by condensation of water vapor, collection between hydrometeors, and finally precipitation of the biggest ones below cloud and their evaporation in subsaturated levels.Implementation of numerical schemes to simulate these processes in multidimensional models raises a few major complications.(i) The supersaturation in warm clouds is on the order of 1%, but the accuracy of the temperature and total water content, that, combined, determine Correspondence to: O. Geoffroy (geoffroy@knmi.nl)the supersaturation, hardly reach 1%.The supersaturation therefore is often parameterized, in particular its peak value that governs CCN activation.(ii) The aerosol response to supersaturation during the activation process depends on the size distribution, chemical composition and mixing state of the aerosol particles.A large number of variables is therefore necessary to fully predict the CCN properties of the aerosol.(iii) Once activated, the diameter of hydrometeors range from 1 µm after activation up to 8 mm for the biggest drops.Accurate description of the hydrometeor size distribution thus calls for an additional large number of model variables, namely from 50 to 200 size bins (Kogan et al., 1994;Liu et al., 2000;Feingold et al., 1994Feingold et al., , 1996Feingold et al., , 1997;;Stevens et al., 1996Stevens et al., , 1998;;Khairoutdinov and Kogan, 2000, referred hereafter to as KK00; Feingold and Kreidenweis, 2002;Jiang et al., 2002).Therefore, microphysical processes are simulated in details using box models with prescribed dynamics (Mordy, 1959;Fitzgerald, 1974;Feingold and Kreidenweis, 2000;Feingold and Chuang, 2002), but their implementation in multidimensional models is numerically very costly, and thus requires drastic simplifications.
Compared to the latter, bulk microphysics schemes offer a very efficient possibility to circumvent these difficulties.In warm clouds, the hydrometeor population is divided in two categories: cloud droplets and precipitating particles.Both are characterized by their mixing ratio (Kessler, 1969), hence reducing to two the number of variables necessary to describe a complete hydrometeor spectrum.The scheme can further be improved by adding two additional variables, describing the number concentration of particles in each category (Berry and Reinhardt, 1974;Beheng, 1994;KK00).With such an improvement, the bulk scheme becomes applicable to studies of the impact of CCN changes, hence the impact of cloud droplet number concentration (CDNC) changes, on precipitation efficiency.The CCN properties can be characterized by the CCN activation spectrum that represents the number concentration of activated CCN as a function of the O. Geoffroy et al.: Parameterization of the mean drizzle rate in stratocumulus peak supersaturation.The parameterization package is finally completed with a diagnostic scheme for the peak supersaturation (Twomey, 1959;Cohard et al., 1998).A few bulk schemes have been developed since the 70s, some being more suited to deep and heavily precipitating convective clouds (Kessler, 1969;Tripoli and Cotton, 1980;Beheng, 1994), and others being more fitted to slightly precipitating extended boundary layer clouds (KK00).The rationale for splitting the hydrometeor size distribution into two categories is to be able to attribute different fall velocities to each category.The definition of the size limit between cloud droplet and precipitating particle is thus crucial.For heavily precipitating clouds, most of the falling mass is contained in millimeter size particles, also referred to as rain drops, and the size limit between the two categories is of the order of 40 µm in radius (Berry and Reinhardt, 1974;Beheng, 1994).In slightly precipitating stratocumulus, most of the falling mass is contained in particles smaller than 50 µm, also referred to as drizzle.The limit is thus set at a smaller value, such as 25 µm in KK00.
In an attempt to simulate precipitation in climate models, bulk schemes were implemented in general circulation models (GCM).Considering that these schemes were initially developed for cloud resolving models (CRM), in which local values of the microphysical fields are resolved, their extension to the GCM coarse spatial resolution is questionable.Indeed, the formation of warm precipitation is a non-linear process.For droplet radii r between 10 and 20 µm, the droplet fall velocity rapidly increases with size (proportional to r 2 ), and the probability of collision increases, leading to droplet coalescence and the formation of precipitation embryos.Once a few embryos have been created, the collection of droplets by falling drops produces precipitation more efficiently.The onset of precipitation is therefore highly sensitive to the size of the biggest cloud droplets, which in turn depends upon the local values of the liquid water content and cloud droplet concentration.In a GCM, the liquid water content is distributed over the cloud fraction of the model grids, i.e. on scale of a few tens of kilometers.Its mean value is thus significantly smaller than the peak values that are simulated by a CRM.CRM bulk schemes have therefore been adapted to the GCM spatial resolution by empirically tuning some key coefficients, e.g. the threshold cloud droplet mean volume radius at which collection starts to be active (Rotstayn, 2000).
The original values of the bulk parameters, such as the threshold radius for collection, are physically based, but the modified values used in GCMs are not.It is therefore not granted that values tuned on present climate and CCN properties still hold in a modified context.The rationale for transposing CRM bulk scheme to the GCM resolution is not obvious either.Indeed, bulk schemes aim at simulating the complete microphysical cycle from CCN activation to precipitation in individual convective cells.At the spatial resolution of a GCM, the precipitation is averaged over a large number of cells at different stages of their cycle.It might therefore be more relevant to develop a parameterization of the precipitation flux averaged over an ensemble of cloud cells, to correctly represent boundary layer clouds in GCMs and particularly to study the aerosol indirect effect.
The horizontal mean precipitation rate of marine stratocumulus clouds has been observed during field experiments and it has been shown that, despite the high spatial and temporal variability and intermittency of the precipitation fields, on average, the precipitation from a large number of convective cells only depends on the liquid water path (LWP) or the cloud geometrical thickness, modulated by the CDNC (Pawlowska and Brenguier, 2003;Comstock, 2004;van Zanten et al., 2005).In this study, a LES model, with bulk parameterizations of the microphysics, is used to replicate these empirical results and provide quantitative assessment of the relationship between the precipitation rate, the LWP and the CDNC.
Bulk parameterizations for CRM are briefly described in the next section and the tuning of the coefficient necessary to extend such parameterizations to GCMs is discussed in Sect.3. Relationships between precipitation rate and cloud properties derived from recent field experiments are summarized in Sect. 4. The LES model used in this study is described in Sect. 5 and simulations of precipitating stratocumulus are compared to the observations in Sect.6.In Sect.7 a parameterization of the precipitation rate is proposed and compared to the ones derived from the observations.

Bulk parameterizations for CRM
In a detailed microphysics model, the collection process between hydrometeors is explicitly solved by calculating the probability for each particle to collide and coalesce with any other particle along its trajectory (Berry, 1967).When the hydrometeors population is split in two categories, the collection process may result in four different scenarios: 1. Cloud droplet self-collection: a cloud droplet collecting a cloud droplet to form a larger cloud droplet.This event does not impact the cloud droplet mixing ratio, but it reduces the cloud droplet number concentration.
2. Cloud droplet autoconversion: a droplet collecting a droplet to form a precipitating particle.This event reduces the cloud droplet mixing ratio and number concentration, and it increases the precipitating particle mixing ratio and number concentration.
3. Accretion: a precipitating particle collecting a cloud droplet to form a larger precipitating particle.This event impacts both mixing ratios as with the previous one, and it reduces the cloud droplet number concentration.
4. Precipitating particle self-collection: a precipitating particle collecting a precipitating particle to form a larger precipitating particle.This event does not impact the precipitating particle mixing ratio, but it reduces the precipitating particle number concentration.
In bulk parameterizations using only two independent variables, the cloud droplet and precipitating particle mixing ratios, only events (ii) and (iii) need to be accounted for, while all of them shall be parameterized if four variables are used (mixing ratios and number concentrations).
Various parameterizations of these collection subprocesses have been developed.Berry and Reinhardt (1974), Tripoli andCotton (1980), andBeheng (1994) use analytical functions to describe the hydrometeor spectra and solve the SCE (Stochastic Collection Equation) numerically to derive parameterized conversion rates.The same technique was applied by Wood (2005) using measured spectra, and by KK00 using spectra simulated with a bin microphysical model implemented in a LES model.Liu and Daum (2004) use analytical representations of the hydrometeor spectra and of the collision kernel to analytically derive functional relationships for the conversion rates.
The main benefit of the analytical approach is to provide functional relationships that can then be used to examine the sensitivity of conversion rate to varying conditions.It is however difficult to fully constrain and some unresolved parameters are tuned empirically, such as the efficiency coefficient in the parameterization of the autoconversion rate (Liu and Daum, 2004).The accuracy thus depends on how realistic are the analytic functions selected for representing the spectra and the collision kernel.In contrast, when using bin spectra, the high accuracy of the explicit microphysical scheme used to build the training data base is transposed to the bulk solutions, but these solutions are limited to the range of variation of the microphysical variables that has been explored in the training data base, either from field data or from numerical simulations.
A summary of the most currently used bulk parameterizations is given in Table 1 of Wood (2005).The main differences are between formulations of the autoconversion process, because of different definitions of the limit between droplet and drop, and because of different ways of parameterizing the rapid increase of the collection efficiency at this limit.

Extension of bulk parameterizations to GCM
In most GCMs, the cloud liquid water content (LWC) is derived from the total water content (a conservative prognostic variable) by saturation adjustment and the precipitating water content is diagnosed (Smith, 1990;Tiedtke, 1993;Lohmann and Roeckner, 1996;Ghan et al., 1997;Rotstayn, 1997;Wilson et Ballard, 1999).Saturation however is assumed to be restricted to a cloud fraction F c of the grid, assuming prescribed distributions of temperature and water content fluctuations at the subgrid scale (Sundqvist, 1978;Tiedtke, 1993;Del Genio et al., 1996;Fowler et al., 1996;Rotstayn, 1997;Boucher et al., 1995;Rasch and Kristjánsson, 1998;Wilson et Ballard, 1999).F c may further be divided into a convective and a stratiform fraction where different parameterizations are applied (Tiedtke, 1993).With a diagnostic treatment of the precipitating water content, long time steps can be used (Ghan and Easter, 1992), but the model is not able to simulate horizontal advection of precipitations.A few models though consider the precipitating water as a prognostic variable (Boucher et al., 1995;Fowler et al., 1996).
Rutledge and Hobbs (1983) use mixing ratios and set the threshold to q ccrit =w ccrit /ρ w =0.7•10 −3 kg kg −1 .Fowler et al. (1996) implemented this parameterization in their GCM, with q ccrit =0.25•10 −3 kg kg −1 , and the same value of the coefficient α as in Kessler (1969).Tripoli and Cotton (1980) also diagnose CDNC and therefore express the threshold for the onset of precipitation in term of a critical cloud droplet mean volume radius, whose value was derived from explicit calculations of cloud droplet collection as r vcrit =10 µm (Manton and Cotton, 1977).This value has been corroborated by observational studies of the onset of precipitation (Gerber, 1996;Boers et al., 1998;Pawlowska and Brenguier, 2003).Rotstayn (1998) however proposed that this value should be reduced to r vcrit =7.5 µm in a GCM for the simulations to better fit satellite climatology.Wilson and Ballard (1999), and Rasch and Kristjánsson (1998) also use smaller critical values of r vcrit =7 µm and 5 µm respectively.Lohmann and Roeckner (1996) use the Beheng (1994) parameterization with no threshold: where d is the cloud droplet spectrum dispersion that is set to a fixed value between 5 and 15.They however increase the autoconversion efficiency α by a factor of 15 to ) 1.75 min: 61 -min: 47 min: 0.02 max: 190 max:275 max: 24.5 better fit satellite and in situ observations.These few examples reveal that the autoconversion rate initially developed for CRMs needs to be intensified in GCMs, by increasing either the autoconversion efficiency or by reducing the threshold for the onset of precipitation.Indeed, the mean LWC over the cloudy fraction of a GCM is obviously smaller than the peak values a CRM is able to simulate in the most active regions of a cloud system.LES simulations with bulk microphysics parameterizations are particularly sensitive to the choice of the autoconversion scheme (Stevens and Seifert, 2008) but this sensitivity is higher in GCM simulations because the autoconversion threshold drives the onset of precipitation over a larger domain corresponding to the cloud fraction of the GCM grid.This feature artificially accentuates the non-linearity of the onset of precipitation: The whole cloud fraction either precipitates or not.In an actual cloud system, when the convective cells approach the depth at which precipitation might start, the most active ones reach the autoconversion threshold and produce a local precipitation, hence a small precipitation rate when averaged over the domain.With the development of convection, the precipitating fraction of the cloudy domain increases progressively.The non-linearity of the onset of precipitation is therefore smoothed out when the precipitation rate is averaged over a large number of convective cells.
Rostayn (2000) carefully examined the impact of the LWC sub-grid distribution by extending the condensation sub-grid scheme to the calculation of the autoconversion rate.The non-linearity of the autoconversion rate with respect to the LWC does not appear to be the main source of discrepancies, which are rather due to the droplet size threshold for the onset of precipitation.For the model to produce the same globally averaged LWP without and with the sub-grid autoconversion scheme, the threshold radius shall be increased from 7.5 to 9.3 µm, a value closer to the original one of 10 µm.Even with the same globally averaged LWP the two schemes however result in different aerosol impacts, especially on the cloud lifetime (the second indirect effect).

ACE-2 Cloudy-Column
The 2nd Aerosol Characterization Experiment (ACE-2) (Raes et al., 2000) took place in the North-East Atlantic, north of the Canary Islands, from 25 June to 14 July 1997.During the ACE-2 Cloudy-Column project, the Merlin-IV instrumented aircraft from Météo-France was sampling stratocumulus clouds over square areas of about 60 km sides and the DLR Do228 was flying the same track above the cloud layer for measurements of the reflected radiation (Brenguier et al., 2000).Microphysical in situ measurements were further analyzed to derive the following parameters: -N act , hereafter referred to as N 1 , the cloud droplet concentration after activation of CCN, represents the cloud system average of the cloud droplet concentration in regions that are not affected by mixing or drizzle scavenging.It was derived by selecting samples in the middle of the cloud layer, with a quasi-adiabatic liquid water content, and no precipitating particles.
-The cloud geometrical thickness H 1 was derived as the 98% percentile of the cumulative frequency distribution of the height above cloud base of the cloudy samples measured during series of ascents and descents through the cloud layer.
-The LWP values were derived by integrating over the cloud depth the values of the LWC measured along the flight path.Two estimations, hereafter referred to as LW P 1rand and LW P 1 max were derived by assuming either random or maximum overlap of the LWC values (Brenguier et al., 2003).An additional estimation of the LWP, hereafter referred to as LW P 1ovid , was obtained by processing the OVID multispectral radiometer radiances, measured concomitantly from above the cloud layer with the DLR Do228, to derive the cloud optical thickness and the LWP (Schüller et al., 2003).
-The cloud system mean precipitation rate R 1 was derived by integrating over the drizzle size distribution measured on a 100 m scale the precipitation rate in each size class, using a formulation of the drizzle fall speed as a function of drop size following Pruppacher and Klett (1997).These values were then averaged over the cloudy fraction of the flight legs (Pawlowska and Brenguier, 2003).The sampling section of the drizzle probe, that was used to derive concentration from particle counting, is difficult to calibrate for small drops (the first three size classes of the instrument, from 15 to 75 µm in diameter).The value used for the ACE-2 data processing was significantly underestimated in the first three size classes that in fact contained most of the drizzle water content.The derived precipitation rates were thus significantly overestimated.The ACE-2 database was further reprocessed using the same sampling sections as during the DYCOMS-II experiment.
Summarizing the observations of 8 cases studies during ACE-2, Pawlowska and Brenguier (2003) established that the precipitation rate scales with the geometrical thickness and the cloud droplet concentration.The best fit derived using the corrected values of the precipitation rate expresses as (Pawlowska and Burnet, personal communication): It predicts significantly lower values of the precipitation rate, compared to the original fit of Pawlowska and Brenguier (2003).
Note also that the flight legs just below cloud base were not long enough nor the drizzle concentration was high enough for deriving a statistically significant estimate of the precipitation rate at cloud base.Therefore, the precipitation rate was averaged over the whole cloud depth.Vertical profiles of the drizzle water content reveal that it is greater in the cloud layer than at the base (Pawlowska and Brenguier, 2003).The cloud depth averaged precipitation rate therefore overestimates the precipitation rate at cloud base.

EPIC
The Eastern Pacific Investigation of Climate Processes in the Coupled Ocean-Atmosphere System (EPIC) took place in the stratocumulus region near Peru in October 2001 (Bretherton et al., 2004).A scanning C-band radar, a vertically pointing millimeter-wavelength cloud radar, a microwave radiometer and radiation flux instruments, aboard the NOAA Ron Brown ship, were used to sample the structure of the drizzling stratocumulus.Additional in situ measurements of cloud microphysics were collected on board the UK Met Office C-130 aircraft.The data were further analyzed to derive, among others, the following parameters (Comstock et al., 2004): -The LWP, hereafter referred to as LW P 2 , from the microwave radiometer (Fairall et al., 1990).
-The mean CDNC <N >, hereafter referred to as N 2 , was derived from the combination of the microwave observed LWP and cloud transmission measurements from the pyranometer (Dong and Mace, 2003).Linear interpolation of the daytime CDNC was used to estimate values during the night.
-The drizzle precipitation rate, hereafter referred to as R 2 , was derived from the millimeter-wavelength cloud radar reflectivity.
Note that these estimates are averaged over a 3 h period, which corresponds to a horizontal scale of about 75 km, similar to the one sampled by the aircraft during ACE-2.
Summarizing about 30 independent cloud samples, Comstock et al. ( 2004) established the following relationship between the drizzle precipitation rate, the LWP and the domain averaged CDNC: R 2 = 24.37 • 10 9 • (

DYCOMS-II
The Dynamics and Chemistry of Marine Stratocumulus cloud experiment took place at the west of the coast of California in July 2001 (Stevens et al., 2003).The NCAR C-130 was instrumented with the UWyo millimeter-wavelength radar and a backscattering lidar, both pointing at the nadir, and a complete suite of instruments for in situ measurements of cloud properties.The data collected during seven flights of the campaign were further analyzed to derive the following parameters (van Zanten et al., 2005): -The mean CDNC <N >, hereafter referred to as N 3 , was directly measured with a PMS (Particle Measuring Systems) FSSP (Forward Scattering Spectrometer Probe) and averaged over flight legs flown above cloud base and just below cloud top, which corresponds to about 2 h of measurements, i.e. about 700 km in cloud.
-The cloud geometrical thickness, hereafter referred to as H 3 , was derived from the average altitude of the cloud top, as measured with the nadir pointing lidar over a 90 min flight leg (500 km), and in situ measurements of the cloud base altitude averaged over the four hours of flight legs just below and above cloud base.
-The drizzle precipitation rate, hereafter referred to as R 3 , was derived from about 150 min of radar sampling from above the cloud layer, i.e. about 900 km, using a relationship between radar reflectivity and drizzle precipitation rate, that was calibrated for each flight with in situ measurements of the drizzle size distributions. (5)

Discussion
These three field experiments on marine stratocumulus tested the same concept, namely that the mean precipitation rate scales with the cloud geometrical thickness or the LWP, and that it is modulated by the CDNC.In ACE-2 and DYCOMS-II, during which measurements were performed with instrumented aircraft, the cloud systems were characterized by their geometrical thickness, derived from in situ measurements in ACE-2 and from lidar (cloud top altitude) and in situ measurements (cloud base altitude) in DYCOMS-II.While lidar measurements in DYCOMS-II provided a continuous monitoring of the cloud top altitude along the horizontal leg above cloud, the cloud top altitude in ACE-2 was measured only once per ascent or descent through the cloud layer (15 to 35 ascents or descents per flight).The ACE-2 estimates of the cloud geometrical thickness are therefore less significant statistically than the DYCOMS-II ones.They are likely to be underestimated because the biggest cells were scarce and the probability of exiting or entering the base and the top during ascent and descent was low.This speculation is supported by the radiation column closure study of Schüller et al. (2003) that revealed that the cloud geometrical thickness derived from radiance measurements and cloud radiative transfer calculations was systematically greater than the in situ characterization.
In EPIC, remote sensing instruments were operated aboard a research vessel and the cloud thickness was represented by the LWP derived from microwave radiometer measurements.The accuracy of these measurements is of the order of 25 g m −2 (Comstock et al., 2004).LWP values were also derived in ACE-2 from remote sensing of cloud radiances with an airborne multispectral radiometer (OVID).They were evaluated against estimates based on in situ measurements of the LWC assuming either random or maximum overlap, and against estimates derived from the measured frequency distribution of cloud geometrical thickness, assuming each cloud column is adiabatic.The discrepancies between these estimates reach values up to 40 g m −2 , with the OVID remote sensing derived LWP always being greater than the adiabatic, the random and the maximum overlap estimates (Schüller et al., 2003).This result corroborates the previous assessment that the cloud geometrical thickness was underestimated during ACE-2.
The cloud droplet concentration was directly measured in situ in ACE-2 and DYCOMS-II, with an accuracy bet-ter than 20%, since the concentrations were low and coincidence effects were limited (Brenguier et al., 1994(Brenguier et al., , 1998)).These measured values were averaged over the cloud system in DYCOMS-II, while an attempt was made in ACE-2 to characterize the initial cloud droplet concentration that results from CCN activation, before it is diluted by mixing and drizzle scavenging.The comparison between N act and <N > reveals that the former is 10 to 40% greater than the latter, on average.In EPIC, the cloud droplet concentration was derived from remote sensing of the cloud during the day and further extrapolated to get nighttime values.The resulting accuracy is difficult to evaluate, but it can reasonably be assumed to be worse than 50%.
In ACE-2, the precipitation rate was not calculated at cloud base because cloud base legs were too few for a statistically significant estimation of this parameter.All legs were combined, including horizontal legs, ascents and descents through the cloud layer.The largest values of drizzle rate were sampled along horizontal legs during which the estimation of the cloud thickness is less accurate than during ascents and descents.In some cases, the contribution of deeper convective cells with a lower cloud base than the stratocumulus cannot be excluded.For these reasons, it can be speculated that the precipitation rate in ACE-2 is significantly overestimated.During DYCOMS-II, a sophisticated airborne Doppler radar technique was developed to estimate the precipitation rate at cloud base with an accuracy varying between 5 to 11%, depending on the flight (van Zanten et al., 2005).During EPIC, the precipitation rate at the cloud base was derived from the millimeter wavelength radar on board the research vessel and its accuracy is of the order of 50% (Comstock et al., 2004).
Table 1 summarizes these observational results.For ACE-2 and DYCOMS-II, the original relationships using the geometrical thickness are expressed as functions of the LWP by assuming an adiabatic profile of liquid water content: where C w ≈2•10 −6 kg m −4 .

Description of the LES model
The non-hydrostatic model Méso-NH, described in detail by (Lafore et al., 1998), was designed to simulate air motions over a broad range of scales, ranging from the synoptic scale to turbulent eddies.The Méso-NH configuration chosen here for LES modeling of marine stratocumulus uses an anelastic system of equations (Lipps and Hemler, 1982) and a 3-D turbulence scheme with a one-and-a-half-order closure, i.e. prognostic turbulent kinetic energy (TKE) and a diagnostic mixing length (Deardorff, 1980).The conservative variables, liquid water potential temperature θ l and total water specific content w t , are advected with a positively definite second order centered scheme.The surface sensible and the latent heat fluxes are proportional to the difference in temperature and specific humidity between the ocean and the air situated just above the surface.The coefficient of proportionality is derived by taking into account the thermodynamic stability above the surface and the roughness length is derived following (Charnock, 1955).

The radiation scheme
The radiative transfer is computed using the ECMWF operational model radiation code (Morcrette, 1991).The radiation package is based on a two-stream formulation that solves separately the longwave (LW) and shortwave (SW) radiative transfers for independent model columns.The radiative fluxes are computed taking into account the absorptionemission of the LW radiation and the reflection, scattering and absorption of solar radiation by the atmosphere and by the earth surface (Morcrette, 1991).In the LW, the radiative transfer is parameterized with a broad band flux emissivity method (Morcrette et al., 1986), while the Delta-Edington approximation is used for the SW (Joseph et al., 1976).
The cloud LW optical properties are computed using the Savijarvi and Raisanen (1998) parameterization.In the SW, the cloud optical thickness τ and the asymmetry factor g are computed following Fouquart (1987).τ is expressed as a function of the LWP and of the cloud droplet effective radius r eff : where ρ w is the density of water.g is set to 0.85 and 0.92 respectively in the visible and the near-infrared ranges of the spectrum.

The microphysical scheme
The model includes a 2-moment bulk microphysical scheme based on the parameterization of KK00, which was specifically designed for LES studies of warm stratocumulus clouds.The limit between cloud droplet and drizzle drop is defined at 25 µm in radius.Four prognostic variables are used for the condensed water: the cloud droplet and drizzle drop concentrations N c and N r , and the cloud droplet and drizzle drop mixing ratios q c and q r .
A fifth prognostic variable is used to account for already activated CCN, following the activation scheme of Cohard and Pinty (2000).This scheme may be seen as an extension of the Twomey (1959) parameterization for more realistic activation spectra (Cohard et al., 1998).The number of CCN, activated at any time step, is equal to the difference between the number of CCN which would activate at the diagnosed pseudo-equilibrium peak supersaturation S max in the grid and the concentration of already activated aerosols N CCN .The aerosols are assumed to be lognormally distributed and the activation spectrum is prescribed as: where N CCN is the concentration of activated aerosol, F (a, b, c, x) is the hypergeometric function, k, µ and β are parameters that can be tuned to represent various aerosol types and C (m −3 ) is the concentration of aerosol (Cohard et al., 2000).S max is diagnosed using vertical velocity and temperature.
The condensation/evaporation rate is derived using the Langlois (1973) saturation adjustment scheme.The cloud droplet sedimentation is computed by considering a Stokes law for the cloud droplet sedimentation velocity and by assuming the cloud droplet size distributions n c (r) fit a generalized gamma law (Cohard and Pinty, 2000): where r is the radius and λ is the slope parameter, which is a function of cloud droplet concentration N c and cloud droplet mixing ratio q c .The parameters α and ν were adjusted using droplet spectra measurements from the ACE-2 database and were set at α=3, ν=2.
The evolution of the prognostic variables by autoconversion, accretion, drizzle precipitation and evaporation is parameterized following KK00.The autoconversion rate and the accretion rate for mixing ratio are expressed respectively as: where N c is given in cm −3 , and q c and q r are in kg kg −1 .
The autoconversion rate and accretion rate of the cloud droplet concentration are defined by assuming that all cloud droplet radii are equal to the mean volume radius of the distribution.The autoconversion rate of the drizzle drop concentration is defined by assuming that all new drizzle drops have a radius equal to 25 µm.
Terminal velocities V N r and V q r , respectively of drizzle drop number concentration and drizzle drop mixing ratio are parameterized as a function of the drizzle drop mean volume radius r vr : where r vr is given in µm and terminal velocities are in m s −1 .

Three dimensional simulations
The Méso-NH model is used here to simulate marine boundary layer clouds over a domain of 10 km×10 km, with a horizontal resolution of 100 m and a vertical resolution of 10 m in the cloud layer and varying from 10 m to 100 m in the free troposphere.Initial profiles and large scale forcing are taken from recent data bases, the 14-17 July 1987 FIRE case (Duynkerke et al., 2004), the ACE-2 26 June and 9 July 1997 cases (Brenguier et al., 2003) and the DYCOMS-II RF02 and RF07 cases, on 11 and 24 July 2001, respectively (Stevens et al., 2003).The simulated cloud fraction is close to unity in all cases.Nighttime and daytime simulations are performed and the initial profiles are slightly modified to generate various values of the LWP.The cloud droplet concentration is also modified by changing the C coefficient in the Cohard et al. (1998) parameterization of the CCN activation spectrum.Each simulation evolves slowly, after a spin-up period of 2.5 h, with successively increasing and decreasing values of the LWP.Samples are taken every 20 min and the cloud fraction averaged values of geometrical thickness, LWP, N act , and drizzle precipitation rate are stored.Overall, the database contains 215 samples with N act values ranging from 45 to 260 cm −3 and LW P values ranging from 20 to 225 g m −2 (Fig. 1).
As an illustration of the simulations, Fig. 2 shows the horizontal distribution of the LWP and of the rainwater path RWP, defined as the vertically integrated drizzle water content.Up to six cloud cells are growing in the domain.The differences between the LWP and RWP horizontal fields reflects the life cycle of the convective cells.Initially, the LWP increases until it reaches a value large enough for the onset of the droplet collection process.Then the RWP increases at the expense of the LWP.In the final stage, the LWP has been depleted, while the RWP remains noticeable until all precipitating particles have reached the ground or evaporated below cloud base.

Comparison with field experiments
In a first step, the results of the simulation are compared to the observations using the same parameters and the same power laws of H or LW P , and N , as proposed by Pawlowska and Brenguier (2003) for ACE2, Comstock et al. (2004) for EPIC, and van Zanten et al. ( 2005) for DYCOMS-II, respectively.
-For ACE-2, the precipitation rate is averaged over the simulated cloud depth.The droplet concentration, represented by N act and the cloud geometrical thickness represented by H g , are derived from the simulated fields using the same conditions (Sect.4.1) as in Pawlowska and Brenguier (2003).
-For EPIC, the precipitation rate is averaged in the model layer below cloud base, the LWP is averaged over the cloudy fraction of the domain, and <N > is averaged over cloudy fraction of the simulation domain.
-For DYCOMS-II, H is the mean simulated cloud thickness and the precipitation flux and <N > are calculated as for EPIC.
The comparison between the simulations (small triangles) and the observations (large grey symbols) is summarized in Fig. 3a, c and e, for ACE-2, EPIC and DYCOMS-II, respectively.The right hand side figures, Fig. 3b, 3d, and 3f, show the same comparison in a log-log scale to emphasize the smallest values.
The ACE-2 observations are well distributed along the proposed power law but they overestimate the simulations by an order of magnitude in precipitation rate.As already discussed in Sect.4, various aspects of the sampling strategy and data analysis contributed to an overestimation of the precipitation rate and an underestimation of the cloud geometrical thickness.
The EPIC data are more scattered than the ACE-2 ones, but they better fit the simulations, with a slight overestimation at low values and underestimation at the largest ones.During EPIC, both the precipitation rate and the LWP were accurately measured with a millimeter wavelength radar and a radiometer respectively.In contrast, the CDNC was poorly characterized.Indeed, day time estimates were derived from remote sensing and then extrapolated to the night time.It can therefore be reasonably speculated that the uncertainty in the characterization of the CDNC is the main source of discrepancies in the comparison.
Overall, the DYCOMS-II observations are well aligned along the proposed power law and closely fit the simulations.Indeed, this experiment benefits of both the accuracy of in situ measurements for characterizing the CDNC and remote sensing for estimating the cloud geometrical thickness and the precipitation rate, hence providing a very accurate and consistent data set.
7 Parameterization of the precipitation rate in stratocumulus clouds

Best fit to the simulations
As mentioned in Sect.3, the extension of CRM bulk microphysics parameterizations to the resolution scale of a GCM is not physically justified.Indeed, bulk parameterizations have been developed for small scale predictions of the LWC and the non-linearity of the onset of precipitation should prevent their application to values averaged over large scales.It is not either efficient because these parameterizations transpose the non-linearity of the internal cloud processes to the GCM scale, while such a non-linearity is smoothed out when averaged over the large number of convective cells that develop in a GCM grid.
Once the LWP, the cloud base and top altitudes have been predicted in a GCM column, instead of tuning a CRM schemes, it would be more efficient to apply an empirical parameterization of the precipitation rate at cloud base that expresses as a function of the large scale properties of the cloud system.During field experiments, it is sometimes more practicable to characterize the cloud geometrical thickness, as in ACE-2 and DYCOMS-II.This parameter however is diffi-cult to precisely predict in a GCM, especially when the vertical resolution is coarse.In contrast, the LWP is directly derived from the conservation of heat and moisture.Similarly, it is easier in field experiments to characterize the mean CDNC value <N > by spatially averaging in situ measurements.However, we rather recommend to use the droplet concentration following CCN activation N act , as defined by Brenguier et al. (2000).Indeed precipitation forms in the core of the convective cells where the CDNC is close to its initial value, while its spatially averaged value can vary significantly depending on the intensity of the mixing processes, without impacting the onset of precipitation.This assessment is corroborated by Fig. 5 in Pawlowska and Brenguier (2000), where the droplet mean volume diameter increases with height above cloud base as predicted with the adiabatic parcel model, using a CDNC value equal to N act .Morevover, parameterizations exist for the prognostic of the concentration of activated droplets from the aerosol properties and a diagnostic of the vertical velocity in convective cores, while there are still no parameterizations of the impacts of the mixing processes and precipitation scavenging on the CDNC to predict its spatially averaged value.
The set of simulations is therefore analyzed to establish a relationship between the precipitation rate at cloud base, the mean LW P and the droplet concentration after CCN activation, N act .The best fit is obtained using the following power law relationship: where R p , LW P , and N act are expressed in kg m −2 s −1 , kg m −2 , and cm −3 , respectively.In Fig. 4a (linear scale) and b (log scale), the simulated rain rates R s are plotted as a function of the parameterized ones R p (triangles).Overall, the above relationship closely fits the simulated cases, especially at small values of the precipitation rate.The larger dispersion at larger values is partly dues to the limited size of the simulated domain.The largest values of the precipitation rate correspond indeed to larger convective cells of a size approaching the size of the domain.When there is only one cell in the domain, the resulting precipitation rate reflects the cycle of precipitation formation.This is illustrated in Fig. 4 by two simulations for which the successive values (every 20 min) have been connected from the start of the simulation (S) to its end (E).Both follow a counter clockwise loop, with sequences of increasing precipitation rate at constant LWP, when precipitation develops, followed by a decreasing LWP at constant precipitation rate, when the cell collapses.Note however, that the mean value over each cycle (the center of the loop) is closer to the proposed fit.In the last step, the proposed relationship is compared to the empirical relationships derived from the three field experiments (Fig. 5).ACE-2 is represented by dotted-dashed lines, EPIC by long-dashed lines and DYCOMS-II by solid lines.Each relationship is plotted for two values of the CDNC, 50 cm −3 (black lines) and 250 cm −3 (red lines) respectively.Each curve is restricted to the range of LWP and precipitation rate values that were documented during each field experiment, as reported in Table 1.
As in Fig. 3, the simulations best agree with the DYCOMS-II fit, even though LW P and N act are used here instead of the geometrical thickness and <N>, as in the original relationship.Discrepancies appear at the limits, for very low or very high precipitation rates.For the largest values, it is difficult to conclude which relationship, the simulated one or the empirical one, is the most suited because the simulation domain used here is too limited.

Conclusions
Since GCM models attempt to simulate the aerosol indirect effect, and more specifically its impact on the life cycle of boundary layer clouds, it is crucial to improve the parameterization of rain formation.Instead of using CRM bulk parameterizations, after tuning their coefficients for larger scales, specific parameterizations should be developed, that represent the mean precipitation production from an ensemble of clouds.Data sets from three field experiments have recently been analyzed to establish empirical relationships between the drizzle precipitation rate on the one hand, and the LWP or the cloud geometrical thickness, and the CDNC on the other hand (Pawlowska and Brenguier, 2003;Comstock et al. 2004;van Zanten et al., 2005).
In this study we use numerical LES simulations to replicate these observations.A large variety of stratocumulus clouds have been simulated, with mean LWP values ranging from 20 to 225 g m −2 and CDNC values ranging from 50 to 250 cm −3 .The resulting drizzle precipitation rates vary from 10 −7 to 4•10 −5 kg m −2 s −1 .The results of the simulation agree with the analysis of the field experiments, especially the DYCOMS-II conclusions and, to a lesser extent, the EPIC ones.They also suggest that the precipitation rates in ACE-2 were significantly overestimated.
In ACE-2 and DYCOMS-II, the clouds were characterized by their geometrical thickness, while the LWP was used in EPIC.For the parameterization of precipitation in GCMs, we recommend to use the LWP instead of the geometrical thickness, because the LWP is directly derived from the thermodynamics, while the prediction of the geometrical thickness requires additional hypotheses about the cloud fraction and the cloud adiabaticity.
In ACE-2, the CDNC was represented by the initial value N act , after CCN activation and before dilution by mixing and precipitation scavenging.In EPIC and DYCOMS-II, the CDNC was represented by its mean value <N > over the cloud system.For the parameterization of precipitation in GCMs, we recommend to use N act , because it can be predicted using a CCN activation parameterization, while the prediction of the cloud system mean CDNC requires a CDNC prognostic scheme that accounts for the impact of entrainment mixing processes and precipitation scavenging.
Such a relationship between R, LW P and N act has been derived from the set of simulations.It closely agree with the relationships derived from EPIC and DYCOMS-II, after replacing the geometrical thickness by the LWP, assuming adiabatic LWC profiles.
The similarity between the results of three field experiments conducted in three different geographical areas, suggest that such relationships are not fortuitous, and the LES simulations reported here demonstrate that they are physically based.

Fig. 1 .
Fig. 1.Parameter space of LW P and N act values simulated in this study.Each data point represents values averaged over the cloudy fraction of each cloud scene.

OFig. 3 .Fig. 5 .
Fig. 3.Scatter plot of the drizzle precipitation rate averaged over the cloud layer as a function of H 4 /N , for the ACE-2 measurements (grey triangles) and for the simulations (black triangles).Linear scale in (a) and log scale in (b).Scatter plot of the averaged drizzle precipitation rate at cloud base as a function of (LW P /N ) 1.75 , for the EPIC measurement (grey squares) and for the simulations (black triangles).Linear scale in (c) and log scale in (d).Scatter plot of the averaged drizzle precipitation rate at cloud base as a function of H 3 /N, for the DYCOMS-II measurement (grey circles) and for the simulations (black triangles).Linear scale in (e) and log scale in (f).

O
. Geoffroy et al.: Parameterization of the mean drizzle rate in stratocumulus 7.2 Comparison with field observations

Table 1 .
Summary of observational results.