Kinetic multi-layer model of gas-particle interactions in aerosols and clouds (KM-GAP)

. We present a novel kinetic multi-layer model for gas-particle interactions in aerosols and clouds (KM-GAP) that treats explicitly all steps of mass transport and chemical reaction of semi-volatile species partitioning between gas phase, particle surface and particle bulk. KM-GAP is based on the PRA model framework (P ¨ oschl-Rudich-Ammann, 2007), and it includes gas phase diffusion, reversible adsorption, surface reactions, bulk diffusion and reaction, as well as condensation, evaporation and heat transfer. The size change of atmospheric particles and the temporal evolution and spatial proﬁle of the concentration of individual chemical species can be modeled along with gas uptake and accommodation coefﬁcients. Depending on the complexity of the investigated system and the computational con-straints,

Atmospheric particles consist of a wide variety of organic and inorganic chemical compounds which can occur as different phase states of liquids, semi-solids or solids (crystalline, amorphous, glassy, ultraviscous, gel-like) depending on their composition and ambient conditions (Marcolli et al., 2004;Zobrist et al., 2008;Mikhailov et al., 2009;Virtanen et al., 2010;Koop et al., 2011). The phase state can influence the gas-particle partitioning, heterogeneous reactions and multiphase processes in atmospheric aerosols (Knopf et al., 2005;Nash et al., 2006;Renbaum and Smith, 2009;Cappa and Wilson, 2011;Pfrang et al., 2011;Shiraiwa et al., 2011a;Vaden et al., 2011). Chemical reactions can proceed in both gas and condensed phase but it is often difficult to discriminate gas, surface and bulk reactions. Moreover, the relative importance of them for secondary organic aerosol formation and aging is poorly understood (Moise and Rudich, 2000;Kalberer et al., 2004;McNeill et al., 2008;Hallquist et al., 2009;George and Abbatt, 2010a;Fry et al., 2011;Salo et al., 2011;Loza et al., 2012). So far resistor model formulations are widely used to describe and investigate heterogeneous reactions and multiphase processes in laboratory, field and model studies of atmospheric chemistry (Hanson, 1997;Finlayson-Pitts and Pitts, 2000;Worsnop et al., 2002;Anttila et al., 2006;King et al., 2009;Xiao and Bertram, 2011). The traditional resistor models, however, are usually based on simplifying assumptions such as steady state conditions, homogeneous mixing, and limited numbers of non-interacting species and processes.
In order to overcome these limitations, Pöschl, Rudich and Ammann  have developed a kinetic model framework (PRA framework) with a doublelayer surface concept and universally applicable rate equations and parameters for mass transport and chemical reactions at the gas-particle interface of aerosols and clouds. Ammann and  provided first examples on how the PRA framework can be applied to describe various physicochemical processes in aerosols and clouds such as reactive gas uptake on solid particles and solubility saturation of liquid droplets under transient or steady-state conditions. Shiraiwa et al. (2009) presented a kinetic double-layer surface model (K2-SURF) for the degradation of a wide range of polycyclic aromatic hydrocarbons (PAHs) by multiple photooxidants (O 3 , NO 2 , OH and NO 3 ). Using K2-SURF Shiraiwa et al. (2011b) showed that long-lived reactive oxygen intermediates are formed in the reaction of ozone with aerosol particles. Springmann et al. (2009) and Kaiser et al. (2011) demonstrated the applicability and usefulness of the PRA framework in an urban plume box model of the degradation of benzo[a]pyrene on soot by ozone and nitrogen dioxide. Pfrang et al. (2010) developed a kinetic double-layer model coupling aerosol surface and bulk chemistry (K2-SUB), in which mass transport and chemical reactions in the bulk are not explicitly resolved but represented by a reacto-diffusive flux analogous to traditional resistormodel formulations. Shiraiwa et al. (2010) developed a kinetic multi-layer model of aerosol surface and bulk chemistry (KM-SUB) that explicitly treats all steps of mass transport and chemical reaction from the gas-particle interface to the particle core, resolving concentration gradients and diffusion throughout the particle bulk. Pfrang et al. (2011) applied and extended KM-SUB to investigate the impact of transformation of diffusivity on the chemical aging of multi-component organic aerosol particles.
Here we present a kinetic multi-layer model of gas-particle interactions in aerosols and clouds (KM-GAP) that builds on KM-SUB and explicitly treats all steps of mass transport and chemical reaction from the gas phase to the particle core, including the evaporation and condensation of semi-volatile species. We demonstrate the applicability of KM-GAP for the condensation of water as well as the oxidation and evaporation of organics. The model results are compared with earlier experimental and theoretical studies.

Model description
As illustrated in Fig. 1, KM-GAP consists of multiple model compartments and layers, respectively: gas phase, nearsurface gas phase, sorption layer, quasi-static surface layer, near-surface bulk, and a number of n bulk layers. The sorption and quasi-static surface layer has a monolayer thickness that corresponds to the (average) effective molecular diameter of semi-volatile species Z i (δ Z i ). The following processes are considered in KM-GAP: gas phase diffusion, gassurface transport (reversible adsorption), surface-bulk transport, surface and bulk reactions, bulk diffusion, vapor pressure change, and heat transfer. The following differential equations can be used to describe the mass or number balance of each molecule Z i for surface and each bulk layers:

Fig. 1.
Kinetic multi-layer model of gas-particle interactions in aerosols and clouds (KM-GAP): (a) model compartments and layers with corresponding distances from the particle center (r), surface areas (A), and volumes (V ). λ Z i is the mean free path of semi-volatile species Z i in the gas phase; δ Z i is the thickness of sorption layer. (b) Transport fluxes (green arrows) and chemical reactions (red arrows). The black arrows indicate that each bulk layer can grow or shrink.
where N Z i ,s , N Z i ,ss , and N Z i ,bk are the absolute number of Z i molecules at surface, at quasi-static surface layer, and in bulk layer k, respectively. The various types of mass transport fluxes (J ) and rates of chemical production and loss (P , L) are defined in the list of symbols (Table A1). A(k) and V (k) are the outer surface area and the volume of the bulk layer k, respectively (Fig. 1). The volume of each layer V (k) can be calculated using N Z i ,bk and molecular volume V Z i assuming ideal mixing (volume additivity).
The ideal mixing assumption is also employed in practically all chemical transport models, and recent experiments are consistent with pseudo-ideal mixing of anthropogenic and biogenic SOA compounds at equilibrium (Hildebrandt et al., 2011). Note that the assumption of ideality can be relaxed in KM-GAP if activity coefficients are available.
For spherical particles, the radius position r(k), the outer surface area A(k) and the layer thickness δ(k) of the bulk layer k and particle diameter d p can be calculated as follows: In this way, each layer can either shrink or grow, in response to mass transport and chemical reactions. KM-GAP can also deal with planar geometry (thin films) assuming A(k) is constant. The surface and bulk number concentrations of Z i can be calculated as follows: The uptake coefficient of gas species Z i can be expressed as a ratio between the net fluxes of Z i from the gas phase to the condensed phase J net,Z i , and the collision flux J coll,Z i .
Note that γ Z i is 0 ≤ γ Z i ≤ 1 when Z i is adsorbing/condensing onto a surface (particle growth), whereas γ Z i can be negative when Z i is desorbing/evaporating from a surface (particle shrinkage). The change of gas phase concentration of Z i can be described using J net,Z i as follows: where N p is the number concentration of aerosol particles and A ss,m is the surface area of the particle m. In case of monodispersed particles the above equation is reduced to: In the appendix we specify the formalisms used to describe and calculate fluxes or rates of gas phase diffusion, gassurface interactions, surface reaction, surface-bulk transport, diffusion and reaction in the particle bulk, and heat transfer (Sects. A1-A5). Note that KM-GAP can treat species of any volatility, without the need to constrain any species to just one phase (aerosol or gas), i.e. it effectively treats each species Z i as semi-volatile by using its actual or estimated vapor pressures (Sect. A4). In the kinetic models previously developed by the authors, the volatile species X and non-volatile species Y were treated separately (Shiraiwa et al., 2009Pfrang et al., 2010), but the distinction between them is removed now in the KM-GAP model presented here.

Model application
To test and demonstrate the applicability of KM-GAP, we applied KM-GAP to three different processes: water condensation (Sect. 3.1), evaporation kinetics of organics (Sect. 3.2), and particle shrinkage due to chemical reactions (Sect. 3.3). In this study the coupled differential equations were solved in Matlab software with an ordinary differential equation solver (ode23tb), which integrates a system of differential equations using second and third order Runge-Kutta formulas. The computational costs for one simulation were less than one minute on a standard desktop computer. In each simulation the KM-GAP results are compared to experimental data available in the literature.

Condensation of water vapor
We simulated the condensation of water vapor in comparison to experimental data from Winkler et al. (2006). In their 2 Figure 2.  Winkler et al. (2006). The blue line is modeled with a low bulk diffusion coefficient (D b,w = 10 −11 cm 2 s −1 ) and α s,0,w = 1. (b) Temporal evolution of surface temperature (T s , black), ambient temperature (T , red) and supersaturation (S, blue). experiments, monodisperse Ag particles with a diameter of 9 nm have been used as condensation nuclei and humidified under initial supersaturation of 37.5 % at 268 K and 737 Torr. The particle number concentration was 4381 cm −3 . The growth of water droplet was observed as shown in Fig. 2a.
We modeled the temporal variation of particle radius and surface as well as ambient temperature by numerically solving the differential equations of mass balance for each model compartment (Eqs. 1-5, 16, A41, A42). The model simulations were performed with n = 50 bulk layers. The thermodynamic and kinetic parameters required for the simulations are summarized in Table 1. The thermal accommodation coefficient of water (α T,w ) is reported to be 1 . The bulk diffusion coefficient of water (D b,w ) of supercooled water is reported to be 2 × 10 −6 cm 2 s −1 at ∼268 K (Debenedetti, 1996). The desorption lifetime of water (τ d,w ) is estimated by the molecular dynamic simulation to be 3.5 × 10 −11 s (Vieceli et al., 2005;Garrett et al., 2006). A sensitivity study on τ d,w showed that the growth curve is insensitive to τ d,w as long as it is less than 10 −7 s. Values of 0.1-1 have been reported for the mass accommodation coefficient of water vapor on pure water (Davidovits et al., 2004;Kolb et al., 2010), and we used this range of α s,0,w to fit the experimental data. For simplicity the kinetic parameters were assumed to be constant throughout each model run.
Note that the molecular dynamic simulations suggest an air/water interface thickness of 0.3-0.6 nm (Garrett et al., 2006), which is equivalent to approximately two water molecules. At the interface each water molecule has on average two hydrogen bonds less than it would have in the bulk (Taylor et al., 1996). KM-GAP captures this interface with monolayers in sorption and quasi-static surface layers.
As shown in Fig. 2b, the simulated temperature of the droplet surface increased to 270 K because of latent heat release by water condensation at the surface. The released heat was then transferred to the ambient air, leading to the gradual increase of the ambient air temperature. The red lines in Fig. 2a illustrate the model results of KM-GAP with different values of α s,0,w . With α s,0,w = 1 (solid line), the simulation was in very good agreement with the experimental observation. In this case, the droplet growth is limited only by the gas phase diffusion of water vapor. Due to the continuous water condensation the ambient supersaturation decreased from 37.5 % to ∼33 %. When α s,0,w values between 0.1-0.5 were used, the simulations substantially underestimated the observed growth. In these cases the growth curves were kinetically limited by the surface accommodation of water vapor.
Our results suggest that the surface accommodation coefficient of water is close to unity, which is consistent with molecular dynamic simulations (Morita et al., 2004;Winkler et al., 2004;Garrett et al., 2006). Indeed,  and (2006) already showed that the accommodation coefficient is close to unity, but the actual value remained under discussion (Kolb et al., 2010). Winkler et al. (2004Winkler et al. ( , 2006) used a continuum model of gas phase diffusion and condensation that does not resolve microscopic information such as the desorption and surface residence times of water molecules. These and other earlier studies did not provide a direct, quantitative link between molecular dynamic simulations and laboratory observations through a single model resolving both the molecular processes at the surface and the macroscopic growth of particles. Thus, our model calculations help to confirm the findings of Winkler et al. (2004Winkler et al. ( , 2006, which is relevant for estimating the number concentration of cloud droplets (Laaksonen et al., 2005).
To investigate the effects of bulk diffusion and surfacebulk transport of water, particle growth was simulated with a smaller D b,w value of 10 −11 cm 2 s −1 , which is a characteristic value for ice and other solids (Huthwelker et al., 2006;Shiraiwa et al., 2011a). As illustrated as blue line in Fig. 2a, the particle growth was retarded substantially being kinetically limited by bulk diffusion and surface-bulk transport of water molecules. The bulk diffusion-limited particle growth may be relevant for condensational growth of amorphous (semi-)solid aerosol particles (Mikhailov et al., 2009;Koop et al., 2011;Tong et al., 2011;Zobrist et al., 2011), which we intend to investigate further in follow-up studies.

Evaporation of dioctyl phthalate
Here we simulated the size-dependent evaporation of singlecomponent organic particles composed of dioctyl phthalate (DOP) and compare the results to experimental data from Vaden et al. (2011). In their experiments, the monodisperse DOP droplets were loaded into a 7 L evaporation chamber containing 1 L of activated charcoal at the bottom of the chamber which continuously removed gas-phase organics. The number concentration of particles in the chamber was ∼150 cm −3 and this value decreased by ∼20 % during the course of the evaporation experiments due to wall losses and dilution during sampling (A. N. Zelenyuk, personal communication, 2011).
We modeled the temporal evolution of particle diameter and gas phase concentration of DOP as a function of evaporation time. The coupled differential equations (Eqs. 1-5, 16) were solved numerically (n = 5). Heat transfer was not considered for simplicity. The required parameters are: vapor pressure of DOP (10 −7 Torr) (Cappa et al., 2008); gas phase diffusion coefficient of DOP (D g,DOP = 4.4 × 10 −2 cm 2 s −1 ) (Ray et al., 1988); bulk diffusion coefficient of DOP (D b,DOP = 3.0 × 10 −8 cm 2 s −1 ) which is estimated from the viscosity of DOP (0.081 Pa s) through the Stokes-Einstein equation (Shiraiwa et al., 2011a), using an effective molecular radius of DOP (δ DOP = 0.86 nm) calculated by density (ρ DOP = 0.986 g cm −3 ) and molar mass (M DOP = 390.56 g mol −1 ). Desorption lifetime of DOP (τ d,DOP ) is assumed to be 10 −6 s which is longer than the one of H 2 O inferred from molecular dynamics simulations (10 −6 s vs. ∼10 −11 s) as DOP molecules are larger and heavier than H 2 O molecules. The surface accommodation coefficient of DOP on free substrate (α s,0,DOP ) is varied to fit to the experimental data. The loss of DOP to the denuder is considered using a first-order wall loss coefficient (k w ) and adding the term −k w [DOP] g to the right hand side of Eq. (16). We assume that the core of the chamber is well-mixed and the loss occurs via molecular diffusion (1st Fick's law) through a layer adjacent to the denuder wall to estimate k w . We varied the thickness of the denuder wall up to the half-height of the chamber (∼3 cm). Particle losses to the wall and coagulation were not considered. Table 1. Thermodynamic and kinetic parameters of air and water used for the simulation of water condensation and droplet growth.

Parameter Description
Value c p,a specific heat capacity of dry air at constant pressure (J g −1 K −1 ) 1.005 a c p,w specific heat capacity of water vapor at constant pressure (J g −1 K −1 ) 1.860 a D g,w gas diffusion coefficient of water vapor (Torr cm 2 s −1 ) thermal conductivity of binary mixtures of air and water vapor (cm 2 s −1 ) K w /(1 + 0.556p a /p w ) + K a /(1 + 1.189p w /p a ) a K a thermal conductivity of dry air (J s Winkler et al. (2006), b Ivanov et al. (2007) The experimentally measured size-and time-dependent evaporation of DOP is modeled well by KM-GAP as shown in Fig. 3, in which the molecular diffusion distance for the denuder loss is assumed to be the half-height of the chamber and α s,0,DOP is taken to be 0.5. The modeled gas phase concentration of DOP is ∼10 9 cm −3 and mainly determined by the vapor pressure of DOP. Due to turbulent mixing in the chamber, the actual molecular diffusion distance might be smaller (Crump and Seinfeld, 1981). For example, the data can also be reproduced with a smaller diffusion distance of 1 mm and α s,0,DOP of ∼0.3, while variations of τ d,DOP in the range of 10 −12 -1 s showed practically no effect. These findings emphasize the importance of mixing effects and chamber geometry. This information and related model sensitivity studies can be used to design and optimize further experiments for efficient probing of specific physical effects (e.g. experiments with variable chamber geometry to probe accommodation vs. mixing).
For the above modeling the particle number concentration was kept constant at the initial concentration for simplicity. The sensitivity studies on particle number concentration showed that a particle loss of 20 % during the course of the experiments does not affect the evaporation behavior significantly. If the number concentration is decreased by a factor of 10, however, the particles shrink faster by ∼5 % because of the lower gas phase concentration of DOP, suggesting that the particle number concentration also plays a role in such evaporation experiments. While traditional evaporation models may suffice to describe simple systems consisting of single components or liquid droplets, KM-GAP may provide further insight into the interplay of mass transport, phase transitions and chemical reactions in complex multicomponent systems of variable composition and phase state such as (secondary) organic aerosols that undergo chemical aging and transformation (Cappa and Wilson, 2011;Pfrang et al., 2011). We intend to investigate the evaporation kinetics of (semi-)solid particles composed of multiple components in follow-up studies. Table 2. Physical and chemical parameters for oleic acid ozonolysis including the assumed mole-based product yields (compare Zahardis and Petrucci, 2007; all minor products are represented by a dimer with the properties given below).

Oxidation and volatilization of oleic acid
Here we simulated the degradation of oleic acid by ozone in comparison to experimental data from Ziemann (2005). The same data set has been used by Pfrang et al. (2010) and Shiraiwa et al. (2010) for simulations with a kinetic doublelayer model (K2-SUB) and a kinetic-multi layer model (KM-SUB) for aerosol surface and bulk chemistry. In those studies the gas-particle partitioning of products was not considered, which is now included in the KM-GAP simulation.
The products of oleic acid oxidation by ozone are mainly 1-nonanal, 9-oxononanoic acid, nonanoic acid, azelaic acid, and peroxidic dimer (e.g. Zahardis and Petrucci, 2007;Vesna et al., 2009). Known and estimated physical properties of oleic acid and its products are summarized in Table 2. The assumed product yield of oleic acid ozonolysis is also given (compare Zahardis and Petrucci, 2007). For simplicity we assume only nonanal is semi-volatile and oleic acid and other products are non-volatile as their vapor pressures are low and all products other than 1-nonanal, 9-oxononanoic acid, nonanoic acid and azelaic acid can be represented by a dimer. We intend to remove these simplifying assumptions in follow-up studies. The gas phase ozone concentration was set to [Z 1 ] g = [Z 1 ] gs = 7.0 × 10 13 cm −3 (corresponding to 2.8 ppm at 1013 hPa and 298 K). The initial surface and bulk concentrations of ozone (Z 1 ) and nonanal were set to [Z i ] s,0 = [Z i ] bk,0 = 0. The initial surface and bulk concentrations of oleic acid (Z 2 ) were set to [Z 2 ] ss,0 = 9.7 × 10 13 cm −2 and [Z 2 ] b,0 = 1.2 × 10 21 cm −3 , respectively , corresponding to 3.15 mol l −1 as reported by Ziemann (2005). Accordingly, the initial value of the total number of oleic acid molecules in a particle with a radius of 0.2 µm was N Y,0 = 4.1 × 10 7 . Note that the bulk concentration of pure oleic acid is actually 1.95 × 10 21 cm −3 based on molecular weight and density, but the lower value of 1.2 × 10 21 cm −3 is used as an initial value because oleic acid was already exposed to ozone until ozone is well-mixed in the chamber (Ziemann, 2005). Therefore, oleic acid is already partly degraded and initial concentration of nonvolatile products are estimated based on the product yield.
We modeled the temporal evolution of particle size, the particle surface and bulk composition, and of the ozone uptake coefficient by numerically solving the differential equations of mass balance in terms of molecular number for each model compartment (Eqs. 1-5). Evolution of gas phase concentration and heat transfer were not resolved for simplicity. The kinetic parameters required for the model simulations are based on parameters of base case 1 given in Shiraiwa et al. (2010) and summarized in Table 2: the surface accommodation coefficient of ozone (α s,0,Z 1 ), the desorption lifetime of ozone (τ d,Z 1 ), the secondorder reaction rate coefficient at surface, quasi-static surface, and bulk (k SLR,Z 1 ,Z 2 = k QSLR,Z 1 ,Z 2 = 6 × 10 −12 cm −2 ; k BR,Z 1 ,Z 2 = 5×10 −17 cm −2 ). The bulk diffusion coefficients (D b ) of oleic acid, nonanal, and nonanoic acid were estimated from viscosity data (Noureddini et al., 1992) using the Stokes-Einstein equation (Shiraiwa et al., 2011a). Note that in our earlier studies D b for oleic acid was assumed to 10 −10 cm 2 s −1 Shiraiwa et al., 2010), but here we provide a better estimate of 1.9×10 −7 cm 2 s −1 . Additional input parameters were the mean thermal velocity of ozone (ω Z 1 = 3.6 × 10 4 cm s −1 ), the Henry's law coefficient of ozone (K sol,cc,Z 1 = 4.8 × 10 −4 mol cm −3 atm −1 ), the vapor pressure of nonanal (∼0.4 Torr) (Verevkin et al., 2003), and desorption lifetime of nonanal which is assumed to be 10 −4 s. For simplicity the kinetic parameters were assumed to be constant throughout each model run. Note that the actual reaction mechanism might be more complex, involving the formation of reactive oxygen intermediates (ROIs) with a different set of kinetic parameters including a shorter ozone desorption lifetime (Shiraiwa et al., 2011a, b).
The model simulations were performed with n = 100 layers to obtain high-resolution results of bulk concentration profiles. To test how the number of model layers in the particle bulk affects the simulation results, we have run the model within the range of n = 1-200. The model results of oleic acid degradation and particle growth were almost identical, demonstrating the robustness of the multi-layer model approach with transport rate coefficients scaled by layer thickness (Eqs. A18, A19, A27) .    Fig. 4a, KM-GAP fits very well to the simulated decay of oleic acid. The particle shrinks with the radius decreasing from 200 to 187 nm in 30 s (∼18 % of original mass is lost), due to the formation and evaporation of the volatile product nonanal. Sage et al. (2009) observed that only 1.5 % of the original oleic acid mass was lost by volatilization of the reaction product nonanal within a reaction time of several hours. The discrepancy might be due to additional secondary chemistry involving the carbon backbone or heterogeneous uptake of nonanal to the particle, which are not considered in this study. These processes can be implemented in KM-GAP when the relevant kinetic parameters of secondary reactions can be determined or estimated, which is a target for followup studies but goes beyond the scope of this study introducing the new model approach.
The simulated ozone uptake coefficient is nearly constant and identical to the surface and bulk accommodation coefficients (γ Z 1 ≈ α s,Z 1 ≈ α ss,Z 1 ≈ α b,Z 1 ≈ 4.2 × 10 −4 ), indicating that the gas uptake is limited by interfacial mass transport, i.e. by the process of bulk accommodation which is in turn limited by surface accommodation. Note that even if evaporation of nonanal is switched off (i.e. no particle shrinkage) the oleic acid degradation is well reproduced by KM-GAP with the same kinetic parameters, indicating the lim-ited impact of evaporation of volatile products on organics degradation in this particular experimental data set.
As shown in Fig. 4b, the surface concentration of ozone at the quasi-static surface layer increases gradually owing to the combination of reversible adsorption, surface reaction, and surface-to-bulk transport driven by the chemical reaction in the bulk. In contrast, the surface concentration of oleic acid decreases due to chemical reaction with ozone and while the concentration of nonvolatile products increases for the same reason. The surface concentration of nonanal stays very low (∼10 6 cm −2 ) at steady-state determined by desorption, chemical production and transport from the bulk. A very similar behavior is observed for the average bulk concentrations as illustrated in Fig. 4c. The bulk concentration of nonanal stays constant at a low level, whereas the concentrations of non-volatile products increase as oleic acid is degraded by bulk reaction. Figure 5 shows the bulk concentration profiles of (a) ozone, (b) oleic acid, (c) nonanal, and (d) other nonvolatile products (i.e. oxononanoic + nonanoic + azelaic acids + dimer). The y-axis displays the radial distance from the particle core. As shown in Fig. 5a, ozone diffuses rapidly into the bulk and a concentration gradient between nearsurface bulk and core is observed, which is determined by the interplay of interfacial transport with bulk diffusion and Atmos. Chem. Phys., 12, 2777-2794  reaction. During the first few seconds, the ozone concentration in the near-surface bulk is almost two orders of magnitude higher than in the particle center, because ozone is depleted due to fast bulk reaction in the near-surface bulk before it can diffuse more deeply into the bulk. Up to ∼30 s, the concentration gradient swiftly relaxes with decreasing abundance of oleic acid. As illustrated in Fig. 5b, the strong concentration gradient of ozone induces almost no gradient for oleic acid because the concentration of oleic acid is ∼8 orders of magnitude higher than that of ozone (Fig. 4c). Thus, oleic acid can effectively be regarded as well-mixed (i.e. concentration differences < 1 %). Figure 5c illustrates the bulk concentration profile of nonanal. It is interesting to note that even though the average bulk concentration of nonanal stays constant at a low concentration level (∼10 14 cm −3 ) as shown in Fig. 4c, a concentration gradient is observed. The concentration at the core is by a factor of ∼5 higher than that in the nearsurface bulk for the first 20 s, when the evaporation rate of nonanal is faster than bulk diffusion, which in turn is faster than the chemical production rate. After ∼20 s, when ozone reaches the particle core, the production rate of nonanal becomes faster than bulk diffusion leading to a steeper concentration gradient. In contrast, the non-volatile products show only a small concentration gradient, but at a high concentration level (∼10 21 cm −3 ).

Summary and conclusions
We present a novel kinetic multi-layer model of gas-particle interactions in aerosols and clouds (KM-GAP) that treats explicitly all steps of mass transport and chemical reaction of semi-volatile species partitioning between gas phase, particle surface and particle bulk. KM-GAP builds on the kinetic multi-layer model for aerosol surface and bulk chemistry which includes gas phase diffusion, reversible adsorption, surface reactions, and bulk diffusion and reaction (KM-SUB; Shiraiwa et al., 2010). The new processes included in KM-GAP are heat flux, evaporation and condensation of semivolatile species. The size change of the particle and the temporal evolution and concentration profiles of semi-volatile species at the gas-particle interface as well as in the particle bulk can be modeled along with surface concentrations and gas uptake coefficients. Depending on the complexity of the investigated system, unlimited numbers of semi-volatile species, chemical reactions, and physical processes can be treated, and the model shall help to bridge gaps in the understanding and quantification of multiphase chemistry and microphysics in atmospheric aerosols and clouds.
To demonstrate the applicability of KM-GAP, the condensation of water vapor on nanoparticles was simulated and compared to experimental data. Resolving the reversible adsorption of water vapor with a water desorption lifetime of picoseconds as suggested by molecular dynamic simulations, we help to confirm the best estimate of the surface accommodation coefficient of water vapor on a water droplet to be unity as reported by Winkler et al. (2006). Sensitivity studies suggest that an artificially slow bulk diffusion can inhibit the particle growth significantly, which may be relevant for condensational growth of (semi-)solid organic aerosol particles.
For further demonstration of the ability of KM-GAP, the size-and time-dependent evaporation of single organic component particles was modeled. Experimental data on the evaporation of dioctyl phthalate particles can be reproduced with an accommodation coefficient of ∼0.5, and the model results emphasize the importance of aerosol particle number concentration and mixing effects in chamber experiments. Moreover, we showed how the formation and evaporation of volatile reaction products like nonanal can cause a decrease in the size of oleic acid particles exposed to ozone.
In conclusion, KM-GAP is a versatile tool for the modeling of multiphase chemical and microphysical processes in atmospheric aerosols and clouds in a self-consistent and realistic manner. The model can resolve the effects of surface and bulk mass transport in particles undergoing condensation and/or evaporation, which may be crucial for the description and understanding of the atmospheric aging of multicomponent mixtures of time-dependent composition and diffusivity. So far the formation and aging of secondary organic aerosol (SOA) has generally been described by thermodynamic models that implicitly assume quasi-instantaneous gas-particle partitioning (Pankow, 1994;Donahue et al., 2006Donahue et al., , 2011Zuend et al., 2010). If, however, the phase state of SOA is not liquid but (semi-)solid, this assumption may break down and the SOA partitioning may be kinetically limited (Virtanen et al., 2010;Cappa and Wilson, 2011;Pfrang et al., 2011;Vaden et al., 2011). We suggest and intend to use KM-GAP for the investigation of such effects.

A1 Gas phase diffusion, interfacial transport, and surface reactions
Based on kinetic theory, the gas kinetic flux of Z i colliding with the surface J coll,Z i can be expressed as where [Z i ] gs is the near-surface gas phase concentration of Z i and ω Z i is the mean thermal velocity given by ω Z i = (8 RT/(π M Z i )) 1/2 , where M Z i is the molar mass of Z i , R is the gas constant, and T is the absolute temperature. The significant net uptake of Z i will lead to local depletion at near-surface gas phase ([Z i ] gs < [Z i ] g ) and gas phase diffusion will influence further uptake. In this case near-surface gas phase concentrations should be corrected using a gas phase diffusion correction factor C g,Z i .
C g,Z i can be described as follows .
γ Z i is the uptake coefficients of Z i , which can be calculated as described in Eq. (14). Kn Z i is the Knudsen number which can be approximated by the gas phase diffusion coefficient D g,Z i and particle diameter d p .
The flux of adsorption of gas molecules on the quasi-static particle surface can be expressed as where α s,Z i is the surface accommodation coefficient and k a,Z i (=α s,Z i ω Z i /4) is a first-order adsorption rate coefficient. Estimation of α s,Z i is based on a Langmuir adsorption model in which all adsorbate species compete for a single type of non-interfering sorption sites on the quasi-static surface . α s,Z i is determined by the surface accommodation coefficient on an adsorbate-free surface α s,0,Z i and the sorption layer coverage θ s , which is given by the sum of the fractional surface coverage of all competing adsorbate species θ s,Z p .
θ s,Z p is the ratio between the actual and the maximum surface concentration value of Z p : θ s,Z p = [Z p ] s /[Z p ] s,max = σ s,Z p [Z p ] s , where σ s,Z p is the effective molecular cross section of Z p . The inverse molecular cross section can be regarded as the overall concentration of non-interfering sorption sites on the quasi-static surface layer : [SS] ss = σ −1 s,Z p . The adsorbed molecules can thermally desorb back to the gas phase. Desorption (the inverse of adsorption) can be described by a first-order rate coefficient k d,Z i , which is assumed to be independent of θ s,Z i . The flux of desorption of gas molecules on the quasi-static particle surface can be expressed using the desorption coefficient k d,Z i as The desorption lifetime τ d,Z i is the mean residence time on the surface in the absence of surface reaction and surfacebulk transport.
The transport of semi-volatile species Z i between sorption layer and quasi-static surface layer (J s,ss,Z i and J ss,s,Z i ) is described as follows: whereas k s,ss,Z i and k ss,s,Z i are the first-order transport rate coefficients. Estimates for k ss,s,Z i can be derived from the corresponding bulk diffusion coefficients D b,Z i based on Fick's first law of diffusion considering that a molecule Z i in the sorption layer needs to travel a distance of δ Z i to move into the quasi-static surface layer: An estimate for k s,ss,Z i can be determined considering mass transport equilibrium conditions. Mass balance implies that k s,ss,Z i [Z i ] s,eq = k ss,s,Z i [Z i ] ss,eq (J s,ss,Z i = J ss,s,Z i ) and k d, [Z i ] s,eq , and [Z i ] ss,eq are the equilibrium or solubility saturation number concentrations of Z i in the gas phase, on the sorption layer, and in the quasi-static surface layer, respectively. It leads to Note that k s,ss,Z i reflects the vapor pressure of Z i (p Z i ) as [Z i ] g,eq is a function of p Z i as described in Sect. 2.4. General rate equations of chemical production and loss by surface layer reactions, which can proceed within the sorption layer (SLR; P s,Z i − L s,Z i ), within the quasi-static surface layer (QSLR; P ss,Z j − L ss,Z j ), and between reactants of the sorption layer and of the quasi-static layer (SQSLR; P s,ss,Z i − L s,ss,Z i ), are given by: Here c SLRv,s,Z i , c QSLRv,ss,Z i and c SQSLRv,ss,Z i stand for the stoichiometric coefficients (negative for starting materials and positive for reaction products) of species Z i in surface reactions (v = 1, ..., v max ) in a system with a total number of v max (photo-)chemical reactions occurring on the surface of the investigated aerosol particles. k SLRv,Z p and k QSLRv,Z p are first-order reaction rate coefficients at sorption and quasistatic surface layer, respectively. k SLRv,Z p ,Z q , k QSLRv,Z p ,Z q and k SQSLRv,Z p ,Z q are second-order reaction rate coefficients. The quasi-static surface accommodation coefficient (α ss,Z i ), i.e. the probability for an individual gas molecule colliding with the surface to enter the quasi-static surface layer is given by: α ss,Z i = α s,Z i J s,ss,Z i J des,Z i + J s,ss,Z i + L s,Z i + L s,ss,Z i (A15)

A2 Surface-bulk transport
The surface-bulk transport of semi-volatile species Z i is defined as exchange between the quasi-static surface layer and near-surface bulk. Based on the PRA framework and following previous studies of kinetic models Pöschl et al., 2007;Pfrang et al., 2010;, the transport flux from quasi-static surface layer to bulk (J ss,b1 ,Z i ) can be described as: k ss,b1,Z i is the first-order transport rate coefficient (s −1 ) of Z i . In the same way, bulk to surface transport (J b,ss,Z i ) can be described as follows: is the first-order transport rate coefficient, which can be regarded as an effective transport velocity. Neglecting surfactant effects, k b1,ss,Z i can be estimated based on Fick's first law of diffusion considering that a molecule Z i in the near-surface bulk layer on average needs to travel a distance of (δ Z i + δ(1))/2 to move into the quasi-static surface layer: Under equilibrium conditions, mass conservation implies k b1,ss,Z i [Z i ] b1 = k ss,b1,Z i [Z i ] ss (J b1,ss,Z i = J ss,b1,Z i ) and for pure Z i the surface and bulk concentrations are given by the inverse of the effective molecular cross section and of the effective molecular volume, respectively: The bulk accommodation coefficient of Z i (α b,Z i ), i.e. the probability for a gas molecule colliding with surface to enter the bulk of the particle, can be derived by considering that a molecule Z i in the quasi-static surface layer needs to travel only over a distance of about one molecular diameter δ Z i to move into the bulk. The corresponding rate coefficient and transport flux can be estimated by Note that the physical parameters k ss,b,Z i , J ss,b1,Z i , and α b,Z i do not depend on the bulk layer thickness, whereas the model parameters k ss,b1,Z i and J ss,b1,Z i vary with the bulk layer thickness.
Z i can enter the bulk after iterative exchange between the sorption layer and the quasi-static surface layer, and α b,Z i can be calculated as follows: where Ψ ss,b and Ψ ss,s are the probabilities for Z i in the quasistatic surface layer to enter the bulk (Ψ ss,b ) or transport back to the sorption layer (Ψ ss,s ); Ψ s,ss is the probability for Z i in the sorption layer to enter the quasi-static surface layer:

A3 Bulk diffusion and reaction
Bulk diffusion is explicitly treated in KM-GAP as the mass transport (J bk,bk±1 ) from one bulk layer (bulk k) to the next (bulk k ± 1). In analogy to surface-bulk mass transport, we describe the mass transport fluxes between different layers of the bulk by first-order rate equations: Estimates for the transport rate coefficients or effective velocities of Z i from layer k to k + 1, k b,b,Z i (k) (cm s −1 ), can be calculated from the corresponding diffusion coefficients based on Fick's first law of diffusion. For this purpose we assume that each layer is homogeneously mixed (no concentration gradient within a layer), but note that this assumption can effectively be overcome by choosing an adequately large number of bulk layers. The lower limit of physically meaningful layer thickness corresponds to the molecular diameter of the involved chemical species, and the upper limit for correct treatment of bulk diffusion depends on the diffusivity. In (semi-)solid matrices with low diffusivity, the mixing times are long and concentration gradients can be steep, so that relatively thin layers should be used (>10 layers) (Shiraiwa et al., 2011a). In liquid matrices with high diffusivity, the mixing times are short and concentration gradients are small, so that relatively thick layers can provide sufficient resolution (<10 layers) (Shiraiwa et al., 2011a). The average travel distance for molecules moving from layer k to k + 1 is (δ(k) + δ(k + 1))/2, and from Fick's first law follows: Note that KM-SUB had used an alternative approach to estimate a transport velocity, but the approach presented here is more straightforward and also used successfully in another model study for similar purposes (Zobrist et al., 2011). This treatment of bulk diffusion yields practically the same results (concentration profiles) as the solving of partial differential equations (Smith et al., 2003;Shiraiwa et al., 2010), but it is more flexible and requires no assumptions about interfacial transport. The influence of changing chemical composition of the particle bulk on diffusion can be taken into account by describing D b,Z i using the obstruction or percolation theory (Pfrang et al., 2011;Shiraiwa et al., 2011a), or from parameterizations of experimental data (Zobrist et al., 2011). Chemical reactions proceeding within the bulk of a particle are defined as bulk reactions (BR). For simplicity, we assume that all relevant bulk reactions proceed via quasielementary steps with straightforward first-or second-order rate dependences on the concentrations within each bulk layer. The following generalized expressions can be used to describe net chemical production (i.e. production minus loss) of bulk species Z i within the bulk layer k.
Here c BRv,Z i stands for the stoichiometric coefficients (negative for starting materials and positive for reaction products) of species Z i in reaction BRv; v = 1, . . . , v max in a system with a total number of v max chemical reactions occurring in the bulk layer k. k BRv,Z p is a first-order reaction rate coefficient and k BRv,Z p ,Z q is a second-order bulk reaction rate coefficient between Z p and Z q in the condensed phase bulk of a system with multiple semi-volatile species which can react with each other. In principle, higher-order reactions might also occur in real systems and could be flexibly included in the model. Moreover, the concentration variables in the rate equations could be replaced by activities or corrected by activity coefficients to account for non-ideal concentration dependencies.

A4 Vapor pressure
The equilibrium (saturation) number concentrations of Z i in the gas phase [Z i ] g,eq can be calculated using the saturation vapor pressure p Z i (d p ) of semi-volatile species Z i , which is a function of the particle diameter as defined by the ideal gas equation:  (Köhler, 1936;Seinfeld and Pandis, 1998): where σ is the solution surface tension, n s is moles of solute in the particle, and p o (T s ) is the vapor pressure of pure Z i over a flat surface. p Z i (d p ) can be also calculated using a single hygroscopicity parameter κ (Petters and Kreidenweis, 2007): where d d is the dry diameter of the particle. Note that is also a function of the surface temperature T s .

A5 Heat transfer
For water condensation the heat flux Q directed away from a single droplet can be expressed as (Vesala et al., 1997;Winkler et al., 2006): where K(T s ) and K(T ) correspond to thermal conductivities of the binary mixture of inert gas and water vapor at the droplet surface temperature (T s ) and the ambient gas temperature (T ), respectively. I is the mass flux directed away from a single droplet: β T is the transitional correction factor for heat transfer given by (Vesala et al., 1997;Winkler et al., 2006) β T = 1 + Kn T 1 + (4/3α T + 0.377)Kn T + 4/3α T Kn 2 T (A34) Kn T is the Knudsen number for heat transfer characterized by particle diameter (d p ) and effective mean free path of the dry air molecules (λ a ): where D T is the thermal diffusivity of air and ω a is the mean thermal velocity of dry air. H v is the specific enthalpy of moist air: where H a and H w are the specific enthalpies of dry air and water vapor, respectively. These can be calculated from the specific heat capacities of dry air (c p,a ) and water vapor (c p,w ), respectively: L w is the latent heat of water vaporization. The humidity ratio x given in units of g g −1 is the ratio between the actual mass of water vapor present in moist air to the mass of the dry air. The droplet surface temperature T s can be obtained from the heat balance equation (Vesala et al., 1997).
H l (T s ) is the specific enthalpy of the liquid at the droplet surface temperature.
Inserting Eqs. (A32) and (A40) into Eq. (A39) yields Parameter definitions and values are given in Table 1. The ambient temperature (T ) can be calculated by the heat balance equation below.
The amount of heat released upon chemical reaction (reaction enthalpy) can be of similar magnitude or larger than the heat release upon condensation. Usually, however, the heat release upon reaction of trace gases with aerosol particles proceeds over much longer time scales than the condensation of cloud droplets, because the mass flux is much smaller. Therefore, the heat released by chemical reactions can be efficiently buffered by the ambient gas and does normally not lead to a substantial increase of particle surface temperature. For example, an oleic acid particle with a radius of 200 nm can be oxidized by ozone on the timescale of a minute. If we insert the corresponding average rate of heat release for L w I in Eq. (A41), the temperature change due to reaction is estimated to be below ∼10 −4 K as long as the reaction enthalpy is less negative than about −2000 kJ mol −1 . Indeed, the reaction enthalpy of oleic acid ozonolysis is estimated to be on the order of −300 to −700 kJ mol −1 based on published enthalpies of formation for the reaction products (Bond, 2007;NIST, 2011), assuming one mole of nonanal, 0.5 mole of nonanoic acid, and 0.5 mole of azelaic acid are formed from one mole of oleic acid. Even if we consider the complete combustion of oleic acid with a reaction enthalpy of −11161 kJ mol −1 (NIST, 2011), this would only lead to a temperature change of less than 10 −3 K. Thus, we ignore heat transfer effects in the modeling of chemical reactions in Sect. 3.3. Table A1. List of symbols.

Symbol
Description Unit α b,Z i bulk accommodation coefficient of Z i α s,Z i surface accommodation coefficient of Z i α ss,Z i quasi-static surface accommodation coefficient of Z i α s,0,Z i surface accommodation coefficient of Z i on an adsorbate-free surface β T transitional correction factor for heat transfer δ(k) thickness of bulk layer k cm δ Z i effective molecular diameter of Z i cm λ Z i mean free path of Z i in the gas phase cm θ s sorption layer surface coverage θ ss quasi-static layer surface coverage σ s,Z i effective molecular cross section of Z i cm 2 γ Z i uptake coefficient of Z i τ d,Z i desorption lifetime of Z i s ω Z i mean thermal velocity of Z i in the gas phase cm s −1 A(k) outer surface area of bulk layer k cm 2 A ss surface area of particle (quasi-static layer) cm 2 C g,Z i gas phase diffusion correction factor for Z i D b,Z i bulk diffusion coefficient of Z i cm 2 s −1 d d particle dry diameter cm d p particle diameter cm H a specific enthalpy of dry air J g −1 H v specific enthalpy of moist air J g −1 H w specific enthalpy of water vapor J g −1 Kn Knudsen number J ads,Z i flux of adsorption of Z i cm −2 s −1 J coll,Z i flux of surface collisions of Z i cm −2 s −1 J des,Z i flux of desorption of Z i cm −2 s −1 J b,ss,Z i flux of bulk-to-surface transport of Z i cm −2 s −1 J b1,ss,Z i flux of bulk layer 1-to-surface transport of Z i cm −2 s −1 J s,ss,Z i flux of sorption-to-quasi-static surface layer transport of Z i cm −2 s −1 J ss,s,Z i flux of quasi-static surface-to-sorption layer transport of Z i cm −2 s −1 J ss,b,Z i flux of quasi-static surface-to-bulk transport of Z i cm −2 s −1 J ss,b1,Z i flux of quasi-static surface-to-bulk layer 1 transport of Z i cm −2 s −1 k a,Z i first-order adsorption rate coefficient of Z i k b,ss,Z i rate coefficient (velocity) of bulk-to-quasi-static surface transport of Z i cm s −1 k b1,ss,Z i rate coefficient (velocity) of bulk layer 1-to-quasi-static surface transport of Z i cm s −1 k BR,Z p ,Z q second-order rate coefficients for bulk reactions of Z p with Z q cm 3 s −1 k b,ss,Z i rate coefficient (velocity) of bulk-to-quasi-static surface transport of Z j cm s −1 k b1,ss,Z i rate coefficient (velocity) of bulk layer 1-to-quasi-static surface transport of Z j cm s −1 k d,Z i first-order desorption rate coefficient of Z i s −1 k s,ss,Z i first-order rate coefficient for sorption-to-quasi-static surface transport of Z i s −1 k ss,b1,Z i rate coefficient (velocity) of surface-to-bulk layer 1 transport of Z i s −1 k ss,s,Z i first-order rate coefficients for quasi-static to sorption layer transport of Z i s −1 k SLR,Z p ,Z q second-order rate coefficients for surface layer reactions of Z p with Z q cm 2 s −1 K sol,cc,Z i gas-particle partitioning coefficient of Z i N Z i ,s absolute number of Z i at surface N Z i ,ss absolute number of Z i at quasi-static surface layer N Z i ,bk absolute number of Z i at bulk layer k P bk,Z i , L bk,Z i production (loss) rate of Z i by reaction in bulk layer k cm −3 s −1 P s,Z i , L s,Z i production (loss) rate of Z i by sorption surface layer reaction cm −2 s −1 P s,ss,Z i , L s,ss,Z i production (loss) rate of Z i by sorption -quasi-static surface layer reaction cm −2 s −1 P ss,Z i , L ss,Z i production (loss) rate of Z i by quasi-static surface layer reaction cm −2 s −1 bulk number concentration of Z i in the bulk layer k cm −3 V b volume of particle bulk cm −3 V (k) volume of bulk layer k cm −3 V Z i molecular volume of Z i cm −3