IMPROVEMENT OF NUMERICAL SOLUTIONS FOR AEROSOL CONDENSATION EQUATIONS
Final Report
Final Report
CRC Project No. A-39
CRC Project No. E551
Prepared by:
Timothy M. Gaydos, Bonyoung Koo, and Spyros Pandis
Department of Chemical Engineering
CarnegieMellonUniversity
July 2002
Prepared for:
Coordinating Research Council, Inc.
3650 Mansell Road, Suite 140 - Alpharetta, Georgia 30022
IMPROVEMENT OF NUMERICAL SOLUTIONS FOR AEROSOL CONDENSATION EQUATIONS
Final Report
Final Report
CRC Project No. A-39
CRC Project No. E551
Prepared by:
Timothy M. Gaydos, Bonyoung Koo, and Spyros Pandis
Department of Chemical Engineering
CarnegieMellonUniversity
5000 Forbes Avenue
Pittsburgh, PA 15213
July 2002
Abstract
Condensation onto existing aerosols can account for over two-thirds of total fine particulate mass. Since aerosols affect human health, visibility, and climate, it is important to be able to model their condensation/evaporation accurately in order to predict how particle mass will change over time. Early models assumed equilibrium between the gas and aerosol phases. Research has shown, however, that the assumption of equilibrium does not hold under many atmospheric conditions, particularly for low temperatures, high relative humidity, and small number concentrations of particles. Consequently, kinetic mass transfer must be accounted for when solving condensation/evaporation. Chock and Winkler (2000) proposed adapting the trajectory-grid method, used for solving the transport equation, to solve the condensation/ evaporation equation. Their preliminary results showed it to be fast and accurate in simple systems (one component, etc.). The approach of Chock and Winkler has been modified for implementation in both one-dimensional trajectory models and the three-dimensional chemical transport model CAMx. The October 1995 PM episode in the SouthCoastAirBasin in CA is analyzed using the three-dimensional host model. The trajectory-grid method is shown to represent a significant advantage in speed over previous methods, while maintaining the same level of accuracy.
Introduction
Atmospheric aerosols negatively affect visibility, human health, and also have an effect on global climate change. Fine and coarse particles can be emitted directly into the atmosphere. A significant portion of aerosol mass, however, is created by gases partitioning into the aerosol phase. Specifically, greater than two-thirds of PM2.5 mass has anthropogenic sources, with much of this mass being created by the condensation of gases onto existing aerosol particles. Air quality standards are based on total PM10 and PM2.5 mass, so it is important to be able to predict how much new material will condense onto existing aerosols. Models must describe this process effectively so they can be used to develop strategies for reducing particulate matter, and the associated harmful effects, in the atmosphere.
In early three-dimensional air quality models, equilibrium was assumed between the gas and particle phase (Pilinis et al., 1987; Russell et al., 1988; Binkowski and Shankar, 1995; Lurmann et al., 1997), rather than accounting for mass transfer limitations. Although this assumption makes these models efficient and stable, several studies have shown that this equilibrium assumption can introduce significant error into aerosol models, particularly for low temperatures, high relative humidity, and larger particles (Wexler and Seinfeld, 1990; Meng and Seinfeld, 1996, Lurmann et al., 1997). To increase the level of accuracy, several models that treat aerosol dynamics and account for mass transfer between phases have been developed (Meng and Seinfeld, 1996; Meng et al., 1998; Jacobson et al., 1996; Jacobson, 1997a, b; Sun and Wexler, 1998a, b; Pilinis et al., 2000; Nguyen and Dabdub, 2002). A hybrid method has also been developed which solves the dynamic equations for the larger particles, which are slower to reach equilibrium, and establishes equilibrium for the smaller particles (Capaldo et al., 2000).
The problem with these methods is that they improve accuracy at the expense of a large increase in computational cost. Models that treat kinetic mass transfer explicitly typically have at least 100 times the CPU requirements of equilibrium models (Seigneur, 2001). These increased computational requirements are caused by large differences in the mass transfer rates of different species and particle sizes, creating a stiff system of differential equations, which need to be solved. More efficient methods of solving this system of equations need to be developed for application to state-of-the-art three-dimensional chemical transport models, which can be used for developing emission control strategies.
Chock and Winkler (2000) proposed adapting the trajectory-grid method, used for solving the transport equation, to solve the condensation/evaporation equation. Their preliminary results in simple systems (one component, etc.) found it to be fast and accurate. In this paper, this method has been extended for use in a three-dimensional chemical transport model. A brief overview of the trajectory-grid method appears in the next section. The method is then applied to two test cases. The first involves application of a one-dimensional trajectory model to an air pollution episode in CA in August 1987. The results with the trajectory-grid method using the hybrid approach of Capaldo et al. (2000) and the fully dynamic case are compared to results using currently used approaches. The second case applies a three-dimensional chemical transport model to a PM episode in the South Coast Air Basin of California. Again, results using the trajectory-grid and current methods are compared, this time only with the hybrid approach.
Trajectory-Grid method
The trajectory-grid method solves the condensation/evaporation equation for a multicomponent aerosol population using the method of characteristics along with equation splitting. The condensation /evaporation equation can be written as
(1)
where pi is the aerosol mass concentration distribution of the ith species such that pid is the mass concentration of i in the size range between and + d, is the logarithm of the aerosol diameter, H = iHi, and p = ipi (Chock and Winkler, 2000). Hi is the condensation rate of species i and is defined as:
(2)
where Di is the diffusion coefficient for species i in air, Mi is its molecular weight, f(Kn,) is the correction due to noncontinuum effects and imperfect accommodation (Dahneke, 1983), p is the density of the particle, Dp is the diameter of the particle, and (ci-ceq,i) is the difference between the bulk and equilibrium vapor pressure for species i. Using the trajectory-grid method, Eq. (1) becomes
(3)
along the characteristic curve ddt = H/3. Finally, it is split into two equations
(4)
(5)
to make it easier to find a solution. A solution for Eq. (4) can be derived (Dhaniyala and Wexler, 1996) assuming that Hi and H are constant for the time-step t:
for H0(6)
for H =0(7)
The solution to Eq. (5) is
(8)
The value of the diameter is updated along the characteristic curve to (t+t) =(t)+(H/3)t.
The approach of Chock and Winkler (2000) was extended in the current work to handle situations that did not arise in the original publication. Previously, the model was only used with a constant condensation rate, i.e. H/ = 0. When H/ 0, it should be calculated explicitly. However, for small time steps and condensation gradients, H/ may be considered negligible and the solution to Eq. (3) is simply Eq. (4).
In order to use the solution to Eq. (4), a small enough step size must be chosen so that Hi and H do not change significantly. A simple algorithm to determine the step size t was developed. The algorithm being used calculates the change in H after each time step. If
(9)
then the nth step is repeated with t = t/2, where is an arbitrary error tolerance. Also, if
(10)
the time steps is increased by 10%, with being an arbitrary lower limit. This allows the time step to increase when the system is under fairly stable conditions, but prevents it from increasing too quickly.
Test Cases
The first test case involves a one-dimensional trajectory framework. Twenty-four trajectories were modeled, each reaching Claremont, CA on a different hour on August 28, 1987. The Secondary Organic Aerosol Model (SOAM v2.0) is the host module and is described in Koo et al. (2002). The MADM framework (Pilinis et al., 2000) is used for solving aerosol dynamics, and results using the trajectory-grid method to solve the condensation/evaporation equation are compared to results using the VODE solver (Brown et al., 1989). Eight size sections are used to represent the aerosol distribution, with 25 organic and inorganic species in the gas and aerosol phase.
The second test case involves a three-dimensional framework. The host model in this case is the chemical transport model CAMx (ENVIRON). The trajectory-grid model was incorporated into the aerosol module of this chemical transport model. The October 1995 PM episode in the South Coast Air Basin of California was studied using both the VODE and trajectory-grid methods and solving the system of equations using the hybrid method. The model inputs are setup for a 325x200 km region with 5x5 km grid resolution and ten vertical layers (65 x 40 grid squares). The aerosol distribution is split into 10 size sections, with 13 species in the aerosol phase and 8 in the gas phase.
Results
For the one-dimensional study, the results of the trajectory-grid and VODE runs were similar, with the trajectory-grid method providing a significant difference in CPU time, including over an order of magnitude decrease for the hybrid case. Fig. 1 shows the results compared to measurement data collected by a filter-based sampler (Fitz et al., 1989). PM10 and PM2.5 mass of various aerosol species are compared over five sampling periods on August 28 (0:00-5:00, 5:00-9:00, 9:00-13:00, 13:00-17:00 and 17:00-24:00 PST).
Little difference is seen between the predicted values using the trajectory-grid method and the VODE solver. Compared to the measured data, only one point has an error greater than 30% (PM2.5 ammonia from 17:00-24:00 PST). The trajectory-grid has a great advantage in efficiency over the VODE solver (Fig. 2). The VODE solver averages 120.6 CPU seconds per simulation hour per cell in the fully dynamic case and 53.9 CPU seconds using the hybrid method. The T-G method averages 48.9 CPU seconds for the dynamic case, beating the hybrid time of the VODE solver. For the hybrid case, the trajectory-grid method requires just 2.3 CPU seconds per simulation hour per cell, over an order of magnitude decrease compared to the other methods.
For the results of the October 1995 episode, the results were again similar between the VODE and trajectory grid methods. The daily average total PM10 and PM2.5 mass on October 17 for each method is shown in Figs. 3 and 4, along with the percent difference between the two methods. The percentage difference between the two methods is less than 1.0% for all locations. The trajectory grid method again has the advantage in speed, taking 0.09 CPU seconds per simulation per cell compared to 0.57 CPU seconds per simulation hour per cell using the VODE solver on October 17. The VODE solver fails during the second day, while the trajectory grid method has times of 0.07 CPUs per simulation hour per cell on October 18 and 19, respectively. For comparison, using the equilibrium assumption takes approximately 0.02 CPUs per simulation per cell on the same machine.
The results using the trajectory grid method are compared to measurements at Anaheim, Los Angeles, Diamond Bar, Fontana, and Riverside. The predicted vs. measured daily average PM10 and PM2.5 mass and nitrate concentrations are seen in Fig. 5. To evaluate the model, the mean normalized bias
(11)
and mean normalized error
(12)
are calculated, where the N represents the measurements at all five sites for October 17, 18, and 19. The results are shown in Table 1. The mean normalized bias is within +/-15% for all four variables, while the mean normalized error ranges from 15.7% to 26.4%. Although no guidelines are available for particulate matter concentrations, these values satisfy the statistical goals of ozone performance given by the ARB and EPA of a mean normalized bias within 15% and a mean normalized error less than 35%, indicating satisfactory performance of the model for this episode.
Conclusions
The trajectory-grid method is shown to be an effective method for solving the condensation/evaporation equation in chemical transport models. It is faster than existing methods of accounting for kinetic mass transfer, while maintaining the same level of accuracy. Combined with the hybrid method, which applies the equilibrium assumption to particles less than 0.625 m, the trajectory-grid method requires just three to five times the computational requirements of applying the equilibrium assumption to all particles in a three-dimensional chemical transport model.
The results obtained with the trajectory-grid method are virtually identical to results obtained with the same framework but using an ODE solver to solve the condensation/evaporation equation, both in one-dimensional and three-dimensional chemical transport models. The performance is satisfactory in both the one-dimensional and three-dimensional models, compared to measured concentrations.
References
Binkowski, F. S., and Shankar, U. (1995). The regional particulate matter model: 1. Model description and preliminary results. J. Geophys. Res. 100:26191-26209.
Brown, P. N., Byrne, G. D., and Hindmarsh, A. C. (1989). VODE: A Variable Coefficient ODE Solver. SIAM J. Sci. Stat. Comput. 10:1038-1051.
Capaldo, K. P., Pilinis, C., and Pandis, S. N. (2000). A computationally efficient hybrid approach for dynamic gas/aerosol transfer in air quality models. Atmos. Environ. 34:3617-3627.
Chock, D. P., and Winkler, S. L. (2000). A trajectory-grid approach for solving the condensation and evaporation equations of aerosols. Atmos. Environ. 34:2957-2973.
Dahneke, B. (1983). Simple Kinetic Theory of Brownian Diffusion in Vapors and Aerosols. In Theory of Dispersed Multiphase Flow, edited by R. E. Meyer. Academic Press, New York, pp. 97-133.
Dhaniyala, S., Wexler, A. S. (1996). Numerical schemes to model condensation and evaporation of aerosols. Atmos. Environ. 30:919-928.
Fitz, D. R., Chan, M., Cass, G. R., Lawson, D. R., and Ashbaugh, L. (1989). A multi-component size-classifying aerosol and gas sampler for ambient air monitoring. Paper No. 89-140.1, presented at the Air and Waste Management Association 82nd Annual meeting.
Jacobson, M. Z. (1997a). Development and application of a new air pollution modeling system- II: Aerosol module structure and design. Atmos. Environ. 31A:131-144.
Jacobson, M. Z. (1997b). Numerical techniques to solve condensational and dissolutional growth equations when growth is coupled to reversible reactions. Aerosol Sci. Technol. 27:491-498.
Jacobson, M. Z., Tabazadeh, A., and Turco, R. P. (1996). Simulating equilibrium within aerosols and nonequilibrium between gases and aerosols. J. Geophys. Res. 101:9079-9091.
Koo, B., Gaydos, T. M., and Pandis, S. N. (2002). Evaluation of the equilibrium, hybrid, and dynamic aerosol modeling approaches. In press.
Lurmann, F. W., Wexler, A. S., Pandis, S. N., Musarra, S., Kumar, N., and Seinfeld, J. H. (1997). Modelling urban and regional aerosols-II. Application to California’s south coast air basin. Atmos. Environ. 31:2695-2715.
Meng, Z., and Seinfeld, J H. (1996). Time scales to achieve atmospheric gas-aerosol equilibrium for volatile species. Atmos Environ. 30:2889-2900.
Meng, Z., Dabdub, D., and Seinfeld, J.H. (1998). Size-resolved and chemically resolved model of atmospheric aerosol dynamics. J. Geophys. Res. 103:3419-3435.
Nguyen, K. and Dabdub, D. (2002). Semi-Lagrangian Flux Scheme for the Solution of the Aerosol Condensation/Evaporation Equation. Aerosol Sci. Technol. 36:407-418.
Pilinis, C., Seinfeld, J. H., and Seigneur, C. (1987). Mathematical modeling of the dynamics of multicomponent atmospheric aerosols. Atmos. Environ. 21:943-955.
Pilinis, C., Capaldo, K. P., Nenes, A., and Pandis, S. N. (2000). MADM-A new multicomponent aerosol dynamics model. Aerosol Sci. Technol. 32:482-502.
Russell, A. G., McCue, K. F., and Cass, G. R. (1988). Mathematical modeling of the formation of nitrogen-containing air pollutants, 1. Evaluation of an Eulerian photochemical model. Environ. Sci. Technol. 22:263-271.
Seigneur, C. (2001). Current Status of Air Quality Models for Particulate Matter. J. Air & Waste Manage. Assoc. 51:1508-1521.
Sun, Q., and Wexler, A. S. (1998a). Modeling urban and regional aerosols near acid neutrality-application to the 24-25 June SCAQS episode. Atmos. Environ. 32:3533-3545.
Sun, Q., and Wexler, A. S. (1998b). Modeling urban and regional aerosols - condensation and evaporation near acid neutrality. Atmos. Environ. 32:3527-3531.
Wexler, A. S., and Seinfeld, J. H. (1990). The distribution of ammonium salts among a size and composition dispersed aerosol. Atmos. Environ. 24A:1231-1246.
Fig. 1. Comparison of predicted and measured concentrations at Claremont, CA on August 28, 1987 for the a) fully dynamic and b) hybrid methods.
Fig. 2. CPU time comparison for one-dimensional case study.
Fig. 3. Daily average PM10 mass (in g/m3) on October 17 for a) VODE and b) trajectory-grid method, and c) the percent difference between the two.
Fig. 4. Daily average PM2.5 mass (in g/m3) on October 17 for a) VODE and b) trajectory-grid method, and c) the percent difference between the two.
Fig. 5. Daily average predictions vs. measured values at five sites for a) PM10 mass b) PM2.5 mass c) PM10 nitrate and d) PM2.5 nitrate.
Table 1.
Statistical comparison of hybrid trajectory-grid predictions with measurements at five locations.
Method / NBIAS / NERRORPM10 mass / 13.34 / 23.36
PM10 NO3 / -11.33 / 24.04
PM2.5 mass / -4.38 / 15.73
PM2.5 NO3 / 13.97 / 26.42
1