Kinetic multi-layer model of aerosol surface and bulk chemistry (KM-SUB): the influence of interfacial transport and bulk diffusion on the oxidation of oleic acid by ozone

Abstract. We present a novel kinetic multi-layer model that explicitly resolves mass transport and chemical reaction at the surface and in the bulk of aerosol particles (KM-SUB). The model is based on the PRA framework of gas-particle interactions (Poschl-Rudich-Ammann, 2007), and it includes reversible adsorption, surface reactions and surface-bulk exchange as well as bulk diffusion and reaction. Unlike earlier models, KM-SUB does not require simplifying assumptions about steady-state conditions and radial mixing. The temporal evolution and concentration profiles of volatile and non-volatile species at the gas-particle interface and in the particle bulk can be modeled along with surface concentrations and gas uptake coefficients. In this study we explore and exemplify the effects of bulk diffusion on the rate of reactive gas uptake for a simple reference system, the ozonolysis of oleic acid particles, in comparison to experimental data and earlier model studies. We demonstrate how KM-SUB can be used to interpret and analyze experimental data from laboratory studies, and how the results can be extrapolated to atmospheric conditions. In particular, we show how interfacial and bulk transport, i.e., surface accommodation, bulk accommodation and bulk diffusion, influence the kinetics of the chemical reaction. Sensitivity studies suggest that in fine air particulate matter oleic acid and compounds with similar reactivity against ozone (carbon-carbon double bonds) can reach chemical lifetimes of many hours only if they are embedded in a (semi-)solid matrix with very low diffusion coefficients (≤10−10 cm2 s−1). Depending on the complexity of the investigated system, unlimited numbers of volatile and non-volatile species and chemical reactions can be flexibly added and treated with KM-SUB. We propose and intend to pursue the application of KM-SUB as a basis for the development of a detailed master mechanism of aerosol chemistry as well as for the derivation of simplified but realistic parameterizations for large-scale atmospheric and climate models.


Introduction
Aerosols are ubiquitous in the atmosphere and have strong effects on climate and public health (Seinfeld and Pandis, 1998;Finlayson-Pitts and Pitts, 2000;Pöschl, 2005). Depending on chemical composition and surface properties, aerosol particles can act as condensation nuclei for cloud droplets and ice crystals, and they can affect the abundance of trace gases through heterogeneous chemical reactions (Ammann et al., 1998;Pöschl, 2005;Fuzzi et al., 2006;Andreae and Rosenfeld, 2008;Hallquist et al., 2009). Gasparticle interactions can also significantly change the physical and chemical properties of aerosols such as toxicity, reactivity, hygroscopicity and radiative properties (Pöschl, 2002;Rudich, 2003;Kanakidou et al., 2005;Rincon et al., 2009;Wiedensohler et al., 2009). Chemical reactions and mass transport lead to continuous transformation and changes in the composition of atmospheric aerosols ("chemical aging") (Schwartz and Freiberg, 1981;Hanson, 1997;Smith et al., 2003;Ammann and Pöschl, 2007).
Atmospheric particles consist of a wide variety of organic and inorganic chemical compounds which can exist in different liquid or (semi-)solid states (crystalline, amorphous, Published by Copernicus Publications on behalf of the European Geosciences Union. glassy, ultraviscous, gel-like) (Seinfeld and Pandis, 1998;Finlayson-Pitts and Pitts, 2000;Mikhailov et al., 2009). Chemical reactions can occur both at the surface and in the bulk of liquid and (semi-)solid particles. It is often difficult to discriminate surface and bulk reactions, and the relative importance of surface and bulk processes is not well understood (e.g., Moise and Rudich, 2000;Hearn et al., 2005;Pfrang et al., 2009). 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., 2008King et al., , 2009; and references therein). 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 double-layer 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. Springmann et al. (2009) 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. Shiraiwa et al. (2009) presented a kinetic double-layer surface model (K2-SURF) and master mechanism for the degradation of a wide range of polycyclic aromatic hydrocarbons (PAHs) by multiple photo-oxidants (O 3 , NO 2 , OH and NO 3 ) through different types of parallel and sequential surface reactions. Pfrang et al. (2009) developed a kinetic double-layer model coupling aerosol surface and bulk chemistry (K2-SUB), in which mass transport and chemical reactions in the particle are not explicitly resolved but represented by a reacto-diffusive flux (Danckwerts, 1951;Hanson, 1997).
Here we present 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. We demonstrate the applicability of KM-SUB for a common model system of organic aerosol chemistry, the ozonolysis of oleic acid droplets (e.g., Moise and Rudich, 2002;Katrib et al., 2004Katrib et al., , 2005King et al., 2004King et al., , 2009Knopf et al., 2005;Ziemann, 2005;Hearn and Smith, 2007;Pfrang et al., 2009;Sage et al., 2009;Vesna et al., 2009), and we compare the results of our numerical simulations with the results of earlier experimental and theoretical studies of this system.

Model description
The kinetic multi-layer model of aerosol surface and bulk chemistry (KM-SUB) builds on the formalism and terminology of the PRA framework Pöschl et al., 2007). A list of symbols is given in Appendix A.
As illustrated in Fig. 1, KM-SUB consists of multiple model compartments and layers, respectively: gas phase, near-surface gas phase, sorption layer, quasi-static surface layer, near-surface bulk, and a number of n bulk layers. The quasi-static surface layer has a monolayer thickness that corresponds to the (average) effective molecular diameter of non-volatile species Y j (δ Y j ). The thickness of bulk layers (δ) follows from the number of layers (n) and particle radius (r p ) : δ = (r p − δ Y j )/n.
The following processes are considered in KM-SUB: gas phase diffusion, gas-surface transport (reversible adsorption), surface layer reactions, surface-bulk transport, bulk diffusion, and bulk reactions. As outlined in the PRA framework, the following differential equations can be used to describe the mass balance of volatile species X i and nonvolatile species Y j for each model layer:

Fig. 1.
Kinetic multi-layer model (KM-SUB): (a) Model compartments and layers with corresponding distances from particle center (r p ± x), surface areas (A) and volumes (V ); λ X i is the mean free path of X i in the gas phase; δ X i and δ Y j are the thicknesses of sorption and quasi-static bulk layers; δ is the bulk layer thickness. (b) Transport fluxes (green arrows) and chemical reactions (red arrows).
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). For spherical particles, A(k) and V (k) can be calculated from the particle bulk radius (r b = r p − δ Y j ), bulk layer thickness (δ) and layer index number k as follows.
For planar geometry (thin films), the surface-to-volume ratio is constant: The surface-to-volume ratios of the layers play an important role in accounting for particle geometry. The absolute number of non-volatile molecules Y j in the particle (N Y j ) can be calculated as follows: where A ss is the surface area of the particle. In the following sub-sections we specify the formalisms used to describe and calculate fluxes or rates of gas-surface interactions, surfacebulk transport, bulk diffusion and bulk reactions.

Gas-surface interactions and surface layer reactions
The principles and details of gas-phase and gas-surface transport (J ads , J des ) and surface layer reactions (L s , L ss P s , P ss ) have been discussed and illustrated in the PRA framework Ammann and Pöschl, 2007) and in the K2-SURF study of aerosol surface chemistry ). Here we just briefly summarize the key aspects and equations. The fluxes of surface collisions, adsorption and desorption of a volatile species X i (J coll,X i , J ads,X i and J des,X i ) are given by: [X i ] gs is the near-surface gas phase concentration of X i , ω X i is the mean thermal velocity, α s,X i is the surface accommodation coefficient, and k d,X i is the desorption rate coefficient or inverse desorption lifetime, respectively (k d,X i = τ −1 d,X i ). α s,0,X i is the surface accommodation coefficient on an adsorbate-free surface, σ s,X i is the effective molecular cross section of X i in the sorption layer, and θ s is the surface coverage by adsorbed species.
The influence of changing chemical composition of the quasi-static particle surface on adsorbate-surface interactions (reversible adsorption) can be taken into account by describing α s,0,X i and τ d,X i as a linear combination of the parameter values α s,0,X i ,Y j and τ d,X i ,Y j that would be observed on pure substrates made up by the different surface components Y j weighted by their fractional surface area θ ss,Y j : For slow gas uptake and in the gas kinetic-regime (small uptake coefficient, large Knudsen number), [X i ] gs is practically the same as the average gas phase concentration of X i in the investigated system, [X i ] g . For rapid gas uptake in the continuum or transition regime (large uptake coefficient, small Knudsen number), the gas phase diffusion correction factor C g,X i can be used to calculate [X i ] gs from . The uptake coefficient (γ X i ) for X i is given by: γ gsr,X i is the reaction probability for elementary gas-surface reactions (potentially relevant for free radicals but assumed to be zero in the remainder of this study). General rate equations of chemical production and loss by surface layer reactions (P s, Pöschl et al., 2007). Different types of surface layer reactions can proceed within the sorption layer (P s,s,X i − L s,s,X i ), within the quasi-static surface layer (P ss,ss,Y j − L ss,ss,Y j ), and between sorption layer and quasi-static layer (P s,ss,X i −L s,ss,X i , P ss,s,Y j − L ss,s,Y j ). The generalized rate equations for reactions between volatile species X i in the sorption layer and non-volatile species Y j in the quasi-static layer are given by: Here c SLRv,s,X i and c SLRv,ss,Y j stand for the stoichiometric coefficients (negative for starting materials and positive for reaction products) of species X i and Y j in reaction SLRv; 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,Xp,Yq is a second-order reaction rate coefficient. In the exemplary simulations of the ozonolysis of oleic acid performed in this study and presented in Sect. 3, we consider only chemical loss and a second-order surface layer reaction between O 3 (X) and oleic acid (Y) with the following rate equation:

Surface-bulk transport
The surface-bulk transport of volatile species X i (J s,b,X i and J b,s,X i ) is defined as exchange between sorption layer and near-surface bulk, and the surface-bulk transport of nonvolatile species Y j (J ss,b,Y j and J b,ss,Y j ) is defined as exchange between quasi-static surface layer and near-surface bulk. Based on the PRA framework, surface to bulk transport can be described by: k s,b,X i and k ss,b,Y j are first-order transport rate coefficients (s −1 ) of X i and Y j , respectively. In the same way, bulk to surface transport (J b,s,X i and J b,ss,Y j ) can be described as follows: k b,s,X i and k b,ss,Y j (cm s −1 ) are transport rate coefficients, which can be regarded as effective transport velocities. Estimates for these effective transport velocities can be derived from the corresponding bulk diffusion coefficients D b,X i and D b,Y j (cm 2 s −1 ) by considering the average distance traveled by molecules diffusing in one direction: x = (4D b t/π) 1/2 (Atkins, 1998;Pfrang et al., 2009), On average, a molecule Y j in the near-surface bulk layer (k = 1) needs to travel a distance of x ≈ (δ+δ Y j )/2 to move into the quasi-static surface layer. The average time required to travel over this distance by diffusion is t = (δ+δ Y j ) 2 π/(16D b,Y j ). By division of travel distance over travel time we obtain the following estimate for the effective transport velocity: Under equilibrium conditions, mass conservation implies ] ss , and for pure Y j the surface and bulk concentrations are given by the inverse of the effective molecular cross section and of the effective molecular volume, respectively: In analogy, the bulk-to-surface transport velocity of X i can be calculated for an average travel distance of x ≈ (δ+δ X i +2δ Y j )/2 from the near-surface bulk layer into the sorption layer: An estimate for k s,b,X i can be determined by matching the rate coefficients for gas-surface and surface-bulk transport Atmos. Chem. Phys., 10, 3673-3691, 2010 www.atmos-chem-phys.net/10/3673/2010/ with a literature value or estimate for the gas-particle equilibrium partitioning coefficient (K sol,cc,X i ) or Henry's law coefficient, respectively : where k d,X i is the first-order rate coefficient of desorption.
From the fluxes and rates of gas-surface transport, surfacebulk transport and surface reaction follows the bulk accomodation coefficient of X i :

Bulk diffusion
Bulk diffusion is explicitly treated in the KM-SUB model 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 X i and Y j , k b,b,X i and k b,b,Y j (cm s −1 ), can be calculated from the corresponding diffusion coefficients. For this purpose we assume that each layer is homogeneously mixed (no concentration gradient within a layer) and that the average travel distance for molecules moving from one layer to the next is the layer thickness δ: This treatment of bulk diffusion yields practically the same results (concentration profiles) as the solving of partial differential equations (Smith et al., 2003), but it is more flexible and requires no assumptions about interfacial transport. As detailed in Appendix B, the ozone concentrations calculated for the near-surface bulk indicate that the assumption of Henry's law equilibrium is not a realistic boundary layer condition for the PDE method when applied to reactive systems. The influence of changing chemical composition of the particle bulk on diffusion can be taken into account by describing D b,X i or D b,Y j as a linear combination of the parameter values that would be observed in pure bulk Y p : The weighting factor b,Y p could be the mole, mass, or volume fraction of Y p in the bulk .

Bulk reaction
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 quasi-elementary 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 X i or Y j within the bulk layer k.
Here c BRv,X i and c BRv,Y j stand for the stoichiometric coefficients (negative for starting materials and positive for reaction products) of species X i and Y j 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,Xp and k BRv,Yp are first-order reaction rate coefficients and k BRv,Xp,Xq , k BRv,Xp,Yr , k BRv,Yp,Yq and k BRv,Xr,Yp are second-order bulk reaction rate coefficients between X p and X q , X p and Y r , Y p and Y q , and X r and Y p respectively, in the condensed phase bulk of a system with multiple 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. In the exemplary simulations of the ozonolysis of oleic acid performed in this study and presented in Sect. 3, we consider only chemical loss and a second-order bulk reaction between O 3 (X) and oleic acid (Y) with the following rate equation:

Model application: oxidation of oleic acid by ozone
To test and demonstrate the applicability of the KM-SUB model, we simulated the oxidation of oleic acid particles by ozone in comparison to experimental data from Ziemann (2005). The same data set has recently been used by Pfrang et al. (2009) for simulations with a kinetic doublelayer model (K2-SUB; Appendix C). The gas phase ozone concentration was set to [X] g = [X] 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 (X) were set to Table 1. Kinetic parameters for the interaction of ozone (X) and oleic acid (Y) in different model scenarios (base cases 1-3).

Parameter
Base case 1 Base case 2 Base case 3 10 −10 10 −10 10 −15 k SLR,X,Y (cm 2 s −1 ) 6.0×10 −12 6.0×10 −12 6.0×10 −12 [X] s,0 = [X] bk,0 = 0, and the initial surface and bulk concentrations of oleic acid (Y) were set to [Y] ss,0 =9.7×10 13 cm −2 and [Y] ss,0 =1.2×10 21 cm −3 , respectively. Accordingly, the initial value of the total number of oleic acid molecules in particles with a radius of 0.2 µm was N Y,0 = 4.1 × 10 7 . We modeled the temporal evolution of the particle surface and bulk composition and of the ozone uptake coefficient by numerically solving the differential equations of mass balance for each model compartment with Matlab software (ode23tb solver with 999 time steps). The rate equations describing the ozone -oleic acid reaction system are listed in Appendix D. The kinetic parameters required for the model simulations are summarized in Table 1: the surface accommodation coefficient of ozone (α s,0,X ), the desorption lifetime of ozone (τ d,X ), the bulk diffusion coefficients of ozone and oleic acid (D b,X , D b,Y ), the secondorder surface layer reaction rate coefficient (k SLR,X,Y ), and the second-order bulk reaction rate coefficient (k BR,X,Y ). Additional input parameters were the mean thermal velocity of ozone (ω X = 3.6 × 10 4 cm s −1 ), the Henry's law coefficient (K sol,cc,X = 4.8 × 10 −4 mol cm −3 atm −1 ), and the molecular diameters of oleic acid (δ Y = 0.8 nm) and ozone (δ X = 0.4 nm) . The model simulations were performed with n = 100 layers corresponding to a layer thickness of δ = 1.99 nm for r p =0.2 µm, unless mentioned otherwise. For simplicity the physico-chemical parameters were assumed to be constant throughout each model run.
The first-generation products of oleic acid oxidation by ozone are mainly 1-nonanal, 9-oxononanoic acid, nonanoic acid, and azelaic acid (Moise and Rudich, 2002;Katrib et al., 2004;Thornberry and Abbatt, 2004;Vesna et al., 2009). 1-Nonanal is highly volatile and likely to evaporate from the particle (Sage et al., 2009). The other products, however, have higher molecular masses and are more polar, less volatile and likely to remain in the condensed phase ). Moreover, they may undergo recombination reactions forming second-generation products of higher molecular mass such as dimers or oligomers in the bulk (Rudich et al., 2007, and references therein). These effects go beyond the scope of the present study, but the gas-particle partitioning of (semi-)volatile species and the effects of chem-ical transformation on particle size and properties shall be incorporated in follow-up studies.
Here we focus on three model scenarios (base cases) characterizing the influence of different parameters and conditions on the uptake of ozone and the decay of oleic acid. The derivation and choice of kinetic parameters were discussed in detail by Pfrang et al. (2009). In view of the uncertainties and limited availability of experimental data we compare three cases covering a range of plausible parameter variations.
In base case 2 (BC2, kinetic limitation by bulk reaction) we assumed slow bulk reaction with a rate coefficient ∼30 times lower than in BC1 (k BR,X,Y = 5×10 −17 cm 3 s −1 ). The adsorption parameters were re-adjusted to match the experimental data (α s,0,X = 8.5 × 10 −4 , τ d,X = 0.001 s), and all other parameters were kept equal to BC1. Note that in BC1 and BC2 the results can be reproduced with different combinations of α s,0,X and τ d,X , that are closer to prediction of molecular dynamic simulations (e.g., α s,0,X ≈ 10 −2 and τ d,X ≈ 10 −9 s) (Vieceli et al., 2005;Shiraiwa et al., 2009). These aspects will be further investigated in follow-up studies.
In base case 3 (BC3, kinetic limitation by bulk diffusion) we assumed slow mass transport in the bulk with diffusion coefficients that are characteristic for amorphous (semi-)solid matrices (Bird et al., 2007;Swallen et al., 2007;Mikhailov et al., 2009) and five orders of magnitude lower than in BC1 and BC2 (D b,X = 10 −10 cm 2 s −1 and D b,Y = 10 −15 cm 2 s −1 ). Figure 2 illustrates the model results of KM-SUB in base case 1 with the kinetic parameters specified in Table 1. As shown in Fig. 2a, the simulated decay of oleic acid is in very good agreement with the experimentally observed decay. Throughout the experimental time scale of 1-15 s, the simulated ozone uptake coefficient is nearly constant and identical to the surface and bulk accommodation coefficients (γ X ≈ α b,X ≈ α s,X ≈ α s,0,X ≈ 4 × 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 Atmos. Chem. Phys., 10, 3673-3691, 2010 www.atmos-chem-phys.net/10/3673/2010/ limited by the process of surface accommodation. After ∼30 s the ozone uptake coefficient rapidly drops off when reaction partner oleic acid is used up by the rapid bulk reaction (γ X /γ X,0 ≈ N Y /N Y,0 < 1% after ∼37 s). As shown in Fig. 2b, the surface concentration of oleic acid also decreases gradually over the first ∼30 s and drops off rapidly afterwards. In contrast, the surface concentration of ozone exhibits a rapid initial increase from zero to a plateau level of [X] s ≈ 10 11 cm −2 , which is reached within less than a second (∼20 ms) and determined by the combination of reversible adsorption, surface reaction, and surfaceto-bulk transport driven by the chemical reaction in the bulk. Over the first 30 s, [X] s gradually increases with the decay of oleic acid and the related decrease of chemical loss. As the chemical loss by reaction with oleic acid rapidly drops off after ∼30 s, [X] s swiftly increases by an order of magnitude to a steady-state level of ∼10 12 cm −2 , which is governed by reversible adsorption.

Base case 1: kinetic limitation by interfacial transport
To test how the number of model layers in the particle bulk affects the simulation results, we have run the model with n = 1, 5 and 100 layers corresponding to layer thicknesses of δ = 199.2, 39.8 and 1.99 nm, respectively. The results obtained with n = 100 and 5 were practically identical, demonstrating the robustness of the multi-layer model approach with transport rate coefficients (velocities) scaled by the layer thickness (Eqs. 28-31). As shown in Fig. 2, the deviations obtained with n = 1 were relatively minor, but they indicate that the particle bulk cannot be regarded as well mixed under the conditions of base case 1. Similar deviations were obtained with the kinetic double-layer model K2-SUB, in which the bulk processes are not explicitly resolved but represented by a reacto-diffusive flux (Appendix C, Pfrang et al., 2009).
The KM-SUB model results for the bulk concentration profiles of ozone and oleic acid were also essentially the same with n = 100 and 5. For high resolution and to avoid congestion of the plot, however, only the profiles obtained 100 layers are shown in Fig. 3. The y-axis indicates the radial distance from the particle center (r) normalized by the particle radius (r p ), ranging from the particle core (δ/r p ≈ 0) to the near-surface bulk ((r p − δ)/r p ≈ 1).
As shown in Fig. 3a, ozone rapidly diffuses into the particle bulk. A gradient between the near-surface bulk and the core is established within less than a second (∼20 ms), and the concentration profile is determined by the interplay of interfacial mass transport (surface and bulk accommodation) with bulk diffusion and chemical reaction. During the first few seconds, the ozone concentration in the particle core is about a factor of ∼50 lower than in the near surface bulk. Up to ∼30 s, the ozone concentration gradient decreases gradually with the decay of oleic acid and the related decrease of chemical loss. As the chemical loss by reaction with oleic acid rapidly drops off after ∼30 s, the ozone concentration gradient swiftly relaxes, and after ∼37 s ozone is well mixed throughout the particle bulk at the concentration level of sol- ubility saturation ([X] b,sat = 8 × 10 14 cm −3 = K sol,CC,X [X] gs with K sol,CC,X = 4.8 × 10 −4 mol cm −3 atm −1 ).
As shown in Fig. 3b, the ozone gradient and the ozonolysis of oleic acid cause a reverse gradient in the bulk concentration of oleic acid. Because the concentration of oleic acid is several orders of magnitude larger than that of ozone, however, the strong gradient of ozone induces only a small gradient in oleic acid. During the first few seconds, the oleic acid concentration in the near-surface bulk is ∼5% lower than in the particle core. As the ozone gradient relaxes after ∼30 s, the small oleic acid gradient also disappears. Thus, oleic acid can effectively be regarded as well-mixed. Figure 4a shows the loss rate of oleic acid and ozone by chemical reaction in the particle bulk 1.1e+21 1e+21 9e+20 8e+20 7e+20 6e+20 5e+20 4e+20 3e+20 2e+20 1e+20 1e+19 (b) . The y-axis indicates the radial distance from the particle center (r) normalized by the particle radius, ranging from the particle core (r/r p = δ/r p ≈ 0) to the near-surface bulk (r/r p = 1). The colored lines are isopleths of bulk concentration with labeled in units of cm −3 (blue=low, red=high). The black solid line represents the reacto-diffusive length for ozone normalized by the particle radius (1 − l rd,X /r p ). The black dotted line represents the 50% isoline for the concentration difference between the particle core and the nearsurface bulk, and the black dashed line shows 63% isoline of L * b .
of KM-SUB with K2-SUB and with the underlying resistor model formulation of the reacto-diffusive flux J b,rd . Due to the overall kinetic limitation by interfacial transport, the volume-average loss rate L b remains near-constant up to ∼30 s. During the first few seconds, the loss rate in the nearsurface bulk (L b1 , r/r p ≈ 1) is more than a factor of 2 higher than the volume-average loss rate, while the loss rates in layers around the particle core (L b40 to L b100 , r/r p < 0.6) are more than a factor of 2 lower. L b1 continually decreases due to progressing consumption of oleic acid, whereas L b40 to L b100 exhibit a pronounced increase up to ∼35 s due to the penetration of ozone to the particle core (see also Fig. 3a).
After ∼35 s all loss rates drop off rapidly and approach zero after ∼37 s as oleic acid is essentially used up. Up to ∼30 s, the strong radial gradient of loss rates L b and the approximate coincidence of L b20 (r/r p ≈ 0.8) with L b in- dicate that the overall rate of conversion or absolute loss rate of oleic acid and ozone molecules in the particle is strongly dominated by the outer bulk layers at r/r p > 0.8. Due to the spherical geometry the outer layers also have much larger volumes, which further enhances their relative importance with regard to the absolute loss rate of molecules in each layer (L * bk = L bk V bk ) and in the entire particle (L The black solid line in Fig. 3a represents the reactodiffusive length for ozone normalized by the particle radius (1 − l rd,X /r p ), which can be regarded as the distance from the surface up to which the chemical reaction proceeds effectively (Finlayson-Pitts and Pitts, 2000;Pöschl et al., 2007). The initial value of l rd,X ≈20 nm is consistent with the values reported by Smith et al. (2003) and Katrib et al. (2004). Over the first ∼30 s of the model run l rd,X remains fairly constant, Atmos. Chem. Phys., 10, 3673-3691, 2010 www.atmos-chem-phys.net/10/3673/2010/ and the reacto-diffusive line (black solid) coincides with the 50% isoline for the concentration difference between the particle core and the near-surface bulk (black dotted), and with the 63% (=1−1/e) isoline of L * b (black dashed). During this time 63% of L * b , i.e., 63% of the overall rate of conversion of oleic acid ozonolysis, occur in the outermost 10 bulk layers, i.e. at r/r p > 0.9 (corresponding to ∼30% of the particle bulk volume). After ∼30 s, l rd,X increases steeply as oleic acid is depleted and ozone can freely penetrate throughout the bulk of the particle.

Base case 2: kinetic limitation by bulk reaction
In base case 2 the second-order bulk reaction rate coefficient was by a factor of ∼30 smaller than in BC1, while the surface accommodation coefficient was by a factor of ∼2 higher (Table 1). As shown in Fig. 5a, the simulation results of BC2 are also in very good agreement with the observed decay of oleic acid. The calculated ozone uptake coefficient behaves very differently than in BC1. The initial value of γ X in BC2 is by a factor of ∼2 higher than in BC1. It reflects rapid uptake into the initially ozone-free particle bulk, which is limited by the kinetics of interfacial transport, i.e., bulk accommodation and surface accommodation, respectively. After less than one second (∼20 ms), however, γ X steeply drops off to a level that is reaction-limited (γ X ≈ 10 −4 < α b,X ≈ 8 × 10 −4 ) and continues to decrease gradually with decreasing abundance of the reaction partner oleic acid.
After a rapid initial increase (∼20 ms), the ozone surface concentration reaches a plateau level which increases only slightly until the end of the simulation, while the oleic acid surface concentration decreases gradually. In contrast to BC1, none of the compositional and kinetic parameters exhibit abrupt changes after ∼30 s. Instead, they all undergo a gradual evolution reflecting the kinetics of the rate-limiting bulk reaction.
The results of multi-layer model KM-SUB with both n = 100 and 1 were practically identical to the results of the double-layer model K2-SUB, indicating that the particle bulk can be regarded as well mixed under the conditions of BC2.
Due to faster interfacial transport and slower chemical reaction, the ozone concentration in the bulk increases much faster than in BC1 (Fig. 6a vs. Fig. 3a). After less than one second (∼20 ms), the ozone concentration in the core of the particle is only ∼30% lower than in the near-surface bulk, and the gradient continues to decrease with decreasing abundance of oleic acid. Due to the slow decay of oleic acid, however, a small ozone gradient persists until the end of the simulation (Fig. 6a) and thus longer than in BC1, where the initially very strong gradient essentially disappears after ∼37 s when practically all oleic acid (>99%) has been consumed by the rapid bulk reaction (Fig. 3a). As illustrated in Fig. 6b, oleic acid can be regarded as well-mixed (concentration differences <1%).  Figure 7a shows the loss rate of oleic acid and ozone by chemical reaction in the particle bulk (L b ) as calculated with KM-SUB for different model layers and with KM-SUB and K2-SUB effectively averaged over the whole bulk volume. The average values of L b calculated by K2-SUB and KM-SUB are almost identical. Due to the kinetic limitation by bulk reaction, the volume-average loss rate L b decreases continuously due to consumption of oleic acid. Initially, L b in the near-surface bulk (L b1 , r/r p ≈ 1) is ∼40% higher than that in the core (L b100 ) corresponding to the concentration gradient of O 3 in the bulk. They become almost same in 20 s as concentration gradient becomes smaller. Figure 7b illustrates L * bk (= L bk V bk ) and L * b (= L b V b ) corresponding to the values of L bk and L b displayed in Fig. 7a (note the linear scaling of Fig. 7a and the logarithmic scaling of Fig. 7b). L * b calculated by KM-SUB and K2-SUB are almost identical.
The black solid line in Fig. 6a represents again the reactodiffusive length for ozone normalized by the particle radius (1 − l rd,X /r p ). Due to the lower reaction rate coefficient, the initial value of l rd,X ≈ 120 nm is by a factor of ∼6 larger . The y-axis indicates the radial distance from the particle center (r) normalized by the particle radius, ranging from the particle core (r/r p = δ/r p ≈ 0) to the near-surface bulk (r/r p = 1). The colored lines are isopleths of bulk concentration with labeled in units of cm −3 (blue=low, red=high). The black solid line represents the reacto-diffusive length for ozone normalized by the particle radius (1 − l rd,X /r p ). The black dotted line represents the 50% isoline for the concentration difference between the particle core and the nearsurface bulk, and the black dashed line shows 63% isoline of L * b .
than in BC1, and after ∼20 s it is already larger than the particle radius (200 nm), indicating that the reaction effectively proceeds throughout the whole bulk. The 50% isoline (black dotted) for the ozone concentration difference between the particle core and the near-surface bulk and the 63% (= 1−1/e) isoline of L * b (black dashed) remain near-constant throughout the simulation. They are from the beginning at the same level reached at the end of BC1, which reflects that both ozone and oleic can be regarded as well mixed in BC2 due to slow reaction and rapid diffusion.
BC1 may be regarded as more realistic, because it uses the only reported value of the bulk reaction rate coefficient of ozone with oleic acid, whereas BC2 assumes a ∼30 times lower value. Nevertheless, both model cases are in good agreement with the available experimental data. Thus, fur- ther experiments covering a wider range of reaction times and conditions are needed to elucidate the actual reaction mechanism (see Pfrang et al., 2009).
As illustrated in Fig. 8a, the decay of oleic acid is very slow in BC3 because the low diffusivity effectively protects oleic acid in the particle bulk from oxidation. Within ∼1 s the calculated ozone uptake coefficient drops from the initial value of ∼4×10 −3 to ∼2×10 −5 , which is due to a rapid decrease of the oleic acid surface concentration from ∼10 14 cm −2 to ∼10 12 cm −2 . After ∼1 s the ozone surface concentration remains practically constant at ∼10 12 cm −2 (adsorption equilibrium). The low levels of ozone uptake Atmos. Chem. Phys., 10, 3673-3691, 2010 www.atmos-chem-phys.net/10/3673/2010/ coefficient and oleic acid surface concentration are maintained by surface bulk exchange and decrease slowly as ozone increases and oleic acid decreases also in the nearsurface bulk. In contrast to BC1, the bulk accommodation coefficient α b , X(≈ 2 × 10 −6 ) is much smaller than α s,X and γ X , reflecting the slow uptake of ozone into the particle bulk, which is limited by the kinetics of surface-bulk exchange and bulk diffusion.
Due to fast reaction and slow diffusion, both ozone and oleic acid exhibit steep concentration gradients near the surface, whereas the inner particle bulk remains nearly unchanged (Fig. 9). During the simulation period ozone penetrates only into the near-surface bulk, and thus the chemical consumption of oleic acid is also restricted to the nearsurface bulk (r/r p > 0.99). Accordingly, the reacto-diffusive length is very small (∼0.1 nm; black solid line in Fig. 9a), and the 50% isoline (black dotted) for the ozone concentration difference between particle core and near-surface bulk as well as the 63% (= 1−1/e) isoline of L * b (black dashed) . The y-axis indicates the radial distance from the particle center (r) normalized by the particle radius, ranging from r/r p = 0.95 to the near-surface bulk (r/r p = 1). The colored lines are isopleths of bulk concentration with labeled in units of cm −3 (blue=low, red=high). The black solid line represents the reacto-diffusive length for ozone normalized by the particle radius (1 − l rd,X /r p ). The black dotted line represents the 50% isoline for the concentration difference between the particle core and the near-surface bulk, and the black dashed line shows 63% isoline of L * b .
remain restricted to the near surface bulk (r/r p >0.99). This is also reflected in Fig. 10, showing that the values of L b and L * b are essentially zero except for k = 1 (near-surface bulk layer).
Compared to KM-SUB, the simulated decay of oleic acid and ozone uptake are faster in K2-SUB. At the end of the model run (40 s), the total amount of oleic acid (N Y ) is lower by ∼10% while the surface concentration of oleic acid ([Y] ss ) and the ozone uptake coefficient (γ X ) are higher by a factor of ∼5. This is because K2-SUB assumes that oleic acid is radially well-mixed and does not resolve its concentration profile. In fact, the diffusion coefficient of oleic acid (D b,Y ) does not appear in the equation for the reactodiffusive flux (Eq. C1), which represents the effects of all bulk processes in the traditional resistor model approach. In case of kinetic limitation by bulk diffusion, however, the reaction rates are influenced by the mass transport of both volatile and non-volatile reactants in the particle bulk.

Chemical half-life of oleic acid and atmospheric implications
In the preceding sections we showed how the KM-SUB model can be used to interpret and analyze experimental data from laboratory studies, which are often carried out with high O 3 concentration to observe measurable effects on short experimental time scales. Here we demonstrate how the results can be extrapolated to dilute atmospheric conditions. We estimated the chemical half-life (t 1/2 ) of pure oleic acid particles, when exposed to ozone at average ambient mixing ratios (<150 ppb). The radii of oleic acid particles (r p ) were set to be 0.2 and 1.0 µm. The chemical half-life was defined as the time when number of oleic acid molecules was degraded to half of its initial value. The solid lines of Fig. 11 show the results of such calculations with BC1. At typical atmospheric O 3 mixing ratios of 30 ppb, t 1/2 was ∼25 min for r p =0.2 µm and ∼130 min for r p =1.0 µm. The chemical half-life depended strongly on particle radius (i.e. number of oleic acid molecules contained in the particle). In the typical polluted air with 100 ppb O 3 , t 1/2 was shortened to ∼8 min for r p = 0.2 µm and ∼40 min for r p = 1.0 µm. The calculations were also made in BC2, which gave almost the same results for r p = 0.2 µm and gave ∼20% smaller value for r p = 1.0 µm compared to BC1. The ozonolysis of oleic acid leads to an increase of oxygen content and hygroscopicity (Asad et al., 2004), and it can convert CCN-inactive particles into CCN (King et al., 2009). Our model results suggest that the transformation of oleic acid particles will occur on a timescale of ∼1 h depending on atmospheric O 3 concentration and particle size. This is in good agreement with the characteristic time of 1.3 h reported by King et al. (2009).
In spite of the short lifetimes found in laboratory experiments, however, oleic acid is detected also in aged atmospheric aerosol particles (Rogge et al., 1991;Morris et al., 2002). Possible explanations include reduced concentrations of O 3 (and other photo-oxidants) in the particle bulk due to the following effects: (a) chemical reaction with other reactive species, self-reaction or catalytic decomposition; (b) competitive co-adsorption and surface reaction of multiple species such as water vapor and nitrogen oxides Springmann et al., 2009); (c) slow mass transport because of low diffusion coefficients in solid or semisolid phases (crystalline, glassy, rubbery, gel-like or ultraviscous: Mikhailov et al., 2009). The formation of semisolid amorphous phases is generally favored by low temperature and low relative humidity (Zobrist et al., 2008), and it can occur in the core of atmospheric particles as well as Atmos. Chem. Phys., 10, 3673-3691, 2010 www.atmos-chem-phys.net/10/3673/2010/ in coatings formed by condensation of secondary particulate matter. For example, oligomerisation reactions may lead to high molecular mass, low hygroscopicity and low diffusivity of secondary organic coatings, effectively shutting off further ozonolysis and oxidation of oleic acid and other organic compounds in the bulk of coated particles. The potential importance of the phase state of aerosol particles (solid, semi-solids, or liquid) for the ozonolysis of oleic acid was already pointed out by Knopf et al. (2005) and Zahardis and Petrucci (2007). To explore and quantify the effects and of kinetic limitation by bulk diffusion, we performed sensitivity studies with low diffusion coefficients that are characteristic for mass transport in highly viscous or semi-solid matrices (Bird et al., 2007;Swallen et al., 2007;Mikhailov et al., 2009).
The dotted lines in Fig. 11 show the effect of reducing the diffusion coefficients of O 3 and oleic acid by three orders of magnitude (D b,X = 10 −8 cm 2 s −1 , D b,Y = 10 −13 cm 2 s −1 ). It increased t 1/2 by up to ∼100%, depending on ozone concentration and particle diameter. The dashed lines shows the effect of reducing the diffusion coefficients to the same level as in the base case 3 presented and discussed above (D b,X = 10 −10 cm 2 s −1 , D b,Y = 10 −15 cm 2 s −1 ). In this case the chemical lifetimes increase by up to a factor of ∼100 and become nearly independent of ozone concentration: ∼3 h for r p = 0.2 µm and ∼3 days for r p = 1.0 µm. The results confirm mass transport can have strong non-linear effects on the chemical aging of atmospheric aerosols, and that (semi-)solid matrices with low diffusion coefficients can effectively shield reactive organic compounds from degradation by atmospheric photo-oxidants.

Conclusions
We presented a novel kinetic multi-layer model (KM-SUB) that explicitly resolves all steps of mass transport and chemical reaction at the surface and in the bulk of aerosol particles. It includes adsorption and desorption, surface reactions and surface-bulk exchange as well as bulk diffusion and reaction.Unlike earlier models of aerosol chemistry (e.g., resistor model), the KM-SUB model approach does not require any simplifying assumptions about steady-state and mixing conditions. The temporal evolution of concentration profiles of volatile and non-volatile species in the particle bulk can be modeled along with surface concentrations and uptake coefficients.
The effects of bulk diffusion on the rate of gas uptake and particle transformation were explored for a popular reference system, the ozonolysis of oleic acid particles, assuming different diameters and different sets of kinetic parameters (surface accommodation coefficients, diffusion coefficients, etc.).
Under conditions where the reaction system was kinetically limited by interfacial transport or chemical reaction, the multi-layer model (KM-SUB) was in good agreement with a double-layer model (K2-SUB) using traditional resistor formulations for bulk processes. In case of kinetic limitation by bulk diffusion, however, the K2-SUB model overestimated the rates of gas uptake and oleic acid degradation, because it does not resolve and account for the concentration profile of oleic acid. Comparison with a partial differential equation (PDE) model confirmed the ability of KM-SUB to resolve diffusion and concentration profiles in the particle bulk.
We demonstrated and discussed how the KM-SUB model can be used to interpret and analyze experimental data from laboratory studies, and how the results can be extrapolated to atmospheric conditions. For liquid oleic acid particles, the chemical half-life of oleic acid molecules with regard to oxidation by ozone was estimated to of the order of one hour. In (semi-)solid particles with low bulk diffusivity, however, the chemical half-life can increase dramatically (up to multiple hours and days). The results of sensitivity studies with different bulk diffusion coefficients confirmed that the phase state of atmospheric particles is highly relevant for chemical transformation and aging (Knopf et al., 2005;McNeill et al., 2007McNeill et al., , 2008Zahardis and Petrucci, 2007;Griffiths et al., 2009;Mikhailov et al., 2009).
Depending on the complexity of the investigated system, unlimited numbers of volatile and non-volatile species and chemical reactions can be added flexibly. Thus, the KM-SUB approach may serve as a basis for the development of a detailed master mechanism of aerosol chemistry and for the derivation of simplified but realistic parameterizations for large-scale atmospheric and climate models.

Symbol
Meaning Unit surface accommodation coefficient of X i on an adsorbatefree surface α s,0,X i ,Y j surface accommodation coefficient of X i on an adsorbatefree surface composed of Y j δ thickness of bulk layers cm A ss surface area of particle (quasi-static layer) c SLRv,s,X i , c SLRv,ss,Y j stoichiometric coefficients of species X i and Y j in surface layer reaction SLRv c BRv,X i , c BRv,Y j stoichiometric coefficients of species X i and Y j in bulk reaction BRv C b,rd,X i reacto-diffusive geometry correction factor for X i C g,X i gas phase diffusion correction factor for X i D b,X i , D b,Y j bulk diffusion coefficient of X i and Y j cm 2 s −1 flux of transport from bulk layer k to bulk layer k-1 for X i and Y j cm −2 s −1 first-order rate coefficients for surface-bulk transport of Y j s −1 k SLR,Xp,Yq second-order rate coefficients for surface layer reactions of X p with Y q cm 2 s −1 K sol,cc,X i gas-particle partitioning coefficient of X i l rd,X i reacto-diffusive length of X i and Y j cm −3 s −1 loss rate of X i and Y j by reaction in the particle bulk cm −3 s −1 L bk,X i , L bk,Y j loss rate of X i and Y j by reaction in bulk layer k cm −3 s −1 L s,X i , L ss,Y j loss rate of X i and Y j by surface layer reaction absolute loss rate of X i and Y j by reaction in the particle bulk absolute loss rate of X i and Y j by reaction in bulk layer k ss,Y j absolute loss rate of X i and Y j by surface layer reaction s −1 n number of bulk layers N Y total number of Y molecules P b,X i , P b,Y j production rate of X i and Y j by reaction in the particle bulk cm −3 s −1 P bk,X i , P bk,Y j production rate of X i and Y j by reaction in bulk layer k cm −3 s −1 P s,X i , P ss,Y j production rate of X i and Y j by surface layer reaction cm −2 s −1 r b particle bulk radius cm r p particle radius cm  Smith et al. (2003) demonstrated that the diffusion and reaction in the bulk can be described by the following partial differential equations: where X=O 3 and Y=oleic acid. The first term describes diffusion based on the Fick's law and the second term describes chemical reaction in the bulk. Note that this partial differential equation (PDE) method can only simulate the processes in the bulk and not coupled with surface processes. In solving these PDE equations, several assumptions are required as boundary conditions (Smith et al., 2003): (1) The concentration of the near-surface bulk is fixed over time, (2) the flux of X at the core is zero (∂[X] b /∂r| r=0 = 0), (3) the flux of Y at both near-surface bulk and core is zero (∂[X] b /∂r| r=r p = ∂[X] b /∂r| r=0 = 0). PDE model simulations were performed by numerically solving the partial differential equations with Matlab software (PDEPE solver). We investigate large particles with a radius of 1.0 µm as investigated by Smith et al. (2003), and we compare the KM-SUB results (n = 200) to the results of the PDE model. The gas phase ozone concentration was set to [X] g = [X] 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 (X) were set to [X] s,0 = [X] bk,0 = 0, and the initial surface and bulk concentrations of oleic acid (Y) were set to [Y] ss,0 = 9.7 × 10 13 cm −2 and [Y] ss,0 = 1.2 × 10 21 cm −3 , respectively. The kinetic parameters used are based on base case 1 summarized in Table 1 with the second-order surface reaction rate coefficient (k SLR,X,Y ) of zero. Unlike KM-SUB, the PDE model depends on the assumption of a fixed ozone concentration in the near surface bulk. Smith et al. (2003) had assumed solubility saturation according to Henry's law. In case of fast reaction, however, chemical loss leads to a substantial decrease and lower value of the actual ozone concentration. For comparability, the concentration of ozone in the near-surface bulk of the PDE model was set to a value of 5.0×10 13 cm −3 , which is similar to the value calculated by KM-SUB. Figure A1 shows the concentration profiles of ozone (top) and oleic acid (bottom) obtained by simulations with KM-SUB (left) and with the PDE model (right). The diffusion of O 3 into the bulk was very slow due to the large particle diameter, so that there was a large concentration gradient of O 3 . The concentration of O 3 at the core remains small (<10 3 cm −3 ) over the reaction time considered. Oleic acid showed also a concentration gradient, which was due to the O 3 concentration gradient and a slow diffusion coefficient of oleic acid (D b,Y = 10 −10 cm 2 s −1 ). If a larger D b,Y of 10 −7 cm 2 s −1 was used in the simulation, the concentration of oleic acid was homogeneous throughout the bulk, which was consistent with the results reported by Smith et al. (2003).
The KM-SUB and PDE methods yield very similar results, but the concentration gradients for both O 3 and oleic acid were slightly steeper in the PDE model. A possible reason for the difference may be that the PDE model requires the assumption of a fixed ozone concentration in the near surface bulk and zero flux values at the surface and at center of the particles as specified above (boundary conditions). Another possible reason for deviations is the estimation of bulk transport velocities (k b,b ). If we assume k b,b = D b /δ instead of Eqs. (35) and (36), the KM-SUB results become nearly identical to the PDE results. Note that varying the number of layers in the KM-SUB model made practically no difference.

Kinetic double-layer model (K2-SUB)
A kinetic double-layer of surface and bulk chemistry (K2-SUB) was developed by Pfrang et al. (2009). K2-SUB uses one layer for bulk, therefore concentration profiles in the bulk is not obtained. The diffusion and reaction in the bulk is considered by the reacto-diffusive flux (J b,rd ) (Danckwerts, 1951;Hanson, 1997;Pöschl et al., 2007). The loss rate of oleic acid in the bulk can be described as follows , which is based on the resistor model formulation (Hanson, 1997): The y-axis indicates the radial distance from the particle center (r) normalized by the particle radius, ranging from the particle core (r/r p ≈ 0) to the near-surface bulk (r/r p = 1). The colored lines are isopleths of bulk concentration with labeled in units of cm −3 (blue = low, red = high).