Aerosols influence cloud and precipitation development in complex ways due to myriad feedbacks at a variety of scales from individual clouds through entire storm systems. This paper describes the implementation, testing, and results of a newly modified bulk microphysical parameterization with explicit cloud droplet nucleation and ice activation by aerosols. Idealized tests and a high-resolution, convection-permitting, continental-scale, 72-h simulation with five sensitivity experiments showed that increased aerosol number concentration results in more numerous cloud droplets of overall smaller size and delays precipitation development. Furthermore, the smaller droplet sizes cause the expected increased cloud albedo effect and more subtle longwave radiation effects. Although increased aerosols generally hindered the warm-rain processes, regions of mixed-phase clouds were impacted in slightly unexpected ways with more precipitation falling north of a synoptic-scale warm front. Aerosol impacts to regions of light precipitation, less than approximately 2.5 mm h−1, were far greater than impacts to regions with higher precipitation rates. Comparisons of model forecasts with five different aerosol states versus surface precipitation measurements revealed that even a large-scale storm system with nearly a thousand observing locations did not indicate which experiment produced a more correct final forecast, indicating a need for far longer-duration simulations due to the magnitude of both model forecast error and observational uncertainty. Last, since aerosols affect cloud and precipitation phase and amount, there are resulting implications to a variety of end-user applications such as surface sensible weather and aircraft icing.
It is well known that aerosols affect cloud microphysics through their role in nucleating cloud and ice particles. An increase in aerosol concentration generally leads to more numerous, but smaller, droplets for a given liquid water content, which results in an increase of the cloud albedo, known as the first indirect effect (Twomey 1974). Further, because of decreases in cloud droplet size, precipitation processes can be delayed and reduced, which is referred to as the second indirect effect (Albrecht 1989). However, numerous feedbacks and interactions with the ice phase and other aspects of cloud dynamics make it difficult to tease out exactly how cloud microphysical changes due to aerosol changes affect the radiative balance, precipitation, and dynamics in a systematic and quantitative way (cf. Levin and Cotton 2009).
The role aerosols play in altering warm-phase clouds has been intensively studied for multiple decades, but, until recently, less attention has been devoted to aerosols affecting mixed-phase clouds. Whereas liquid-only clouds tend to be somewhat simpler to examine, aerosols may impact mixed-phase clouds by changing the overall population and/or size of droplets that potentially alter freezing (Bigg 1953) and riming (Saleeby et al. 2009) processes as well as the vertical profile of latent heat release (cf. Khain et al. 2008). Modeling studies that focused on single-convective cloud systems or simulations performed for short time periods have found precipitation differences from a few to several hundred percent due to aerosols [e.g., review by Tao et al. (2012)] including either/both increases or decreases in precipitation. To complicate matters, some precipitation sign differences may be responding to differing environmental conditions. For example, dry continental versus moist maritime convective clouds respond differently to changes in aerosols (Khain et al. 2008, 2009; Teller and Levin 2006; Cui et al. 2011) and even the environmental wind shear was found to play a role in how aerosols affect convective clouds and resulting precipitation (Fan et al. 2009; Lee et al. 2012; Lebo and Morrison 2013).
Recently, large-scale, high-resolution, and long-duration modeling studies have been conducted (Seifert et al. 2012; van den Heever et al. 2011; Grabowski and Morrison 2011) and found that aerosol impacts to cloud systems interplay with the dynamics in a “naturally buffered” system. Even when relatively large changes in aerosols were simulated, the resulting surface precipitation differences were only a few percent overall; although larger impacts may occur locally.
Other mixed-phase cloud types such as orographic clouds may reveal systematic precipitation impacts in varying aerosol conditions. For example, Saleeby et al. (2009) found a shift in the location of the precipitation, with a reduction on the windward slope and increase on the leeward side of a mountain barrier in wintertime orographic clouds. Although mountain range total snow amount remained mostly unchanged, the distribution over the crest of a mountain range potentially impacts specific water basins. Similar findings were reported by Igel et al. (2013) in relation to aerosol impacts on precipitation in the vicinity of a warm front: precipitation reduced near the front but increased farther northward as aerosol concentration was increased. For the orographic cloud system, Saleeby et al. (2009) argued that the increase in precipitation on the leeward side was attributed to reduced riming on ice crystals due to reduced water droplet size in the more polluted conditions, which allowed the crystals to be carried farther downwind before reaching the ground. The warm frontal study also showed a less efficient snow riming process, but the precipitation increase distant from the warm front was attributed to increased accretion of cloud droplets by rain as aerosols increased the droplet number and liquid water content (LWC) but decreased overall droplet size.
Changes to cloud properties by aerosols are not only important to radiation, precipitation, and dynamics but also to any weather applications in which the phase and amount of water and ice content may be highly susceptible to small changes. For example, the amount of aircraft icing is directly dependent on the LWC and size of supercooled water droplets. Rosenfeld et al. (2013) attributed frequent incidences of aircraft icing near the U.S. West Coast to clean maritime air with low concentrations of cloud condensation and ice nuclei and stressed the importance of including aerosols when modeling aircraft icing.
To address a complex and uncertain problem that affects storms from convective to synoptic scales, the Thompson et al. (2008) bulk microphysics scheme was recently updated to incorporate aerosols explicitly in a simple and cost-effective manner. The scheme nucleates water and ice from their dominant respective nuclei and tracks and predicts the number of available aerosols. Using the Weather Research and Forecasting (WRF) Model (Skamarock and Klemp 2008), the scheme was tested in a high-resolution (4-km grid spacing) simulation of a 3-day winter storm event over the entire contiguous United States. The previous version of the scheme is widely used and well tested in WRF for quantitative precipitation forecast (QPF) applications (Rasmussen et al. 2011; Ikeda et al. 2010; Molthan and Colle 2012) because it consistently compares well to precipitation measurements (Liu et al. 2011) and was developed with inflight aircraft and ground icing applications in mind (Kringlebotn Nygaard et al. 2011; Podolskiy et al. 2012). Therefore, we believe it is well suited to address potential connections of aerosol impacts to cloud properties that subsequently affect radiation, precipitation amount and type, and aircraft icing.
This paper is organized as follows: A description of the numerical model is found in the next section along with more detailed descriptions of the activation of water and ice by two aerosol species. Results of the newly coupled aerosol–cloud physics parameterization tested under idealized flow conditions are presented in section 3. Next, a synoptic-scale, multiday winter cyclone is presented in section 4 along with results from a suite of high-resolution, continental-scale model simulations including sensitivity experiments using different aerosol number concentrations. The final section contains a summary and conclusions.
2. Numerical model
The simulations in this study were performed using the WRF Model, version 3.4.1, with modifications discussed below. The WRF Model includes many choices for various physical parameterizations of radiation, boundary layer, microphysics, and land surface interactions, but we avoided the use of a cumulus parameterization by applying a high-resolution grid as discussed in section 4. Most pertinent to this study, we used the Thompson et al. (2008) bulk microphysics scheme that treats five separate water species: cloud water, cloud ice, rain, snow, and a hybrid graupel–hail category. To minimize computational cost, prior versions of this scheme utilized one-moment prediction of mass mixing ratio for some species (cloud water, snow, and graupel) mixed with two-moment prediction (addition of number concentration) of cloud ice and rain. The number concentrations of single-moment species could be diagnosed from mixing ratio and various diagnostic relations between size distribution shape and other parameters found in the scheme. Cloud water was assumed to follow a generalized gamma distribution with a diagnostic but variable “shape parameter” based on a predetermined value of droplet number concentration set in the code, which was always intended to be changed by users for specific cases.
The scheme has now been updated to incorporate the activation of aerosols as cloud condensation (CCN) and ice nuclei (IN) and, therefore, explicitly predicts the droplet number concentration of cloud water as well as the number concentrations of the two new aerosol variables, one each for CCN and IN. Rather than determine a priori the specific aerosol types and chemical composition of multiple aerosol categories, which can lead to high computational burden and significant complexity, we simply refer to the hygroscopic aerosol as a “water friendly” aerosol (Nwfa) and the nonhygroscopic ice-nucleating aerosol as “ice friendly” (Nifa), although the latter is primarily considered to be dust. Consistent across all forms of water species (vapor, liquid, or solid), each species mass mixing ratio or number concentration follows the same governing conservation equation:
where Φ is mass mixing ratio or number concentration of any water species, t is time, ρ is the air density, U is the 3D wind vector, z is height, is the appropriately weighted fall speed of Φ, δΦ represents the subgrid-scale mixing operator, and represents the various microphysical process rate terms. Descriptions of the numerous process rate terms for previously existing species are found in Thompson et al. (2004, 2008), while the terms for newly predicted variables of cloud droplet number concentration Nc and the number of each aerosol species Nwfa and Nifa are provided in Eqs. (2)–(4) below, along with more detailed descriptions of specific terms found in subsequent paragraphs:
As compared to the prior Thompson et al. (2008) scheme with eight microphysics species to advect–predict, the new scheme with its three additional variables increases computational cost by approximately 16%. The most significant increase in computing time is due to the advection of new species, not the additional coding of various source–sink terms. In contrast, the simplest of Weather Research and Forecasting Model with Chemistry (WRF-Chem) options available at the time of writing increases the number of variables by over a factor of 2, which would massively impact computer memory and time. The subsections below describe the aerosol activation methods and input aerosol dataset used in simulations discussed in subsequent sections.
a. Cloud droplet nucleation
Cloud droplets nucleate from explicit aerosol number concentration (Nwfa) using a lookup table of activated fraction determined by the model’s predicted temperature, vertical velocity, number of available aerosols, and predetermined values of hygroscopicity parameter (0.4 in experiments performed in this research) and aerosol mean radius (0.04 μm). The lookup table was created by explicit treatment of Köhler activation theory using these five variables within a parcel model by Feingold and Heymsfield (1992) with additional changes by Eidhammer et al. (2009) to use the hygroscopicity parameter (Petters and Kreidenweis 2007). This implementation follows the same basic structure used by the Regional Atmospheric Modeling System (RAMS) as described by Saleeby and Cotton (2004, 2008) in which an assumed lognormal distribution with different values of aerosol mean radius and a constant geometric standard deviation (1.8) were preset as parameters when creating the table. The activation of aerosols as droplets in the new scheme is done at cloud base as well as anywhere inside a cloud where the lookup table value is greater than the existing droplet concentration. Nucleation of new droplets is prevented when existing ice hydrometeors would otherwise grow by water vapor deposition in a single time step that causes sufficient depletion of vapor to result in water subsaturated conditions; however, this rarely occurs in most updrafts because ice growth processes are relatively slow. Upon nucleation, the participating aerosols are removed from the population [third term in Eq. (3)], though they can be restored as regenerated aerosols, to the parent category via hydrometeor evaporation in which one aerosol is returned to Nwfa for each cloud or rain drop evaporated [represented by the fourth term in Eq. (3)]. Furthermore, aerosols are removed by precipitation scavenging [first term in Eq. (3)] using aerosol collision efficiencies computed following Wang et al. (2010) using a standard geometric sweep-out such as performed by other parts of the microphysics such as rain collecting cloud water.
For this study, any effects of subgrid turbulence on vertical velocity and nucleation of water droplets or ice were neglected; however, the newly added variables of aerosol number concentrations were mixed consistently with heat, moisture, and momentum fluxes produced by the boundary layer parameterization [represented by δΦ in Eq. (1)]. The simulations discussed below in section 4 used relatively high-resolution grid spacing of 4 km and primarily included well-organized clouds forced by large-scale ascent that suffice to exclude subgrid-scale vertical motions; however, simulations with coarser resolution should consider potential contributions by a distribution of vertical velocities within a model single grid box. Possible alternatives to relate a distribution of vertical velocities coupled with model-predicted turbulent kinetic energy (TKE) or eddy diffusivity variables to nucleate droplets (Ghan et al. 1997, 2011; Morrison and Pinto 2005; Morrison and Gettelman 2008; Morales and Nenes 2010) were bypassed to keep the new version consistent with the simpler one-moment cloud water scheme. Likewise, the RAMS simulations of Saleeby and Cotton (2004, 2008) and Igel et al. (2013) activated aerosols as droplets using grid-scale velocities only. A future version will likely incorporate the subgrid scales using guidance from large-eddy simulations (LES) to parameterize CCN activation due to turbulence, but we avoided this complication at this early stage. However, since the model may have a small downward vertical velocity and yet be fully saturated, although likely to be brief, the CCN activation by the lookup table assumes a minimum upward velocity of 1 cm s−1.
The water-friendly aerosol category was designed to be a combination of sulfates, sea salts, and organic matter because these aerosols represent a significant fraction of known CCN and are found in abundance in clouds worldwide. At this time, black carbon was ignored. More sophisticated aerosol treatments could be incorporated into future versions, but the competition for water vapor to nucleate cloud droplets with many aerosol constituents of unknown chemical composition is poorly understood and subject to further research before incorporation into a mesoscale numerical weather prediction model. Additionally, several studies have concluded that chemical properties are not nearly as important as the assumed aerosol number/size distribution (e.g., Dusek et al. 2006; Ward et al. 2010). The scheme is currently capable of representing different aerosol populations by altering the hygroscopicity and aerosol mean radius variables, although, for this study, these variables were held constant throughout. Additionally, the simulations presented here are intended to be sensitivity experiments using first-order approximately representative aerosol number concentrations, mean size, and hygroscopicity, while we do not claim to be forecasting precise aerosol amounts–composition.
b. Cloud ice activation
Cloud ice activates based on the number concentration of mineral dust aerosols since this species is considered to be highly active and most abundant naturally occurring ice nuclei in the atmosphere (DeMott et al. 2003; Cziczo et al. 2004; Richardson et al. 2007; Hoose et al. 2010; Murray et al. 2012). While other constituents may act as ice nuclei, the best direct correlation of activated ice crystals and aerosols acting as nuclei appears to be dust. Similar to CCN activation, the addition of more aerosol species acting as IN leads to unnecessary complications as multiple species compete for the water vapor in complex ways. A future version of the scheme may incorporate other ice nuclei when future research clearly indicates such a requirement. The number of dust particles that nucleate into ice crystals is determined following the parameterization of DeMott et al. (2010) when above water saturation to account for condensation and immersion freezing and by the parameterization of Phillips et al. (2008) when less than water saturated to account for deposition nucleation. In addition, the freezing of deliquesced aerosols using the hygroscopic aerosol concentration is parameterized following Koop et al. (2000), shown as the second term in Eq. (3).
The freezing of existing water droplets continues to follow the Bigg (1953) volume and temperature parameterization as previously used in Thompson et al. (2008), except that the dust aerosol concentration alters the “effective temperature” to freeze more (or fewer) water drops when more (or fewer) dust particles are present. This connection intends to increase nucleation by considering the higher likelihood of contact nucleation by Brownian motion causing a dust particle to come into contact with a supercooled water droplet. As currently implemented based on an inspection of typical background dust concentration of 0.1 particles per liter of air, there is no alteration of the freezing of water droplets due to dust when compared against the prior Bigg’s freezing implementation in the Thompson et al. (2008) scheme. However, for each order of magnitude increased (decreased) number concentration, the effective temperature for freezing of droplets is decreased (increased) by 1°C. Quantitatively, we have no basis for instituting this ad hoc method of altering the water drops freezing point by the presence of dust, only qualitative belief that some effective increase in freezing water drops due to the presence of high dust concentration is likely due to an increased likelihood of freezing by contact nucleation or immersion via an embedded dust nuclei inside of water drops. All of the ice nucleation mechanisms by dust can be readily switched off in favor of using the previous ice nucleation scheme as found in Thompson et al. (2008). Separate ice nucleation sensitivity experiments are beyond the scope of this study and will be reported in the future since this study focuses exclusively on aerosols acting as CCN, except a single test of the old versus new ice nucleation techniques was performed for an idealized test discussed in section 3. However, since the freezing of water drops contains explicit dependence on their size (volume), there are implicit links to aerosol sensitivities found in the mixed-phase region discussed in detail below (section 4), even without altering the ice nucleation methods explicitly.
c. Aerosol input data
The aerosol number concentrations in the winter storm simulations in section 4 were derived from multiyear (2001–07) global model simulations (Colarco et al. 2010) in which particles and their precursors are emitted by natural and anthropogenic sources and are explicitly modeled with multiple size bins for multiple species of aerosols by the Goddard Chemistry Aerosol Radiation and Transport (GOCART) model (Ginoux et al. 2001). The aerosol input data we used included mass mixing ratios of sulfates, sea salts, organic carbon, dust, and black carbon from the 7-yr simulation with 0.5° longitude by 1.25° latitude spacing. We transformed these data into our simplified aerosol treatment by accumulating dust mass larger than 0.5 μm into the ice-nucleating, nonhygroscopic, mineral dust mode Nifa and combining all other species besides black carbon as an internally mixed cloud droplet–nucleating, hygroscopic, CCN mode Nwfa. Input mass mixing ratio data were converted to final number concentrations by assuming lognormal distributions with characteristic diameters and geometric standard deviations taken from Chin et al. (2002, their Table 2).
For simplicity, we implemented a variable lower boundary condition that represents aerosol emissions based on the starting near-surface aerosol concentration and a simple mean surface wind to calculate a flux (constant through time) using the following relation applied only to the model lowest level, represented by the last term in parenthesis in Eq. (3): dNwfa/dt = 10[log(Nwfa)−3.698 97] which results in 0.01 × 106 kg−1 s−1 for Nwfa = 50 cm−3, 0.1 × 106 kg−1 s−1 for Nwfa = 500 cm−3, and 1.0 × 106 kg−1 s−1 for Nwfa = 5000 cm−3, for example. A 3-day averaging test revealed that the aerosol number concentration remained very close to the climatological condition over most of the domain, revealing that this simple assumption is more advanced than holding initial aerosol concentration constant as other studies have done (Igel et al. 2013). An earlier test that held only the lowest model level constant in time everywhere led to an obvious and unrealistic domainwide increase in aerosols over the 3-day simulation. Our oversimplification can be remedied in future versions using more explicit aerosol emissions inventories or coupling with a full chemistry model such as WRF-Chem (Grell et al. 2005) or WRF–Community Multiscale Air Quality (CMAQ; Wong et al. 2012).
In this application, the aerosols represent a monthly climatology sufficient for running a series of sensitivity experiments. It was not our intent to produce a proper simulation of the aerosol conditions of a particular event since such measurements are not widely available in space and time over a scale of the simulations in section 4. Samples of the climatological aerosol dataset are shown in Fig. 1 and were interpolated to the WRF Model horizontal and vertical points for initial and lateral boundary condition data.
3. Idealized tests
For fundamental testing, WRF was configured using simple two-dimensional flow over a hill similar to tests of prior versions of the scheme described in Thompson et al. (2004, 2008). We used a domain with 600 points spaced 1 km apart and 72 vertical levels spaced from 50 m near the surface to 250 m near the midtroposphere, and the model was integrated for 6 h. To test the warm-rain process in complete isolation from potential complications of the ice physics, the initial temperature profile was warmed to avoid cloud tops reaching glaciation temperatures. Horizontal wind linearly increased from 0 m s−1 at the surface to 10 m s−1 at 1 km and above and impinged on a 1-km high and 25-km half-width mountain barrier that produced a maximum updraft velocity of 0.2 m s−1 (see Fig. 2). Other sensitivity tests with a steeper 10-km half-width mountain increased the maximum updraft velocity to 0.5 m s−1 in order to test higher aerosol activation rates. Two initial aerosol conditions with exponentially decreasing profile of concentration from the surface to 2 km and constant amount above were used to test aerosol sensitivity. In one experiment, a surface aerosol concentration of 250 cm−3 decreased to 50 cm−3 aloft, which might be typical of a clean maritime air mass, whereas a second experiment started with 1000 cm−3 near the surface and decreased to 250 cm−3 aloft (see Fig. 2), which might be more typical of a continental air mass. Additional experiments including the ice phase were performed in which the thermodynamic profile was cooled to match the same temperature and moisture profile used in Thompson et al. (2008); otherwise, the conditions shown in Fig. 2 were maintained, but only the steeper mountain profile was used. To refer to the sensitivity experiments with abbreviated reference names, we use the following nomenclature: warm or cold describes the simulations excluding and including ice phase, respectively; 25 and 10 refer to the mountain half-width; and Mar and Con refer to the maritime and continental aerosol concentrations, respectively, as shown in Table 1.
As expected, all simulations with higher aerosol concentration caused a corresponding increase in cloud droplet number concentration since the updraft strength and attendant LWC remained nearly constant. Figure 3 shows that the low-aerosol-concentration experiments had about twice as many grid points with droplet concentrations below 50 cm−3 compared to simulations with a higher number of aerosols. Also, the continental experiments produced a flat range of droplet concentrations from 25 to 200 cm−3 because the updrafts were relatively weak, while the maritime experiments produced no droplet concentration exceeding 100 cm−3. Table 1 shows that the computed mean cloud droplet concentration remained quite low, on average, at only 54 to 72 cm−3 for the continental experiments, which was 2–3 times larger than the 25 cm−3 of the maritime experiments.
Since the LWC was nearly constant regardless of aerosol concentration, the larger droplet concentration in the continental experiment must contain smaller overall droplet sizes that greatly affect the ability of the cloud to form rain from the collision–coalescence process. This is readily confirmed in Table 1 that shows the median size of cloud droplets was only 9.5 μm in the continental experiments versus 13.5 μm in the maritime experiments. Furthermore, we see that the maritime experiment with the steeper mountain slope was the first to produce rain (16:40), taking roughly half the time that was needed in the higher aerosol number concentration of experiment continental (27:55). The broader, 25-km, half-width mountain required nearly twice the time in each experiment to produce rain due to its weaker updraft. In other words, the strong updrafts associated with the steeper terrain simply supply condensing water at a more rapid pace that enhances droplet growth to rain sizes far more effectively than impacts due to changing aerosol concentrations, droplet number, or size combined with the weaker updrafts.
When the temperature profile was cooled to match the sounding used in Thompson et al. (2004), the simulation produced a cloud with a temperature of −13°C at its top and ice formed in the simulations. However, rather than initiating ice solely from a temperature-dependent function following Cooper (1986), the mineral dust aerosol was responsible for ice initiation as described earlier. The inclusion of ice roughly halved the time to produce precipitation from 16 min 40 s to 8 min 50 s in the maritime experiments or from 27 min 55 s to 17 min 55 s in the continental experiments. The smaller overall number concentration and mean size of cloud droplets was due to riming of droplets onto snow. Another test (not shown) in which the scheme was changed back to the original Cooper (1986) ice nucleation did not alter the precipitation timing or amounts noticeably. A final experiment (not shown) in which the dust aerosol was increased by a factor of 3 between 3 and 6 km, as shown in Fig. 2, also had a negligible effect on precipitation. These additional tests do not reveal significant impacts of ice initiation sensitivity because the cloud was simply too shallow and warm to contain significant ice. Analysis of ice sensitivities remains as future work.
4. Winter cyclone simulations
Between 31 January and 2 February 2011, a large extratropical winter cyclone developed in the central United States and moved eastward across the Appalachian Mountains. With this storm came a variety of surface weather including near-record low temperatures and light snow in the northern high plains; regions of lake-effect snow in the Great Lakes; moderate to heavy snowfall in the central plains; a mixture of snow, ice pellets, and freezing rain from the Ohio River Valley to New England; and moderately strong convection in the Southeast in advance of the cold front. In addition, a weak upper-air low pressure system moved slowly eastward across the desert Southwest region producing mostly light snowfall in the southern Rocky Mountains. In the central part of the country, the storm has been called the “Groundhog Day Blizzard” [National Climatic Data Center (NCDC); http://www.ncdc.noaa.gov/billions/events.pdf], since 1–2 ft (30.5–61 cm) of snow paralyzed the city of Chicago and caused an estimated $1.8 billion worth of damage along with 36 deaths in a multistate region.
For each 6-h period within the 3-day storm, between 350 and 500 surface weather reporting sites recorded a trace or more of precipitation and nearly 1000 sites across the whole country recorded precipitation in the 3-day period. Because of the widespread impact of the storm, and myriad of cloud and precipitation forcing mechanisms, including synoptic, mesoscale, and orographic, we believe the event is well suited for extensive modeling sensitivity experiments.
a. Model configuration
WRF was configured with a single convection-permitting grid of 4-km horizontal spacing with 1200 × 825 grid points covering the entire contiguous United States, also using reasonably high vertical resolution with 72 levels up to model top at 73 hPa with stretched spacing from 50 m near the surface to 750 m near the tropopause. Input and lateral boundary condition atmospheric data were supplied by the Rapid Update Cycle (RUC) model analyses every 3 h. The simulations utilized the Noah land surface model (Barlage et al. 2010), Yonsei University planetary boundary layer scheme (Hong et al. 2006), and the Rapid Radiative Transfer Model (RRTM-G; Iacono et al. 2000) radiation scheme and no convective parameterization, since the grid was sufficiently high resolution to predict most clouds explicitly.
At the time of writing, no existing radiation scheme in publicly available WRF code utilizes fully coupled effective radii of all hydrometeor species as known within the microphysics scheme into the radiative computations involving clouds, which is insufficient for performing aerosol–cloud sensitivity experiments. Therefore, this missing link was remedied by explicitly computing the effective radii of cloud water (cf. Slingo 1989), ice, and snow (cf. Stephens et al. 1990) directly in the microphysics scheme and passing those values to the RRTM-G scheme to calculate the cloud optical depth parameter. At present, the sulfates, sea salt, organic carbon, and dust aerosols used by the microphysics scheme to activate water droplets and ice crystals do not scatter or absorb radiation directly, and only the typical background amounts of gases and aerosols present within the RRTM-G scheme were considered for scattering–absorption–emission of direct radiation in this study.
Figure 4 shows results of the WRF Control simulation at 42 h, valid at 1800 UTC 1 February 2001 and reveals broad regions of snow (blue color fill: 1-h snow amount) and rain (green color fill: 1-h rain amount) with an overlapping region of both plus graupel (red color fill: 1-h graupel amount). The gray-shaded regions represent the accumulated total precipitation amount thus far in the simulation, and the various red- and blue-shaded dots represent the difference between observed and WRF-simulated 6-h precipitation. The storm obviously impacted a very large portion of the United States and includes very typical extratropical cyclone characteristics of a synoptic-scale warm and cold front as well as less obvious orographic forcing, lake-effect snow, and convection.
b. Sensitivity experiments
A suite of sensitivity experiments was run to test in a robust and comprehensive manner the physics of the new aerosol-aware microphysics in contrast to the non-aerosol scheme as well as the impact of changing the amount of aerosols on cloud and precipitation development. First, to create a set of benchmark tests, the non-aerosol scheme with original one-moment cloud water (Thompson et al. 2008) was run with constant and extremely low droplet concentrations of 50 cm−3, followed by a moderately high value of 750 cm−3. These simulations represent very clean versus moderately polluted conditions and provide reasonable bounds for the simulations using explicit aerosols. We will refer to these as Nc50 and Nc750, respectively. Second, the aerosol-aware scheme was run with the input and lateral boundary condition data as described in section 2c with the first simulation representing aerosol conditions that should be representative of conditions present in the current era. This simulation will be referred to as Control. Next, a simulation that reduced at all model grid points the aerosol number concentration to one-tenth the Control concentrations (Clean) and a final simulation with 10 times the Control concentrations (Polluted) were performed. Changes to aerosol characteristics such as chemical composition, hygroscopicity, or mean radii were not tested for these simulations, solely the aerosol number concentrations. Furthermore, there were no changes made to the nonhygroscopic aerosol (dust) number concentration in order to minimize any changes due to ice nucleation in these tests.
Although the benchmark tests used single values of droplet number concentration that were constant in space and time, the computed radiative effective size was fully coupled into the radiation scheme as described previously. If the aerosol-aware scheme produced results that varied wildly in comparison to the benchmark experiments with low and high droplet concentrations, then almost certainly an error in coding would be indicated. Furthermore, the benchmark experiments provide bounds to cloud, precipitation, and radiation properties and impacts for simulations where aerosols were explicitly introduced.
c. Cloud property impacts
Consistent with the results of the 2D idealized tests, the fully 3D simulation showed the expected result that number concentration of cloud droplets increased with increasing aerosol concentration. Along with the increase in droplet concentration, a very prominent decrease in the mean size of droplets was noted, and, since the radiation scheme properly accounted for the radiative effective radius, there was an absolute indication of the first aerosol indirect effect: the “cloud albedo” effect (Twomey 1974). Figure 5 shows mostly positive differences of cloud droplet concentration (top panel; warm colors), mostly negative differences of mean effective radius (middle panel; cool colors), and mostly increased outgoing shortwave radiation (bottom panel; warm colors) when subtracting the less-polluted Clean simulation from the higher-aerosol-concentration Control simulation. Numerically, the average difference of reflected shortwave radiation in these two simulations was a 5.4% increase in cloud albedo due to higher aerosol concentration of the Control versus Clean experiments when computed from 6-h intervals during daylight hours. This behavior was entirely consistent when any of the experiments with lower aerosol or droplet concentration was subtracted from a corresponding experiment with more aerosols. Likewise, consistent behavior was found in the difference of longwave radiation reaching the ground below clouds (average 0.47% increase) as well as top of the atmosphere outgoing longwave radiation (average 0.11% decrease), although not shown.
d. Water droplet distribution changes
It is difficult to encapsulate all of the changes to water droplet sizes and amounts for a series of simulations with millions of spatial grid points over 72 h, but we believe the next set of three figures best illustrates the changes to water droplet distributions as aerosol number concentrations were changed. In Figs. 6–8, we plotted a random sampling of points containing any liquid water, either cloud droplet or rain, in terms of their median volume diameter (MVD) versus LWC. On the left portion of each figure, cloud droplets are shown with a linear MVD scale, while points with rain are shown using a logarithmic scale on the x axis. Each dot is color coded by temperature with gray dots for any temperature value above 0°C, then red, orange, green, and blue dots for each 10° increment below 0°C. This color coding provides insights into possible changes to size as well as frequency of finding water drops in specific temperature ranges in the mixed-phase region as aerosols were changed in the various experiments. In addition, the solid black lines on the left portion of Figs. 6–8 represent the results of the benchmark simulations, Nc50 and Nc750, while we omit the rain drops because they are redundant with those found in the other two figures. Nc50 and Nc750 collapse to a single line because a constant number concentration of droplets gives only one value of MVD for any particular LWC using the simple mass–diameter power-law relation of m(D) = aDb.
Note in Fig. 6, created from the Control simulation, that nearly all points containing cloud water lie within the bounds of the benchmark simulations, and note how the highest LWC and largest MVD correspond to the highest temperatures. Also note how the number of grid points of cloud water decreases sharply with decreasing temperature, which we would expect since liquid water more likely freezes as temperature decreases and larger drops freeze before smaller drops in general (Bigg 1953). Where the MVD of rain is relatively small, the corresponding LWC is small, and the preponderance of these points were produced via the collision–coalescence process and subsequent accretion of other cloud droplets in the warm-rain process, while the narrowing diagonal region into higher LWC and larger MVD dominantly represents grid points of rain produced from melting ice, which would be expected to be larger.
Then, to see alterations to water distributions with the decreased aerosol number concentration in the Clean experiment, refer to Fig. 7 and note how the distribution of cloud droplets significantly shifts to the right side of the Nc50 line, indicating a notable increase of MVD and corresponding decrease in droplet number. Also note the upper extent of LWC as the larger mean size of water leads to more rapid rain production by enhancing the warm-rain processes, which is easily confirmed by the indicated number of rain points at all temperatures. In fact, a factor of 10 more grid points with rain between −20° and −30°C (green dots) appears along with a factor-of-4 increase between −10° and −20°C (orange dots) when reducing aerosols by a factor of 10 between Control and Clean. Note the larger y-axis vertical extent of LWC (rain) by colored dots between 200 and 400 μm in Fig. 7 compared to Fig. 6. A more subtle feature appears in the narrow diagonal region of rain with higher LWC and larger MVD as relatively fewer grid points appear in this region in the Clean experiment as compared to the Control experiment. We will refer to this narrowing region toward the upper right of these graphics as the “flame tip” and provide a physical connection for differences seen between Figs. 6 and 8 in the next subsection.
In the final sensitivity experiment (Polluted), in which the aerosols were increased by a factor of 10 more than Control, the most notable change of Fig. 8 is the dramatic shift of grid points with cloud water toward much smaller MVD and slightly higher LWC. This makes physical sense since the increased aerosol concentration is leading to smaller overall mean size of cloud water that subsequently hinders the warm-rain processes. There are also more cloud droplets surviving to lower temperatures due to their lower likelihood of freezing as their mean size decreases.
e. Precipitation impacts
The changes to water droplet populations by changing aerosols definitely resulted in changes to surface precipitation, but not in entirely obvious ways. Figure 9 shows the individual differences of rain, snow, and graupel amounts for the second day of the simulation between the Control and Clean experiments. Table 2 also contains precipitation amounts by type from all the experiments along with various differences and percentage change between high- and low-aerosol-concentration experiments. Other time periods (not shown) confirmed similar patterns. Overall, there are mixed signals of both increased and reduced rain and snow amounts due to evolution and location differences of narrow precipitation bands; however, the primary signals were a reduction of rain in the southeast portion and an increase of snow in the northern portion as aerosols were increased. The reduction of rain seems logical since the warm-rain processes were hindered by overall smaller droplets (Albrecht 1989), but the very widespread and obvious increase of snow with higher aerosol concentration was not expected.
We believe that the increase in snow was due to the generally reduced warm-rain processes in the southern United States permitting many more cloud droplets, albeit smaller, to be transported northward (and possibly lofted higher) into the snow-producing clouds found to the north. While the overall mean size of droplets may have been smaller when aerosols were more numerous, the geometric sweep-out of those droplets increases because there were so many more droplets to intercept as well as larger LWC, even though there was a general decrease in collision efficiency (Hindman et al. 1992) between snow and cloud water. This was confirmed by calculating the horizontal flux of cloud water crossing through four parallel WRF x–z (west–east) grid planes during four 6-h time periods on day 2. Note in Fig. 10 how the flux of water was largest through each plane and for each 6-h interval in the simulations with the highest aerosol concentration, and the flux was percentagewise larger in the x–z planes to the south and smaller to the north. An additional contribution to the increase of snow in the north was also possible from an enhanced Wegener–Bergeron–Findeisen process as some of the cloud droplets could have evaporated to vapor that subsequently migrated to the ice and snow; however, individual process rates were not captured during the model simulations to confirm this hypothesis.
Further evidence and confirmation that rain and graupel generally decreased while snow increased when aerosols were increased is provided in Fig. 11. Differences of individual rain, snow, and graupel precipitation amounts between experiments with higher aerosol concentration minus experiments with lower aerosol concentrations are shown for each day as well as the sum of all 3 days. The largest decrease in rain and corresponding increase in snow occurs between the experiments with the greatest difference of aerosol concentration (Polluted minus Clean). Comparisons between experiments with less drastic aerosol change produced less drastic precipitation differences, showing consistent and robust behavior of the aerosol effects. Furthermore, the decrease in rain amount exceeded the increase in snow and differences of graupel were quite small, but also consistently less graupel with increasing aerosols. We speculate that this was due to the overall reduction in number of points with rain and overall smaller droplet size that hindered freezing of rain drops into graupel particles in this scheme.
While the amount of rain reaching the surface decreased with higher aerosol concentrations, the most common reductions occurred primarily in association with extremely light precipitation. Figure 12 shows distributions of rain, snow, and graupel in precipitation bins of varying amounts for each hour of the 72-h simulation. Whereas the count of grid points with hourly rain went down as aerosols increased, there was hardly any noticeable change in counts of hourly amounts greater than 2.5 mm over the entire synoptic storm-scale regions. Similarly, decreases in the amount and frequency of light rain but not heavier rain was noted in relation to significant increases of aerosol concentration in an observational and modeling study over eastern China by Qian et al. (2009). Also, Sorooshian et al. (2010) found much greater impact of aerosols to light precipitation in contrast to heavier precipitation and attributed it to cloud thickness property since deep clouds offer plenty of opportunity for rain to accrete cloud droplets over a large cloud depth as compared to relatively thin clouds.
As mentioned in the previous subsection, there was another change to the mixed-phase precipitation region worthy of mention, although more subtle than the preceding noted effects. The flame tip region shown in Figs. 6–8 shows a decrease in grid points with relatively high LWC and large MVD in Clean (Fig. 7) as compared to Control (Fig. 6) with lesser differences seen between Control and Polluted (Fig. 8). The scattering of points oriented vertically along MVD = 300 μm was dominantly produced by warm-rain processes, whereas gray dots (T > 0°C) toward the flame tip were dominantly produced by melting snow/graupel. Consistent with the snow increase due to aerosol increase found in the region north of the warm front (Fig. 9), which was dominated by glaciated clouds filled with snow, it appears that increasing aerosols increased the overall size/mass of snow aloft that subsequently melted into rain before reaching the surface; however, Fig. 11 showed that the additional melted ice does not compensate for the loss of rain by warm-rain processes.
Another interesting aerosol effect in regions of mixed-phase surface precipitation is noted in Fig. 13. For any model grid point containing a mixture of rain and snow–graupel during an hour, we computed the fraction of liquid precipitation as rain/(rain + snow + graupel) and counted each 10% bin. After normalizing by the number of grid points with any precipitation, we found that as aerosols increased, there was a relatively higher fraction of liquid precipitation. One potential hypothesis for this effect is corollary to the increased snow to the north of the warm front, which is that the less efficient rain production in the south allowed more cloud droplets to transport northward into the zone of mixed-phase region near the warm front where rain accreted more cloud droplets simply due to a higher number of them, albeit smaller size, and resulted in a disproportionate increase in rain reaching the ground compared to graupel–snow. This hypothesis is supported by similar results seen in Igel et al. (2013), where they attributed slightly higher surface precipitation amounts approximately 150 km north of the warm front to higher rates of rain accreting cloud droplets.
A final aspect of precipitation was analyzed to determine if using a simpler microphysics scheme without aerosols and constant cloud droplet concentration or the new scheme with low, moderate, or high aerosol concentration produced any improvement as compared to observations. Unfortunately, errors in precipitation observations (Rasmussen et al. 2012) and errors in the model forecast at single sites (even 1000 sites with precipitation during a large-scale winter storm) far outweigh the scale or magnitude of changes seen in our five sensitivity experiments. Figure 14 illustrates that model forecast errors were rather large and extremely variable and each experiment produced very similar error statistics. In fact, our results indicate no statistically significant differences among the five experiments as evidenced by the overlapping means and confidence intervals shown in Fig. 14. And, since the fidelity of observed snow water equivalent data in automated precipitation measurements, especially during moderate to high winds, lacks credibility, we excluded most snow reports from the data used to create Fig. 14. As examples of the measurement problem, Quincy, Illinois, reported 559 mm of snow, yet only 2-mm liquid equivalent; Moline, Illinois, reported 467 mm of snow with 4-mm liquid equivalent; and Chicago–Midway airport reported 457 mm of snow with only 5 mm of liquid, yet Chicago–O’Hare airport reported 508 mm of snow and 41-mm liquid equivalent. Massive errors such as these are rampant in automated reporting stations during snowstorms in recent decades, and evaluators of model forecasts should remember to question observational data quality when assessing model performance. A massive number of the deep blue dots in Fig. 4 representing serious model overforecasts of precipitation are likely to be far lower error than it superficially appears. Regardless, when we exclude measurements that likely coincided with snow, we did find that our WRF simulation produced a noticeable bias of underforecasting the highest precipitation amounts, indicating frequently missed convective events combined with near-0 mean bias of light precipitation (<13 mm over 72 h) with a slight model overforecast problem for the light amounts.
To emphasize a main point about aerosols affecting precipitation amounts, even though aerosols changed the water size distributions as dramatically as seen in Figs. 6–8, which subsequently affected at least six microphysical processes including autoconversion, collection of cloud water by rain, snow, and graupel, and freezing of cloud water and rain, the accumulation of all these processes remained negligible as compared to combined errors in observations and model precipitation forecasts. Perhaps the only way to know for certain if the more complex physics with more realistic spatial and vertical distributions of aerosols improves forecasts of precipitation is to perform far longer integrations over months, seasons, or years.
To address a complex and uncertain research problem that affects storms from convective to synoptic scales, the Thompson et al. (2008) bulk microphysics scheme was updated to incorporate aerosols explicitly. The scheme nucleates water and ice from their dominant respective nuclei and fully tracks and predicts the number of available aerosols. Using the Weather Research and Forecasting (WRF) Model, the scheme was tested in a high-resolution (4-km spacing) simulation of a 3-day winter storm event over the entire contiguous United States. A Control simulation with climatological aerosol conditions and two sensitivity experiments with Clean (one-tenth) and Polluted (10 times) conditions were used to evaluate the magnitude of various aerosol–cloud–precipitation interactions. Additional experiments that ignored aerosols and used the older, one-moment cloud water prediction combined with constant (in space and time) high and low droplet number concentration revealed entirely consistent behavior with the aerosol-included experiments and gives credence to the robustness of results and physics of the scheme.
There were numerous notable and fundamental changes to water droplet size distributions and subsequent precipitation and radiation impacts from varying aerosol number concentration that were consistent with expected aerosol indirect effects. Increasing aerosol concentration produced consistently more droplets of overall smaller size that hindered the warm-rain processes (Albrecht 1989) and increased cloud albedo (Twomey 1974). When comparing the Control versus Clean experiments, the cloud albedo increased by 5.4% in the experiment with the higher aerosol concentration. Differences in longwave radiation reaching the ground due to cloud property changes were more subdued, as expected, increasing only 0.47% while outgoing longwave radiation to space decreased even less, −0.11%, due to cloud opacity changes by the different aerosol concentrations affecting droplets and ice crystal sizes.
Space- and time-integrated surface precipitation differences between experiments with more or fewer aerosols revealed rather modest effects overall (3%–8% reduction of rain, 2%–5% increase of snow) for this 72-h winter cyclone simulation. Findings in section 4e were consistent with a study by Igel et al. (2013) in which the precipitation amount in the immediate vicinity of a synoptic-scale warm front decreased slightly, whereas amounts north of the warm front increased. This was due to higher cloud droplet number concentration and LWC being transported northward as aerosol concentration increased and subsequent capture by falling snow and rain increased due to higher available LWC even though collisions efficiencies reduced due to smaller overall droplet size. This may have broad and important implications for overall water transport being affected by aerosols and provide shifts in precipitation patterns on a continental scale.
Although it is clear from Fig. 9 that very specific locations may have changed precipitation amounts more drastically, most of the changed rain regions involved shifts in location while the amounts nearly offset, especially in moderate to high precipitation bands, since Fig. 12 showed that only the lightest amounts of precipitation showed high susceptibility to aerosol changes. Therefore, we speculate, that if we simulated an entire season’s worth [similar to Seifert et al. (2012)] of real weather across an entire continent, most of the location shifts in precipitation due to different aerosol conditions would be likely to smear out with successive storms due to changing wind directions, convergence features, and dynamical interactions. The basic behavior of dominantly less rain and slightly more snow is a plausible outcome for numerous extratropical winter cyclones such as the one studied here, but we would expect only modest changes to surface precipitation from changing aerosol concentration when using reasonable estimates of typical aerosol concentrations and integrated over an entire season and a large region, especially considering our experiments used factors of 10 above/below the typical values.
Clear from Fig. 14 is that there were no statistically significant differences in the model’s surface precipitation forecasts when using different aerosol conditions and comparing to observations. We point out the following difficulties in verifying model forecasts of surface precipitation to validate sensitivity experiments such as ours:
errors in model cloud forecast timing and location may greatly outweigh differences among sensitivity experiments;
observational uncertainty can be massive, particularly with liquid equivalent snow measurements in blowing- snow conditions; and
sensitivity of aerosols to resulting precipitation is potentially weaker or stronger in models than what is truly found in nature but determining such a bias is exceedingly difficult.
Perhaps more important to validation efforts are the changes to water droplet distributions such as those illustrated by Figs. 6–8, although insufficient aircraft data exist to perform an objective analysis. However, the general trend of cloud droplet concentrations shown in the Control simulation (Fig. 6) as compared to the Clean (Fig. 7) or Polluted (Fig. 8) experiments gives at least subjective positive comparison to previously published aircraft data (e.g., Cober and Isaac 2012; Sand et al. 1984; Politovich and Bernstein 2002).
While extensive research continues to focus on aerosol effects on surface precipitation, this study also shows explicitly how aerosols affect the water droplet size distribution aloft. This is an important consideration for any inflight aircraft-icing applications because the liquid water content and size of drops are critical to the accumulation of ice on airplane control surfaces (Arenberg 1943). Therefore, using the data from these experiments, we performed relatively simple ice accretion calculations intended to predict aerosol effects on a final application to aircraft icing. The equations used for ice accretion on a standard cylinder followed Makkonen (2000) where the change in mass with time is a product of collision efficiency, LWC, velocity, and cross-sectional area of the cylinder (details found in the appendix). Using the WRF simulation data, we calculated a dM/dt value for any model grid point with either cloud water or rain at temperatures below 0°C each hour from 6 to 72 h from all five sensitivity experiments. Next, we calculated the frequency of occurrence of each order of magnitude bin of ice accretion rate, shown in Fig. 15. The figure shows that as aerosols were increased, there was generally an increase in ice accretion by cloud water (left panel) up until the largest ice accretion rate when the trend reversed direction. In contrast to the smaller cloud droplets, the effect of increasing aerosols generally reduced the ice accretion from larger rain drops (right panel) except in the highest category of ice accretion. These appear to be logical because an increase in aerosol concentration led to more numerous (but smaller) cloud water droplet number concentrations with higher LWC because of the hindered warm-rain production. The increase in LWC overcompensated for the decrease in cloud droplet size since collision efficiency of droplets decreases as size decreases since the smallest droplets pass around a moving object and follow the airstream rather than impinge on the surface of the cylinder–wing. The general decrease in frequency of ice accretion due to rain as aerosols increase follows from the decrease in grid points with rain as aerosols increased due to the reduced warm-rain processes. The increased frequency of small droplet icing may imply that more “rime” icing may occur as aerosols are increased, while the frequency of “clear–glaze” icing may decrease if aerosols are increased.
A natural extension of this study would be to run a series of similar sensitivity experiments for multiple months, seasons, and years to capture the breadth of precipitation systems across most of a continent and study the resulting changes in regional precipitation, especially in water-sensitive areas of the western United States. Also this study did not break down various mesoscale forcing mechanisms such as orographic forced snow (or rain), lake effect, sea-breeze areas, or strong convective regions to investigate aerosol effects on more localized precipitation, but the foundation for these tests was demonstrated. Additionally, we believe this scheme is well suited to simulate long-duration convective events including typical nonsevere shallow convection along with deep convective squall lines, supercells, and mesoscale convective systems (MCSs), since all inherently linked dynamics and feedbacks are present in this type of configuration using a well-established convection-permitting model (WRF). Such simulations could be used to validate many claims of aerosol invigoration of shallow and deep convection (cf. Li et al. 2011; Tao et al. 2012) and perhaps reveal if aerosols’ effects are causing specific responses in convection or are simply correlated with various convective weather situations (Morrison and Grabowski 2011, 2013).
The authors wish to thank David Gill and Jimy Dudhia for their advice and support of various WRF code modifications. Alison D. Nugent, Roy Rasmussen, and Hugh Morrison are gratefully acknowledged for their discussions that ultimately improved this work, and Kyoko Ikeda is thanked for helping to create Fig. 14. This research is in response to requirements and funding by the Federal Aviation Administration. The views expressed are those of the authors and do not necessarily represent the official policy or position of the FAA.
Ice Accretion Rates
Equation (A1), used for ice accretion, followed Makkonen (2000) where the change in mass with time is a product of a collision efficiency α1 [computed using Eqs. (A2)–(A8) following Finstad et al. (1988)], LWC, velocity υ, and cross-sectional area of the cylinder A. For simplicity, we used standard values of a 76.2-mm (3 in.)- diameter cylinder assumed to be moving at 89.1 m s−1 (200 mph), consistent with values used for decades by the aircraft-icing research community (Jeck 2001):
and K is the Stokes number, Re is Reynolds number, ϕ is Langmuir’s parameter, μ is dynamic viscosity of air, ρw is the density of water, and ρa is air density.