Current State of Atmospheric Aerosol Thermodynamics and Mass Transfer Modeling: A Review

Abstract

A useful aerosol model must be able to adequately resolve the chemical complexity and phase state of the wide particle size range arising from the many different secondary aerosol growth processes to assess their environmental and health impacts. Over the past two decades, significant advances in understanding of gas-aerosol partitioning have occurred, particularly with respect to the role of organic compounds, yet aerosol representations have changed little in air quality and climate models since the late 1990s and early 2000s. The gas-aerosol partitioning models which are still commonly used in air quality models are separate inorganics-only thermodynamics and secondary organic aerosol (SOA) formation based on absorptive partitioning theory with an assumption of well-mixed liquid-like particles that continuously maintain equilibrium with the gas phase. These widely used approaches in air quality models for secondary aerosol composition and growth based on separated inorganic and organic processes are inadequate. This review summarizes some of the important developments during the past two decades in understanding of gas aerosol mass transfer processes. Substantial increases in computer performance in the last decade justify increasing the process detail in aerosol models. Organics play a central role during post-nucleation growth into the accumulation mode and change the hygroscopic properties of sulfate aerosol. At present, combined inorganic-organic aerosol thermodynamics models are too computationally expensive to be used online in 3-D simulations without high levels of aggregation of organics into a small number of functional surrogates. However, there has been progress in simplified modeling of liquid-liquid phase separation (LLPS) and distinct chemical regimes within organic-rich and inorganic-rich phases. Additional limitations of commonly used thermodynamics models are related to lack of surface tension data for various aerosol compositions in the small size limit, and lack of a comprehensive representation of surface interaction terms such as disjoining pressure in the Gibbs free energy which become significant in the small size limit and which affect both chemical composition and particle growth. As a result, there are significant errors in modeling of hygroscopic growth and phase transitions for particles in the nucleation and Aitken modes. There is also increasing evidence of reduced bulk diffusivity in viscous organic particles and, therefore, traditional secondary organic aerosol models, which are typically based on the assumption of instantaneous equilibrium gas-particle partitioning and neglect the kinetic effects, are no longer tenable.

1. Introduction

Aerosols have substantial impacts on the radiative balance of the atmosphere both directly and indirectly via their influence on cloud composition and evolution [1]. They are a central element in the feedback of chemical and microphysical processes on atmospheric dynamics and transport. Aerosols affect the chemical composition of the atmosphere and thus its radiative balance through heterogeneous chemistry [2,3]. As particulate pollution they also adversely impact human health [4,5] and both agricultural and silvicultural productivity through the formation of haze and by acid deposition [6,7,8]. Due to their importance, representation of aerosol processes has been a major component of air quality (AQ) and climate modeling. Aerosol submodels in climate models are typically of reduced complexity compared to AQ models, but there is a trend towards more comprehensive treatment of aerosols in order to better capture direct and indirect aerosol effects (e.g., [9]). Although important, due to computational expense, coupling between prognostic aerosol submodel and radiative transfer and cloud microphysics in weather prediction models has been very limited until recently [10].

Gas-aerosol mass exchange is a kinetic process that involves the diffusivity of constituents in the gas phase, at the gas-particle interface and in the condensed phase, thermodynamics (reversible phase changes and hydration), and chemical reactions on aerosol surfaces and in the condensed phase [11]. Ambient weather conditions and physicochemical properties of chemical constituents in the gas and aerosol phases determine the driving force for gas-aerosol partitioning. Thus, accurate simulation of gas-aerosol mass transfer in models requires proper representation of aerosol physicochemical characteristics. This involves at least three degrees of freedom: size, mixing state (representing the chemical composition) and morphology (representing the physical distribution of the chemical components), in line with the recent review by Riemer et al. [12]. These properties are important for the proper representation of the growth of aerosols as well as their radiative and hygroscopic properties. Secondary aerosols which are generated via nucleation in the atmosphere are important since they account for about 50% of cloud condensation nuclei (CCN) through activation, with the remaining fraction arising from primary aerosol emissions [13]. Smaller particles are more effective at scattering electromagnetic (EM) radiation and it is accumulation mode aerosol particles that act as CCN and impact cloud droplet number concentration and albedo [14]. In the case of rural secondary aerosol which is dominated by the organic fraction, and likely for other aerosol types, the details of the growth process has a direct impact on CCN concentrations [15]. In the case of cloud ice condensation nuclei (IN) it is primary emissions of insoluble aerosols that are important but aged secondary soluble aerosols play a role as well [16].

Temperature and relative humidity (RH) conditions affect gas-particle mass transfer rates for different chemical components such as nitrate and oxidized non-methane volatile organic compound (VOC)s (e.g., [17,18]). For low temperatures and/or RH conditions, the aged organic component has amorphous solid characteristics with a strong impact on organic particle viscosity which can substantially affect gas-aerosol mass exchange [19,20,21,22,23].

Particle composition (or mixing state) affects the exchange of constituents, including water, between the gas and aerosol phase. Observational studies have demonstrated that aerosols have complex mixing states [24,25,26]. Li et al. [27] and Dong et al. [28] find that submicron aerosols in the Tibetan Plateau region are speciated by different mixtures of primary aerosol (BC and dust), secondary inorganic aerosol (SIA) and SOA which occur in distinct phases with SOA often forming a coating on an SIA core. It is possible for larger SIA and SOA particles (less affected by coagulation) to coexist in a given air mass due to mixing of air parcels with different emissions and aerosol growth histories. Long range transport of aerosol plumes from biomass burning was not found to result in pure internally mixed state but rather to maintain a complex external mixture by Moroni et al. [29]. There is indication of competitive growth behavior between SIA and SOA (e.g., [30]). The aerosol mixing state cannot be simply characterized as internally mixed even if a significant amount of internal mixing occurs downstream from emissions sources (e.g., [31]). Healy et al. [32] used an aerosol time-of-flight mass spectrometer (ATOFMS) and find a range of primary and secondary organic aerosol mixing states of carbonaceous particles in Paris. Measurements in the Arctic indicate a significantly externally aerosol mixed state since important sources are typically remote (e.g., [33,34]).

Distinct species of aerosols can combine to form composite particles consisting of distinct particle domains and exhibit complex non-spherical structure depending on the degree of hydration. Particle morphology also reflects ambient meteorological conditions and composition associated processes such as LLPS and the physical state of constituents (whether they are in a solid, semi-solid phase or liquid solution phase) (e.g., [35,36]). For mixed organic-inorganic atmospheric aerosols LLPS involves the separation into an electrolyte-rich subdomain and an organics-rich subdomain (e.g., [37]). Song et al. [38] investigated the morphological behavior of liquid mixed organic-inorganic aerosols in a laboratory and find a range of behavior depending on the organic species involved. Phase-separated particles can have core-shell and partially engulfed (lens) morphologies with distinct growth kinetics as well as optical properties (e.g., [39]). However, the occurrence of complex quasi-solid phases and LLPS is not universal and they can be in a liquid state under common atmospheric conditions, for example, in the southeastern USA, in the Amazon and elsewhere due to ambient RH levels [40,41]. Sulfate, nitrate, organics and water can accumulate on black or elemental carbon (BC) and substantially change BC radiative properties (e.g., [42]). So inclusions of BC in sulfate aerosols impact the thermal behavior, resulting in more rapid dehydration [43].

Over the last two decades there has been significant improvement in laboratory and field data on aerosol properties and chemistry ([44], and references therein). It has become increasingly evident that organic chemical species play a central role in secondary aerosol formation and evolution accounting from 18% to 90% by mass of accumulation mode aerosols (sizes typically less than 1 micron) [45,46]. The implication of field studies is that pure inorganic secondary aerosols are not typical in the troposphere [47]. In non-urban regions organic compounds are key to nucleation given low background levels of H${}_2$SO${}_4$ and for facilitating initial growth from the molecular cluster scale to the continuum droplet scale (e.g., [48,49]). However, even large aerosol particles over remote forested regions have substantial organic fractions [50]. SOA also forms a substantial fraction of secondary aerosol in polluted urban environments [51]. Organic compounds affect water uptake and the deliquescence and efflorescence of inorganic compounds [52,53]. Organic compounds can have different impacts on aerosol water-uptake depending on the type of inorganic salt involved. For example, the impact on water absorption by NaCl solutions is different from that on ammonium sulfate solutions. In addition, some organic species have little impact on deliquescence relative humidity (DRH) and crystallization relative humidity (CRH), also termed efflorescence relative humidity (ERH), while others introduce significant changes to the inorganic salt hysteresis curve [50,54,55].

Chemical interactions of gas-phase species with the surface and bulk interior of liquid droplets and primarily the surface of solid particles of aerosols or clouds have been termed “multiphase chemistry” and “heterogeneous chemistry”, respectively, by Ravishankara [56]. However, in the general atmospheric chemistry literature the two terms are not distinguished and “heterogeneous chemistry” is used to refer to both liquid and solid particle chemical processes both on the surface and in the bulk [57], we adopt this terminology in this review. Heterogeneous chemistry influences new particle formation and subsequent evolution in composition and mass through the formation of new compounds with different volatility characteristics and different affinity for water [40,58,59]. A large body of research demonstrates that heterogeneous chemistry in aerosols needs to be taken into account as it affects the growth of inorganic aerosol (e.g., [3,60,61,62,63,64]) and SOA [65,66,67] in urban and remote regions. In addition, organic chemistry both in the gas and aerosol phases can produce hydrophobic and hydrophilic compounds and result in heterogeneous aerosol composition (multiphase morphologies) [37,68,69,70,71].

The growth of particles in the nucleation, Aitken and small diameter tail of the accumulation mode is sensitive to the Kelvin (curvature) effect [11] and physicochemical properties of particles such as surface tension (e.g., [72]). Owing to the Kelvin effect, nanoparticles have high saturation vapor pressures at the gas-particle interface and only very low volatility compounds can partition into the condensed phase for freshly nucleated particles. Small particles are generally supersaturated at the gas-particle interface making physicochemical measurements difficult (e.g., experimental probes can trigger phase transitions) and only recently has progress been made in filling the data gap for particles in the Aitken and nucleation modes [73]. Another feature of nanoparticles is that they have much higher number densities compared to larger aerosol particles and are more surface active. Combined with the large surface area density of the nanoparticle population, heterogeneous surface chemistry is particularly important in this size range.

Ansari and Pandis [74] and Zhang et al. [75] reviewed and compared inorganic thermodynamics models for 3-D large-scale models available at the time. To the best of our knowledge no reviews have been published on this topic since then. Here we present highlights of both organic and inorganic thermodynamic model development after the year 2000. The purpose of this review is to assess the current state of aerosol process modeling which pertains to gas-aerosol partitioning and consider which aerosol processes can be reasonably improved or included in 3-D models, which are not yet commonly implemented based on current progress. A second objective is to highlight the gap between current understanding of gas-aerosol partitioning from measurements, comprehensive box models and parameterizations for large-scale modeling. Throughout the review, we extensively refer to up-to-date observational studies and discuss the limitations of presently available models of related processes. Since the focus of the review is gas-aerosol partitioning for air quality and climate models, detailed presentation is restricted to large-scale models and more detailed process-based box models are mentioned for the sake of discussion. We also touch on the issue of computational performance of process models or parameterizations.

This review is organized as follows. Section 2 presents the theoretical basis and modeling status of inorganic and organic thermodynamic models including parameterizations for water uptake and LLPS processes in mixed organic-inorganic aerosol and a discussion on the limitations and computational performance of these models. Section 3 describes modeling of gas-aerosol mass transfer kinetics and introduces available condensation-evaporation models for inorganic aerosol and condensed phase models for organic aerosol. The physicochemical properties of aerosols such as pure liquid vapor pressure, surface tension, particle diffusivity and mass accommodation coefficient are important sources of uncertainties in thermodynamics models. Current approaches for these parameters are also presented in Section 2 and Section 3. In addition, nanosize effects on aerosol thermodynamics and morphology are addressed in Section 2. In Section 4, the role of organic and inorganic heterogeneous chemistry on aerosol growth and its current implementations in models are briefly presented. For detailed reviews of inorganic and organic heterogeneous chemistry, the reader is referred to George et al. [43] and Shrivastava et al. [67], respectively. Section 5 concludes the review with a summary of the status, recommendations based on existing capacity and research gaps and priorities for future development for the gas-aerosol partitioning modeling.

2. Thermodynamic Equilibrium

Particle composition affects the exchange of constituents, including water, between the gas and aerosol phase. Based on the concentrations of various chemical constituents in the gas and aerosol phases, there is a repartitioning of these constituents between the two phases. Gas-aerosol mass exchange or partitioning is a kinetic process that involves the diffusivity of constituents in the gas and intra-aerosol or condensed phases, as well as the influence of either reversible or irreversible chemical reactions in the condensed phase. Aerosol composition also determines if aerosol particles undergo phase separation and the physical state of constituents (whether they are in a solid, semi-solid phase or liquid solution phase). These composition-associated process are treated within the core of thermodynamics model framework which simulates the effect of chemical mixtures on the saturation vapor pressures of aerosol constituents at the gas-aerosol interface and their phase state in the particle [11,37,76]. It should be noted that phase separation is not exclusive to liquid aerosols and desiccated aerosol samples can exhibit substantial separation of the constituents (e.g., [77]).

The traditional treatment of aerosols in 3-D AQ and climate models has been guided by considerations of the overall mass distribution which is weighted to particles with diameters in excess of 1 micron. The fact that larger secondary aerosols are dominated by inorganic components (e.g., [78]) combined with the lack of observational data on the role of organic compounds in secondary aerosol formation and growth has led to the development of pure inorganic secondary aerosol submodels which are the core of the schemes in current 3-D AQ and climate models (e.g., [9,79]). If SOA submodels are included, then they are treated as an independent process from inorganic aerosols and simplified schemes, that assume dilute solutions rather than complex organic mixtures that undergo phase separation and form amorphous solids at low RH, are employed to model them (e.g., [80,81]). A substantial research base has accumulated over the last two decades which makes this approximation no longer tenable, especially for accumulation mode and smaller size aerosol particles. For example, organic compounds reduce the uptake of ammonia by sulfate aerosols and in the process increase aerosol acidity and deviation from equilibrium conditions [82].

The implementation of sophisticated SOA schemes based on thermodynamics models for mixed inorganic-organic systems is limited due the significant uncertainties in the physics and chemistry and the large computational expense for such models ([83,84], and references therein). There are hundreds of organic compounds that participate in SOA formation and they are subject to significant ageing inside aerosols (via oxidation and oligomerization) and in the gas phase (via oxidation) [85]. Ageing changes the physical properties of organic species including their volatility. Some fraction of oxidized products becomes less volatile, but part can become more volatile and escape back to the gas phase or negligibly partition into the condensed phase. The chemistry of most of these organic species is uncertain and not all of the reaction rates and physical properties of the chemical species are known [83]. In addition, the activity and uptake coefficients of most organic species associated with SOA formation are not known (e.g., [84]). However, understanding of organic chemistry is sufficiently advanced that comprehensive offline mixed SIA and SOA thermodynamics models have been developed that use generalized properties of organic molecules based on their functional group structure [86,87,88,89].

2.1. Theoretical Basis

Aerosol thermodynamics models are employed to determine the physical state of aerosols as controlled by the amount of water and the chemical composition of the solutes [76]. Thermodynamics principles are used to describe chemical mixtures and hence the name for these models [11]. The uptake/loss of water by aerosols is affected by the presence of inorganic and organic chemical species via the Raoult effect. The Raoult effect or the solute effect involves the reduction of both solute and solvent saturation vapor pressure at the surface of liquid aerosol particles. The Kelvin effect becomes more important as aerosol particle diameter decreases due to the increasing curvature of the aerosol surface which impacts the saturation vapor pressure of aerosol constituents at the surface. In particular, for a given gas phase vapor pressure of a constituent, smaller aerosols with a larger Kelvin effect undergo more evaporation and less condensation compared to larger aerosols. Via the Raoult effect, solutes act to reduce the evaporation of the liquid aerosol particle. Together the opposing Kelvin and the Raoult effects give rise to the Köhler equation which governs the activation of aerosols and their equilibrium size for a given RH. A major component of thermodynamics models is the model of activity coefficients which are necessary to account for deviations from ideal behavior of chemical mixtures [76,90,91].

Here we present a brief summary of the theory on which thermodynamics models are built. A more detailed treatment can be found in Seinfeld and Pandis [11] and Jacobson [92]. Gas phase species that dissolve into aerosols constitute solutes in a solvent which can be water, sulfuric acid or a liquid mixture of organics. Inorganic solutes tend to dissociate into ions in aqueous solutions and most of the inorganic constituents of aerosols are strong electrolytes. Most organic acids and bases are weak electrolytes and undergo only partial dissolution. Dissolution of inorganic and organic compounds can occur in polar organic solvents, but to a lesser extent compared to water. Dissolution of solids occurs in solvents as well and typically inorganic acids and NH${}_3$ from the gas phase form compounds inside aqueous aerosols that are solids under atmospheric conditions (e.g., ammonium nitrate, ammonium sulfate).

These soluble solid products, or salts, which can also be composed of purely organic or a combination of inorganic and organic precursors, can crystallize (effloresce) or dissolve (deliquesce) depending on the solvent content of the aerosol particle. In the case of aqueous aerosols the liquid water content is a function of RH. Efflorescence and deliquescence typically do not occur at the same RH and exhibit hysteresis as RH is varied. A dry particle consisting of several soluble constituents that is subjected to increasing RH starts to absorb water at the mutual deliquescence relative humidity (MDRH) and soluble solids are fully dissolved at the separation relative humidity (SRH). A fully deliquesced particle subjected to decreasing RH remains in a metastable supersaturated liquid state until spontaneous efflorescence occurs at RH below MDRH, the CRH. The values of RH for which these phase transitions occur are influenced by the presence of other solutes. Prediction of CRH is challenging and there is currently a lack of any comprehensive parameterization [76]. Under ideal conditions crystallization can be treated as a homogeneous nucleation process in a classical nucleation theory (CNT) framework (e.g., [93]) but can occur via other pathways such as heterogeneous nucleation (e.g., [94]) and contact efflorescence (e.g., [95]). In the case of the latter two pathways, efflorescence can occur at much higher RH. Aqueous mixed organic and inorganic solutions formed at high RH can undergo LLPS as RH decreases into an aqueous electrolyte domain and organics mixture domain. This process also reflects the solute composition and pH can have a substantial effect by increasing the solubility of the organic component [96].

Liquid-solid transition hysteresis and LLPS can occur in organic solvent systems as well [97,98]. However, under atmospheric conditions there is no dominant organic compound that acts as a solvent substitute for water. Instead, even pure organic aerosols from primary emissions undergo chemical and physical ageing that involves incorporation of water and inorganic constituents. Similarly, SOA that may have nucleated as pure organic particles undergoes ageing and uptake of water and incorporation of inorganic components such as sulfate and nitrate via coagulation and condensation. The presence of hydrophobic compounds and compounds with mixed hydrophobic-hydrophillic properties results in complex mixtures with a heterogeneous, phase-separated distribution of constituents.

Both the chemical and physical processes described above, except for the metastable branch of the hydration hysteresis curve, can be represented by reversible reactions and the system tends towards an equilibrium state under a given set of conditions (pressure, temperature and RH):

$${\nu }_{M}M+{\nu }_{N}N+\cdots \rightleftharpoons {\nu }_{A}A+{\nu }_{B}B+\cdots$$

(1)

where ${\nu }_X$ are stoichiometric coefficients, M and N are reactants and A and B are products. The amount of reactants and products is determined by the equilibrium coefficient which is given by

$$K_{eq}\left(T\right)=\frac{{\left\{A\right\}}^{{\nu }_A}{\left\{B\right\}}^{{\nu }_B}\cdots }{{\left\{M\right\}}^{{\nu }_M}{\left\{N\right\}}^{{\nu }_N}\cdots }=\prod _i{\left\{a_i\right\}}^{k_i{\nu }_i}$$

(2)

where $\left\{X\right\}$ represents the thermodynamic activity of constituent species X, $\left\{a_i\right\}$ is the thermodynamic activity of species i, $k_i=+1$ for products, and $k_i=-1$ for reactants. This relation is derived by minimizing the Gibbs free energy of the system at fixed pressure and temperature, which takes the form

$$\mathrm{d}G=V\mathrm{d}p-S\mathrm{d}T+\sum _i{k}_i{\mu }_i\mathrm{d}{n}_i$$

(3)

where $\mathrm{d}{n}_i$ is the number of moles of species i, ${\mu }_i$ is the chemical potential of the species, $k_i$ are as defined previously, S is the entropy and V is the volume of the system. The chemical potential of species i can be expressed as:

$${\mu }_i={\left(\frac{\partial G}{\partial n_i}\right)}_{T,p}={\mu }_i^{\circ }\left(T\right)+RTln\left\{a_i\right\}$$

(4)

where the first term on the right represents the standard chemical potential and the second term represents the effect of the activity of the species. The activity of a chemical constituent represents deviation from ideal solution behavior and is defined as one for an infinitely dilute (ideal) solution. It can be less than or greater than one for concentrated solutions. Activity coefficients are functions of temperature.

Under equilibrium conditions $\mathrm{d}G=0$ so that for fixed temperature and pressure

$$\sum _i{k}_i{\mu }_i\mathrm{d}{n}_i=0$$

(5)

The equilibrium reaction Equation (1) imposes a constraint between the reactants and products such that a change $\mathrm{d}n$ in the molar amount of any species determines the changes in all the other species:

$$\frac{\mathrm{d}{n}_M}{-{\nu }_M}=\frac{\mathrm{d}{n}_N}{-{\nu }_N}=\cdots =\frac{\mathrm{d}{n}_A}{{\nu }_A}=\frac{\mathrm{d}{n}_B}{{\nu }_B}=\cdots$$

(6)

Thus, dividing Equation (5) by the smallest of the $\mathrm{d}{n}_i$ gives

$$\sum _i{k}_i{\nu }_i{\mu }_i=0$$

(7)

For a fixed temperature $T_0$ substituting for ${\mu }_i$ we have

$$\begin{array}{c}\hfill \sum _i{k}_i{\nu }_i{\mu }_i^{\circ }\left(T_0\right)+R{T}_0\sum _i{k}_i{\nu }_{i}ln\left\{a_i\right\}=\sum _i{k}_i{\nu }_i{\mathsf{\Delta }}_f{G}_i^{\circ }+R{T}_{0}ln\prod _i{\left\{a_i\right\}}^{k_i{\nu }_i}=0\end{array}$$

(8)

where ${\mathsf{\Delta }}_f{G}_i^{\circ }$ is the molal Gibbs free energy of formation of compound i. Thus

$$\begin{array}{c}\hfill K_{eq}\left(T_0\right)=\prod _i{\left\{a_i\right\}}^{k_i{\nu }_i}=exp\left(-\frac{1}{R{T}_0}\sum _i{k}_i{\nu }_i{\mathsf{\Delta }}_f{G}_i^{\circ }\right)\end{array}$$

(9)

In the derivations above, only bulk solutions are considered. For small particles, this approximation breaks down as the contribution of surface terms to the total Gibbs free energy of the system become more important. In particular, terms for the surface tension between various phases must be included. Surfaces can also interact since they are effectively thin in the small particle limit. A water layer on a solid core (e.g., dissolving salt particle) exhibits a disjoining pressure which will contribute an additional non-negligible term to the Gibbs free energy (e.g., [99]). These size dependent effects are not included in existing thermodynamics models which reduces their accuracy for particles in the nucleation and Aitken modes [73]. There is more discussion of size effects in Section 2.8.

Activity coefficients can be derived from laboratory analysis or theoretical models and reflect both long range and short range interactions between solutes. Since single ions cannot be isolated from solutions laboratory studies typically measure mean binary activity coefficients. A commonly used parameterization for mean binary activity coefficients in aqueous solutions is the method of Pitzer [100,101]. However, it is only applicable for the high RH range when solute molalities are low. The method of Pitzer–Simonson–Clegg (PSC) [102,103] applies to higher solute molalities and thus a broader range of RH values. For mixtures of many electrolytes in an aqueous solution, such as found in aqueous aerosol with dissolved sulfuric acid, nitric acid and ammonia, a more general activity coefficient model is required. There are empirical mixing rules which use mean binary activity coefficient expressions. The methods of Bromley [104] and Kusik and Meissner [105] have been used in aerosol thermodynamics models. Some models include their own mixing rule formulations such as the Equilibrium Simplified Aerosol Model (EQSAM) [106] and the Multicomponent Equilibrium Solver for Aerosols (MESA) [107].

For organic mixtures activity coefficients are usually estimated based on the functional groups from which the organic compounds are composed. The immense diversity of organic compounds and lack of data on the physical properties of every organic species necessitates such an abstracted approach. A widely used model is UNIQUAC Functional-group Activity Coefficients (UNIFAC) derived from the Universal Quasichemical (UNIQUAC) model [108], which is able to characterize mixtures of most organics relevant for the atmosphere [109,110]. LIFAC [111] is an extension of UNIFAC to handle mixtures of organics and inorganics via an approach similar to that of Pitzer [100]. However, this model was not developed for atmospheric science applications but for chemical engineering purposes. Zuend et al. [89] introduced the Aerosol Inorganic-Organic Mixtures Functional Groups Activity Coefficients (AIOMFAC) model, a modified version of LIFAC, that is applicable to the atmosphere. In particular, it can handle high ionic strengths associated with low RH conditions and extends the list of electrolyte species to include sulfuric acid and ammonium bisulfate. This model is able to calculate vapor-liquid, liquid-liquid and solid-liquid equilibria. In contrast to the inorganic activity coefficient parameterizations, these mixed organic-inorganic models are much more computationally expensive and cannot typically be used in 3-D atmospheric models without substantial simplification of the organic chemistry representation into a small number of organic surrogates. There is a trade-off for using such complex models online. For example, the Extended AIM (E-AIM) [112] employs a chemical surrogate approach that attempts to preserve accuracy and is not used online to the best of our knowledge.

The water content of aqueous solutions reflects the hydration property of solutes such as sulfuric acid. Aqueous aerosol droplets can exist at RH below 100% without evaporating due to the reduction of saturation vapor pressure at the droplet surface (Raoult’s law or solute effect). A commonly used method of estimating the aerosol liquid water content is the Zdanovskii–Stokes–Robinson (ZSR) method [113]. ZSR is not restricted to inorganic electrolytes. However, ZSR does not consider organic-inorganic interactions in mixed systems. Zuend et al. [89] compared ZSR to the more comprehensive treatment in AIOMFAC and found a reasonable agreement.

2.2. Inorganic Thermodynamics Models

Inorganic thermodynamics models deal with a relatively small number of species which have the most atmospheric significance aside from organic species. These include H${}_2$SO${}_4$, HNO${}_3$ and NH${}_3$ from gas phase inputs, NaCl from sea salt, and CaCO${}_3$, K, and Mg from lithogenic (dust) sources. There are many other species that can dissolve into aerosols but they are not included for reasons of numerical tractability and due to the fact that they are generally secondary in terms of concentration. This is true even for many species with significant gas phase concentrations, e.g., O${}_3$, since the Henry’s law constants for these species are small. However, this simplification breaks down if there are fast reactions on dissolution. For O${}_3$ there are such reactions with organic molecules such as chlorophyll (e.g., [114]). Reactive uptake can significantly increase the dissolution of gases into the aerosol phase (e.g., [11]).

Sulfuric acid partitioning is typically treated as a pure condensation process due to its very low volatility [115]. To save computational expense the explicit kinetic partitioning can be replaced with an equilibrium assumption [116,117] but this can degrade the results. In particular, the equilibrium approach lacks the Kelvin effect, which is affects aerosol growth in the nucleation mode, Aitken mode and the small-radius tail of the accumulation mode (e.g., [118]). The non-equilibrium approach has been implemented other models (e.g., [76]).

Many thermodynamics models for inorganics were developed from the 1980s onward (e.g., [75]). Some of the most commonly used or advanced variants are considered here. From the 1990s period we have EQUISOLV II [91] and ISORROPIA [116,119] (updated to version II by Fountoukis and Nenes [117]). From the 2000s we have UHAERO [120] and MESA [107,121], which attempt to improve on their predecessors in terms of realism and numerical efficiency.

A primary function of thermodynamics models is to enable accurate gas-aerosol mass transfer calculations for certain chemical species and for water. These mass transfer calculations are a distinct model component that is typically packaged with the thermodynamics model. In this section we describe the models listed in

Table 1. Inorganic Models and Components.

Thermodynamics Model	Acitivity Coefficient Scheme	Condensed Phase Equilibrium Solver	Gas-Aerosol Partitioning Scheme	Host Model Availability
EQUISOLV II [90,91]	Mean mixed activity coefficients based on empirical mixing rules of Bromley [104]. Necessary mean binary activity coefficients from various sources Parker [122], Hamer and Wu [123]. Temperature dependent sulfate and bisulfate activity coefficients obtained using the method of Stelson et al. [124] to the mixed coefficients from Clegg and Brimblecombe [125]. Zdanovskii-Stokes-Robinson (ZSR) used for liquid water content.	Analytical equilibrium iteration. Modification of the mass flux iteration scheme where individual equilibrium reactions are solved or approximated analytically. Approximate treatment of deliquescence relative humidity of mixtures (MDRH). NH${}_4$-Na-Ca-Mg-K-NO${}_3$-SO${}_4$-Cl-CO${}_3$	Kinetic approach with solution of mass transfer equations. Kelvin effect is included. Noniterative, unconditionally stable solution schemes for both dissolutional and condensational growth [126].	GATOR-GCMM [127]
				[email protected]
ISORROPIA II [116,117]	Mean mixed activity coefficients based on Bromley [104]. Binary activity coefficients obtained using Kusik and Meissner [105] following evaluation of Kim et al. [128]. ZSR for liquid water content.	The full multicomponent system is approximated by several subdomains identified on the basis of composition. Only a subset of component mixtures and reactions is solved. Analytical solutions are used where possible and iteration with the bisection method is used for the rest. MDRH obtained via fast custom scheme. MDRH temperature dependence determined by a formula similar to Wexler and Seinfeld [129] for deliquescence relative humidity (DRH). NH${}_4$-Na-Ca-Mg-K-NO${}_3$-SO${}_4$-Cl	Components are assumed to be under instantaneous equilibrium based on the observation that transient partitioning effects occur mostly for coarser aerosol sizes. No Kelvin effect included.	WRF-Chem [130], GEOS-Chem [131]
				https://www.epfl.ch/labs/lapi/isorropia/
EQSAM [132]	Uses an empirical binary-type relation which takes into account the Kelvin effect.	Structured similar to ISORROPIA. Simplification into composition subdomains and approximate MDRH approach. NH${}_4$-Na-Ca-Mg-K-NO${}_3$-SO${}_4$-Cl	Instantaneous equilibrium without Kelvin effect.	Can be used online.
				[email protected]
MOSAIC [76]	Mixing rules based on the Pitzer, Simonson and Clegg (PSC) model optimized using the Taylor’s series approach of Wagner [133]. Superior to other models since there is no limitation on saturation and valid over the full RH range. ZSR for liquid water content using PSC parameters.	Uses pseudo-transient continuation (PTC) methods to solve the nonlinear algebraic reaction equations of the full multicomponent system. General solver that can model equilibrium salt formation characteristics which enables an accurate parameterization of MDRH. NH${}_4$-Na-Ca-NO${}_3$-SO${}_4$-HSO${}_4$-Cl	Kinetic mass transfer approach using the adaptive step time-split Euler method. Includes the Kelvin effect.	WRF-Chem
				[email protected]
UHAERO [120]	Mixing rules based on choice of PSC or the extended UNIQUAC model [134]. Water content is based on Clegg et al. [135,136].	Solves full multicomponent system via Gibbs free energy minimization using a primal-dual active-set algorithm for speed [137,138]. NH${}_4$-Na-NO${}_3$-SO${}_4$-Cl	Instantaneous equilibrium without Kelvin effect.	CMAQ (experimental)
				Co-authors
PD-FiTE [139]	Activity coefficients taken from ADDEM [86,140] and represented by a partial derivative fitted Taylor expansion. A more general approach than MOSAIC since the fitting is done to the whole composition space. Binary activity coefficients determined using PSC.	Not included. NH${}_4$-Na-NO${}_3$-SO${}_4$-HSO${}_4$-Cl	Not included.	Experimental
				[email protected]

.

2.2.1. EQUISOLV

This model was introduced by Jacobson et al. [90] and updated by Jacobson [91] for computational speed improvement. It has two modes of operation: (1) solution of equilibrium equations between gas phase and multiple size bins of aerosol with diffusion limited mass transfer ignored, (2) solution of internal aerosol equilibrium to provide saturation vapor pressure terms for diffusion limited mass transfer equations between the gas phase and multiple aerosol size bins.

Mean mixed activity coefficients for solutes are calculated using the empirical mixing rule of Bromley [104]. Mean binary activity coefficients are required by this method. The liquid water content is determined with the ZSR method (e.g., [141]).

Solids can form in two ways: (1) RH increasing, solid phase is permitted but not required if the RH is below the DRH and no solid phase is allowed if RH exceeds the DRH even for multicomponent mixtures. (2) RH decreasing, below the DRH solutions are assumed to stay supersaturated until CRH. CRH is temperature dependent.

The DRH treatment in EQUISOLV is not fully accurate as in general the DRH of a mixture is less than the minimum of all the individual DRH values for each component. The MDRH is not unique and depends on the mixture.

This model uses two analytical equilibrium iteration (AEI) methods in lieu of mass flux iteration (MFI) to solve the system of equilibrium equations. The first AEI method pertains to reactions of the type D ⇔ A, D ⇔ A + B, D + E ⇔ A + B, the second AEI method handles reactions of the type D ⇔ 2A + B and D ⇔ A + 2B. The single reaction of the type D ⇔ 3A + B + C is solved by the MFI method.

EQUISOLV includes the Kelvin effect in the gas-particle transfer equations for every species and bin size (see Section 3 in [90]). The solution approach does not make use of extensive approximations but there are multiple nested iterations which make this scheme numerically expensive compared to ISORROPIA and MESA.

2.2.2. ISORROPIA

This model was introduced by Nenes et al. [116] and updated in Fountoukis and Nenes [117] for more species (e.g., Ca, K, Mg and CO${}_3$). ISORROPIA lacks the Kelvin effect based on the assumption of equilibrium motivated by the consideration that most atmospheric aerosol mass is in particles with diameters greater than 100 nm and that smaller particles equilibrate rapidly. With the neglect the of the Kelvin effect the water activity is equal to the RH.

ZSR is used to calculate the water content of aerosols. Multicomponent activity coefficients are calculated following Bromley [104]. The necessary binary activity coefficients are calculated following Kusik and Meissner [105] based on the evaluation of Kim et al. [128,142]. ISORROPIA uses a custom and fast approach to calculate the MDRH instead of the full and expensive solution approach of Potukuchi and Wexler [143,144] but this entails approximations. The DRH temperature dependence is based on Wexler and Seinfeld [129]. The MDRH temperature dependence is calculated with a similar formula. It was not possible to find MDRH for all mixtures and these cases were approximated by other mixtures.

ISORROPIA assumes plenty of available water so the ambient RH is not affected by the deliquescence or efflorescence of aerosol particles. This assumption may not be valid for low RH conditions. EQUISOLV does not make such an assumption but does impose a (small) minimum aerosol water to simplify the treatment. The Adaptive Step Time-Split Euler Method (ASTEM), used in the Model for Simulating Aerosol Interactions and Chemistry (MOSAIC), makes no assumptions about ambient water and can solve for mass transfer for solid and mixed-phase aerosol [76].

The solution of the system involves finding as many analytical solutions as possible which reduces the total number of equations requiring a numerical solution. Most cases in ISORROPIA require only one iteration with the bisection method. The binary activity coefficients are not calculated at run time but are precalculated.

ISORROPIA reduces the computational burden by decomposing the multicomponent system into several subdomains based on relative composition and restricting the number of components and associated reactions in each subdomain. By contrast, EQUISOLV takes a more general and rigorous approach which does not make such simplifying assumptions. But EQUISOLV is more numerically expensive as it solves one equilibrium equation at a time (analytically or iteratively) and all equations are iterated many times to get convergence. There can be up to three levels of nested iteration loops.

The neglect of the Kelvin curvature effect imposes a limitation which can reduce the applicability of this model, since much of the aerosol mass is secondary and originates from nucleated particles growing to be larger than 100 nm. The physicochemical evolution of aerosols can be affected by this approximation since the Kelvin effect restricts the condensation of gas phase constituents based on their volatility. Debry and Sportisse [145] evaluate the hybrid methods employed by ISORROPIA and other models to handle the stiff nature of the dynamical system they are solving (e.g., the bulk equilibrium method of Pandis et al. [146]) and find that the neglect of the Kelvin effect is consistent with the approximation of separation of timescales where fast components equilibrate instantaneously. However, for size-resolved equilibrium approaches (e.g., [90]) this neglect of the Kelvin effect is not an obviously valid assumption. Debry and Sportisse [145] find that the quasi-Newton Broyden–Fletcher–Goldfarb–Shanno (BFGS) minimization algorithm [147] is best suited to solving the size resolved equilibrium system.

2.2.3. EQSAM

Metzger et al. [106] introduced a computationally efficient gas-aerosol partitioning model based on a simplified representation of the activity coefficients. EQSAM produces results similar to ISORROPIA partly due to the fact that it uses a similar concentration subdomain and MDRH approach. EQSAM avoids numerically expensive iteration for the activity coefficients by using an empirical expression:

$$\gamma =[\frac{{\mathrm{RH}}^N}{(\frac{1000}{N}(1-\mathrm{RH})+N)}{]}^{\frac{1}{\xi }}$$

(10)

where N represents the number of water molecules that participate in dissolution of a salt solute in a binary solution, and $\xi$ represents the charge and stoichiometric coefficient of the dissolved electrolyte. These are concentration domain dependent parameters chosen to fit the results of ISORROPIA. Mixed solutions satisfy the same relation with an appropriate choice of N and $\xi$. A limitation of this relation is that it is strictly valid only for particles larger than 100 nm since the Kelvin effect is ignored by assuming the water activity equals RH. This is the same approximation that is made in ISORROPIA and MESA.

The analytical approach was extended further in the EQSAM Version 3 (EQSAM3) [148]. Non-ideal solution properties such as activity coefficients, aerosol water, and DRH and CRH of binary and multi-component mixtures are only analytical functions of the solubility of aerosol constituents at a given RH and temperature. The introduction of a stoichiometric coefficient for water (see Equation (19) in [148]) proved to be controversial [149]. But Metzger et al. [132] effectively address the criticism with a general solubility based approach for calculating the water activity which includes the Kelvin effect. In this framework the water activity takes the form:

$$\begin{array}{ccc}\hfill a_{\mathrm{w}}& =& (A+{\mu }_s^o·M_{\mathrm{w}}·{\nu }_i[\frac{{\mu }_s}{{\mu }_s^o}+B{]}^{{\nu }_i}{)}^{-1}\hfill \end{array}$$

(11)

$$\begin{array}{ccc}\hfill A& =& (1+{\nu }_i·{\mu }_s·M_{\mathrm{w}})·exp\left(-M_{\mathrm{w}}·{\mu }_s^o·{\nu }_i·{\left(\frac{{\mu }_s}{{\mu }_s^o}\right)}^{{\nu }_i}\right)\hfill \end{array}$$

(12)

$$\begin{array}{ccc}\hfill B& =& \left(1+\frac{1}{{\nu }_i·{\mu }_s·M_{\mathrm{w}}}\right)·{\left({\nu }_i·\frac{{\mu }_s}{{\mu }_s^o}\right)}^{-\frac{1}{{\nu }_i}}\hfill \end{array}$$

(13)

where ${\mu }_s$ (mol(solute) kg${}^{-1}$ (H${}_2$O)) is the solute molality, ${\mu }_s^o=1$ (mol kg${}^{-1}$) is a reference/scaling concentration, $M_{\mathrm{w}}$ (kg mol${}^{-1}$) is the molar mass of water, and ${\nu }_i$ is a solute specific constant that accounts for non ideal behavior. The functional form of $a_{\mathrm{w}}$ is motivated by empirical expressions derived in previous research as summarized by Rose et al. [150].

The single solute constant ${\nu }_i$ can be solved by using laboratory data for DRH and the saturation molality, ${\mu }_s^{\mathrm{sat}}$, which is related to the mass fraction solubility, $w_s$, and the fact that

$$\begin{array}{ccc}\hfill a_{\mathrm{w}}& =& \mathrm{RH}exp\left(-\frac{4M_{\mathrm{w}}{\sigma }_{\mathrm{sol}}}{RT{\rho }_{\mathrm{w}}{D}_{\mathrm{wet}}}\right)\hfill \end{array}$$

(14)

where the second term on the right hand side is the Kelvin term which depends on ${\sigma }_{\mathrm{sol}}$ (J m${}^{-2}$), the surface tension of the solution, ${\rho }_{\mathrm{w}}$, the density of water, $D_{\mathrm{wet}}$ the diameter of the aerosol solution droplet, R (J mol${}^{-1}$ K${}^{-1}$), the ideal gas constant, and the temperature, T (K). This system is solved for ${\nu }_i$ by assuming that ${\sigma }_{\mathrm{sol}}$ and ${\rho }_{\mathrm{w}}$ are pure water constants (Rose et al. [150]) and using a root-finding method such as bisection.

This formulation is evaluated by Metzger et al. [151] using the EQSAM Version 4 (EQSAM4) which is compared to EQUISOLV II, ISORROPIA II and a reference model, the Aerosol Inorganic Model (AIM) [152], for various cases. ISORROPIA II was used to fit mixed solution cases based on the single solute framework introduced by Metzger et al. [132]. In terms of aerosol water content and solid–liquid separation EQSAM4 compares well to EQUISOLV II in general and in many cases outperforms ISORROPIA II. There are cases where EQSAM4 deviates from both ISORROPIA II and EQUISOLV II, but not to an extreme degree. EQSAM4 behaves similarly to ISORROPIA II when compared to AIM. This behavior reflects the ISORROPIA-like framework on which EQSAM4 is based, namely the projection onto composition subdomains and simplified treatment of MDRH. EQSAM4 performs well relative to ISORROPIA II when compared to observations (see Figure 8 in [151]) from the Mediterranean Intensive Oxidant Study (MINOS) [153] in terms of fine and coarse mode particulate matter and water. In particular it has significantly less particulate matter outlier events where a peak or trough occurs that does not agree with the observational time series data. However, ISORROPIA II is able to capture some peaks which EQSAM4 cannot. Considering the computational performance gain offered by EQSAM4, it is a worthwhile replacement for ISORROPIA II in 3-D models.

2.2.4. MOSAIC: MTEM, MESA and ASTEM

MOSAIC, introduced by Zaveri et al. [76], has two submodels to solve for the activity coefficients and condensed phase equilibrium reactions and a submodel for non-equilibrium gas-aerosol partitioning which are described below.

The Multicomponent Taylor Expansion Method (MTEM) was developed by Zaveri and Easter [121] to efficiently compute activity coefficients for mixed systems. The Bromley and Kusik and Meissner mixing rules are easy to implement and computationally inexpensive but they are reasonably accurate only for subsaturated solutions. Zaveri and Easter [121] note that for saturated and supersaturated conditions multicomponent solutions are limited by the maximum ionic strengths that mean binary activity coefficient parameterizations are able to handle and this is a problem for low RH. In addition, since the binary activity coefficients are functions of ionic strength, the mixing rules require repeated evaluations in an iteration as the ionic strengths change in the course of a numerical solution of the equilibrium equations.

The PSC model [103,154] is superior to the other mixing rules since it has no limitations on saturation and is valid over the whole RH range, but is numerically expensive. The method of Zaveri and Easter [121] has the accuracy of PSC but is much faster. It uses an extension of the Taylor’s series approach of Wagner [133] for dilute alloy systems which is applicable to aqueous electrolytes at any concentration.

PSC model parameters are not known for all combinations of species so there are gaps for some component combinations (e.g., H${}^{+}$, Ca${}^{2+}$, NO${}_3^{-}$, H${}_2$O). Expressions for binary activity coefficients are obtained by polynomial fits using PSC model output over the designated RH range (assuming activity of water equals RH).

MTEM assumes the water activity equals RH which is only true if the Kelvin effect is neglected. This is a common approximation valid in the bulk limit and as noted above, thermodynamics models do not consider small particle effects such as the surface tension terms in the Gibbs free energy of the system. However, as noted above the computation of the equilibrium water content in MESA (see below) in liquid particles is size resolved and takes the Kelvin effect into account [76].

On introduction, MESA [107] treated only the system H${}^{+}$–NH${}_4^{+}$–Na${}^{+}$–Ca${}^{2+}$–SO${}_4^{2-}$–HSO${}_4^{-}$–NO${}_3^{-}$–Cl${}^{-}$ without lithogenic (dust) aerosol contributions. Later versions found in WRF-Chem [130] include calcium carbonate from dust. The MESA solver is general enough to model equilibrium salt formation characteristics and this allows a parameterization of MDRH including its temperature dependence using cubic polynomials.

MESA deals with the same nonlinear algebraic equilibrium reaction equations as ISORROPIA and EQUISOLV I but instead of multiple nested iteration like EQUISOLV II, it uses pseudo-transient continuation methods [155] which are much faster. MOSAIC uses ASTEM to dynamically integrate the mass transfer equations [76]. H${}_2$SO${}_4$ and MSA (CH${}_3$OS${}_3$H) are irreversibly transferred to pre-existing particles using a time-dependent uptake. HNO${}_3$, HCl and NH${}_3$ are semivolatile and the mass exchange is largely reversible. A series of complex conditions is used to model the partitioning which we do not restate here. The Kelvin effect is included for all species.

2.2.5. UHAERO

Amundson et al. [120] introduced a model based on a computationally efficient minimization of Gibbs free energy and was designed to be incorporated into 3-D AQ models. The use of Gibbs free energy minimization puts the model in the superior class with the least approximations out of the models being considered and in the same class as expensive offline thermodynamics models such as AIM [112] which are considered reference models. The aim of the model is to capture both deliquescence and crystallization based purely on thermodynamics without approximations. This model does not include the Kelvin effect, but the authors state that the effect could be easily included.

This model handles the crystallization of metastable solutions with a sophisticated approach but only homogeneous crystal nucleation is considered due to the complexity of heterogeneous nucleation.

In terms of performance the model averages less than 3.5 Newton iterations per grid point for liquid-solid equilibrium and gas-aerosol equilibrium. This is for so-called warm starts where some initialization state exists as would be the case in an AQ model. This model is faster than ISORROPIA in some aspects (depending on convergence criteria imposed).

The gas phase to aqueous phase exchange is treated as a Gibbs free energy minimization problem just as the aerosol liquid phase to solid phase equilibrium. There are no peer-reviewed references that could be found comparing the model to other models. However some work was done to compare UHAERO with ISORROPIA in the CMAQ model [156] in a conference presentation by Cheng et al. [157].

2.2.6. PD-FiTE

Topping et al. [139] introduced a more sophisticated approach for activity coefficient calculation compared to MTEM. The Partial Derivative Fitted Taylor Expansion (PD-FiTE) scheme uses the comprehensive thermodynamic model ADDEM [140] to optimize the Taylor expansion fit for activity coefficients. It was found empirically that the logarithm of any binary activity coefficient is nearly linear over the entire RH range. This justifies a linear expansion of the logarithm of an activity coefficient for species A, ${\gamma }_A$, as

$$\begin{array}{ccc}\hfill \ln {\gamma }_A(x_B,x_C,x_D,\dots ,\mathrm{RH})& =& \ln {{\gamma }^0}_A\left(\mathrm{RH}\right)+\sum _{B\ne A}^N\left(\frac{\partial \ln {\gamma }_A}{\partial x_B}\right)\left(\mathrm{RH}\right)x_B\hfill \end{array}$$

(15)

$$\begin{array}{ccc}\hfill x_A& =& \frac{{\nu }_A{N}_A}{{\sum }_E{\nu }_E{N}_E}\hfill \end{array}$$

(16)

where $\ln {{\gamma }^0}_A\left(\mathrm{RH}\right)$ is the contribution at a given RH with no other solutes present. The variables $x_E$ for a set of species E are ionic mole fractions as defined in (16) where ${\nu }_E$ is the stoichiometric coefficient and $N_E$ is the number of moles of solute. Topping et al. [139] do not make additional simplifications to the above formulation like Zaveri and Easter [121] for MTEM. The fitting is done for the whole composition space of H${}^{+}$–NH${}_4^{+}$–Na${}^{+}$–SO${}_4^{2-}$–HSO${}_4^{-}$–NO${}_3^{-}$–Cl${}^{-}$.

The PD-FiTE approach does not involve consideration of sulfate-poor and sulfate-rich regimes like in MTEM since the full spectrum is covered by ADDEM and captured in the fit. The binary activity coefficient terms, $\ln {{\gamma }^0}_A\left(\mathrm{RH}\right)$, are calculated via the PSC activity model and expressed as polynomials of RH. The treatment of ammonia is more explicit and faster compared to MTEM.

PD-FiTE does not currently deal with solids based on the assumption that multicomponent aerosols in the troposphere are liquid. This includes organic aerosol. But solid formation can be included due to the general nature of the formulation. At present, no thermodynamics model based on PD-FiTE is available.

2.2.7. Performance Comparison of Inorganic Thermodynamics Models

Here we summarize the published information on the relative performance, both accuracy and computational efficiency, of the activity coefficient and equilibrium phase state solvers presented above. Evaluations are not available for every model and we did not perform our own evaluations of these schemes. As noted previously, Ansari and Pandis [74] and Zhang et al. [75] did evaluate various models available at the time, including ISORROPIA I. However, there have been no followup inter-comparisons since that time even though there have been new and improved models introduced.

Zaveri et al. [107] evaluated MESA (with the MTEM activity model) against EQUISOLV II using the comprehensive offline AIM Version III (AIM III) [152] as a reference. Test cases were restricted to mixed solid-liquid conditions for all aerosol salt components. Only single salt precipitation was evaluated and double salts and salt hydrates were not included. Two types of accuracy evaluations were conducted, for convergence criteria and thermodynamics deviations associated with uncertainties in the thermodynamics parameters used by both models. The convergence error was evaluated using the mass growth factor (MGF). The thermodynamics deviations were evaluated using MDRH and MGF. The MDRH limitations of EQUISOLV II were taken into account when evaluating MGF to produce a fair comparison. Overall, MESA is more accurate in both accuracy aspects, i.e., exhibits lower errors compared to the reference model, than EQUISOLV II. However, neither model had excessive disagreement with AIM III.

Over the entire mixed-phase RH range MESA had an average root mean square (RMS) relative convergence error of 0.14% while EQUISOLV II had an error of 3.67%. Excluding data points for RH with discrepancies in MDRH the error for EQUISOLV II falls to 0.38%. The MDRH weakness of EQUISOLV II unsurprisingly manifested itself in the thermodynamics evaluation. For one test the absolute error was 15%. The average error for MESA was 1.17% and for EQUISOLV II it was 4.42%. The RMS error in MGF was 2.31% for MESA and 4.43% for EQUISOLV II. In terms of total water content the errors were 2.46% and 3.75%, respectively.

In terms of computational performance MESA is substantially faster than EQUISOLV II. EQUISOLV II was tested in both scalar mode and vector mode. MESA lacks a vector mode. Over the entire RH range the average speed ratio of MESA relative to EQUISOLV II is 2.7 (vector) and 5.8 (scalar) for the test cases used to evaluate accuracy. Timing tests were also conducted accounting for the MDRH weakness of EQUISOLV II by restricting RH to mixed-phase ranges. In this case the average speed ratio drops to 1.4 (vector) and 2.8 (scalar) but is still a significant improvement.

The scalar single grid tests conducted by Zaveri et al. [107] have to be taken in context as the performance of the scalar mode EQUISOLV II increases with the grid size due to the optimal structuring of inner loops [75]. In addition, the convergence criteria in EQUISOLV II can be relaxed to yield performance in terms of accuracy and speed similar to ISORROPIA [158]. At the same time, the accuracy differences between thermodynamics models become obscured in the context of 3-D AQ models since other factors affect aerosol distributions and characteristics [158].

Zaveri et al. [76] evaluated MOSAIC against ISORROPIA. It should be noted that ISORROPIA II has the same computational performance as ISORROPIA used in these tests [117]. Box model simulations were conducted for 14 idealized test cases covering monodisperse particles which were purely solid, purely liquid, and mixed phase, and polydisperse particles in two RH limits. Accuracy was tested by using AIM III as a benchmark model. Three different Fortan compilers on two different computer platforms were used to determine CPU times for a single timestep and single bin (MOSAIC) and single equilibrium calculation (ISORROPIA) and the results averaged. MOSAIC was found to be comparable overall in performance to the less accurate ISORROPIA. For some tests ISORROPIA was faster but for other tests MOSAIC was faster and several tests had similar timing results (see Table 5 in [76]). The box model timings were validated with a 3-D Eulerian model applied to trace gas and aerosol evolution in the Los Angeles basin for a period spanning three days. The average CPU time per timestep per bin was found to be 125 $\mathsf{\mu }$s on the same platform using the same compiler where the box model timings spanned 13 to 274 $\mathsf{\mu }$s for MOSAIC and 3.5 to 593 $\mathsf{\mu }$s for ISORROPIA. It appears that these tests were conducted with the default convergence settings in ISORROPIA which correspond to the “low” setting in Makar et al. [159]. The execution time of ISORROPIA increases rapidly as convergence criteria are made tighter [120,159].

MOSAIC was found to predict equilibrium concentrations and pH similar to the benchmark model AIM, whereas ISORROPIA showed substantial differences. In three of the cases ISORROPIA predicted a pH near 8 when AIM and MOSAIC predicted it to be below 3. In general ISORROPIA exhibits non-negligible errors compared to AIM and MOSAIC at low and moderate RH. These deviations were attributed to errors in the activity coefficients using the Bromley mixing rule at high ionic strengths, and approximation errors in the treatment of water content and solid-liquid equilibria in the MDRH region. However, a recent study of pH for fine particles during Beijing haze events by Song et al. [160] found coding errors in ISORROPIA for its forward solution regime that result in unrealistically high (over 7) pH for several sub-cases where the model fails to properly account for partitioning of ammonia between the gas and aerosol phase. A revised version of ISORROPIA gave substantially lower pH values more consistent with AIM and would likely address the issues found by Zaveri et al. [76].

The performance of UHAERO was evaluated against ISORROPIA by Amundson et al. [120]. Solid-liquid phase equilibria, with supersaturated salts precipitated, were calculated in ammonium-RH space spanned by a 100 point grid for a fixed sulfate fraction. UHAERO with the Extended Universal Quasichemical (ExUNIQUAC) activity model [134] was found to take about twice the time of ISORROPIA with the default “low” convergence criteria settings. ISORROPIA with the “high” convergence criteria [159] was about nine times slower relative to the “low” settings. UHAERO using the PSC activity model was about four times slower than ISORROPIA with the default convergence criteria. Considering the accuracy of UHAERO with either activity model, the performance relative to ISORROPIA is notable. Given the various approximations in ISORROPIA it cannot reach the accuracy of UHAERO even with the “high” convergence criteria so the value of running ISORROPIA with these settings is limited.

Pye et al. [161] have used E-AIM, MOSAIC, ISORROPIA II and EQUISOLV II to evaluate predictions of pH against idealized simulations and observations. E-AIM accounts for the single-ion activity coefficient of H${}^{+}$ but the other models use mean molal ion activities for their computations and do not compute single-ion activities. For a system composed of water, ammonium bisulfate, sulfuric acid and ammonia, the predicted pH of MOSAIC was found to deviate from E-AIM with decreasing RH, from 0.03 pH units at 99% to 0.46 pH units at 40%. ISORROPIA II and EQUISOLV II exhibit more substantial deviation from E-AIM at low RH with the former being too low and the latter being too high. A more detailed evaluation of various aspects can be found in their study.

2.3. Organic Thermodynamics Models or Components

The SOA schemes in AQ models deal with very simple representations of the aerosol physical state. Organic aerosol (OA) are vastly more complex than inorganic aerosols in terms of physical characteristics (e.g., [19]). Depending on the organic species profile of OA, they may form amorphous instead of crystalline phases upon dehydration. The characteristics of these amorphous phases can be those of glass, rubber, gel or viscous liquids and are sensitive to temperature. The shape and porosity of these amorphous phases reflect the history of the dehydration process and the chemical composition. Multiple separated phases with differing morphologies can exist in OA and organics can significantly impact the behavior of inorganic fractions including the total suppression of efflorescence of inorganic salts [162]. These microphysical characteristics have a direct impact on the heterogeneous chemistry associated with OA (e.g., [163]). They also impact the hygroscopic growth and activation of OA [19] which has implications for OA as CCN and IN and hence their role in cloud formation.

The details of the partitioning into multiple liquid and solid phases of organic and mixed organic-inorganic aerosol particles impact their evolution. Phase partitioning can occur thermodynamically (e.g., LLPS) and through aerosol dynamics resulting in composite aerosol characteristics (e.g., [88,164]). The latter case is possible in the atmosphere since aerosol particles reflect multiple emissions sources such as lithogenic, sea salt, combustion soot, and other types. In this section we present some examples of thermodynamics models for organics. Currently, only offline comprehensive coupled organic-inorganic models are available. However, some attempts have been made towards reducing computational cost through the use of chemical surrogates [112,165,166] or developing simplified models which could be incorporated into 3-D AQ models [167,168].

Table 2. Organic Models and Components.

Thermodynamics Model	Acitivity Coefficient Scheme	Condensed Phase Evolution	Gas-Aerosol Partitioning Scheme	Host Model Availability
MPMPO [169,170]	Uses the UNIFAC [109] model for activity coefficients. MPMPO is for SOA in used in conjunction with an inorganic thermodynamics scheme (e.g., SCAPE2).	Aerosol particles are assumed to be of two types: single phase (mixed inorganic and organic) or separated into an aqueous phase and an organic phase.	Pankow [171] equilibrium absorption model is used for partitioning into the organic phase. Activity-corrected effective Henry’s law scheme is used for partitioning into the aqueous phase. An iterative solution approach is used for both phases since activity coefficients are functions of composition and to conserve mass.	UCI-CIT [170]
				[email protected]
UHAERO -organics [87]	Hybrid scheme based on the Clegg-Seinfeld-Brimblecombe (CSB) model [172,173,174] where the Pitzer-Simonson-Clegg (PSC) method is used for inorganics and UNIFAC for organics. The Clegg et al. [135,136] solvent activity model is used to determine water content instead of Zdanovskii-Stokes-Robinson (ZSR).	Includes organic compounds which are not water soluble and multiple liquid phases are allowed to form. Uses an efficient Gibbs free energy minimization approach with a hybrid solver based on the primal-dual active-set algorithm for the gas-aerosol equilibrium and the primal-dual interior-point algorithm for the liquid phase equilibrium.	Equilibrium	Offline
				Co-authors
PD-FiTE-organics [88]	Follows the Taylor expansion approach of the inorganics version but uses UNIFAC for organic activity coefficients. The model is for organics only and not for mixed inorganic-organic systems. ZSR is used for water content.	Not included.	Not included.	Can be used online.
				[email protected]
UNIPAR [175]	Activity coefficients are assumed to be unity in the organics particle phase. Semi-empirical approximations are used for activity coefficients in the inorganic phase.	Uses the Bertram et al. [176] scheme for liquid-liquid phase separation. Predicts aerosol water content. Organic chemistry is cast as oligomerization in the organics phase and acid-catalyzed reactions in the inorganic phase. Treats formation of organosulfates.	Uses a lumped gas-phase chemistry based on the Master Chemical Mechanism [177] with 30 groups to represent volatile organic compound precursor oxidation products. Pankow [171] equilibrium partitioning between the gas and both condensed phases is assumed.	Offline.
				[email protected]
SOAP [165]	Uses UNIFAC for short-range interactions and AIOMFAC [89] for medium-long range interactions between electrolytes and organics. This model needs to be coupled to an inorganic thermodynamics model to provide aerosol water content, pH, concentrations of inorganic species and the ionic strength.	Organic compounds are represented by surrogates to reduce complexity. Phase separation based on Gibbs free energy minimization. With dynamic option, a multilayer scheme for aerosol composition is used. Morphology factors are used to account for aerosol physical heterogeneity which impacts diffusion of organics.	Solves full reaction-diffusion system and includes the Kelvin effect. Two solution options: equilibrium or dynamic. Capable of treating organic aerosol high viscosity effects. Uses a thin interface layer without diffusion to avoid artificially limiting absorption (dynamic option).	Can be used online with the equilibrium option.
				http://cerea.enpc.fr/soap/index.html
Pye et al. [167]	Based on binary mixtures of individual organic constituents and water using a one-constant Margules equation. Coefficient scaling was required to get realistic water uptake.	Calculates the phase separation relative humidity based on You et al. [178]. Water uptake to the organic phase is predicted using $\kappa$-Köhler theory [179]. Uses ISORROPIA II to calculate water uptake by the inorganic phase.	Modified Pankow [171] equilibrium partitioning into the organic phase. The volatility constant is defined to be proportional to the total moles of the absorbing medium including water and the inorganics following Zuend et al. [180].	CMAQ v5.1 [181]
				[email protected]
BAT [168]	Based on binary mixtures of individual organic constituents and water based on the Duhem-Margules equation. Model parameters were optimized using AIOMFAC.	Solves for liquid-liquid phase separation. Calculates equilibrium water uptake given the relative humidity. Uses detailed or averaged physicochemical properties of organics.	Coupled to a non-ideal volatility basis set scheme. Uses an equilibrium solver for co-condensation of organics.	Can be used online.
				https://github.com/Gorkowski/Binary_Activity_Thermodynamics_Model

summarizes the models presented in this section.

2.3.1. MPMPO

The Model to Predict the Multiphase Partitioning of Organics (MPMPO) was developed by Griffin et al. [169] and is a fully coupled hydrophobic-hydrophillic organic gas-aerosol partitioning model. This model allows the inorganic and organic components to interact. Equilibrium partition coefficients are calculated for each precursor organic that is absorbed into the organic aerosol phase similar to the approach of Pankow [182]. The Myrdal and Yalkowsky [183] method is used to estimate the vapor pressure of organic compounds. UNIFAC is used to calculate activity coefficients. For the aqueous aerosol phase species specific Henry’s laws and activity coefficients from UNIFAC are used.

Aerosol liquid water content (LWC) affects both the organic and aqueous phase partitioning of organic compounds. For low LWC partitioning into the organic phase is favored. As LWC increases partitioning into the aqueous phase increases as well. In MPMPO pH influences partitioning as well. Inorganics are the primary controls on pH in the aerosol. As pH increases so does the aqueous phase concentration of dissociating organics. LWC and pH act in concert to control which aerosol phase accepts gas phase organic compounds. The aqueous phase becomes the primary recipient at high LWC and pH, and conversely the organic phase dominates at low LWC and pH. However, total condensed organic concentrations in the aerosol show a weak dependence on pH consistent with observations.

Primary organic aerosol impacts SOA formation by affecting the activity coefficients. For primary organics without functional groups the activity coefficients in the organic phase are increased which leads to the reduction of concentrations of gas phase precursors in the organic aerosol phase. For oxygenated primary organics the effect is to reduce the activity coefficient in the organic aerosol phase which leads to an increase in organic phase concentrations of secondary organics.

Compared to decoupled approaches such as that of Pun et al. [35], the coupled approach in MPMPO increases the SOA production by about 10%. MPMPO interacts with the inorganic thermodynamics model in two stage fashion via updates to LWC. The inorganic model provides the initial state for MPMPO but is called again to determine the effect of organic anions and additional water on the inorganic components of the aerosol.

2.3.2. UHAERO For Mixed Organic-Inorganic Systems

Amundson et al. [87] extended UHAERO to include organic species and LLPS. Activity coefficients are determined using a hybrid scheme based on the Clegg–Seinfeld–Brimblecombe (CSB) model [172,173,174] where PSC is used for inorganics and UNIFAC is used for organics. Instead of ZSR the more accurate Clegg [135,136] solvent activity model is used to determine aerosol water content.

As with the inorganic version, gas-aerosol exchange is assumed to be in equilibrium. The condensed phase state of the aerosols is modeled by a hybrid Gibbs free energy minimization scheme. The gas-aerosol equilibrium is determined with the primal-dual active-set algorithm but the liquid phase equilibrium is determined with the primal-dual interior-point algorithm.

This version of UHAERO, coupled with a non-equilibrium gas-aerosol mass transfer model such as that of Zaveri et al. [84] could form the basis of a comprehensive thermodynamics framework. However, the numerical cost will be substantially higher than for pure inorganic formulations.

2.3.3. PD-FiTE for Organics

Topping et al. [88] introduced an organic version of their inorganic activity coefficient scheme. The fitting procedure is applied to UNIFAC instead of ADDEM. This model is for organics only and not for mixed inorganic-organic systems. It uses ZSR to calculate the water content of organic aerosols.

2.3.4. UNIPAR

Im et al. [175] have developed an intermediate complexity model that takes into account the physical characteristics of SOA such as LLPS and treats chemistry in an explicit fashion. The Unified Partitioning–Aerosol Phase Reaction (UNIPAR) model uses the Pankow equilibrium absorption model to calculate partitioning of organics into both the the organic and inorganic aerosol phases. UNIPAR uses the Bertram et al. [176] semi-empirical scheme for LLPS, efflorescence and deliquescence in mixed organic ammonium sulfate systems. Organic and inorganic aerosol phases are assumed to be fully separated if there is phase separation. Activity coefficients for organics are set to unity in the organic phase. The activity coefficients of organics are semi-empirically approximated in the inorganic phase. It is assumed that organic matter (OM) formed in the aerosol phase is non-volatile and chemically irreversible.

The explicit gas phase chemistry from the Master Chemical Mechanism (MCM) [177] is used to derive a lumped mechanism. Oxidation products of precursor VOCs are lumped into 30 groups based on their vapor pressures and aerosol phase reactivity levels (fast, medium, slow, nonreactive and multi-alcohol). Partitioning between the gas-phase and either of the condensed phases is based on the Pankow absorption equilibrium approach. Organic chemistry in the aerosol phase is partitioned into oligomerization in the organic phase and acid-catalyzed reactions with inorganics in the inorganic phase. UNIPAR also includes organosulfate species (OS) formation. Aerosol acidity is tracked dynamically based on ammonia titration and OS concentrations.

Beardsley and Jang [184] have successfully applied UNIPAR to model SOA formation from isoprene photooxidation both with and without inorganic seed. Isoprene is the largest single natural source of OA globally. SOA formation from isoprene is sensitive to aerosol acidity and LWC. This requires a model that can keep track of aerosol phase chemistry and water amount. Simple SOA models such as the two-product scheme of Odum et al. [185], SORGAM [80] and gas phase ageing parameterizations such as those based on the volatility basis set (VBS) approach [186] are not able to do this.

2.3.5. SOAP

Introduced by Couvidat and Sartelet [165], SOAP is in the class of more comprehensive thermodynamics models for organics such as ADDEM. It also treats the full reaction-diffusion system for gas-aerosol mass transfer, as in Zaveri et al. [84]. Only the effect of inorganic compounds on organic compound formation is considered but in general there is a two way interaction which affects the aerosol size and composition. SOAP does not treat inorganic aerosol formation and if such a component was added it could handle the two way interaction. SOAP has two options to model condensation and evaporation of organic aerosols: equilibrium and dynamic.

SOAP uses a multilayer approach to model the particle-phase diffusion of organics. By contrast, inorganic thermodynamics models make the implicit assumption that the condensed phase is well mixed. Pure inorganic aerosol particles do not exhibit the large viscosities that can be found in SOA (e.g., [21]). The multi-layer approach involves tracking whether organics are hydrophobic, hydrophillic or hybrid. It is an approximation to the complex multiphase behavior of OA which are not diffuse mixtures of electrolytes as in the inorganic case but rather compound entities composed of phase-separated organic and water mixtures.

Couvidat and Sartelet [165] evaluated the performance of SOAP in a 3-D AQ model and found the computational cost to be too high for the dynamic option. For the equilibrium option, SOAP appears to be tractable in 3-D AQ models which is supported by its use (with the equilibrium option) in a subsequent study using CHIMERE [166]. It should be noted that this is made possible by the substantial degree of simplification of the organic species spectrum in SOAP through the use of lumped surrogates which reduces the numerical cost of activity coefficient calculations by UNIFAC and AIOMFAC employed in SOAP. Given that various VOC precursors result in thousands of daughter products through oxidation and photolysis in the real atmosphere (e.g., [187]), and that these daughter products have a wide spectrum of functional group composition, lumping into surrogates necessarily loses information about the physicochemical properties of the system.

2.3.6. Updated Pankow Scheme

Pankow et al. [188] revisit the application of the equilibrium theory of Pankow [182] but with refined treatment of key parameters. Instead of a simplistic approximation of the activity coefficient for each species, use is made of UNIFAC. In contrast to previous approaches of Odum et al. [185] and Donahue et al. [186] water is included as a constituent of SOA. The addition of water affects the partition coefficient for every organic species by changing the mean molecular weight of the SOA particle. Inclusion of water and inorganic constituents was recommended by Zuend et al. [180] for simplified SOA mass transfer models. It is likely that SOA are not pure organic particles and inorganics will affect SOA partition coefficients through the activity coefficients and mean molecular weight of the SOA particle. Unlike with previous simplifications there is at least partial accounting for organic polarity effects via the activity coefficient. So the composition of the SOA particle can hinder or enhance the partitioning of gas phase organics depending on their polarity. Simulations conducted using CMAQ with comprehensive SOA parameterizations indicate a substantial contribution of water to SOA mass in the southeastern USA, so this aspect of SOA formation cannot be ignored.

2.3.7. Pye CMAQ Extensions

Pye et al. [167] updated the CMAQ model [189] to study the influence of aerosol water on OA growth. This is the first 3-D AQ model study to consider the role of water in OA growth which includes the mixing of organic and inorganic phases. Like water soluble inorganics, hydrophillic organic compounds increase the uptake of water by OA (including SOA and mixed aerosol species). This in turn promotes additional aerosol growth through uptake of organics including loss of gas phase compounds through aqueous reactions (e.g., glyoxal).

LLPS occurs below the SRH. SRH is primarily a function of the inorganic salt present in the aerosol particle with the highest values associated with ammonium sulfate. Pye et al. [167] use the experimental results of You et al. [178] for SRH in the presence of ammonium sulfate. In terms of the ratio of OM to organic carbon ($\mathrm{OM}/\mathrm{OC}$) it is given by:

$$\mathrm{SRH}={\left[1+exp\left(7.7\frac{\mathrm{OM}}{\mathrm{OC}}-15.8\right)\right]}^{-1}$$

(17)

This choice was motivated by the dominance of ammonium sulfate and ammonium bisulfate for the observational campaign used to evaluate the model performance, the Southern Oxidant and Aerosol Study (SOAS) 2013 [190]. In general, there will be spatial and temporal variation of the dominant inorganic salt so a more general expression for SRH would be needed.

Water uptake to the organic phase is predicted using $\kappa$-Köhler theory [179] and solving for the volume equivalent diameter, D:

$$\mathrm{RH}-\frac{D^3-D_{core}^3}{D^3-D_{core}^3(1-\kappa )}exp\left(\frac{4{\sigma }_w{M}_w}{RT{\rho }_{w}D}\right)=0$$

(18)

where $D_{core}$ is the volume equivalent accumulation mode diameter excluding water associated with organic species, ${\sigma }_w$ is the surface tension of water, $M_w$ is the molecular weight of water, ${\rho }_w$ is the water density, R is the universal gas constant, and T is the temperature. The mass of particle liquid water associated with organic compounds per volume of air (W${}_o$) was obtained using:

$${\mathrm{W}}_o=\frac{\pi N_p{\rho }_w}{6}(D^3-D_{core}^3)$$

(19)

where $N_p$ is the number of aerosol particles per volume of air. The mass of water associated with inorganics per volume of air (W${}_i$) was calculated using ISORROPIA v2.2 [117]. Total aerosol water is given by W${}_o$ + W${}_i$.

The hygroscopicity parameter, $\kappa$, is determined as a volume weighted sum of the individual constituent parameters, ${\kappa }_i$:

$$\kappa =\sum _{i=1}^n\frac{{\kappa }_i{V}_{core,i}}{V_{core}}$$

(20)

The Lambe et al. [191] CCN-based ${\kappa }_i$ values were used even for subsaturated conditions for technical reasons.

Instead of the Pankow [171,182] definition of volatility, this work uses:

$$C_i^{*}=\frac{G_i{M}_{i}N}{A_i}$$

(21)

where $G_i$ is the gas phase concentration of constituent i, $A_i$ is its aerosol phase concentration, $M_i$ is its molecular weight and the total moles of the absorbing medium, N, is given by:

$$N=N_{other}+\sum _i\frac{A_i}{M_i}$$

(22)

As suggested by Zuend et al. [180] simplified models should include water and inorganic constituents in the absorbing phase to better capture more detailed calculations. $N_{other}$ represents this inorganic fraction.

Activity coefficients for organic species were derived from a one-constant Margules equation [192]. The temperature dependent constant in this equation corresponds to an infinitely dilute mixture and was thus estimated using a combination of Henry’s and Raoult’s laws. The resulting expression for the activity coefficients, $a_i$, is:

$$ln\left(a_i\right)=x_w^{2}ln\left(a_{\infty }\right)=x_w^2[ln\left(M_i{\rho }_w\right)-ln\left(H_i{C}_{0,i}^{*}RT{M}_w\right)]$$

(23)

where $x_w$ is the mole fraction of water in the absorbing/partitioning medium, $H_i$ is the Henry’s law constant for species i and $C_{0,i}^{*}$ is the pure species saturation concentration at T. It can be shown that the saturation concentration as a function of water is given by:

$$C_i^{*}=C_{0,i}^{*}{\left(a_{\infty }\right)}^{{\left(\frac{N_w}{N}\right)}^2}$$

(24)

where $N_w$ is the aerosol water in moles. This equation is applicable for any aerosol water fraction.

The activity coefficient formulation adopted by Pye et al. [167] results in a competitive response in SOA growth from addition of aerosol water. More water increases the amount of partitioning medium which act to increase SOA mass, but at the same time there is an increase in the volatility of organic constituents as seen from the above expression which acts to reduce SOA mass.

2.3.8. BAT

A reduced complexity organics thermodynamics model has been developed by Gorkowski et al. [168]. The Binary Activity Thermodynamics (BAT) model can use explicit (composition-resolved) or aggregated (limited chemical composition information) organic species characteristics such as the oxygen to carbon ratio, hydrogen to carbon ratio, molar mass and vapor pressure. BAT captures LLPS into aqueous-rich and organics-rich phases. The primary approximation to reduce numerical cost in BAT is that only non-ideal interactions with water are treated for each organic constituent (or lumped proxy) instead of non-ideal interactions between all constituents.

Gas-particle partitioning is simulated by coupling BAT to a non-ideal VBS model. Non-ideal VBS is distinguished from regular VBS in that $C^{*}$ values are no longer constant due to dependence on activity coefficients and mole fraction of organics and water. The non-ideal VBS is consequently more numerically expensive. Gorkowski et al. [168] reduce the computational burden by using a neural network to generate an optimal initial guess for the iterative solver for the gas-liquid partitioning coefficients.

A neural network approach is also used to simplify the solution of the water activity, $a_w$, for a given mole fraction of organic species, $x_{org,i}$, where i is the species index. Typically the former is known for a given RH and the latter is not. To avoid using an expensive iterative solver, a neural network is trained to guess $x_{org,i}$. There is close agreement for $a_w$ less than 0.95, but an iterative refinement is required for higher values. Overall, the neural network approach makes the numerical cost insignificant compared to the iterative solution approach.

Comparison of $\alpha$-pinene and isoprene SOA systems modeled by AIOMFAC and VBS + BAT in terms of co-condensation of organic matter and hygroscopic growth indicates a high level of agreement as long as the properties of organic compounds are represented with sufficient accuracy. In the isoprene SOA case, it was found that an overly-simplified description of constituents as multifunctional hydroperoxides underestimated the hygroscopic properties of isoprene-derived epoxydiol (IEPOX) oxidation products. So the VBS + BAT model system requires a more detailed view of the SOA chemistry to be properly tuned.

2.4. Hygroscopic Behavior of Mixed Organic-Inorganic Aerosols and Liquid-Liquid Phase Separation Effects

A primary deficiency of the pure inorganic aerosol formulation common in AQ models is that it fails to address the actuality that a substantial fraction of aerosols are mixtures of organic and inorganic constituents. In the case of aerosols in the nucleation, Aitken and accumulation modes, the organic fraction is generally significant and thus these aerosols have distinct hygroscopic characteristics compared to pure inorganic idealizations. Organic constituents increase particle viscosity and under low RH can result in amorphous or semi-solid characteristics [81,193,194]. A significant organic fraction in sub-micron aerosol particles generally acts to reduce their hygroscopic growth rate for subsaturated RH conditions [23]. Mixed organic-inorganic aerosols also exhibit distinct physical characteristics in the form of heterogeneous internal mixing states including LLPS. The occurrence of phase separation has been established to occur in atmospheric aerosols by You et al. [164].

Phase separation affects the uptake of organics by mixed aerosols since hydrophobic organic species can partition into the organics-rich phase much more readily than into the mixed aqueous phase that occurs at higher RH [180,195] and for aqueous aerosol particles without phase separation. The increase in organics uptake can be as high as 50%. AQ model idealizations where SOA is modeled as a distinct aerosol species maintained by an equilibrium partitioning of organics cannot account for this growth pathway.

Organics change the thermodynamic behavior of inorganic constituents. They generally act to suppress efflorescence [196,197,198]. Under low RH conditions, hydrophillic organics delay evaporation and the onset of efflorescence of inorganic salts and enable more water to be retained for a given RH level compared to pure inorganic aerosols. This effect will occur even in the case of LLPS since the aqueous rich phase will retain a fraction of hydrophillic organics. Bodsworth et al. [199] find that less organic content in mixed sulfate-organic aerosol is need to inhibit efflorescence of ammonium sulfate at colder temperatures. Given the relatively high fraction of organics in the sulfate aerosols it is possible that a larger fraction of sulfate aerosol particles are liquid in the upper troposphere than commonly assumed. So the contribution of solid ammonium sulfate to heterogeneous ice nucleation is likely to be smaller. The organic fraction of aqueous IN is also believed to be responsible for the observed supersaturation of water inside cirrus clouds in the cold tropical tropopause layer (TTL) [200]. Organic-rich aqueous aerosols become amorphous under cold temperatures and are effective for heterogeneous ice nucleation. This ice nucleation pathway can explain the observed low number density of ice crystals and supersaturation RH. Heterogeneous ice nucleation occurs at a much lower initial RH than homogeneous ice nucleation and is facilitated by the very low temperatures in the TTL. However, heterogeneous nucleation allows RH to increase after the onset of freezing whereas the RH is reduced after homogeneous nucleation. Ultimately, less ice crystals are formed through heterogeneous nucleation and they grow to a larger size and significantly affect the cloud radiative properties.

As organics influence phase transitions for inorganic aerosol constituents, the reverse is also true. Baustian et al. [162] find that the glass transition temperature for organic aerosol particles was substantially reduced when mixed with ammonium sulfate in a 1:1 ratio. The degree of reduction depends on the organic species involved. But there is always a reduction which results in mixed particles having a liquid state over a wider range of temperature and RH. The effect of inorganic salt solutions is to reduce the bulk particle viscosity. It was found that particles of the organic species considered could act as effective heterogeneous ice nuclei at temperatures of 220–234 K. However, mixtures with ammonium sulfate only exhibited freezing at temperatures of 200 K and lower.

Baustian et al. [162] also investigated differences in crystalline and amorphous deliquescence of organics. Compared to crystalline organic particles, deliquescence of amorphous particles initiated at a lower RH and spanned a significant range of RH. There was little change in particle volume over the transition RH range. This reflects effects of high particle viscosity and it is even possible for particles to experience reduction in volume due to micro-structural rearrangement [19].

Hodas et al. [201] demonstrate that no simplified representation of the hygroscopic behavior of mixed phase aerosols can be made that avoids having to deal with the details of aerosol phase composition. In particular, the errors in hygroscopic growth factors from simplifying assumptions are substantial. Depending on the the organic system being considered the peak error ranged from 10% to 26%. Peak error was found for RH around 80% which corresponds to a phase transition as determined by AIOMFAC. Systems that underwent LLPS showed the most deviation. As a consequence the particle optical and radiative properties deviate significantly as well since the particle size and refractive index are sensitive to growth and composition [202,203,204]. Recent work indicates that there can be regional impacts of organics on CCN activation for example in the boreal forest regions where aerosols have substantial biogenic organic fractions [205]. Sensitivity tests based on adjustment of the hygroscopicity factor ($\kappa$-Köhler framework of [179]) in two climate models showed a non-negligible impact on top-of-the-atmosphere radiative flux with differences between −1 and −0.25 W${\mathrm{m}}^{-2}$ [206]. The range of variation of the hygroscopicity factor was based on field measurements. In addition to the hygroscopic growth, aerosol composition affects the evolution of surface tension during activation [207]. As noted by Davies et al. [207], due to the approximate cancellation effects of bulk to surface partitioning of surface-active constituents, which increases particle water activity, and dilution due to rapid water vapor uptake, the composition effect on surface tension has been neglected in commonly used activation parameterizations. However, Davies et al. [207] demonstrate that the simplistic one-parameter formulation of these parameterizations cannot capture the full surface tension effect and makes its implementation in such schemes difficult.

It is likely that biases in aerosol optical depth found by the Atmospheric Chemistry and Climate Model Intercomparison Project (ACCMIP) [208] reflect in part the simplistic process representations related to the hygroscopic growth and activation of aerosols. When aspects pertaining to hygroscopic growth and cloud wet scavenging are improved there is improvement in model radiative biases [203,209]. Latimer and Martin [210] find that adjusting hygroscopic growth of SIA and OA below 65% RH along with the size distribution in a physically plausible way can remove substantial biases (40–80%) in the aerosol mass scattering efficiency in the GEOS-Chem model. However, there are many other model issues pertaining to emissions accuracy, lack of aerosol chemistry and other simplifications that make attribution of improvements difficult. Given advances in our understanding in recent years, there is need and opportunity for studies which quantify errors associated with processes representations in climate and AQ models.

The influence of the organic content of submicron aerosols manifests as discontinuous hygroscopic growth behavior for subsaturated and supersaturated conditions as established by field and laboratory measurements ([211], and references therein). Following Hodas et al. [211] we list potentially relevant processes that result in this behavior. For subsaturated conditions organics act to reduce hygroscopic growth but for supersaturated conditions there is enhanced CCN activity. The origin of this discontinuity is likely to involve several effects. For semi-solid organic-rich particles water uptake occurs via adsorption for subsaturated conditions [23]. This reflects high particle viscosity at low RH. As RH increases so does aerosol water content and it is possible for a transition to low viscosity conditions to occur near 100% RH [193,212]. LLPS is likely to occur at higher RH due to the nonlinear uptake of water by organic-rich aerosols and LLPS will increase CCN activity of aerosol particles [213]. Above 95% RH, water uptake becomes sensitive to the surface tension of liquid particles [214]. Surface active organics act to lower the surface tension and may significantly contribute to the observed discontinuous water uptake for high RH [215]. Supersaturated conditions are often associated with decreasing temperatures which allow condensation of higher volatility organic compounds. Thus it is possible for enhancement of water uptake to result from the partitioning of water soluble organic compounds into the aerosol phase [88,216].

Hodas et al. [211] investigated hygroscopic growth and CCN activity of aerosols composed of pure polyethylene glycol (PEG) and solutions of polyethylene glycol and ammonium sulfate (PEG-AS). PEG was chosen as a proxy for atmospheric organic compounds for several reasons. PEG is available in a wide range of molecular weights which enables study of volatility, solubility and viscosity effects. PEG-AS mixtures can exhibit LLPS at certain RH levels. The pure component surface tension of PEG compounds is much lower than water. So LLPS and surface tension effects can be studied in this system. Both pure PEG and peg-as mixtures were found to have lower hygroscopic growth factors (HGF) compared to pure ammonium sulfate for molecular weights of PEG spanning the range 190 g/mol to 11,500 g/mol (PEG-200, PEG-1000 and PEG-10000). PEG-AS mixtures exhibited higher HGF at RH over 75%. In the case of CCN activity (for 0.8% supersaturation), aerosol particles involving heavier molecular weight PEG compounds exhibited higher activation fractions for dry diameters less than 120 nm. Pure PEG particles were resistant to activation for dry diameters less than 60 nm. By contrast, PEG-AS mixtures experience rapid activation for particles smaller than 60 nm dry diameter. However, smaller molecular weight PEG compounds were associated with smaller activation fractions compared to pure ammonium sulfate particles for particles with a dry diameter of less than 120 nm. Both pure PEG and PEG-AS mixed particles converged to the CCN activation fraction of pure ammonium sulfate as the dry particle diameter increased toward 200 nm. In the case of PEG with molecular weight in the 10,000 g/mol range the PEG-AS mixture exhibited enhanced activation over pure ammonium sulfate. The most likely explanation is that a fraction of the PEG compounds in this size range are surface active.

Hodas et al. [211] find that non-ideal behavior as measured by activity coefficients using AIOMFAC, is a function of molecular weight and exerts a substantial effect on HGF for subsaturated conditions. In the saturation limit the non-ideal behavior disappears. LLPS was found to introduce a kink in the hygroscopicity parameter ($\kappa$-Köhler theory of [179]) as a function of RH for PEG-AS systems involving PEG compounds with a molecular weight range near 1000 g/mol. This kink corresponds to an effective increase in hygroscopicity. In general, $\kappa$ declines with increasing RH for pure PEG and PEG-AS mixed systems. Similar decline in hygroscopicity with increasing RH is found for aerosols involving oxidation products of $\alpha$-pinene [23,213].

Phase separation will have different effects on aerosol water depending on the LLPS morphology. A hydrophobic organic coating on an aqueous core will inhibit both evaporation and condensation of water. If a lens forms instead of an organic shell there will be an impact on the evaporation rate due to the reduced effective surface area of the aerosol particle. However, Davies et al. [217] conclude based on laboratory studies that surface shell coatings for certain organics are likely to be unstable and will break up, thereby allowing rapid exchange of water with the aerosol particle.

The impact of organics on the surface tension and the consequent impact on aerosol CCN activity was found by Ovadnevaite et al. [218] to depend on LLPS. Previous research indicated that surfactants would reduce both the surface tension and the solute effect by reducing organic solutes in the aerosol particle interior. So that both the Kelvin effect and the Raoult effect were reduced, thereby counteracting any enhancement of CCN activity. However, mixed organic-inorganic aerosols tend to have a core-shell structure (e.g., [38]) and Ovadnevaite et al. [218] find that even for partially engulfed particles the impact on surface tension is large without an associated negative effect on the Raoult effect. Taking this into account results in CCN number concentrations an order of magnitude larger than currently predicted by climate models. Davies et al. [207] find that surface tension evolution during aerosol growth affects activation. Even if the surface tension in the system reaches that of pure water at activation, the history of the surface tension at prior stages is reflected in the particle diameter and critical supersaturation.

A comprehensive model approach based on, for example, AIOMFAC which accounts for LLPS does a reasonable job reproducing the observed behavior [37]. However, thermodynamics models that treat both organics and inorganics, for example ADDEM [86] and AIOMFAC [89], are numerically expensive which generally inhibits their use in 3-D AQ models. In this light, Zuend and Seinfeld [37] considered approximations. Of the approximate methods considered in Zuend and Seinfeld [37], only one had the best fit to the full solution, namely the assumption of no LLPS but with non-ideal interactions taken into account through activity coefficients together with solid-liquid phase equilibrium permitted for ammonium sulfate. However, this approximation is the most numerically expensive as it requires the AIOMFAC model to be used to determine the activity coefficients. The ZSR approach with the assumption of full LLPS was second best but failed for systems which did not exhibit phase separation.

Bertram et al. [176] introduced a model for LLPS of mixed ammonium sulfate and organic aerosols. The simplified approach adopted was based on the mass ratio of organic material to sulfate and the elemental ratio of oxygen to carbon of the organic material. A training data set was used to parameterize LLPS, RH, (SRH), ERH and DRH. The SRH parameterization predictions were found to be within 15% for 88% of the measurements. The organic material modeled was based on the oxidation products of isoprene, $\alpha$-pinene and $\beta$-caryophyllene. A core-shell morphology for LLPS is assumed in this model.

Kwamena et al. [219] introduced a thermodynamic model that predicts the equilibrium morphology of mixed phase aerosol. This model applies to two liquid phases only but accounts for lens formation in addition to the core-shell regime. At present no comprehensive model system exists that captures the full morphological and compositional inorganic-organic aerosol formation and evolution process. A relatively unexplored path is to model the impact of the detailed aerosol particle state on its water content and size without explicitly resolving the microphysics. It may be possible to use machine learning methods applied to extensive measurement data to obtain a parameterization based on variables such as RH, temperature and VOC precursor class (e.g., $\alpha$-pinene or anthropogenic hydrocarbons).

2.5. Heterogeneous Efflorescence

One of the most uncertain aspects of aerosol thermodynamics is the onset of crystallization. The liquid to solid phase transition is complex and occurs through several pathways [220] similar to ice nucleation [221]. This aspect of aerosol models has not been subjected to a detailed parameterization approach (e.g., Zaveri et al. [76]). For inorganic aerosols, homogeneous nucleation described by CNT is assumed to be the dominant process for solid phase formation. However, the solid phase transition in sulfate aerosols can be triggered by seeds such as mineral dust particles that are incorporated into their interior during coagulation and other growth pathways [94]. A feature of heterogeneous efflorescence is that it occurs at much higher RH compared to homogeneous efflorescence and the phase transition RH depends on the physicochemical characteristics of the seed particle. Efflorescence on the metastable branch of the phase transition between liquid and solid states can also be triggered by the contact with an external solid particle. Davis et al. [95] find that contact with soluble inorganic crystalline particles triggered efflorescence of aqueous ammonium sulfate (at 80% RH), sodium chloride (at 75% RH), and ammonium nitrate (at 62% RH). The phase change RH under contact efflorescence is close to the DRH and in the case of ammonium nitrate there is solid phase formation even though no such formation is observed to occur via the homogeneous pathway to RH as low as 1%.

The above results imply that the inclusion of mineral solids into aqueous inorganic or mixed organic-inorganic particles needs to be tracked to properly capture particle phase transitions. This means differentiating primary dust emissions into relevant mineral fractions. The thermodynamic state of aerosol particles must also be taken into account during aerosol dynamics calculations (e.g., coagulation) since collisions can trigger efflorescence. Contact efflorescence has been accounted for by Zaveri et al. [222] in an aerosol dynamics model with a comprehensive mixing state representation, PartMC [223,224]. Although, the highly resolved mixing state modeled by PartMC is not tractable in AQ and climate models, a similar contact efflorescence approach can be adopted. However, since pure inorganic aerosols are not characteristic of the atmosphere, the effects of organics need to be considered.

Heterogeneous efflorescence is not excluded for mixed organic-inorganic particles since solid transitions can be triggered in the inorganic aqueous phase in particles subject to LLPS. Phase separation tends to offset the effect of organics on ERH and DRH of inorganic constituents [225]. However, phase separation can be suppressed as the organic volume fraction is increased. For the system consisting of ammonium sulfate and toluene oxidation products, phase separation was restricted to organic volume fractions under 20%. For higher fractions the effect on ERH and DRH increased substantially due to the onset of a single mixed state. By contrast, systems based on oxidation products of biogenic compounds such as $\alpha$-pinene change the ERH and DRH by less than 4% for organic volume fractions up to 0.96 [176,226,227]. This indicates that biogenic VOC oxidation products are less miscible with inorganic aqueous phases than anthropogenic VOC oxidation products. Realistic simulation of efflorescence in the presence of organics requires identification of the nature of the organic species, the organic volume fraction and the potential presence of LLPS.

2.6. Models for the Pure Liquid Vapor Pressure of Organics

In order to determine the partitioning of an organic vapor between the gas and condensed phase it is necessary to determine its surface vapor pressure. This not only requires the activity coefficient but its vapor pressure when condensed into a liquid, $P^{\circ }\left(T\right)$. The Raoult expression for the surface vapor pressure, $P_c$, of an organic solute is given by (e.g., [180]):

$$P_c=x_c{a}_c{P}^{\circ }$$

(25)

where $x_c$ is the mole fraction of the organic compound in the condensed phase and $a_c$ is its activity coefficient expressed on a mole fraction basis. For particles smaller than 100 nm in diameter, surface tension related effects become important so the above expression needs to be corrected as described in Section 2.8.1.

However, pure component vapor pressure, $P^{\circ }$, laboratory data are lacking for the majority of condensable multifunctional organic compounds predicted to be formed from gas-phase VOC oxidation. This has prompted the development of models that use the functional characteristics of organic compounds to estimate this property. Mishra and Yalkowsky [228] developed a formulation dependent on chemical structure, namely the rotational symmetry and the conformational flexibility number, and did not require regression fitting for each compound. Myrdal and Yalkowsky [183] extended this formulation to account for hydrogen bonding compounds via the entropy of boiling term and used updated data to modify the heat capacity change upon boiling [229,230]. Chickos et al. [230] developed a group contribution approach to estimating heat capacities for organic liquids and solids. The Mishra and Yalkowsky [228] scheme and its update Myrdal and Yalkowsky [183] requires the boiling point temperature of organic compounds, which in most cases must be estimated. Stein and Brown [231] have developed a group contribution based scheme to estimate boiling point temperatures for a large number of organics with an average error of less than 5%. Nannoolal et al. [232] improved upon the previous boiling point temperature modeling using the group contribution approach such that the percentage of estimates with a deviation over 20 K was 3% compared to 16% for the Stein and Brown [231] scheme and the average error was under 2%.

Pankow and Asher [233] introduced a simple group contribution based liquid vapor pressure scheme for organic compounds. The expression for the liquid surface vapor pressure $P^{\circ }$ (atm) is:

$${log}_{10}{P}^{\circ }=\sum _k{\nu }_k{b}_k\left(T\right)$$

(26)

where ${\nu }_k$ is the number of groups of type k and $b_k$ is the associated group contribution. The temperature dependence of $b_k$ is given by

$$b\left(T\right)=B_1{T}^{-1}+B_2+B_{3}T+B_{4}ln\left(T\right)$$

(27)

The number of structural groups considered was 30 and the values of the temperature coefficients, B, were obtained by fitting 272 compounds for which $P^{\circ }$ versus temperature data was available.

Nannoolal et al. [234] have developed and explicit group contribution based model for the vapor pressure of organics. The role of groups and group interactions is to shift the slope of the curve of the logarithm of the liquid surface vapor pressure $P^{\circ }$ (atm):

$$\begin{array}{ccc}\hfill {log}_{10}{P}^{\circ }& =& (4.1012+dB)\left[\frac{T-T_b}{T-0.125T_b}\right]\hfill \end{array}$$

(28)

$$\begin{array}{ccc}\hfill dB& =& \left(\sum N_i{C}_i+GI\right)-0.176055\hfill \end{array}$$

(29)

$$\begin{array}{ccc}\hfill GI& =& \frac{1}{n}\sum _{i=1}^m\sum _{j=1}^m\frac{C_{i-j}}{m-1}\hfill \end{array}$$

(30)

where $N_i$ is the number of groups of type i, $C_i$ is the group contribution of group i, $GI$ is the total group interaction contribution. $C_{i-j}$ is the group interaction contribution between groups i and j ($C_{i-j}=C_{j-i}$ and $C_{i-i}=0$), n is the number of atoms excluding hydrogen in the molecule, and m is the total number of interaction groups in the molecule. The total number of groups considered was 207 which included ones making secondary contributions. These secondary group interactions are necessitated by the fact that not all groups are additive, specifically strongly associating groups exhibiting hydrogen-bonding interactions which are associated with multifunctional compounds. A total of 1663 compounds were used to perform the regression fitting for this model.

Moller et al. [235] modified the Nannoolal et al. [234] model to improve estimates for aliphatic alcohols and carboxylic acids. The updated vapor pressure expression is:

$$\begin{array}{ccc}\hfill {log}_{10}{P}^{\circ }& =& B^{\prime }\left[\frac{T-T_b}{T-C\left(T_b\right)}\right]+D^{\prime }{log}_{10}\left(\frac{T}{T_b}\right)\hfill \end{array}$$

(31)

$$\begin{array}{ccc}\hfill C\left(T_b\right)& =& -2.65+\frac{T_b^{1.485}}{135}\hfill \end{array}$$

(32)

where $B^{\prime }$ is a group and group interaction contribution expression as in Nannoolal et al. [234], $D^{\prime }$ is a group contribution expression dependent on heavy atoms that is set to zero if no aliphatic alcohols or carboxylic acids are involved.

For the atmospheric case polyfunctional organic molecules such as the products of VOC oxidation are of interest and recent experimental data relevant for ambient conditions has shown existing vapor pressure models to be inadequate [236,237]. This has led to the development of a new model by Compernolle et al. [238]. The vapor pressure expression is

$${log}_{10}{P}^{\circ }=A+\frac{B}{T^{1.3}}$$

(33)

where

$$\begin{array}{ccc}\hfill A& =& \left\{\begin{array}{cc}\sum _k{c}_k{a}_k,\hfill & \mathrm{for}\text{zero-}\mathrm{and}\text{mono-functional}\mathrm{compounds}\hfill \\ \sum _l{c}_l{a}_l+\sum _m\frac{c_m{a}_m}{\sqrt{N_{CL}}}+\sum _n\frac{c_n{a}_n}{\sqrt{N_{HB}}}\hfill & \mathrm{for}\mathrm{multi-functional}\mathrm{compounds}\hfill \end{array}\right.\hfill \end{array}$$

(34)

$$\begin{array}{ccc}\hfill B& =& \sum _k{c}_k{a}_k\hfill \end{array}$$

(35)

where different subscripts are assumed to represent different subsets of groups, $N_{CL}$ is the total number of carbonyl, ester and peroxy acyl nitrate groups, $N_{HB}$ is the total number of hydrogen bonding (hydroxyl, acid, hydroperoxide, and peracid) groups. The CL and HB terms for multifunctional compound reflect non-additivity of group contributions and are associated with the diminishing incremental impact on molecular order of successive functional group addition.

2.7. Surface Tension Models for Aqueous Aerosols

Surface tension of aqueous solutions is a key element for aerosol growth, surface heterogeneous chemistry and phase transitions. Laboratory data for all possible mixtures is lacking and the full range of temperature and supersaturation found in the atmosphere has not been measured (e.g., [239]). In fact, surface tension data is insufficient to properly formulate aerosol thermodynamics models for particles smaller than 100 nm [73]. Atmospheric aerosols exhibit substantial compositional variation which evolves with time. Thus, there is a need for surface tension models that can fill the measurement gaps and can be used in aerosol models.

Dutcher et al. [239] have developed a semi-empirical model for aqueous and non-aqueous inorganic electrolyte solutions that applies to ternary mixtures. This model spans the full concentration range for electrolytes and the full range of atmospheric temperatures. However, pure inorganic mixtures are not sufficient to represent the characteristics of atmospheric aerosols, which have a non-negligible organic component. Topping et al. [240] considered the model options for mixed organic-inorganic aqueous solutions. It was found that for higher order organic-inorganic mixtures that using binary mixture data with a mixing rule [241] produced reasonable fits to laboratory data. There is a need for a composition-dependent surface tension formulation to properly calculate the critical supersaturation using the Köhler equation which is important for CCN evolution. The critical point is also increasingly affected with decreasing aerosol particle size by the relative partitioning of solutes between the interior and the surface layer [242]. For spherical particles the surface area grows as the inverse of the radius relative to the volume. This effect is important for CCN and for Aitken and nucleation mode aerosol particle growth.

Boyer et al. [243] have developed a statistical mechanical model for the surface tension of ternary aqueous solutions of organics and inorganics. The model is able to capture laboratory measurements for aqueous mixtures involving two electrolytes, two organics or a combination of an electrolyte and an organic. There are indications that this model can be extended to more components.

2.8. Nanosize Effects on Aerosol Thermodynamics and Morphology

At the nanoparticle scale, aerosol thermodynamics becomes very sensitive to the Kelvin curvature effect and size-dependence of intrinsic particle properties such as surface tension, dry-particle density, and the shape factor accounting for non-sphericity of particles. The impact of observational uncertainties in these properties on thermodynamics is also the largest for nanoparticles. For example, the DRH and ERH of 6 nm NaCl particles increase to 87% and 53%, respectively, from 75% and 45% in the case of particles larger than 40 nm [244]. In addition, hygroscopic growth factors for nanoparticles decrease due to the Kelvin effect. The hygroscopic behavior of nanoparticles is affected by the increased contribution of the surface to the total Gibbs’ free energy of the system. This term increases substantially for particles smaller than 100 nm [245,246].

2.8.1. Chemical Potential in the Small Size Limit

Existing aerosol thermodynamics models use bulk formulations for the chemical potential of solutes (see Section 2.1). In particular, they neglect the contribution of interfacial energy of solid phase solutes in aqueous aerosol particles as well as the energy of the gas-particle interface. These approximations are reasonable for large aerosol particles where the bulk dominates surface contributions to the Gibbs free energy of the system. However, for nanoparticles the contribution of phase interface surfaces increases owing to geometry. As the particle radius decreases an increasing fraction of molecules reside in the surface layer versus those in the bulk. Molecules at phase interfaces are subject to higher energy and entropy compared to the bulk (e.g., [247]) and make a positive contribution to the chemical potential.

Interface curvature has a major impact on the vapor-liquid phase equilibrium [248]. Liquid drop nanoparticles exhibit a large difference between the ambient gas phase pressure and the internal pressure according to the Young-Laplace equation [249]. This surface tension related term modifies the chemical potential of the solute relative to its vapor phase:

$$\begin{array}{ccc}\hfill {\mu }_{g,i}(T,P_g)& =& {\mu }_{g,i}(T,P_{\mathrm{sat},i})+RTlog\left(y_i\frac{P_g}{P_{\mathrm{sat},i}}\right)\hfill \end{array}$$

(36)

$$\begin{array}{ccc}\hfill {\mu }_{l,i}(T,P_l)& =& {\mu }_{l,i}(T,P_{\mathrm{sat},i})+V_{l,i}(P_l-P_{\mathrm{sat},i})+RTlog\left(x_i{a}_i\right)\hfill \end{array}$$

(37)

where ${\mu }_{g,i}$ and ${\mu }_{l,i}$ are the chemical potential functions for chemical constituent i in the gas and liquid phases, respectively. The mole fractions of i in the gas and liquid phases are given by $y_i$ and $x_i$, respectively. $P_{\mathrm{sat},i}$ is the bulk pure component vapor pressure of i, $P_g$ is the gas phase pressure and $P_l$ is the liquid phase pressure. The molar volume of i in the liquid phase is denoted by $V_{l,i}$ and the activity coefficient of i is denoted by $a_i$. The second term in the liquid phase chemical potential becomes negligible in the flat bulk limit. For a spherical droplet the difference between the exterior gas pressure and the interior droplet pressure is given by the Young-Laplace equation:

$$P_l-P_g=\frac{2\sigma }{r}$$

(38)

where $\sigma$ is the surface tension and r is the droplet radius. This equation breaks down in the molecular cluster limit, but for continuum droplets under 10 nm in diameter the internal pressure can reach ten thousand atmospheres.

As a result of the additional surface tension related term in the liquid phase, the Raoult equation needs to be modified so that:

$$y_i{P}_g=x_i{a}_i{P}_{\mathrm{sat},i}exp\left(\frac{V_{l,i}(P_l-P_{\mathrm{sat},i})}{RT}\right)$$

(39)

Shardt and Elliott [248] investigate the behavior of two-component fluids (methanol-ethanol and ethanol-water) using the two-constant Margules equation for the activity coefficients. The surface pressure of solutes is found to increase relatively uniformly for all solute mole fractions ($x_i$). The amount of deviation from the bulk limit increases rapidly for droplets smaller than 50 nm in diameter. The exponential term in Equation (39) is the Kelvin curvature effect. Omission of the curvature term in the Raoult equation cannot be justified by steady state arguments since the Raoult equation already represents a thermodynamic equilibrium. Some thermodynamics models include the Kelvin effect in gas-aerosol partitioning (e.g., [76,250]). However, these models do not take fully into account the size effect on solubility [73]. In particular, the size effect on liquid-solid phase equilibrium which is critical for ERH and DRH is not commonly considered.

A number of other theoretical studies have considered the detailed thermodynamic behavior of aqueous inorganic salt systems in the small size limit [99,245,246,251,252,253]. But not all of them accounted for the effect of the formation of a surface liquid layer. One study advanced the concept of prompt deliquescence [245] based on the assumption that solid particles remain dry until the DRH point, which is not consistent with the observed behavior of solid NaCl particles under increasing RH. NaCl crystals undergo the formation of a liquid surface layer at RH lower than DRH [254,255] and experience gradual dissolution at the DRH rather than instantaneous dissolution as expected from a first order phase transition. The formation of a liquid surface layer for RH below DRH can be expected for other inorganic salts as well [254].

McGraw and Lewis [252] considered a nanoparticle inorganic salt aerosol system that accounted for surface tension contributions to the system Gibbs free energy within the classical capillary drop model. McGraw and Lewis [252] invoke a thin layer criterion to model the liquid layer surrounding a dissolving salt crystal and separating it from the gas phase. With this criterion the chemical potential of the dissolved salt and the salt in the core are equal as are the chemical potentials of the water in the liquid shell and the gas phase. Unlike the bulk limit, it is found that deliquescence and efflorescence are very similar energetically. For nanoparticles deliquescence can be triggered by a nucleation process in the metastable regime just below DRH or like a barrier free spinodal transition near the threshold for instability at DRH. The results of McGraw and Lewis [252] are informative but their approach did not include measurable thin film effects such as the disjoining pressure. The disjoining pressure cannot be neglected in the small size limit since all surface films become effectively thin.

Djikaev et al. [246], Russell and Ming [251], Shchekin et al. [99] and Hellmuth and Shchekin [253] included the surface liquid layer in their thermodynamics models and the associated disjoining pressure effect [256]. The disjoining pressure, ${\mathrm{\Pi }}_D$, reflects the molecular interaction of the two interface boundaries, which give rise to surface tension of a thin film and introduces a substantial contribution to the Gibbs free energy of the system. The disjoining pressure is the deviation from bulk pressure in a thin film due to the overlap of the interface boundaries and decreases rapidly with increasing layer thickness. If the liquid film engulfs a solid core, i.e., wets it, then ${\mathrm{\Pi }}_D>0$. If ${\mathrm{\Pi }}_D<0$, then the liquid will not fully coat the core and will collect in a lens. The former case will be true for the soluble inorganic salts common in atmospheric aerosols.

Djikaev et al. [246] identified the disjoining pressure as the key element controlling the steepness of the deliquescence curve. A higher disjoining pressure leads to a more gradual dissolution of the solid phase. The coupling between the interface layers that results in the disjoining pressure is due to long-range molecular interactions. As a result there is some effect even for relatively thick (versus the core diameter) surface liquid coatings. This is particularly true for small particles where even relatively thick liquid coatings are thin in absolute terms and so the impact of the disjoining pressure on the thermodynamics increases with decreasing particle size.

Hellmuth and Shchekin [253] have produced a very detailed thermodynamic model of a single-solute solid core particle with an aqueous surface layer. They applied this framework to sodium chloride nanoparticles and found that the humidification of these particles at the initial hydration stages, associated with increasing RH, was very sensitive to the disjoining pressure. The disjoining pressure acts to enhance the wetability of the particle surface and results in an earlier onset of hygroscopic growth. This study did not consider the effect of solution concentration on the surface tension but found that for nanoparticles, in the 10 nm dry diameter size range, that it was necessary to consider the non-ideal behavior of the solution in the aqueous surface layer. Tuning of other parameters within uncertainty bounds in the model failed to produce reasonable results if ideal limit solute and solvent activities were used.

2.8.2. Effect of Particle Size on Phase Transitions and Hygroscopic Growth

The thermodynamics behavior of supersaturated inorganic aerosols such as the ammonium sulfate salt exhibits significant deviation from model predictions of phase transitions and hygroscopic growth factors for aerosol particles smaller than 100 nm [73]. This is true even for advanced models such as AIM [112]. In the case of DRH, the error is very large for particles smaller than 20 nm. Without accounting for the size effect, DRH is 100% for ammonium sulfate and salt particles with diameters less than 10 nm. But observations indicate that the DRH is below 85% for particles in this size range. In the case of ammonium sulfate, the AIM model makes substantial errors in predicting the growth factor over the full RH range for particles with a diameter of 6 nm.

The model-measurement disagreement partly reflects the lack of laboratory data required to characterize nanoscale droplets as outlined by Cheng et al. [73]. Small droplets can reach higher levels of supersaturation compared to bulk solution samples. This makes measurements difficult with existing laboratory techniques which involve the use of probes since probes can trigger heterogeneous nucleation of supersaturated droplet solutions.

Another reason for the disagreement is that thermodynamics models have not been formulated to account for nanoparticle effects such as the large contribution of surface tension at interfaces and coupling of interfaces. Due to the quadratic dependence of the particle surface area on radius and the cubic dependece on radius of the particle volume, small particle liquid surface boundary layers increase in relative thickness versus the bulk interior in the small radius limit [99]. Cheng et al. [73] introduce small size corrections by using the Ostwald-Freundlich equation [257] to determine the deliquescence concentration of the solute.

Cheng et al. [73] introduce the differential Köhler analysis (DKA) method to enable determination of the water activity ($a_w$) and the surface tension (${\sigma }_{lv}$). The Köhler equation for the equilibrium saturation over a spherical aerosol particle is

$$\begin{array}{c}\hfill s_w=a_{w}exp\left(\frac{4{\sigma }_{lv}{\nu }_w}{RT{D}_{sol}}\right)=a_{w}exp\left(\frac{4{\sigma }_{lv}{\nu }_w}{RT{D}_d{g}_f}\right)\end{array}$$

(40)

where ${\nu }_w$ is the partial molar volume of water, R is the universal gas constant, T (K) is the temperature, $D_{sol}$ (m) is the wet diameter of the aerosol particle, $D_d$ (m) is the dry diameter and $g_f=\frac{D_{sol}}{D_d}$ is the hygroscopic growth factor. The above equation is an approximation that neglects the size dependence of ${\sigma }_{lv}$ but Cheng et al. [73] find that this error is very small (<1%) for particles between 6 nm and 100 nm in diameter. High precision measurements using a hygroscopic tandem differential mobility analyzer (HTDMA) can be used to determine $s_w$, $D_d$, and $g_f$. So there are two unknowns left, namely $a_w$ and ${\sigma }_{lv}$. Taking the natural logarithm of both sides in the above equation gives

$$\begin{array}{c}\hfill ln{s}_w=ln{a}_w+\frac{A}{D_d}{\sigma }_{lv},A=\frac{4{\nu }_w}{RT{g}_f}\end{array}$$

(41)

Measurements of $s_w$ for two dry diameters, $D_{d1}$ and $D_{d2}$, but for the same value of $g_f$ are made which allows determination of $a_w$ and ${\sigma }_{lv}$ via

$$\begin{array}{c}\hfill a_w=\frac{s_{w1}^{c1}}{s_{w2}^{c2}},c1=\frac{D_{d1}}{D_{d1}-D_{d2}},c2=\frac{D_{d2}}{D_{d1}-D_{d2}}\end{array}$$

(42)

and

$$\begin{array}{c}\hfill {\sigma }_{lv}={\left(\frac{A}{D_{d1}}-\frac{A}{D_{d2}}\right)}^{-1}(ln{s}_{w1}-ln{s}_{w2})\end{array}$$

(43)

Both $a_w$ and ${\sigma }_{lv}$ are functions of $g_f$ and so are fixed with only $s_w$ and $D_d$ varying during this procedure.

Using the above procedure, Cheng et al. [73] find that the surface tension for ammonium sulfate exhibits asymptotic behavior for high molality, approaching a constant value. All of the established parameterizations for the surface tension deviate substantially from these new findings and from each other in the supersaturated regime. Cheng et al. [73] also find that water activity is underestimated by available parameterizations in the supersaturated regime.

The experimental study of Tikkanen et al. [258] finds that the reference thermodynamic model AIM cannot reproduce the observed hygroscopic behavior of DMA-SO${}_4$-NH${}_4$-H${}_2$O particles much smaller than 80 nm. The authors consider size-dependent composition effects such as base evaporation but cannot account for the discrepancy. Zieger et al. [259] also find that thermodynamic models, including ISORROPIA with Kelvin curvature corrections, underpredict the hygroscopic growth of sea salt particles smaller than 100 nm with a progressive divergence as the size is reduced. Since the deviation cannot be explained by Kelvin curvature effects alone, the authors posit that other size dependent effects are occurring such as changes in solubility of sea salt components and reduction of surface tension via calcium enrichment [260]. These thermodynamic models, including AIM, do not account for important size effects in the chemical potential formulation, especially thin liquid layer effects (see Section 2.8.1).

Biskos et al. [261] find that NaCl particle hygroscopic growth factors decrease with decreasing particle size with the largest change in the high RH limit. This behavior was reproduced by a model which accounted for the Kelvin effect and included an additional size-dependent shape factor that was based on theoretical arguments rather than tuning. The shape factor is associated with the change in particle geometry from cubic to spherical upon deliquescence. The shape factor increases from 1.08 in the continuum regime (Kn < 0.1) to 1.24 in the free-molecule regime (Kn > 10).

Bahadur and Russell [262] noted that the deviation of NaCl DRH from predictions is in part due to uncertainties in the solid-liquid and solid-vapor surface tensions and the curvature dependence of the surface tension. They conducted molecular dynamics (MD) simulations to determine the surface tensions (NaCl-air, NaCl-solution and solution-air) of particles in the 2 to 10 nm size range and used the Tolman expression [263] to relate the calculated surface tension ($\sigma$) to the zero curvature limit typically used:

$$\begin{array}{c}\hfill \sigma ={\sigma }_0\left(1-\frac{2\delta }{R}\right)\end{array}$$

(44)

where $\delta$ is the Tolman length, R is the particle radius and ${\sigma }_0$ is the zero curvature limit surface tension. Two parameter fitting of the MD data with both $\sigma$ and ${\sigma }_0$ allowed to vary is used to determine a size dependent expression for the surface tension. The surface tension remains close to ${\sigma }_0$ except for particles smaller than 2 nm. Two MD simulation approaches are used, the energy difference method and the test area method. The two methods do not yield identical surface tension values and are used as uncertainty bounds.

To determine the impact of the size-dependent surface tension on the DRH an iterative procedure with the Russell and Ming [251] wetted particle thermodynamic model was applied to NaCl particles with sizes in the 5 to 150 nm range. In accordance with laboratory measurements, it was found that size dependence of the surface tension can reduce the DRH for nanoparticles but the actual impact depends on surface tension for all three phase interface surfaces in this system. Even with the MD simulations there is a high degree of uncertainty in the surface tension values and the laboratory measurements of Biskos et al. [261] serve as an important empirical constraint.

2.8.3. Size Impact on LLPS in Mixed Organic-Inorganic Aerosol

Cheng et al. [73] only consider ammonium sulfate and salt solutions but their method should be valid for other solutions including organics as long as the necessary measurements can be made. They find that the critical diameter for aerosol particles to be liquid and well mixed has a very similar linear inverse dependence on the bulk phase transition temperature for salt and ammonium sulfate solutions. The glass transition of short chain polystyrene, for which data is available, falls on the same line. This indicates that the relationship between interfacial energy (surface tension) and the enthalpy of phase transition in the aerosol bulk follows a similar relationship for inorganic salts and at least some organic compounds. The results suggest that SIA and SOA particles with diameters under 100 nm should be in a liquid and internally mixed state under atmospheric conditions.

These results are consistent with the findings of Veghte et al. [264] that mixtures of ammonium sulfate and succinic acid do not exhibit phase separation for sizes smaller 170 nm. If succinic acid is replaced by pimelic acid then the limit for phase separation grows to particles larger than 270 nm. This indicates that the size dependence of the morphology reflects the characteristics of the organic species. Altaf et al. [265] find that certain mixtures of ammonium sulfate (AS) and poly(ethylene-glycol) 400 (PEG-400 having molecular weight in the 400 g/mol range) fail to show any size dependence and particles of all sizes are phase separated. This occurs for 50/50 mixtures. By contrast 20% AS and 80% PEG-400 mixtures show lack of phase separation for particles smaller than 86 nm. There is an overlap region between 86 nm and 182 nm where phase separated and homogeneous particles can occur. Altaf et al. [265] conclude that size dependent particle morphology reflects LLPS governed by nucleation with an associated activation barrier. The size-independent LLPS of the 50/50 mixture indicates a spinodal decomposition which lacks an activation barrier.

In a follow-up study, Altaf and Freedman [266] find that the transition region from phase separation to no phase separation depends on the drying rate of the particle. The study conducted in Altaf et al. [265] subjected particles to drying rates (i.e., ≈99.7% RH/s) orders of magnitude greater than occur in the ambient atmosphere. By contrast, during the SOAS 2013 campaign the range of drying rates for particles collected over a one-day sampling cycle was between 6.5 × 10${}^{-4}$ and 1.4 × 10${}^{-3}$ RH/s [267]. As the drying rate is reduced into the atmospheric limit (3.5 × 10${}^{-3}$ RH/s) it is found that particles remain homogeneous only for sizes smaller than 40 nm. Phase separation is found to occur for all sizes larger than this. In addition, the width of the transition region separating the two regimes of phase behavior attenuates substantially from around 121 nm in the high drying rate limit to 4 nm in the low limit. Altaf and Freedman [266] conclude that both kinetic and size dependent thermodynamic effects are responsible. The kinetic effect consists of the impact of particle viscosity (as determined by the drying rate) on phase separation and the increase in the timescale for the transition with increasing viscosity. However, the laboratory data exhibits an asymptotic behavior in the small size limit which indicates a decoupling from the drying rate. This implies another size dependent effect manifesting itself likely associated with the thermodynamics of this system. In our view, it is the surface contributions to the Gibbs free energy discussed in Section 2.8.1 that account for this observation.

Mixed organic-inorganic particles can exhibit complex behavior that affects their surface tension and maximum solubility of organic compounds. Werner et al. [268] find that in aqueous succinic acid and ammonium sulfate or sodium chloride aerosols there is a redistribution of the succinic acid from the interior to the outer layers. For saturated succinic acid solutions the surface concentration is about 10 times larger than in the interior. Inorganic salts enhance this effect. As a result the organic fraction of the particle increases with decreasing particle diameter. The effect is apparent for particles smaller than 200 nm and the deviation factor from bulk behavior exceeds 20 for 10 nm particles with dilute succinic acid solutions and an inorganic mole fraction of about 0.05. These results are consistent with the conclusion of Cheng et al. [73] that small particles can remain liquid and exceed the bulk solubility limit. However, these results indicate that there is some degree of phase separation for the succinic acid and ammonium sulfate aerosols for all sizes below 170 nm, in contrast to the findings of Veghte et al. [264]. The disagreement appears to be related to the lack of an interface between the organic-rich and inorganic-rich phases. It is evident that this aqueous organic-inorganic system cannot be described as well mixed. Other mixtures involving relatively low solubility organic compounds should behave similarly.

3. Gas-Aerosol Mass Transfer Kinetics

Thermodynamic instability between the gas and aerosol phases induces molecular mass transport of chemical species across the interface that divides them. The surface gas-aerosol partitioning is a kinetic process and can exhibit non-negligible time dependence reflecting the aerosol size, composition, and phase, as well as the nature of the partitioning gas. In many cases the timescale for gas-aerosol mass transfer is fast and an instantaneous equilibrium approach can be employed. However, the instantaneous equilibrium approximation is not generally valid under atmospheric conditions [269,270,271].

There are two types of gas-aerosol mass exchange models used in aerosol modeling. The most commonly used in AQ models with various thermodynamics submodels consider only the surface mass exchange (condensation and evaporation) between the gas and aerosol phase. This approach assumes that constituent redistribution inside the condensed phase is rapid and that a uniform condensed phase concentration can be defined. To simulate the nonequilibrium mass transfer regime requires a multilayer model of the condensed phase to resolve the diffusion of surface concentrations into the interior (e.g., [272]). The delayed equilibration of the interior affects the surface vapor pressure of constituents and can potentially be just as important as non-ideal solution activity effects.

Although the surface mass exchange approach is sufficient for pure inorganic aerosols, in reality, aerosol particles (e.g., ones containing substantial SOA fractions) can have very high viscosity and the condensed phase concentration evolves towards equilibrium slowly. These non-ideal properties that vary with RH and temperature can substantially impact the particle mass transfer including the interaction between particles via condensation and evaporation. Ye et al. [273] have conducted a study using isotopic tracking that investigates the interaction of particle populations formed from different biogenic and anthropogenic SOA precursors. It is found that isoprene and $\alpha$-pinene derived SOA particles undergo rapid exchange of mass via semi-volatile constituents even for RH less than 10%. But for $\beta$-caryophyllene derived SOA there was limited exchange even for RH much higher than 30%. The authors point to possible effect of low diffusivity of the latter or its inability to form ideal solutions with toluene and $\alpha$-pinene SOA constituents as possible explanations.

3.1. Surface Condensation-Evaporation Models for Inorganic Aerosol

An early dynamic gas-aerosol portioning formulation was introduced by Wexler and Seinfeld [129]. It includes a correction factor in the molecular transition regime for the Maxwellian flux [274]. The mass transfer equations take the following form (e.g., [275]):

$$J_c=4\pi D_i{R}_p{f}_i({Kn}_i,{\alpha }_{s,i})(C_{i,\infty }-C_{i,s})$$

(45)

$$f_i({Kn}_i,{\alpha }_{s,i})=\frac{0.75{\alpha }_{s,i}(1-{Kn}_i)}{{Kn}_i^2+{Kn}_i+0.283{Kn}_i{\alpha }_{s.i}+0.75{\alpha }_{s,i}}$$

(46)

where $J_c$ is the particle molar growth rate contribution of species i, $D_i$ (m${}^2$ s${}^{-1}$) is the gas-phase diffusivity of the species, $C_{i,\infty }$ (mol m${}^{-3}$) is the species concentration far from the particle surface, $C_{i,s}$ (mol m${}^{-3}$) is the species saturation gas-phase concentration at the particle surface, ${Kn}_i$ is the dimensionless species Knudsen number, ${\alpha }_{s,i}$ is the species surface-bulk accommodation coefficient, and $R_p$ (m) is the particle radius. The function $f_i({Kn}_i,{\alpha }_{s,i})$ is a correction factor for the molecular transition regime. Time-dependence and nonlinearity is introduced through the surface saturation concentration ($C_{i,s}$) which depends on the composition of the aerosol particle including non-ideal solution effects and particle curvature effects.

As described by Zaveri et al. [76], due to the dependence of equilibrium time scale on particle size, the system of ordinary differential equations (ODEs) describing this dynamic partitioning is very stiff and its solution has posed a challenge [126,250,276,277]. In the case of H${}_2$SO${}_4$ the solution is easier since the surface concentration on liquid aerosol can be approximated by zero. For more volatile species (e.g., HNO${}_3$, HCl and NH${}_3$) this is not the case and requires a more complex approach for dynamic gas-aerosol transfer.

In the case of acid-base systems, Wexler and Seinfeld [270] introduced a generalization of the surface mass transfer formulation above, called coupled transport:

$$J_c=4\pi \frac{D_i{R}_p}{f_i}\overline{C}\left(1-\sqrt{1-4f_i^2\frac{C_{b,\infty }{C}_{a,\infty }-K_{ab}}{{\overline{C}}^2}}\right)$$

(47)

where subscript a represents the acid and b the base species, and

$$\begin{array}{ccc}\hfill \overline{D}& =& \sqrt{D_b{D}_a}\hfill \end{array}$$

(48)

$$\begin{array}{ccc}\hfill \overline{C}& =& \frac{D_b{C}_b+D_a{C}_a}{\overline{D}}\hfill \end{array}$$

(49)

$$\begin{array}{ccc}\hfill K_{ab}& =& C_{b,s}{C}_{a,s}\hfill \end{array}$$

(50)

However, not all thermodynamics models solve the time-dependent surface gas-aerosol partitioning problem. ISORROPIA (I and II) is an equilibrium model and approximates the gas-aerosol transfer as an instantaneous process and neglects the Kelvin effect. The equilibrium approach is also adopted by EQSAM [106]. Metzger et al. [106] invoke the results of Meng and Seinfeld [278] that aerosol size determines the time scale for gas-aerosol equilibration. The reasoning goes that small particles have a larger surface area versus volume and so reach equilibrium faster than large particles. However, small aerosol particulates are also subject to coagulation processes and the Kelvin effect so that transience can be imposed by the size redistribution of the aerosol population. Ambient gas and vapor concentrations cannot be assumed to be under equilibrium. Rapid transport processes and gas phase composition changes can invalidate the assumption of separation of timescales between gas-aerosol partitioning and gas phase evolution.

3.1.1. Formulations in EQUISOLV I and II

The dynamic approach to surface gas-aerosol partitioning is adopted by EQUISOLV I and II [90,91,279]. However, an effort was made to minimize the computational cost by finding analytical solutions where possible. Jacobson [126] presents a method to reduce the computation burden of solving the dissolutional and condensational surface mass transfer using a quasi-analytical approach. The approach consists of the following steps: (a) the ordinary differential equation for dissolutional and condensational growth can be integrated to give an analytical expression for the concentration in the gas phase of the species under consideration, (b) a mass balance equation is imposed on all species in the gas and condensed phase, (c) use of the mass balance equation and the analytically derived concentration of each species determines the combined concentration in the gas phase which upon substitution into the expression from (a) gives the condensed phase concentration. These two surface mass transfer schemes are called analytical predictor of dissolution (APD) and analytical predictor of condensation (APC).

Jacobson [279] considered the numerical problem of simultaneous nonequilibrium gas-aerosol mass transfer of multiple dissociating acids and the ammonia base for long time steps. Previous formulations resulted in oscillations for time steps over ∼15 s. This resulted from inadequate time for the pH to reach equilibrium during the time step leading to oscillations between excessive and inadequate growth over successive timesteps. Smaller timesteps involve smaller concentration changes which allows for faster equilibration of pH. The stiffness issue also manifested itself in the equilibration of solids for low LWC (<0.01 $\mathsf{\mu }$g m${}^{-3}$). The amount of gas permitted to transfer into solid was too small for long time steps.

Various approximations were attempted to deal with the stiffness problem [250,276,280] but introduced their own limitations. The predictor of nonequilibrium growth (PNG) scheme detailed in Jacobson [279] avoids the problems of previous schemes while being stable and accurate for long time steps (over 300 s) for all RH and over the full size distribution. This is achieved for semi-volatile acids by applying dissolutional growth at high LWC and condensational growth at low LWC, for non-volatile acids by applying condensational growth at all LWC, with the subsequent equilibration of both ammonia and pH simultaneously between gas and condensed phases at all LWC based on updated acid species content. An operator split equilibrium calculation using EQUISOLV II is performed to obtain the condensed phase composition, updated pH and LWC.

Hu et al. [281] evaluate different schemes for condensation-evaporation on the aerosol surface. They find that kinetic schemes are the most accurate and that of these the kinetic APC scheme [282] is the best in terms of computational efficiency.

3.1.2. ASTEM Module in MOSAIC

Zaveri et al. [76] describe the ASTEM scheme for solving the surface mass transfer. ASTEM operates in two stages spanning $h_{ASTEM}$ seconds each. During the initial stage the nonvolatile species H${}_2$SO${}_4$ and MSA are partitioned pseudo-analytically. During the final stage adaptive time steps, h, are used to cover the h${}_{ASTEM}$ interval. The Jacobson [279] approach is used to discretize the liquid phase cases. Solid phase cases are discretized using the explicit Euler method. At the end of each adaptive timestep the updated concentrations in the liquid and solid phases are summed to obtain the particle-phase concentrations and fed into the thermodynamics subroutines to update the intra-particle phase equilibrium. $h_{ASTEM}$ is chosen to divide the overall chemistry time step (which in many models is the transport timestep) such that ASTEM is called more than once during each transport timestep.

The various cases reflecting different inorganic species and aerosol phase regimes are detailed in Zaveri et al. [76].

3.1.3. Hybrid Dissolution Solver HyDiS-1.0

Benduhn et al. [283] introduced a hybrid scheme for surface mass transfer of semi-volatile species (HNO${}_3$, NH${}_3$ and HCl). The hybrid approach attempts to avoid the numerical cost of explicitly solving stiff equations. This involves assuming instantaneous equilibrium partitioning for certain particle sizes. As noted by the authors, hybrid schemes suffer from limitations related to the equilibrium assumption itself. Under certain conditions, small particles may take much longer to reach equilibrium due to competitive partitioning effects influenced by particle chemistry [271]. Hybrid models typically assume that small particles are most quick in reaching equilibrium. The HyDiS scheme considers four regimes of stiffness based on concentration gradients and particle chemistry (specifically the pH). Equilibrium timescales vary depending on conditions and require a species-specific approach to determine if the equilibrium or dynamical solver is applied.

HyDiS is found to produce solutions with very good agreement with a fully dynamical reference solver. Under highly polluted conditions results start to degrade when HCl is included as a third semi-volatile species. Compared to ASTEM, the hybrid solver may be up to three times faster but a proper quantification is currently lacking.

3.2. Mass-Transfer Formulations for Organic Aerosols

The most common approach to determining SOA mass is to use the Pankow [182] formulation popularized by Odum et al. [185]. This is a thermodynamic equilibrium approach which has recently come into question due to observations of kinetically limited SOA growth [15,21]. However, some models apply the Fuchs and Sutugin [274] kinetic theory that is used for condensation of inorganics [271]. For example, GLOMAP-mode [284] treats secondary organics as nonvolatile like H${}_2$SO${}_4$ and irreversibly condenses them by assuming a zero aerosol surface concentration. Even though this approximation is valid for extremely-low volatility VOC (ELVOC) ($C^{*}\sim {10}^{-5}$ $\mathsf{\mu }$g m${}^{-3}$) and low volatility VOC (LVOC) ($C^{*}\sim {10}^{-2}$ $\mathsf{\mu }$g m${}^{-3}$) this is not necessarily the case for semi-volatile VOC (SVOC) ($C^{*}\sim 10$ $\mathsf{\mu }$g m${}^{-3}$) and intermediate volatility VOC (IVOC) ($C^{*}\sim {10}^4$ $\mathsf{\mu }$g m${}^{-3}$) even for larger particle sizes with reduced Kelvin effect [272]. However, Shiraiwa and Seinfeld [272] find that SVOC can rapidly mass-transfer in the reversible equilibrium manner envisioned by Pankow [182] even in the case of semi-solid SOA particles.

Surface mass transfer models make significant assumptions about particle mixing properties which do not apply in general. This is due to the assumption that the limiting step is the surface mass transfer and not transport and chemical sinks inside aerosol particles. Particles with substantial viscosity have a delayed uptake of constituents into the particle bulk. The impact of aerosol-phase chemistry is highlighted by the study of Shen et al. [285]. They find that both the Pankow absorptive partitioning theory and effective Henry’s law both grossly underestimate the particle-phase concentration of carbonyl compounds by over 10${}^5$ and 10${}^3$, respectively. Field-derived partitioning coefficients for formaldehyde, acetaldehyde, acetone, propionaldehyde, glyoxal and methylglyoxal were found to increase with decreasing RH. It is likely that oligomer formation in the particle-phase exerts a substantial impact on the mass transfer balance when ambient RH is not high. Water dilution at high RH appears to limit this chemical effect likely due to reduced low-volatility oligomer formation rates and hence reduction of this chemical sink inside the particle.

The numerical and laboratory study by Shiraiwa et al. [66] is consistent with the conclusions of Shen et al. [285]. An IVOC (dodecane) precursor system is shown to form oligomers in the aerosol phase which enhance SOA particle growth. Oxidative ageing of VOC products outside the particle phase is not sufficient to account for the observed growth.

3.2.1. Condensed Phase Reactive Mass Transfer Models

To deal with mass transfer in aerosols an additional model is needed to simulate slow diffusive transport in the aerosol particle interior and condensed phase chemistry. Such a model requires spatial discretization of the aerosol particle which typically involves concentric shells below the surface layer and the surface layer itself is represented by a transient sorption layer and a quasi-static surface layer which together can capture surface saturation and competitive co-adsorption [194]. Pöschl et al. [286] developed a kinetic model framework that enables detailed representation of aerosol mass transfer and chemical processes. This so-called PRA framework has served as the basis for detailed aerosol mass-transfer model development (Pöschl [194]) including KM-SUB by Shiraiwa et al. [287] and its successor KM-GAP [288] which included condensation, evaporation and heat transfer.

This class of offline aerosol models with detailed chemistry and physical composition representation includes the ADCHAM model of Roldin et al. [275]. The mass-transfer component of this model consists of a surface condensation-evaporation model based on Fuchs and Sutugin [274] theory and employs the APC and PNG schemes [279], and a kinetic multilayer model similar to KM-GAP [288]. ADCHAM includes a detailed gas phase chemistry and particle phase chemistry model and employs AIOMFAC [89] for activity coefficients.

Such offline models have been shown to capture diffusion and reaction effects in solid and semi-solid aerosols with varying physical properties [194]. However, these detailed kinetic models are too computationally expensive to be run online in 3-D AQ models.

Zaveri et al. [84] introduced a simplified framework designed for use in 3-D models for SOA mass transfer that includes condensed phase diffusion and chemistry effects. This framework avoids the complexity of offline models such as KM-GAP by applying a quasi-steady-state formulation of the reactive diffusion system in the condensed phase for fast-reacting organics and a two-film approximation for slowly reacting and non-reacting organics. Tests of competitive growth between the Aitken and accumulation modes using this framework in the MOSAIC aerosol thermodynamics and mass transfer model show that the size distribution dynamics is sensitive to the volatility of the organic species, condensed phase diffusivity, and condensed phase reactivity. In the case of an initial value problem involving condensation of a plasticising component (e.g., water) on viscous organic particles there is an attenuation of the particle number peak and broadening of the distribution at later times due to increased equilibration times for decreased particle diffusivity.

Zaveri et al. [289] conducted experiments in a Teflon chamber to investigate the impact of SOA viscosity on particle growth and decay dynamics. The system investigated consisted of isoprene SOA formed under low RH conditions in the presence of a bimodal particle background consisting of Aitken mode ammonium sulfate and accumulation mode isoprene and $\alpha$-pinene SOA. The observed growth of the Aitken mode SOA cannot be explained if SOA particles have liquid-like bulk diffusivities. But semi-solid SOA can explain the observations since there is a limitation on semi-volatile organic compound uptake by larger particles compared to smaller particles which have shorter diffusion time scales.

The above work highlights the importance of aerosol phase physicochemical processes for the size evolution of aerosol populations, which in turn affects CCN formation and gas phase chemical composition due to heterogeneous chemistry impacts from aerosol surface area density. Therefore, it is necessary to know the details of the condensed phase chemistry and physical properties of the various constituents to be able to model the SOA mass transfer problem. At the present time, such details are sparse reflecting the partial understanding of SOA chemistry and specifics of the organic species involved. But there is active research in this area aimed at filling these gaps.

3.2.2. Organic Aerosol Diffusivity Models

The emerging consensus from the laboratory and field studies of organic aerosols is that they are characterized by high viscosity which is sensitive to RH [290]. SOA can exhibit amorphous or semi-solid behavior under low RH conditions but transition to low viscosity liquids at high RH. Relatively high viscosity SOA respond more rapidly to changes in ambient RH compared to initial estimates due to the plasticizing effect of water. In addition, changes in gas phase composition can affect aerosol composition on short timescales. However, high viscosity aerosols do have smaller diffusivities and take longer to respond to changes in ambient conditions than assumed in widely used dilute SOA models. This affects the chemical and size evolution of SOA. Recent work attempts to model the diffusion coefficient for different organics likely to be present in SOA.

O’Meara et al. [291] review three models with explicit composition-dependent diffusion to evaluate model dependence of the results. The review also evaluated the equilibration timescales for water soluble particles and potential limitation of hygroscopic growth under atmospheric conditions. The three models considered were ETH [292], KM-GAP [288] and Fi-PaD [293]. All models assume that Fick’s second law of diffusion applies and that the aerosol particle is spherical and this domain is discretized into concentric shells. A thin surface shell is subjected to instantaneous equilibrium with the gas phase. The three models were subjected to a test system consisting of a non-volatile and semi-volatile component representing sucrose and water, respectively. The diffusion for each component was assumed to be independent of location inside the particle. Three cases for diffusion dependence on composition were investigated: (a) independent of the semi-volatile component and therefore constant during the simulation, (b) with a logarithmic dependence on the semi-volatile mole fraction and (c) with a sigmoidal dependence on the semi-volatile mole fraction. The system was evolved to equilibrium from an initial viscous state subject to an instantaneous change in the gas-phase saturation ratio of the semi-volatile component from 1% to 90%, 60% to 80%, 10% to 20% and vice versa. The ETH model was found to be the fastest and reached the e-fold fraction of the equilibrium time about 100 times faster compared to the other two models. All three models produced very similar results for various ranges of diffusion of the non-volatile and semi-volatile components. A key finding by this study is that when semi-volatile component saturation ratio transitions from 60% to 80% the equilibrium state was reached almost instantaneously relative to the residence time of the particle.

The aerosol-phase models considered by O’Meara et al. [291] are too numerically expensive to include in 3-D air quality models. Going from the constant diffusion case (a) to the two variable diffusion cases involves a 100 fold increase in computing time. A deficiency of existing air quality models is that they do not account for composition dependence of aerosol-phase diffusion for organics. A simplified analytical approach for composition-dependent diffusion is considered in O’Meara et al. [294]. It is found that the analytical solution compares well to the numerical solution for cases when the partitioning species has a constant gas-phase saturation ratio and the partitioning rate is limited by aerosol-phase diffusion. When the gas-phase saturation varies, no general analytical solution is found.

A limitation of the above diffusion models is that they assume Fickian diffusion. However, glassy polymers experience additional physical changes such as swelling upon uptake of a constituent such as water that cannot be modeled by concentration gradients with fixed boundary conditions. Such non-Fickian behavior reflects the finite time required for glassy polymers to adjust their structure to compositional change [295]. The non-Fickian regime becomes dominant when the temperature approaches the glass transition temperature ($T_g$) and when the glassy polymer is not diluted sufficiently by a solvent such as water.

Evoy et al. [296] consider the application of the Stokes-Einstein relation (SER) for predicting the diffusivity of SOA organics. SOA was represented by proxies consisting of citric acid, sorbitol and a sucrose—citric acid mixture. The SER links viscosity to diffusivity (D) via the equation:

$$D=\frac{kT}{6\pi \eta R_H}$$

(51)

where k is the Boltzmann constant, T is the temperature in Kelvin, $R_H$ is the hydrodynamic radius of the diffusing species, and $\eta$ is the bulk solution medium viscosity based on laboratory measurements. It is found that the SER does not fit the SOA proxy data. Based on previous work [297] with sucrose-water systems the fractional SER is found to give a much better fit:

$$D=C\frac{kT}{6\pi {\eta }^t{R}_H}$$

(52)

where C and t are empirical fit parameters. Evoy et al. [296] find a best fit to their and literature data for $C=1.66$ and $t=0.93$. The SER is known to break down for glassy polymer systems in the limit of $T_g$ or high viscosity (e.g., [298,299]).

When applied to atmospheric conditions in the planetary boundary layer and above, it is found that the diffusive mixing time based on the regular SER overestimates the result based on the fractional SER by a factor of three and higher in regions where mixing time is an hour or more. The mixing time controls the uptake of gas phase constituents into SOA and their rate of growth and decay. The Evoy et al. [296] findings indicate that the off-line aerosol-phase models described above overestimate the mixing time, especially at lower temperatures or lower RH. Thus, even for relatively low RH, SOA particles should have mass-transfer equilibrium times substantially smaller than the worst cases expected for semi-solids.

3.2.3. Organic Aerosol Mass Accommodation Coefficients

An important physical parameter that limits the mass transfer from the gas phase to the aerosol phase is the mass accommodation coefficient ($\alpha$). For many atmospherically relevant inorganic constituents, the $\alpha$ for uptake on water droplets is around 0.1 [300]. Only some strong acids such as H${}_2$SO${}_4$ and MSA have $\alpha =0.7$ or higher [301]. By contrast, SOA models have assumed that organic constituents have $\alpha =1$ based mostly on liquid evaporation measurements and lack of data for specific organics [302,303].

Julin et al. [302] use molecular dynamics simulations to investigate the uptake of several organic species with atmospherically relevant structural characteristics by media consisting of the same species representing SOA surface conditions ranging from liquid to high viscosity solid and semi-solid. Both surface and bulk mass accommodation coefficients were considered. The species considered were adipic acid, succinic acid, naphthalene and nonane. An expansion chamber experiment was conducted to measure nonane accommodation coefficients for uptake on liquid nonane particles. The Kelvin curvature effect was considered by simulating adipic acid uptake by a 400 molecule (5 nm) adipic acid cluster. No significant difference was found compared to a planar surface. It was found that surface accommodation coefficients for the organics on the range of surfaces considered fell between 0.96 and 1. In the case of the bulk mass accommodation coefficients the upper bound was 1 but the lower bounds were much lower specifically for the cases of adsorption onto solid and semi-solid surfaces. This reflects the transport delay from the surface to the bulk and depends on the boundary between the surface and the bulk. Julin et al. [302] argue that the lower bounds are very conservative. The transition from surface to bulk makes the bulk mass accommodation coefficient difficult to compare to other experimental studies but highlights the impact of sub-surface processes on gas to aerosol phase mass transfer.

Krechmer et al. [303] find a volatility dependence of $\alpha$ for organics using a chamber study. Photochemical oxidation products of a precursor system consisting of alkanol VOCs (hexanol, octanol, nonanol, decanol and dodecanol) together with CH${}_3$ONO and NO were seeded onto dioctyl sebacate (DOS) particles. DOS was selected as a seed particle since it remains a fluid under experimental conditions due to its low vapor pressure and is non-reactive with the oxidation products being considered. Detailed box modeling is done to account for chamber wall losses but with sufficient surface area of DOS seeds it is found that oxidation products partition at much higher rates into the seed particles. The box model was also used to fit the data to infer the accommodation coefficient for different volatility class compounds. For the chemical system considered, it was found that $\alpha$ was close to 0.9 for effective concentration ($c^{*}$) values less than 2 $\mathsf{\mu }\mathrm{g}{\mathrm{m}}^{-3}$, with a transition to 0.4 between 10 and 100 $\mathsf{\mu }\mathrm{g}{\mathrm{m}}^{-3}$. Previous experimental work has pointed to $\alpha$ being close to 1 for $c^{*}$ less than 1000 $\mathsf{\mu }\mathrm{g}{\mathrm{m}}^{-3}$. But the results of this study show this is not the case for $c^{*}$ exceeding 10 $\mathsf{\mu }\mathrm{g}{\mathrm{m}}^{-3}$. So volatility of organics needs to be considered when dealing with SOA mass transfer.

3.2.4. Size Dependence of Mass Accommodation Coefficients

The mass accommodation coefficient is typically treated as a constant in aerosol modeling. However, it can be a function of temperature depending on the molecule being considered (e.g., [304,305,306]) and also exhibits substantial size dependence in the molecular cluster limit [307]. Barclay and Lukes [307] use molecular dynamics simulations to investigate the mass accommodation coefficient, $\alpha$, on curved surfaces for different temperatures and particle sizes. They find a universal particle radius dependence of $\alpha$ for different compounds after scaling:

$$\alpha \left(r\right)={\alpha }_{\infty }\left(1-\frac{r_0}{r}\right)$$

(53)

where ${\alpha }_{\infty }$ is the asymptotic limit of $\alpha$ corresponding to zero curvature, which is typically used in aerosol models, r is the radius (nm) and $r_0$ (nm) is a parameter with a positive linear dependence on temperature and and liquid-vapor interfacial width. For the compounds they considered (Lennard-Jones fluids, Buckingham fluids, benzene, methane, butane, methanol and water) $r_0$ exceeded 0.82 nm for all compounds except the Lennard-Jones fluid with $r^{-12}$ repulsive forces at the lowest temperature (74 K). It is possible to view the 1–2 nm radius limit as the breakdown of the continuum regime and the emergence of the molecular cluster regime. The mass accommodation coefficient becomes zero at this transition.

The size dependence of $\alpha$ has implications for recent efforts to include more explicit treatment of post nucleation growth of secondary aerosols (e.g., [308]). It is tempting to apply the Fuchs and Sutugin [274] theory to drive pure surface condensation growth via extremely low volatility compounds such as sulfuric acid and ELVOC. However, only particles greater than 2–4 nm in diameter can likely be considered. To get within 70% of ${\alpha }_{\infty }$ requires particles with diameters in excess of 6–12 nm.

3.2.5. Organic Aerosol Mass Transfer Summary

For organics kinetic growth aspects can be more critical than for inorganics, especially for small particle sizes. At low RH the organics fraction of aerosols forms very high viscosity amorphous phases which have no resemblance to the low viscosity SOA assumed in the Pankow [182] equilibrium partitioning theory. This is true even if very small particles remain liquid as they will retain high viscosity [73]. It is also possible that even for moderate RH, under conditions of LLPS the organics rich phase has a viscosity high enough to require kinetic mass transfer modeling. So the Meng and Seinfeld [278] size-dependent equilibrium assumption can break down for organic and mixed organic-inorganic aerosols. That SOA growth is to a significant part kinetic, i.e., being proportional to the surface area and not the volume, was established by Riipinen et al. [15] through field studies. Scott et al. [309] evaluated kinetic versus thermodynamic growth using GLOMAP-mode in the TOMCAT CTM [310]. It was found that a kinetic formulation of condensation promotes new particle formation (NPF) as it distributes mass to the finer modes with higher surface area density. By contrast, thermodynamic equilibrium growth distributes more mass to particles with larger volumes. Given the results of Mikhailov et al. [19], Perraud et al. [21] and Abramson et al. [81] SOA typically has significant viscosity and exhibits substantial growth transients and irreversible effects. So it is unlikely that size determined by thermodynamic equilibrium between the condensed and gas phases is a realistic approximation for SOA and mixed organic-inorganic aerosols in the nucleation, Aitken and accumulation modes. This has also been noted by Pöschl [194] and Shiraiwa and Seinfeld [272].

For liquid aerosol particles with relatively low viscosity the mass transfer behavior depends on the volatility of the organic compound being considered. ELVOC, which is critical for post-nucleation growth, exhibits kinetically limited, essentially irreversible, condensation like sulfuric acid and can be modeled using the Fuchs and Sutugin [274] theory assuming zero surface concentration. As aerosol particles grow in size more volatile compounds can partition into the aerosol phase due to the attenuation of the Kelvin effect. So condensation of organic compounds becomes more compatible with thermodynamic equilibrium partitioning at larger particle sizes. Shiraiwa and Seinfeld [272] demonstrate this behavior using detailed model simulations. They also investigate the impact of surface accommodation and particle bulk diffusivity on mass transfer. As the bulk diffusivity transitions into the semi-solid phase regime there is a rapid growth in the timescale required to reach equilibrium for any oxidized VOC. The increase in equilibrium timescale is linear with bulk diffusivity decrease in this diffusivity limited regime. In the liquid and semi-solid bulk diffusivity range the equilibrium timescale has a logarithmic dependence on the accommodation coefficient. In this accommodation coefficient limited regime there is basically no dependence of the equilibrium timescale on the bulk diffusion. So the bulk diffusivity needs to reach a critical point, dependent on the value of the accommodation coefficient, where the mass transfer behavior changes. For SVOC with a volatility (C${}^{*}$) of 10 $\mathsf{\mu }$g m${}^{-3}$ the transition is at a bulk viscosity of about 3 × 10${}^{-13}$ cm${}^2$ s${}^{-1}$, which is within the semi-solid phase and indicates that the thermodynamic equilibrium regime holds even for non-liquid particles. But for very low RH and larger particle sizes this is not likely to be the case.

Li and Shiraiwa [311] evaluate the gas-particle mass equilibrium timescales for SOA under a large range of temperature and RH conditions found in the troposphere using the KM-GAP model. They find that for semi-volatile compounds the equilibrium timescale is on the order of seconds or minutes under conditions found in the planetary boundary layer (PBL), excluding low temperature environments. If particles form glassy or amorphous states due to low RH, then the equilibrium timescale can exceed one hour. In the lower troposphere above the PBL the timescale can exceed days due to low temperatures. These results are consistent with the findings of Evoy et al. [296] in the PBL and to some extent at higher altitudes, but deviations from Stokes-Einstein diffusion still need to be accounted for in the free troposphere.

Although the Shiraiwa and Seinfeld [272] and Li and Shiraiwa [311] model results indicate equilibrium mass transfer even for semi-solid SOA in the case of SVOC, the picture is incomplete. Zhou et al. [312] find that the heterogeneous reaction of benzoapyrene (BaP) with ozone in SOA coated ammonium sulfate particles is significantly affected by the viscosity of the SOA shell, which depends on RH. For RH under 50% the exchange of reactants through the viscous SOA shell was inhibited. This RH range corresponds to a bulk diffusivity for $\alpha$-pinene SOA of under 8 × 10${}^{-14}$ cm${}^2$ s${}^{-1}$. This raises the aspect of particle morphology. Although Shiraiwa and Seinfeld [272] did not detail the distribution of SVOC in the semi-solid condensed phase, it is possible that most of the SVOC accumulated near the surface and established a liquid phase which was subject to the Pankow [182] equilibrium mass-transfer regime. Li and Shiraiwa [311] do discuss this transient near surface bulk effect due to kinetic limitations and concentration gradients. They point to evaporation experiments [313] that indicate formation of a local thermodynamic equilibrium with the near-surface bulk that is distinct from full equilibrium with the entire bulk. The Shiraiwa and Seinfeld [272] results also indicate that there is a transition as the SOA core becomes more viscous and solid-like to a state which takes a very long time to reach mass-transfer equilibrium. This is a transition from a volumetric absorption to a slow kinetically limited surface absorption. If the SVOC dilutes the viscosity of the surface layers of the semi-solid phase, then it enhances the volumetric uptake of SVOC and at the same time experiences reduced rates of evaporation due to a lower mean volatility of these layers on account of pre-existing low volatility organics.

Much like the inhibition of mass transfer of reactants found by Zhou et al. [312], viscous SOA (for RH less than 50%) should inhibit the dehydration of sulfate aerosol if it forms a contiguous surface layer of sufficient depth. Whether such particles form naturally is not so clear. Davies et al. [217] find that formation of coatings with long-chain alcohols can be disrupted by water condensation. For some organics such as isoprene oxidation products no LLPS occurs at any RH (Smith et al. [314]). But for $\alpha$-pinene oxidation products LLPS occurs and the organic shell is semi-solid [315]. In addition, there is the impact of size dependence of morphology discussed in Section 2.8 and phase separated particles may not exist at sizes less than 100 nm depending on the types of organic species involved and their amount relative to inorganic components.

4. Aerosol Heterogeneous Chemistry

In this section we give a brief overview of inorganic and organic heterogeneous aerosol chemistry and present recent work on development of chemistry parameterizations for 3-D models. Heterogeneous or multiphase chemistry [56,57] includes the processing that occurs on the surface and inside of both activated aerosol particles (cloud droplets) and unactivated aerosol particles. For a comprehensive review of heterogeneous photochemistry we refer the reader to George et al. [43] and Shrivastava et al. [67] for SOA processes including organic heterogeneous chemistry and particle growth behavior. In-depth reviews of cloud phase chemistry from measurements and modeling are available in Herrmann et al. [316] and Ervens [317], respectively.

AQ models generally do not include heterogeneous chemistry for unactivated aerosol particles. However, they do include in-cloud aqueous inorganic chemistry since it is an important source of sulfate. Partition of soluble species from the gas phase into the aqueous phase is typically dealt with via Henry’s law constants without detailed calculations for activity coefficients which applies to dilute solutions. Model cloud aqueous chemistry affects aerosol composition via the impact on the primary soluble species (sulfate and nitrate). One of the main sources of sulfate aerosol mass is via in-cloud growth from aqueous phase SO${}_2$ oxidation, and the subsequent survival upon cloud and rain evaporation [11].

In addition to the thermodynamics processes described in Section 2, aerosol growth, hygroscopicity and phase state are affected by both gas phase and multiphase irreversible chemical processes. To model aqueous chemistry associated with liquid aerosols requires accounting for non-dilute solutions, which involves determination of activity coefficients, and reactions in a higher acidity environment which is the case for accumulation mode and smaller sulfate particles.

In the case of SOA both gas phase and condensed phase chemistry results in the formation of reduced volatility products including SVOC, LVOC and ELVOC that serve to enable post-nucleation growth and subsequent growth into the accumulation mode [85]. Chemistry also modifies the uptake of even relatively volatile species. Conversion to lower volatility products in the condensed phase establishes a mass sink and concentration gradient that promotes more partitioning from the gas phase. Such chemical transformation in the particle phase also accounts for the increased resistance to evaporation of SOA particles relative to the simplistic VBS estimates of the volatility of the compounds contributing to SOA growth [318]. At the same time, reactions that result in the break-down of organics into more volatile compounds can facilitate the evaporation of aerosol particles. Heterogeneous chemistry also acts to age hydrophobic primary emissions such as black carbon in soot particles (e.g., [319]). This facilitates formation of aqueous and organic coatings that changes the growth and coagulation properties of these particles. Photochemical ageing is also found to substantially enhance the lensing effect on radiative absorption by black carbon through formation of highly oxidized SOA coatings [320].

There are synergistic effects in coupled organic-inorganic aerosol systems. An interesting finding of Wang et al. [60] is that SOA plays an important role for fine particle modes during heavy smog events by forming aqueous seeds for sulfate formation. Concurrently, rapid hygroscopic growth drives hydration and oligomerization of glyoxal and methylglyoxal associated with emissions of aromatic hydrocarbons which forms SOA [321,322].

4.1. Brief Theoretical Overview

There are three distinct chemical domains in the atmosphere: the gas phase, particle surfaces and particle interiors or bulks. Each domain offers its own set of chemical pathways which is relevant for aerosol evolution. The morphology of the aerosol bulk varies from all-liquid to all-solid and includes mixed states such as a liquid layer fully or partially encasing a solid core and an organic coating either forming a shell or a lens on an aqueous core. Mass transfer links all three domains and is itself affected by chemistry in each domain which determines concentrations of constituents. Uptake of constituents into the bulk is constrained by particle viscosity. In the case of solid or highly viscous particles, such as organics-rich aerosols at low RH, chemistry is confined to either the surface or the outer layer.

The uptake of gas phase constituents is enhanced by reversible and irreversible reactions in the bulk phase since the concentration of the partitioning compound is effectively reduced in the condensed phase; which maintains a concentration gradient and promotes mass transfer. As outlined in Section 2, the dissolution reactions are of the reversible equilibrium type and this includes the primary inorganic salt reactions. However, there are irreversible reactions as well for inorganics and especially organics which experience photolytic and oxidative ageing in the liquid phase as well as in the gas phase.

Surface chemistry is distinct from bulk chemistry due to the effect of unbalanced van der Waals forces (e.g., [323,324,325]) and electric moments of solution electrolytes. As discussed in the cited references, the surface energy effects are comparable to the molecular kinetic energy which is central to reaction kinetics. The chemical potential at the particle surface is thus substantially modified and chemical pathways can occur that are not found in the bulk. For example, gas phase hydrophobic organics can adsorb and react on the surface of aqueous aerosols (e.g., [326]).

Surface reactions on both solids and liquids can exhibit a saturation type behavior following a Langmuir-Hinshelwood (LH) reaction mechanism [326,327]. Using the examples from Valsaraj [326], the pseudo first order reaction rates of naphthalene and anthracene with gaseous ozone on aqueous surfaces take the form

$$k_{1,\sigma }=k_{\max }\left(\frac{C_{ox}\left(g\right)}{C_{ox,1/2}+C_{ox}\left(g\right)}\right)$$

(54)

where $C_{ox}\left(g\right)$ is the gas phase concentration of ozone, $C_{ox,1/2}$ reflects the binding energy of ozone with the water surface and $k_{\max }$ is the bimolecular reaction rate constant which varies depending on surface composition and organic species involved. The LH reaction rate form has more general applicability and is valid for PAH reactions with ozone on surfaces which involves the formation reactive oxygen intermediates which have a much longer residence time on the surface [328].

The heterogeneous reaction of various atmospheric gases on mineral dust aerosol was reviewed by Tang et al. [329]. The reaction of ozone with dust metal oxides follows an Eley-Rideal mechanism (e.g., [330]) since ozone is never observed to desorb from such surfaces implying that the reaction proceeds via a collision with surface chemical constituents. Other cases such as HCHO reactions with dust TiO${}_2$ follow the LH mechanism.

Bulk reactions involving liquid atmospheric aerosols are typically aqueous which involves dissolution of inorganic and organic constituents into strong and generally weak electrolytes, respectively. Inorganic constituents which contribute to aerosol formation such as the acids H${}_2$SO${}_4$ and HNO${}_3$, and the base NH${}_3$ form strong electrolytes and form stable salt compounds through acid-base reactions. The aqueous environment promotes rapid proton transfer (e.g., [331]) which facilitates acid-base reactions compared to the gas phase. In the atmosphere water plays a role at all stages of secondary aerosol formation, and nucleation is likely enhanced by proton transfer in sulfuric acid-water clusters [332,333,334] through the increase in cluster polarity and cluster stability.

Oxidation in the aqueous phase is distinct from the gas phase as well. For example, Munoz-Santiburcio et al. [335] describe the effect of solvation on oxidation of methanol by gold/titania nanocatalysts. Water molecules participate in the catalytic oxidation process and change its character. In the gas phase it is necessary for O${}_2$ and methanol to be co-adsorbed on proximate surface sites but in the aqueous phase this is no longer necessary. Proton transfer, or Grotthuss diffusion, can transfer the acidic hydroxy proton to an atomic oxygen from the dissociated O${}_2$. This requires only for the methanol to come into contact with the catalyst surface in a way to release a proton from its methyl group. Another difference between the gas-phase and the aqueous-phase oxidation process is that charge transfer can occur in stages in the latter compared to the former. The enhanced charge of the solvated catalyst enables O${}_2$ activation without having to first accept a H${}^{-}$ from the methanol. Non-catalytic oxidation can similarly be affected by the presence of water molecules and the associated degree of freedom offered by proton transfer.

As to be expected, organic aqueous reactions follow distinct chemical pathways compared to the gas phase. Photo-oxidative ageing of VOC-derived products in aqueous aerosol yields different products [336,337,338]. Secondary atmospheric aerosols typically have a complex composition including inorganic acids and organic compounds. In addition, these can be concentrated with respect to water for RH less than 80% and exhibit high pH. Thus the chemical environment inside liquid atmospheric aerosols is fundamentally different from that of the gas phase and offers distinctive chemical pathways such as oligomerization and organo-sulfate formation.

4.2. Inorganic Aerosol Chemistry

The AQ model thermodynamics schemes discussed in Section 2 include reversible aerosol phase equilibrium chemistry needed to determine the composition of inorganic salts based on the dominant set of inorganic species (sulfuric acid, nitric acid, ammonia, etc.). However, even with a restricted set of species excluding organic chemistry not all inorganic chemical processes relevant for aerosol particle evolution are represented in such models. Heterogeneous reactions impact the aerosol evolution as well.

Inorganic heterogeneous chemistry acts to age hydrophobic BC particles to hydrophillic (e.g., [319]). Treatment of chemical aging of carbonaceous aerosols (BC) into hydrophillic aerosol facilitated by the surface reaction of ozone with surface adsorbed material is included in some global chemical transport models (e.g., [339,340,341]). Kaiser et al. [342] use the Particle Monte Carlo model coupled to the MOdel for Simulating Aerosol Interactions and Chemistry (PartMC-MOSAIC) system to describe the surface oxidation reactions of generic polycyclic aromatic hydrocarbons (PAHs) with O${}_3$, H${}_2$O, NO${}_2$, OH and NO${}_3$.

Chang et al. [343] provide a review of measurements and model parameterizations of heterogeneous chemistry of N${}_2$O${}_5$ and its impacts on atmospheric chemistry. The significant impact of N${}_2$O${}_5$ hydrolysis on nitrogen oxides and the formation of particulate nitrate at night is shown in several modeling studies (e.g., [344,345,346]). Organic coatings are found to suppress N${}_2$O${}_5$ hydrolysis [224,347].

Sulfate formation by aqueous-phase S(IV) oxidation by O${}_3$ at higher pH in sea-salt and mineral dust aerosols is shown to impact the global tropospheric sulfur distribution (e.g., [348,349,350]). Transition metal-catalyzed heterogeneous oxidation of SO${}_2$ inside or on the surface of particles is likely to be important [351]. Experimental studies suggest strong oxidation of SO${}_2$ by NO${}_2$ in the presence of NH${}_3$ and alkaline dust within haze particles during heavy haze episodes in China [60,352]. There is also evidence for heterogeneous SO${}_2$ oxidation on existing low pH aerosol surfaces and RH-dependent parameterizations have been applied to reproduce the observed high concentrations of sulfate [39,353]. Reactive uptake of SO${}_2$ onto mineral dust particles dependent on NH${}_3$, NO${}_2$ and RH has been modeled by Zhang et al. [354]. Shao et al. [355] parameterize aqueous-phase S(IV) oxidation in wintertime haze aerosols in Beijing, China, by various oxidants such as H${}_2$O${}_2$, O${}_3$, NO${}_2$ and O${}_2$ catalyzed by transition metal ions and conclude that transition metal ion chemistry is important.

Dimethyl sulfide (DMS) oxidation is a dominant source of sulfate aerosols in the remote marine boundary layer and global upper troposphere but multiphase oxidation of DMS can reduce the formation of sulfate particles [356]. Heterogeneous bromine chemistry in sea-salt aerosol is a major source of gaseous reactive bromine in the lower troposphere currently included in a few models using empirical approaches [357,358,359]. Reactive bromine acts to oxidize SO${}_2$ into sulfate in the aqueous phase [360].

4.3. Organic Aerosol Chemistry

The complexity of organic chemistry has resulted in only a small part of the organic species believed to be present in the atmosphere being identified and physically characterized by laboratory studies. This applies both to the gas phase VOC reactions and aerosol phase reactions. Lack of a clear view of the aerosol related chemistry processes has prompted the development of stand-alone SOA models which are idealizations that do not take into account the coupling with inorganic aerosol constituents.

However, there has been progress in characterizing the basic processes involved in SOA formation. Starting from various emitted VOCs multiple stages of oxidation of products occur both in the gas and aerosol phases with repartitioning between these phases at each stage [85]. Even though most of the species and reactions are uncertain, there are chemical pathways leading to higher molecular weight compounds with lower volatility which are eventually trapped in the aerosol phase. At the same time, a fraction of the reactions produces lower molecular weight products that have high volatility which escape the aerosol phase. The VOC classes that are most likely to lead to the formation of less volatile organic oxygenates are the cyclic ones. Cyclic compounds have the feature that fragmentation products associated with ring-opening reactions have a similar number of carbon atoms as the parent compound [85]. Compounds with fewer carbon atoms tend to be more volatile and less likely to contribute to SOA formation and growth. For cycloalkenes, aromatic hydrocarbons and terpenes oxidation occurs primarily via addition of polar functional groups, which reduces volatility [85].

4.3.1. Gas Phase Processes

The gas phase oxidation or ageing in successive sages described above is the foundation of the VBS scheme [186]. A limitation of the standard VBS framework is the lack of any fragmentation pathway which would increase volatility. The impact of fragmentation reactions has been taken into account by Koo et al. [361] and Shrivastava et al. [362] with updated VBS frameworks. Another limitation of VBS schemes is that the link to actual chemical processes is obscured since the chemistry is parameterized and simplistic. It would be preferable to have an explicit lumped chemistry directly tied to the oxidative ageing volatility evolution.

Recent laboratory and field work has revealed that the progression of volatility reduction envisaged in the VBS formulation is not always valid. In the case of $\alpha$-pinene oxidation, a major source of biogenic SOA, the lowest volatility products, ELVOC, which are key to overcoming the Kelvin effect growth barrier after nucleation, are formed in the initial oxidation stage [363]. This first stage oxidation ELVOC decomposes into more volatile fractions in the condensed phase. For $\alpha$-pinene it is clear that the progression of oxidation state is not sequential. The rapid autoxidation formation of ELVOC can occur from other terpenes but not for isoprene [364].

Generalizations of VBS which allow for evolution in both volatility and oxidation state (C:O ratio) have been developed [365,366,367] which substantially improve the representation of the chemical processes compared to the original VBS scheme. This 2-D VBS approach allows more generalized trajectories in oxidation evolution, including fragmentation and increase in volatility of products. Zhao et al. [367] used chamber experiments to calibrate the 2-D VBS model. They use a three level approach, one level models anthropogenic VOC oxidation while another level models biogenic VOC oxidation, and an additional level is dedicated to photo-oxidation of primary organic aerosol (POA) and IVOC from all sources. The need for this arises from the differences in underlying chemistry which cannot be generally parameterized even within a 2-D VBS framework. As indicated by their results with CMAQ this approach manages to overcome the serious SOA yield deficit problem.

Stadtler et al. [368] have developed a semi-explicit gas phase isoprene oxidation chemistry model that they find to substantially improve the global isoprene-derived SOA amounts in a 3-D chemistry climate model. A group contribution scheme was used to obtain the saturation vapor pressure for the compounds considered. Reactive uptake of IEPOX was also included. By tracking the volatility of each precursor in every size section it was possible to identify ISOP(OOH)${}_2$ as the primary SOA formation channel. IEPOX was found to be the secondary SOA formation channel. This model did not include particle phase chemistry but highlights the need for more detailed gas phase chemistry formulations compared to VBS.

The improved results in the absence of particle phase chemistry seen by Stadtler et al. [368] suggest a cancellation of missing processes. As in the case of the $\alpha$-pinene SOA system, there may be a decomposition of gas phase low volatility products in the condensed phase [363]. In addition, certain compounds have short thermal breakdown lifetimes such as is the case for peroxide-containing highly oxygenated molecules formed from terpene precursors [369]. Carbonyl compounds formed during terpene oxidation are also sensitive to photolysis which affects the composition of SOA and acts to reduce its mass [369]. In general, particle phase chemistry cannot be ignored for SOA evolution.

4.3.2. SOA Formation Due to Aerosol Phase and Cloud Chemistry

Shiraiwa et al. [370] identify molecular corridors and kinetic regimes in the heterogeneous chemistry of SOA. There is a rather tight dependence between molecular weight and volatility in the species generated by organic reaction chains in the aerosol phase. Different classes of reactions and their products fall within a confined region of the 2-D molecular weight and volatility space. This work suggests a VBS-like approach for the condensed phase and offers guidance to the development of more explicit chemistry schemes for SOA formation.

The aerosol phase chemistry is a critical component of the evolution of SOA but it has been omitted in SOA models including those based on the VBS approach. Shiraiwa et al. [66] demonstrate that condensed phase chemistry plays an active role in the dynamics of the size distribution of organic aerosol. They investigated the dodecane SOA system and found that heterogeneous surface reactions drive SOA mass increase. It is likely that other SOA systems are subject to similar condensed phase chemistry growth enhancement. This growth mechanism for SOA is consistent with the field study of Riipinen et al. [15] that identifies kinetic or surface area dependent growth as accounting for about 50% of SOA size evolution for ultrafine aerosol (<100 nm diameter). However, this kinetically limited growth behavior can also reflect the physical characteristics of atmospheric SOA particles for medium and low RH which can have very high viscosity. These semi-solid and solid particles are subject to a slow diffusion regime and even concentration front formation and take a long time to equilibrate to changes in gas phase composition.

One of the pathways for SOA formation is via aqueous chemistry in clouds [316,371]. Isoprene is the most abundant biogenic VOC emitted into the atmosphere [372] and thus can make a significant contribution to SOA formation even if it contributes lower yields than other VOC precursors such as $\alpha$-pinene [373]. Glyoxal is a key intermediate product of isoprene oxidation in the gas phase. It effectively partitions into the aqueous phase due to its high solubility (Henry’s law constant of 3.5 × 10${}^5$ M/atm) and undergoes further oxidation via reaction with OH. For relatively dilute cloud droplets the main product is oxalic acid which readily forms the ammonium oxalate salt in the presence of ammonium which greatly reduces its volatility and accounts for oxalate being present primarily in the particle phase in the atmosphere [374]. Upon droplet evaporation any residual glyoxal can self-oligomerize and contribute to SOA mass [375].

Oxalic acid interaction with inorganic salts was found to form products with low volatility and hygroscopicity by Drozd et al. [376] in submicron particles. Lee et al. [377] find that glyoxal facilitates formation of highly oxidized, low volatility products from biogenic precursors (e.g., $\alpha$-pinene) in cloud droplets.

For submicron wet particles (not necessarily associated with clouds), Gaston et al. [378] find a strong dependence on particle acidity and reactive uptake of IEPOX. This basically reflects that the aqueous phase IEPOX reaction mechanism is that of acid-catalyzed epoxide ring opening. However, the acidity sensitivity is greatest for very acidic aerosols which would be relevant for urban emissions sources. Gaston et al. [378] also identify a self-limitation effect for IEPOX uptake if it is associated with the formation of organics that act to produce a coating on the aerosol particle. Pre-existing aerosol may include an organic coating which will limit IEPOX uptake. The surface-active species that from coatings are typically comprised of both hydrophillic and hydrophobic moieties [379]. The IEPOX uptake inferred by Gaston et al. [378] is much higher than that specified by Pye et al. [380] who investigate the impact of IEPOX. Higher uptake would contribute to improved modeling given the evaluation of Karambelas et al. [381]. However, the results of Gaston et al. [378] should be regarded as being valid only over short timescales.

Liu et al. [382] find that ambient SO${}_2$ enhances the formation of SOA in the OH ageing of eugenol via conversion into sulfate. Eugenol is used as a proxy for biomass burning products dominated by methoxyphenols. The effect is at least in part due to condensed phase acid driven oligomerization or other low-volatility product formation. But sulfate is also hygroscopic and acts to increase the SOA particle size. The increased particle surface area and volume can facilitate additional heterogeneous chemistry and organic species uptake. As with isoprene-derived SOA, phase separation and the formation of a mass-transfer limiting organic shell appears not to have occurred [314]. But there was a decrease of the SOA enhancement for higher OH concentrations while the sulfate amount increased linearly. This could indicate some limitation on the uptake of organics, but could also reflect the relative increase in fragmentation reactions for high levels of OH.

Pye et al. [383] investigate the impact of organic and inorganic aerosol phase interaction on gas-particle partitioning. They find that organics tend to increase the pH. At the same time, the uptake into the aerosol phase of highly oxygenated compounds (O:C ≥ 0.6) was increased and the uptake of low oxygenation compounds was reduced. Pye et al. [383] also consider the model-measurement discrepancy for ammonium found by Silvern et al. [384] and conclude that it reflects variations in data sets rather than a failure of thermodynamics models.

Silvern et al. [384] find that thermodynamics models cannot predict the ratios of ammonium relative to sulfate in eastern US aerosol measurements. This ratio should approach a value of 2 when ammonia is in excess but is found to be 1.04 ± 0.21 (mol mol${}^{-1}$) under these conditions. In spite of SO${}_2$ emission reductions between 2003 and 2013, the acidity of aerosol is increasing while ammonia emissions remain unchanged. The authors point to the increasing OA concentrations during this period from 1.1 to 2.4 (g g${}^{-1}$) as the explanation. Organic coatings on sulfate aerosol would act to inhibit uptake of ammonia and such coatings are becoming more substantial with time.

The Pye et al. [383] results suggest that phase separation is sufficiently variable that the ammonia uptake by the aerosol phase is not systematic. This is consistent with the results of Smith et al. [314] that isoprene SOA does not exhibit LLPS and isoprene is a major precursor for SOA in the south-eastern USA [385].

Jo et al. [386] have developed a parameterized IEPOX-SOA model which is based on a combination of analytical and fitting approximations of a full gas phase chemistry model. Simulations with the global chemistry and transport model GEOS-Chem [131] show this scheme to be superior to both a simple fixed yield formulation and the VBS scheme. In addition to the parameterized IEPOX-SOA chemistry, the authors included organic coating induced limitation of reactive IEPOX uptake by sulfate aerosol. The total available surface area of mixed inorganic-organic aerosol was used. The addition of organic coatings increased the available surface area which counteracted the reduction in uptake.

5. Summary

5.1. State of Thermodynamics Models

Thermodynamics models for inorganic chemical species that can be implemented in 3-D models have improved in terms of accuracy and numerical performance [76]. The MOSAIC model by Zaveri et al. [76] does not assume mass transfer equilibrium which agrees better with observations and it removes many other simplifying assumptions applied in the popular ISORROPIA model [116]. The computational expense of MOSAIC is close to that of ISORROPIA and it is significantly faster than the vectorized EQUISOLV II [91]. UHAERO [120] is another model that makes fewer approximations, but it appears not to be used in current AQ models.

There have also been significant advances in aerosol thermodynamics modeling which includes both inorganic and organic components and accounts for LLPS [36,37]. However, these models are computationally expensive and cannot be incorporated into AQ or climate models. Currently no mixed inorganic-organic thermodynamics submodel is available for 3-D modeling. Nevertheless, simplified approaches have been implemented which parameterize LLPS and apply distinct gas-particle mass exchange and thermodynamics formulations to the organics-rich and electrolyte phases, and are able to model a more realistic hygroscopic behavior for mixed organic-inorganic aerosols such as those by Pye et al. [167] in CMAQ and UNIPAR [175]. The non-ideal VBS + BAT model of Gorkowski et al. [168] appears to currently be the best choice for 3-D transport models (for organics-only thermodynamics) in terms of the representation of co-condensation of organics and LLPS. It can be coupled to inorganic constituents similarly to UNIPAR.

The gas-aerosol partitioning problem for inorganics is handled with a kinetic approach in MOSAIC and EQUISOLV II. These models correctly take into account the Kelvin effect. In contrast, ISORROPIA assumes steady state equilibrium partitioning and does not include the Kelvin effect based on the argument that most aerosol mass is distributed in particles larger than 100 nm in diameter [116]. Considering the importance of aerosol with diameters less than 100 nm for CCN and for post-nucleation growth this is a major shortcoming. Aside from the lack of the Kelvin effect, it is not clear that the steady state approximation is valid for a rapidly evolving aerosol size distribution even if the typical diameters are small. The separation of timescales between growth by coagulation and condensation and compositional equilibrium inside aerosol particles could be too short. At the same time, the mass exchange between the gas and aerosol phases reflects the instantaneous aerosol composition so that simulation of partitioning requires a dynamical approach ([76], and references therein).

5.2. Need for Improvement in the Small Particle Limit

Ultimately, surface tension terms that are typically neglected need to be included in the Raoult equation which substantially impacts the thermodynamic equilibrium between the gas and liquid aerosol phases for particles much smaller than 100 nm. At the present time both the online inorganic models and offline mixed organic-inorganic models lose accuracy for particles in the nucleation and Aitken modes. These models are formulated based on bulk limit Gibbs free energy for the aerosol particle chemical system which is not valid in the small particle limit due to the lack of explicit surface tension and disjoining pressure terms in the chemical potentials. There is also a lack of surface tension and growth factor data for small supersaturated particles [73]. So even reference models such as AIM employ extrapolated data and exhibit significant deviation from observations for particles smaller than 100 nm.

The deviation of the chemical potential from the bulk limit affects both the particle mass transfer and phase transition behavior. This includes LLPS and solid-liquid phase transitions. Efflorescence and deliquescence can occur via a nucleation type process that can be externally triggered (e.g., [94]) or as spinodal type transitions without any energy barrier (e.g., [252,265]). The plethora of organic constituents makes characterization of phase transition behavior very case specific and challenging to formulate in a model.

The breakdown of thermodynamics models in the small particle limit appears to not to be common knowledge. For example, MOSAIC has been applied to aerosol particles as small as 1 nm in diameter (e.g., [308,387,388]). The improved data gathering method for particles smaller than 100 nm introduced by Cheng et al. [73] has not yet been used to update existing inorganic thermodynamics models. Size effects on solubility also have not been incorporated. As noted by Cheng et al. [73] there is a lack of surface tension data for nanoparticles due to the difficulty of directly measuring such supersaturated particles with probes and lack of alternative unintrusive techniques. Taking into account the size effect on solubility, particles of ammonium sulfate with diameters less than 10 nm can be predicted to deliquesce at less than 85% RH, in contrast to standard thermodynamics formulations which predict deliquescence only at 100% RH. The Cheng et al. [73] work is also relevant for the thermodynamics of organic aerosols.

The nucleation mode is an important stage in secondary aerosol growth. However, additional corrections are required for particles smaller than 3 nm aside from those mentioned above pertaining to the lack of surface tension effects on solubility in existing thermodynamics models. The surface tension begins to exhibit rapid variation with particle radius in the molecular cluster limit. Terms associated with rapid variation of the surface tension with particle radius are generally neglected (see the supplementary material of Cheng et al. [73] for detailed examples). The surface tension exhibits a very large deviation from its bulk flat surface value in the molecular cluster limit as the surface layer gets progressively thicker relative to the particle size and the distinction between the particle surface and interior breaks down. So in order to properly model liquid aerosol growth from nucleation it is necessary to have a surface tension model. Some mixed composition surface tension models have been developed but this remains an evolving field [240,243].

5.3. Need for Improved Representation of Organics in 3-D Model Thermodynamics

Until recently, the treatment of organics in thermodynamics and heterogeneous chemistry has been non-existent or simplistic. For SOA the most common formulation of gas-aerosol partitioning is that of Pankow [182] which is an equilibrium absorption model that does not capture the effects of high viscosity and the Kelvin effect (e.g., [389]). Some studies use a Fuchs and Sutugin [274] type condensation formulation for organics which is a kinetic approach that includes molecular regime effects and takes into account the Kelvin effect (e.g., [390]). A general kinetic framework for SOA gas-aerosol exchange has been developed by Zaveri et al. [84]. However, due to lack of relevant data this framework has not been deployed in models. Pure SOA or mixed organic-inorganic thermodynamic modeling is not currently done in 3-D models due to the numerical expense.

Riipinen et al. [15] show that exclusively thermodynamic partitioning of organics following the Pankow [182] model is not able to explain observations for ultrafine aerosol (under 100 nm) in a boreal environment, which is critical for CCN formation. Around 50% of the partitioning of organics onto fine mode aerosol is kinetic in nature. Their use of the terms “kinetic” and “thermodynamic” refers to the growth dependence for LVOC and SVOC, respectively. LVOC partitioning depends on the surface area of the pre-existing SOA and proceeds as a kinetic condensation process. In contrast, SVOC partitioning depends on the volume of pre-existing SOA and the assumption that mixing is rapid into the SOA bulk. This suggests that LVOC is forming coatings or lenses on SOA and SVOC is partitioning deeper into the SOA interior or core. This is consistent with LVOC having higher viscosity than SVOC and exhibiting slower diffusion into the bulk phase of aerosol particles. Shiraiwa and Seinfeld [272] find that SVOC can follow the Pankow equilibrium absorption regime even for aerosol particles with very high viscosity and exhibiting semi-solid characteristics. By contrast LVOC, and specifically ELVOC, undergoes irreversible condensation. The non-equilibrium, kinetic growth mode can also reflect heterogeneous chemistry in the aerosol phase [66].

Another major implication of the differences in SOA growth from LVOC and SVOC is that the latter partitions preferentially into larger SOA particles with greater volumes, whereas the former partitions into smaller SOA particles which account for most of the SOA surface area and number density. Details of the mass transfer affect the number distribution since a pure thermodynamic growth regime would distribute mass into larger particles which have smaller number densities. A surface area driven process increases the size of numerous small particles in the nucleation and Aitken modes and leads to higher CCN number densities, even after accounting for coagulation.

Since aerosol indirect effects are associated with CCN production in the ultrafine mode, it is necessary to understand and properly account for aerosol growth processes of nanoparticles, including non-equilibrium condensation and heterogeneous chemistry. Submicron secondary aerosol particles have substantial organic fractions so pure inorganic models are limited in their ability to represent real-world aerosol growth and properties. In the case of SOA, the low viscosity liquid state may only occur at very high RH and for more common lower RH conditions semi-solid and solid states are more likely. In this case diffusive transport limitations into the particle bulk need to be accounted for as well. In the list of organic thermodynamics models covered in this paper there are options for explicit modeling of kinetic partitioning from the gas to aerosol phase, including reactive uptake, (e.g., [84]) and LLPS) (e.g., AIOMFAC and KM-GAP). However, these models are numerically expensive or lack sufficient observational data and thus are not practical for incorporation into AQ and climate models at the present time.

LLPS is observed to be a common feature of mixed organic-inorganic particles [38,391] with an organic shell surrounding an aqueous core. Given the high viscosity nature of SOA (e.g., [81]), it is reasonable to ask if this organic shell inhibits exchange of gas phase constituents with the aqueous core. Zawadowicz et al. [391] find that there is little inhibition of mass exchange between the core and the gas phase in the presence of an organic shell. This implies that the organic constituents that they investigated did not form a high viscosity phase. This likely reflects the hygroscopic properties of these organics and it is possible that the organic shell has a high water fraction at RH values which do not involve efflorescence of the inorganic ammonium nitrate core. It is not clear if their experiments subjected phase-separated particles to one or several dehydration cycles which would have changed the physical properties of the organic shell.

SOA with a low inorganic fraction will likely exhibit distinct behavior compared to more common mixed particles since exposure to low RH would entail the formation of an amorphous phase with a very high viscosity in the particle bulk. The timescale for this phase to fully dissolve at higher RH can be expected to be long [272]. However, such particles can reach quasi-equilibrium conditions with the ambient environment due to surface adsorption of water and semi-volatile organics [272]. Sulfate aerosol particles with higher inorganic fractions and an organic shell subjected to dehydration could also involve transition to high viscosity in the shell that would act as a barrier to mass exchange between the gas phase and the inorganic core. But for the core of such particles to undergo deliquescence, a sufficient dissolution of the organic shell is required so as to render this case academic. The dissolution of a viscous organic shell would likely take much less time than for the interior bulk of dehydrated SOA particles of similar dry diameter.

An important aspect highlighted by the experiments of Zawadowicz et al. [391] is that small (e.g., under 100 nm) mixed organic-inorganic aerosol particles do not undergo LLPS regardless of organic species, consistent with other studies (e.g, [264,265]). A more recent study by Altaf and Freedman [266] indicates that the transition occurs below 40 nm. This means that organics have maximal effect on the ERH and DRH of nucleation and Aitken mode particles and potential inhibition of mass exchange with the gas phase is suppressed. Due to lack of LLPS, mixed organic-inorganic particles with sufficient organic content in the nucleation and Aitken modes will have a propensity to remain in the liquid phase over a wide range of RH due to suppression of efflorescence [73,196,197,198,199]. Coupled with the fact that small particles are surface active, this implies that bulk phase aqueous chemistry in these particles is maintained when it shuts down in larger particles due to efflorescence.

Current Options for 3-D Models

Based on the current understanding and model limitations the following suggestions can be made for 3-D AQ and climate models. Pure inorganic aerosols are an idealization for submicron secondary aerosol particles. Instead, secondary aerosols in this size range are mixtures of organic and inorganic constituents with substantial organic fractions. For particles larger than 100 nm, LLPS is likely to be a common feature and parameterized approaches such as those of Pye et al. [167] and Im et al. [175] offer a numerically tractable aerosol submodel approach. For particles smaller than 100 nm, size effects are important and there is limited applicability of existing thermodynamics models. Aerosol particles exist in a liquid state at these sizes under regular atmospheric conditions and LLPS is unlikely. The need for mixed organic-inorganic thermodynamics cannot be avoided. Any thermodynamics model would need an empirically constrained surface tension model to explicitly account for surface terms and thin layer effects in the Gibbs free energy of the system.

The Pye et al. [167] scheme can be extended further into the mixed organic-inorganic regime to some extent. The $\kappa$-Köhler theory can be used to predict water uptake by the mixed aerosol. A modified Pankow [182] type mass partitioning between the gas and aerosol phase that accounts for the Kelvin effect can be used for relatively volatile compounds but for LVOC compounds a surface condensation scheme (e.g., [274]) should be adopted. A neural-network trained two-parameter or four-parameter Margules equation [168,392] can be used to obtain more accurate binary organic activity coefficients. The Margules equation is based on the observation that $log\left(a_i\right)$ is approximately a linear function of the mole fraction of constituent i in the liquid phase. This is the basis of the PD-FiTE approach as well [88,139]. A surface tension model is needed for mixed organic-inorganic aerosols that is valid in the small size limit. At present there are no established models that could serve this purpose and laboratory derived data is lacking [73].

5.4. Need for Aerosol Heterogeneous Chemistry

Heterogeneous chemistry is a critical part of the general mass transfer problem for aerosols. However, the currently available models in this field are incomplete especially when organic constituents are involved. Fully coupled organic-inorganic heterogeneous chemistry models that include all VOC precursors or even a representative selection (e.g., most important biogenic and anthropogenic VOCs and their gas phase oxidation products) that can be incorporated into 3-D models are not available. Existing simplified chemistry models that improve on the SOA yield problem indicate that more sophisticated gas phase chemistry representations are necessary. This can include semi-explicit mechanisms or 2-D VBS tuned for particular VOC precursor classes. Similar conclusions are likely to apply for condensed phase chemistry as well. The fact that organics are a substantial component of submicron particles makes the need for organic chemistry modeling unavoidable. The approach common in most AQ models of today which involves a complete process splitting of organics from inorganics is not adequate. Organics have a leading order impact on the thermodynamics and mass transfer of submicron aerosol particles and hence on their size evolution.

5.5. Alternative Modeling Approaches to Deal With Aerosol Physicochemical Complexity

Given the complex and varied nature of aerosol composition and dynamics in the presence of organics, it is unlikely that traditional parameterization approaches are viable to advance aerosol process representation in AQ and climate models. Thus, there is a need for application of machine learning [393] and other meta-modeling (e.g., [394]) methods. As noted by Rothenberg and Wang [394] the meta-model framework captures physical process complexity where standard “physically based” approximations common to parameterizations cannot. There is active research using these methods in many fields and this includes to some extent aerosol science. For example, work has been done to process observational data of aerosol composition (e.g., Christopoulos et al. [395]) and to produce a mixing state scheme usable by climate models based on the processing of the output of a numerically expensive and detailed off line model [396].

However, to the best of our knowledge these methods have not been applied in the case of thermodynamics and other numerically expensive parts of aerosol models until recently. It is reasonable to assume that applying meta-modeling methods to complex offline mixed organic-inorganic thermodynamics models can produce a reduced numerical cost scheme usable in AQ and climate models as long as there is some coherent representation of the process with a reduced description of the aerosol composition. Machine learning could be applied to determine the optimal reduced aerosol composition that is required. Ramírez-Beltrán et al. [392] trained a neural network combined with a group-contribution approach to obtain parameters for the four-parameter Margules equation. The activity coefficient predictions were slightly better than those of the UNIFAC model. Gorkowski et al. [168] trained two neural networks to solve for aerosol water content and mass transfer of semi-volatile species, respectively, to develop their reduced complexity organics thermodynamics model. Smirnov et al. [397] have applied a tensor train decomposition for solving the multicomponent Smoluchowski equation for coagulation. The numerical expense of coagulation increases rapidly with the number of sectional bins and aerosol constituents so efforts to reduce numerical costs are justified [398]. This work has yet to propagate into AQ and climate models.